This is a read only copy of the old FEniCS QA forum. Please visit the new QA forum to ask questions

Consistent IO of Mesh and Function in parallel

+3 votes

Is there a working way how to store Mesh and Function in parallel for reading in future? Combination of XDMF for Mesh and XML for Function does not work because some entity indices gets messed. This is demonstrated by following snippet when run on 2 proccesses.

from dolfin import *

# output
mesh = UnitSquareMesh(5, 5)
V = FunctionSpace(mesh, 'CG', 1)
u = project(Expression('x[0]'), V)
File('mesh.xdmf') << mesh
File('u.xml') << u

# input
mesh = Mesh('mesh.xdmf')
V = FunctionSpace(mesh, 'CG', 1)
u = Function(V, 'u.xml')
plot(u, interactive=True)
asked Jun 22, 2013 by Jan Blechta FEniCS Expert (51,420 points)

2 Answers

+1 vote
 
Best answer

This code snippet works as expected with recent master - probably because of commit

Garth Wells  3541c5e  Merge branch 'chris/hdf5-cell-indices'  2013-06-21
answered Jun 22, 2013 by Jan Blechta FEniCS Expert (51,420 points)

The latest master now has some code for I/O of Functions in HDF5, which should allow you to write on "m" and read on "n" processes, m != n.
It would be good if it gets tested out...

Complete output code:

from dolfin import *

mesh = UnitSquareMesh(5, 5)
V = FunctionSpace(mesh, 'CG', 1)
u = project(Expression('x[0]'), V)

f = HDF5File('u.h5', 'w')
f.write(mesh, 'mesh')
f.write(u, 'u')

and input code:

from dolfin import *

f = HDF5File('u.h5', 'r')

mesh = Mesh()
f.read(mesh, 'mesh')

V = FunctionSpace(mesh, 'CG', 1)
u = Function(V)
f.read(u, 'u')

plot(u, interactive=True)

u_e = Expression('x[0]')
print 'errornorm=', assemble((u-u_e)*(u-u_e)*dx)

Seems to work ok even when output and input is run on different number of processes. Maybe this could be incorporated into test suite if you tweak dolfin/test/unit/test.py to prepend mpirun -n 3 to output code and mpirun -n 2 to input code and force it somehow to run it in necessary order. Maybe by manually invoking

mpirun -n 3 python output_code.py && mpirun -n 2 python input_code.py

Do you consider adding filename constructor for Mesh, MeshFunction and Function?

Some of that exists for XDMF, e.g.

mesh = Mesh("mfile.xdmf")

should read the first (only) Mesh in mfile.h5. The same should work for MeshFunction.

The problem with Function is that the output is fundamentally different for XDMF than HDF5, as it needs to be a visualisation format. HDF5 stores all the DOFs, so higher-order Functions can be reconstructed. In any case, a Function needs to define its Mesh and FunctionSpace before reading.

+3 votes

HDF5 support for Function IO will land in the development version of DOLFIN in the next few days. It works in parallel.

answered Jun 22, 2013 by Garth N. Wells FEniCS Expert (35,930 points)

can anyone point me to the correct dorsal commit for this functionality?

I don't know if HDF5 has been added to Dorsal, but it's easy to install on most Linux distributions, and on OSX with MacPorts or Homebrew. Just install the HDF5 development package (with MPI support) via a package manager.

...