Hi all,
I am trying to load a mesh and a scalar function from a .xdmf/.h5 pair in 3D, assign the function data to a new Function, and then save it again for visualization in ParaView. However, the result looks visually very different from the original. I suspect it might be due to incorrect DOF ordering, missing scatter_forward(), or mismatch between the original and reconstructed mesh.
Here is a minimal example of what I’m doing:
from dolfinx.io import XDMFFile
from dolfinx.fem import Function, functionspace
from basix.ufl import element
from mpi4py import MPI
import h5py
# Load the mesh
with XDMFFile(MPI.COMM_WORLD, "3d_example.xdmf", "r") as xdmf:
mesh = xdmf.read_mesh()
P0 = element("Lagrange", mesh.topology.cell_name(), 1)
V0 = functionspace(mesh, P0)
# Explore the HDF5 file structure
with h5py.File("3d_example.h5", "r") as h5file:
def print_hdf5_structure(name, obj):
print(name)
h5file.visititems(print_hdf5_structure)
# Load the function data manually
with h5py.File("3d_example.h5", "r") as datafile:
data = datafile["Function/u"][:]
print("Data shape:", data.shape)
data = data[:, 0] # Assuming shape (N, 1)
# Assign to a new Function
u = Function(V0)
u.x.array[:] = data
u.x.scatter_forward()
# Save again
with XDMFFile(MPI.COMM_WORLD, "3d_example_new.xdmf", "w") as file1:
file1.write_mesh(mesh)
file1.write_function(u)
When I visualize 3d_example_new.xdmf in ParaView, the function looks completely different from the original output (even though I believe the mesh and data match). In the following, 3d_example.xdmf on the left and 3d_example_new.xdmf on the right.
In the original simulation, the mesh was created using:
from dolfinx.mesh import create_box, CellType
import numpy as np
mesh = create_box(MPI.COMM_WORLD, [np.array([0, 0, 0]), np.array([L, W, H])], [nx, ny, nz], CellType.tetrahedron)
So it is a structured tetrahedral mesh.
Has anyone encountered similar issues?
Could this be due to DOF reordering between the original FunctionSpace and the reconstructed one?
Is there a better or safer way to reconstruct a Function from .xdmf/.h5 if read_function() is not available in FEniCSx?
Any help or insight would be greatly appreciated.
Many thanks!