two different fields into the same XDMF file


471
views
2
6 months ago by
Dear all,

I have the 2017.2 version, that supports xdmf i/o. However I cannot save two fields (say, displacement and pressure) on the same xdmf file.

from dolfin import *

mesh = UnitSquareMesh(10,10)
fileO = XDMFFile(mesh.mpi_comm(), "fileO.xdmf")

Vh = VectorElement("CG", mesh.ufl_cell(), 2)
Qh = FiniteElement("CG", mesh.ufl_cell(), 1)
Hh = FunctionSpace(mesh,MixedElement([Vh,Qh]))

(u, p) = TrialFunctions(Hh)
(v, q) = TestFunctions(Hh)

GammaU = CompiledSubDomain("near(x[1],1) || near(x[0],0)")

bc = DirichletBC(Hh.sub(0), Constant((0.0,0.0)), GammaU)

f = Constant((0.0,1.0))

a = 2*inner(nabla_grad(u),nabla_grad(v)) * dx \
    -div(v) * p * dx - div(u) * q * dx
b = dot(f,v) * dx

sol = Function(Hh)
solve(a==b,sol,bc)
u_h,p_h = sol.split()
u_h.rename("u","u"); fileO.write(u_h,0.0)
p_h.rename("p","p"); fileO.write(p_h,0.0)
​

Of course, I can generate two separate files and save each field on it, replacing the relevant lines in the snippet above by

fileU = XDMFFile(mesh.mpi_comm(), "fileU.xdmf")
fileP = XDMFFile(mesh.mpi_comm(), "fileP.xdmf")​

and

u_h.rename("u","u"); fileU.write(u_h,0.0)
p_h.rename("p","p"); fileP.write(p_h,0.0)​

However one needs the mentioned feature to observe the value of the scalar field plotted on the deformed configuration. I could not do that in Paraview (at least not in versions 4.4 and 5.0.1) from separate xdmf files (and the same happens with pvd files).

Any help would be welcome.
Community: FEniCS Project

1 Answer


4
6 months ago by
In order to open multiple functions in one xdmf file in ParaView, try opening with option xdmf3 Reader (!!! without Top Level Partition). This worked for me with the latest ParaView.

Another option, much more safe I would say, is to really save one function per one xdmf file. Read then both of them into ParaView and apply Append attributes filter, https://www.paraview.org/Wiki/ParaView/Users_Guide/List_of_filters#Append_Attributes
1
Just for the record: I always do what OP was trying to do, and it works fine using ParaView's "Xdmf Reader". So at least in my experience, "xdmf3" is not required. I've been doing this since ParaView 5.3.0. Currently I use ParaView 5.4.1.

This makes me still wonder what the bug is for OP's case; but of course if OP is happy, then I'm happy.

For the sake of completeness, this is how I write my solution file (though it's exactly what OP was doing):

def write_solution_to_xdmf(solution, time, file):
        """Write the solution to a file.
        
        Parameters
        ----------
        solution: fenics.Function
        time: float
        file : fenics.XDMFFile
        
            write_solution should have been called from within the context of this open fenics.XDMFFile.
        """
        pressure, velocity, temperature = solution.leaf_node().split()
    
        for var in [pressure, velocity, temperature]:
        
            file.write(var, time)​
written 6 months ago by Alexander G. Zimmerman  
This works for me, however when I try to look at already computed timesteps while the simulation is still running paraview crashes with
"Number of dimensions in light data description in Xdmf does not match number of dimensions in hdf5 file."

After the simulation is finished the problem vanishes..
Is this the same in your case?
written 6 months ago by Lukas O.  
1
@Lukas O. The flush_output option is what you need, see https://fenicsproject.org/qa/3997/previewing-xdmf-file/
written 6 months ago by Michal Habera  
Thank you Michal! Together with the reload feature in the newer paraview versions I really like this.
written 6 months ago by Lukas O.  
Yes. Perhaps if one closes the file after writing at each timestep the problem disappears?
written 6 months ago by Ever Ardo RB  
I never trust ParaView to read files that I'm still writing to. Maybe there's a safe way to do it; but I just copy the files somewhere else.
written 6 months ago by Alexander G. Zimmerman  
1
Some more experience for future reference. I found a problem using Paraview (5.4.1)'s 'plot over line' filter wherein only the last function written to the file was shown. I fixed this by setting the parameters:

file0.parameters['rewrite_function_mesh']=False
file0.parameters["functions_share_mesh"] = True

which reuses the mesh for each time step and also for each function (also reducing file size).
written 6 months ago by Mike Welland  
1
Brilliant! This is precisely what I needed. Thank you!
written 6 months ago by Ever Ardo RB  
1
Thanks for this! I recently did some refactoring and ParaView would no longer draw streamlines. Setting file0.parameters["functions_share_mesh"] = True fixed the problem.

On the other hand, I use the adaptive solver, and setting file0.parameters['rewrite_function_mesh']=False seems to then only write the coarse mesh. Thankfully, I just skip this part, and "Plot Over Line" still works fine for me.
written 6 months ago by Alexander G. Zimmerman  
That makes sense that it only writes the first mesh I suppose. Are you using the adaptive mesh in parallel? What is your experience? Last time I tried it it didn't parallelise very well (but I haven't tried it recently)...
written 6 months ago by Mike Welland  
I suppose this is getting off topic, but shortly: Last I checked (a few months ago), the adaptive solver was not implemented for distributed-memory parallelism.
written 6 months ago by Alexander G. Zimmerman  
Great! The append attributes filter does the trick. Thank you so much!
written 6 months ago by Ever Ardo RB  
You are welcome. Please select the answer as the correct, it could help the others ;)
written 6 months ago by Michal Habera  
I might add that the 'save/load state' feature in Paraview is useful for setting up a pipeline for appending attributes and reusing it. (Set it up once, save state, then load it whenever you change the mesh or the files).

You might also check out cbcpost:
https://bitbucket.org/simula_cbc/cbcpost
written 6 months ago by Mike Welland  
Please login to add an answer/comment or follow this question.

Similar posts:
Search »