Importing mesh gives error

dolfinx
1

Hi,
I am trying to import mesh file created by using third-party software. And I am simply trying to reproduce the example posted here: Hyperelasticity — FEniCSx tutorial by using externally (importing) created mesh.

# Scaled variable
L = 1
W = 0.2
mu = 1
rho = 1
delta = W/L
gamma = 0.4*delta**2
beta = 1.25
lambda_ = beta
g = gamma


import dolfinx
import numpy as np
from mpi4py import MPI
from dolfinx.cpp.mesh import CellType
from dolfinx.io import (XDMFFile, extract_gmsh_geometry,
                        extract_gmsh_topology_and_markers, ufl_mesh_from_gmsh)


print(f"DOLFINx version: {dolfinx.__version__} is installed")


#mesh = Mesh("solid-mec-00.xdmf");                                             

import dolfinx
import numpy as np
from mpi4py import MPI
from dolfinx.cpp.mesh import CellType


filename="solid-mec-00.xdmf"
with dolfinx.io.XDMFFile(MPI.COMM_WORLD, "solid-mec-00.xdmf", "r") as xdmf:
       mesh = xdmf.read_mesh(name="Grid")

def clamped_boundary(x):
    return np.isclose(x[0], 0)

fdim = mesh.topology.dim - 1
boundary_facets = dolfinx.mesh.locate_entities_boundary(mesh, fdim, clamped_boundary)

u_D = dolfinx.Function(V)
with u_D.vector.localForm() as loc:
    loc.set(0)
bc = dolfinx.DirichletBC(u_D, dolfinx.fem.locate_dofs_topological(V, fdim, boundary_facets))

T = dolfinx.Constant(mesh, (0, 0, 0))

import ufl
ds = ufl.Measure("ds", domain=mesh)

def epsilon(u):
    return ufl.sym(ufl.grad(u)) # Equivalent to 0.5*(ufl.nabla_grad(u) + ufl.nabla_grad(u).T)
def sigma(u):
    return lambda_ * ufl.nabla_div(u) * ufl.Identity(u.geometric_dimension()) + 2*mu*epsilon(u)

u = ufl.TrialFunction(V)
v = ufl.TestFunction(V)
f = dolfinx.Constant(mesh, (0, 0, -rho*g))
a = ufl.inner(sigma(u), epsilon(v)) * ufl.dx
L = ufl.dot(f, v) * ufl.dx + ufl.dot(T, v) * ds

problem = dolfinx.fem.LinearProblem(a, L, bcs=[bc], petsc_options={"ksp_type": "preonly", "pc_type": "lu"})
uh = problem.solve()

and I am getting the following error:

DOLFINx version: 0.3.1.0 is installed
Traceback (most recent call last):
  File "LE.py", line 39, in <module>
    mesh = xdmf.read_mesh(name="Grid")
  File "/python-venv/lib/python3.8/site-packages/dolfinx/io.py", line 50, in read_mesh
    cell_shape, cell_degree = super().read_cell_type(name, xpath)
RuntimeError: Cannot recognise cell type. Unknown value: mixed

Could you please tell me where I am making a mistake?
Thanks!

2

Dolfinx does not support xdmf files with multiple cell types in one mesh. See for instance section 2 of: Mesh generation and conversion with GMSH and PYGMSH | Jørgen S. Dokken

3

Hi,
Thanks! During the mesh creation, I do not name any surfaces or domains (I do not do the labeling)
As you say, for example, would the following work for tetrahedral elements?

import meshio
mesh_from_file = meshio.read("solid-mec-00.msh")

import numpy
def create_mesh(mesh, cell_type, prune_z=False):
    cells = mesh.get_cells_type(cell_type)
    cell_data = mesh.get_cell_data("gmsh:physical", cell_type)
    points = mesh.points[:,:2] 
    return out_mesh

triangle_mesh = create_mesh(mesh_from_file, "tetrahedron")
meshio.write("solid-mec-00.xdmf", tetrahedron_mesh)

I have just ignored the z-values in the mesh coordinates, it is a3D mesh.

Have you put out an example, how to import the external mesh, which is converted .xdmf from .msh.

I have used:

meshio-convert example.msh example.xdmf

And I would like to import the example.xdmf file to FEniCS-X.

4

See: Defining subdomains for different materials — FEniCSx tutorial
Note that it is using a slightly outdated version of meshio (v.4.4.6) which has the prune_z_0 command.
As long as you are converting 3D meshes (tetrahedral meshes) this should not affect you.

5

Thank you!

I have .cgns mesh (only tetrahedron elements), then I convert them into .vtk, from there I use the gmsh to get .msh file.
Finally I use meshio-convert .msh .xdmf

In that case, the following should work, but it is not?

import dolfinx.io
with dolfinx.io.XDMFFile(MPI.COMM_WORLD, "solid-mec-00.xdmf", "r") as xdmf:
    mesh = xdmf.read_mesh(name="Grid")

Or am I still missing some concept here?

6

I depends on what the msh conversion to change it.
Note that you can read vtk files directly with meshio and create the xdmf file. As you can see in the code above. The only thing that needs to change is how you fetch physical tags:

import meshio
def create_mesh(mesh, cell_type, prune_z=False):
    cells = mesh.get_cells_type(cell_type)
    # This line needs to get cell data from what ever the markers have been called in the vtk file.
    # Expore this by printing mesh.cell_data
    # cell_data = mesh.get_cell_data("gmsh:physical", cell_type)
    if prune_z:
         points = mesh.points[:,:2] 
    else:
         points = mesh.points
    out_mesh = meshio.Mesh(points=points, cells={cell_type: cells})
    #, cell_data={"name_to_read":[cell_data]})
    return out_mesh

