I require to read mesh and initial conditions from a file which is generated by another code (outside of FEniCSX). If I write them in a xdmf format, to read mesh is trivial, however, reading function to initialize fields is challenge. If I correctly understood the challenge is attributed to reading in parallel. If so, my question is that how I can read mesh and function in serial, and then FENICSX takes care of domain decomposition and parallel run (any minimal script or pseudo script is helpful)?
I would suggest to read:
You have not specified what the input field is tied to though, is it tied to the mesh nodes or the cell?
Hi dokken,
Thanks for your suggestion. The data is defined on vertexes of simple CG Hex elements (8-DoFs per element).
The target is initialize a function with a distance function that is computed using VTK library (distance field is computed from surface triangulated geometry; STL file).
I am curios, whether it is possible to do this job without using adios4dolfinx module?
You would have to implement the precise code that is provided by adios4dolfinx for communicating the data from an input process to the one owning the degree of freedom.
Thanks for your feedback. I just followed your first recipe and use the same scripte you refereed, i.e. this one. I did not modify the code but change input stuff based on my data. However, running the script I face the following error. I did not realize if there are something specific to a element types in this code to be modified. Could you please help me, whether I missed something etc (I put my python script together with data in a zip archive on this link if I can help to realize the problem):
Traceback (most recent call last):
File “/home/rtavakoli/Downloads/tesio/test_init.py”, line 59, in
uh = read_mesh_data(infile, mesh, data)
File “/home/rtavakoli/Downloads/tesio/test_init.py”, line 32, in read_mesh_data
element = basix.ufl.element(
File “/root/miniconda3/envs/dolfinx-checkpoint/lib/python3.10/site-packages/basix/ufl.py”, line 1671, in element
e = _basix.create_element(family, cell, degree, lagrange_variant, dpc_variant, discontinuous)
TypeError: create_element(): incompatible function arguments. The following argument types are supported:
1. (family_name: basix._basixcpp.ElementFamily, cell_name: basix._basixcpp.CellType, degree: int, lagrange_variant: basix._basixcpp.LagrangeVariant = <LagrangeVariant.unset: -1>, dpc_variant: basix._basixcpp.DPCVariant = <DPCVariant.unset: -1>, discontinuous: bool = False, dof_ordering: List[int] = ) → basix._basixcpp.FiniteElement
Invoked with: <ElementFamily.P: 1>, <CellType.hexahedron: 5>, 1, 1, <DPCVariant.unset: -1>, False
Traceback (most recent call last):
File “h5py/_objects.pyx”, line 201, in h5py._objects.ObjectID.dealloc
RuntimeError: Can’t decrement id ref count (can’t close file, there are objects still open)
Exception ignored in: ‘h5py._objects.ObjectID.dealloc’
Traceback (most recent call last):
File “h5py/_objects.pyx”, line 201, in h5py._objects.ObjectID.dealloc
RuntimeError: Can’t decrement id ref count (can’t close file, there are objects still open)
Many thanks for your consideration
I would need to know what version of dolfinx you are running. From
it looks like you are running through conda? so you are probably using v0.7.x?
I would appreciate if you could put your code here, instead of in a zip file (for the sake of transparency and accessibility in the future).
What I do assume in the link I referenced is that the data is associated with the vertices of the mesh, and that one has to map these to an appropriate 1st order Lagrange space.
Yes I used conda (conda install -n dolfinx-checkpoint -c conda-forge fenics-dolfinx mpich) and installed ADIOS4DOLFINx v0.7.3
And data are attributed to vertices, in fact number of data are identical to number of vertices, is it an issue?
The code is here (identical to your former comment but changing input data; however, it is not possible to attach data to this forum; I put the contents of xdmf file here too)
from pathlib import Path
import adios4dolfinx
import dolfinx
import numpy as np
from mpi4py import MPI
filename = "scaffold.xdmf"
encoding=dolfinx.io.XDMFFile.Encoding.HDF5
with dolfinx.io.XDMFFile(MPI.COMM_WORLD, filename, "r", encoding=encoding) as file:
mesh = file.read_mesh(name="Grid")
import basix
import h5py
def read_mesh_data(file_path:Path, mesh: dolfinx.mesh.Mesh, data_path:str):
assert file_path.is_file(), f"File {file_path} does not exist"
infile = h5py.File(file_path, "r", driver="mpio", comm=mesh.comm)
num_nodes_global = mesh.geometry.index_map().size_global
assert data_path in infile.keys(), f"Data {data_path} does not exist"
dataset = infile[data_path]
shape = dataset.shape
assert shape[0] == num_nodes_global, f"Got data of shape {shape}, expected {num_nodes_global, shape[1]}"
dtype = dataset.dtype
# Read data locally on each process
local_input_range = adios4dolfinx.comm_helpers.compute_local_range(mesh.comm, num_nodes_global)
local_input_data = dataset[local_input_range[0]:local_input_range[1]]
# Create appropriate function space (based on coordinate map)
assert len(mesh.geometry.cmaps) == 1, "Mixed cell-type meshes not supported"
element = basix.ufl.element(
basix.ElementFamily.P,
mesh.topology.cell_name(),
mesh.geometry.cmaps[0].degree,
mesh.geometry.cmaps[0].variant,
shape=(shape[1],),
gdim=mesh.geometry.dim)
# Assumption: Same doflayout for geometry and function space, cannot test in python
V = dolfinx.fem.FunctionSpace(mesh, element)
uh = dolfinx.fem.Function(V, name=data_path)
# Assume that mesh is first order for now
assert mesh.geometry.cmaps[0].degree == 1, "Only linear meshes supported"
x_dofmap = mesh.geometry.dofmap
igi = np.array(mesh.geometry.input_global_indices, dtype=np.int64)
global_geom_input = igi[x_dofmap]
global_geom_owner = adios4dolfinx.utils.index_owner(mesh.comm, global_geom_input.reshape(-1), num_nodes_global)
for i in range(shape[1]):
arr_i = adios4dolfinx.comm_helpers.send_dofs_and_recv_values(global_geom_input.reshape(-1), global_geom_owner, mesh.comm, local_input_data[:,i], local_input_range[0])
dof_pos = x_dofmap.reshape(-1)*shape[1]+i
uh.x.array[dof_pos] = arr_i
infile.close()
return uh
infile = Path("scaffold.h5")
for data in ["data2", "data3"]:
uh = read_mesh_data(infile, mesh, data)
with dolfinx.io.XDMFFile(mesh.comm, f"{data}.xdmf", "w") as xdmf:
xdmf.write_mesh(mesh)
xdmf.write_function(uh)
scaffold.xdmf file:
<Xdmf Version="3.0"><Domain><Grid Name="Grid"><Geometry GeometryType="XYZ"><DataItem DataType="Float" Dimensions="238328 3" Format="HDF" Precision="8">scaffold.h5:/data0</DataItem></Geometry><Topology TopologyType="Hexahedron" NumberOfElements="226981" NodesPerElement="8"><DataItem DataType="Int" Dimensions="226981 8" Format="HDF" Precision="4">scaffold.h5:/data1</DataItem></Topology><Attribute Name="signed-dist" AttributeType="Scalar" Center="Node"><DataItem DataType="Float" Dimensions="238328 1" Format="HDF" Precision="8">scaffold.h5:/data2</DataItem></Attribute><Attribute Name="level-set" AttributeType="Scalar" Center="Node"><DataItem DataType="Float" Dimensions="238328 1" Format="HDF" Precision="8">scaffold.h5:/data3</DataItem></Attribute></Grid></Domain></Xdmf>
Thanks
The same question. I use dolfinx v0.7.3. And I save the mesh and its function (tied to cells) to “.xdmf” and “.h5" file through the method dolfinx.io.XDMFFile.write_mesh and dolfinx.io.XDMFFile.write_function. Now I want to read the mesh data and the function to initialize a dolfinx.fem.Function object. What should I do with it？
You can’t you need to use ADIOS4dolfinx to read and write checkpoints. See:
http://jsdokken.com/adios4dolfinx/README.html
I cannot download the zip file any longer, and thus I cannot get your mesh.
Could you generate some dummy data that could be put somewhere accessible?
It also seems that you are having issues defining the finite element. For me to known how to solve this I would need:
- A mesh file
- Information about what dolfinx version you are running.