in_mesh = meshio.read("mesh.vtk")
out_mesh = create_mesh(in_mesh, "tetra")
meshio.write("mesh.xdmf", out_mesh)
7

Thank you!

I have followed up with your suggestion. But, I am still getting some errors:

DOLFINx version: 0.3.1.0 is installed
WARNING:py.warnings:/fenicsx-iris-master-r23/python-venv/lib/python3.8/site-packages/h5py/__init__.py:36: UserWarning: h5py is running against HDF5 1.10.7 when it was built against 1.12.1, this may cause problems
  _warn(("h5py is running against HDF5 {0} when it was built against {1}, "

Warning! ***HDF5 library version mismatched error***
The HDF5 header files used to compile this application do not match
the version used by the HDF5 library to which this application is linked.
Data corruption or segmentation faults may occur if the application continues.
This can happen when an application was compiled by one version of HDF5 but
linked with a different version of static or shared HDF5 library.
You should recompile the application or check your shared library related
settings such as 'LD_LIBRARY_PATH'.
You can, at your own risk, disable this warning by setting the environment
variable 'HDF5_DISABLE_VERSION_CHECK' to a value of '1'.
Setting it to 2 or higher will suppress the warning messages totally.
Headers are 1.12.1, library is 1.10.7

.
.
.
Features:
---------
                   Parallel HDF5: yes
Parallel Filtered Dataset Writes: yes
              Large Parallel I/O: yes
              High-level library: yes
                Build HDF5 Tests: yes
                Build HDF5 Tools: yes
                    Threadsafety: no
             Default API mapping: v110
  With deprecated public symbols: yes
          I/O filters (external): deflate(zlib),szip(encoder)
                             MPE: 
                      Direct VFD: no
                      Mirror VFD: no
              (Read-Only) S3 VFD: no
            (Read-Only) HDFS VFD: no
                         dmalloc: no
  Packages w/ extra debug output: none
                     API tracing: no
            Using memory checker: no
 Memory allocation sanity checks: no
          Function stack tracing: no
                Use file locking: best-effort
       Strict file format checks: no
    Optimization instrumentation: no
Bye...
Aborted

import dolfinx
import numpy as np
from mpi4py import MPI
from dolfinx.cpp.mesh import CellType
from dolfinx.io import (XDMFFile, extract_gmsh_geometry,
                        extract_gmsh_topology_and_markers, ufl_mesh_from_gmsh)

print(f"DOLFINx version: {dolfinx.__version__} is installed")

import dolfinx
import numpy as np
from mpi4py import MPI
from dolfinx.cpp.mesh import CellType

import meshio
def create_mesh(mesh, cell_type, prune_z=False):
    cells = mesh.get_cells_type(cell_type)
    # This line needs to get cell data from what ever the markers have been called in the vtk file.
    # Expore this by printing mesh.cell_data
    # cell_data = mesh.get_cell_data("gmsh:physical", cell_type)
    if prune_z:
         points = mesh.points[:,:2] 
    else:
         points = mesh.points
    out_mesh = meshio.Mesh(points=points, cells={cell_type: cells})
    #, cell_data={"name_to_read":[cell_data]})
    return out_mesh

in_mesh = meshio.read("solid-mec-00.vtk")
out_mesh = create_mesh(in_mesh, "tetrahedron")
meshio.write("solid-mec-00.xdmf", out_mesh)

print(f"DOLFINx version: {dolfinx.__version__} is installed")
8

You should install h5py and meshio with the following options:

pip3 install --no-binary=h5py h5py meshio

as DOLFINx already has a version of HDF5 that it is working against.

See for instance: Problem with trying to get dolfin, meshio and GMSH python interface to coexist - #4 by BillS

9

Thanks! my pip list shows:

$$ pip list
.
h5py               3.5.0
.

And I have tried your suggestion:

$$ pip3 install --no-binary=h5py h5py meshio
Requirement already satisfied: h5py in /python-venv/lib/python3.8/site-packages (3.5.0)
Requirement already satisfied: meshio in /python-venv/lib/python3.8/site-packages (4.4.6)
Requirement already satisfied: numpy>=1.17.5 in /python-venv/lib/python3.8/site-packages (from h5py) (1.20.3)

Still, I am not able to read/convert the mesh.

10

You need to uninstall the current version of h5py first, ie. pip3 uninstall h5py

11

Thanks! sorry, my bad!

Now it is working with newly installed h5py. But however, I am getting the error in the mesh_out option.

Successfully built h5py
Installing collected packages: h5py
Successfully installed h5py-3.6.0
(python-venv) $ python3 mesh-read.py
DOLFINx version: 0.3.1.0 is installed
Traceback (most recent call last):
  File "mesh-read.py", line 42, in <module>
    meshio.write("solid-mec-00.xdmf", out_mesh)
  File "/python-venv/lib/python3.8/site-packages/meshio/_helpers.py", line 141, in write
    return writer(filename, mesh, **kwargs)
  File "/python-venv/lib/python3.8/site-packages/meshio/xdmf/main.py", line 546, in write
    XdmfWriter(*args, **kwargs)
  File "/python-venv/lib/python3.8/site-packages/meshio/xdmf/main.py", line 369, in __init__
    self.write_cells(mesh.cells, grid)
  File "/python-venv/lib/python3.8/site-packages/meshio/xdmf/main.py", line 429, in write_cells
    xdmf_type = meshio_to_xdmf_type[meshio_type][0]
KeyError: 'tetrahedron'
12

Try "tetra" instead of tetrahedron.

13

Thank you very much! it solved the issue.