Can dolfinx read a GMSH file containing 1D mesh of an ellipse?

dolfinx
1

I have created a 1D mesh of an ellipse, i.e., a 1D mesh immersed in 2D space, using the following .geo file in Gmsh.

//+
SetFactory("OpenCASCADE");
Ellipse(1) = {0, 0, 0, 5, 2.5, 0, 2*Pi};

Then I make a 1D mesh and export it to a MSH file. I use meshio to convert the MSH file to XDMF as follows

import meshio
m = meshio.read('ellipse.msh')
m.write('ellipse.xdmf')

I try to read this file in dolfinx as follows:

from mpi4py import MPI
from dolfinx import io

# Define parameters
meshfile = "ellipse.xdmf"

# Read the mesh
with io.XDMFFile(MPI.COMM_WORLD, meshfile, "r") as file:
    domain = file.read_mesh(name='Grid')

I get an error

Traceback (most recent call last):
  File "/home/fenics/galvanotaxis_ellipse.py", line 24, in <module>
    domain = file.read_mesh(name='Grid')
  File "/usr/local/dolfinx-real/lib/python3.10/dist-packages/dolfinx/io/utils.py", line 154, in read_mesh
    cell_shape, cell_degree = super().read_cell_type(name, xpath)
RuntimeError: Cannot recognise cell type. Unknown value: mixed

How can I read this mesh in dolfinx?

2

Noet that you should use Physical Tags when using GMSH, i.e.

SetFactory("OpenCASCADE");
Ellipse(1) = {0, 0, 0, 5, 2.5, 0, 2*Pi};
Physical Curve(1) = {1};

Then, in turn, when using meshio, you should use the script from
https://jorgensd.github.io/dolfinx-tutorial/chapter3/subdomains.html?highlight=meshio#read-in-msh-files-with-dolfinx
and then you have the full pipeline:

from dolfinx import io
from mpi4py import MPI
import meshio


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] if prune_z else mesh.points
    out_mesh = meshio.Mesh(points=points, cells={cell_type: cells}, cell_data={
                           "name_to_read": [cell_data]})
    return out_mesh


mesh = meshio.read("ellipse.msh")
meshio.write("ellipse.xdmf", create_mesh(mesh, "line", True))


# Define parameters
meshfile = "ellipse.xdmf"

# Read the mesh
with io.XDMFFile(MPI.COMM_WORLD, meshfile, "r") as file:
    domain = file.read_mesh(name='Grid')

print(domain.topology.cell_type, domain.topology.index_map(
    domain.topology.dim).size_local, domain.geometry.dim)

which returns

CellType.interval 22 2
1 Like
3

How can we also load physical tags of the nodes and subdomains?

4

See Defining subdomains for different materials — FEniCSx tutorial

Same procedure with any sub entity (i.e, same strategy for facet tags as for nodes).

5

The following MWE;

import numpy as np
from mpi4py import MPI
from dolfinx.io import gmshio
import gmsh

gmsh.initialize()
gdim = 1
mesh_comm = MPI.COMM_WORLD
model_rank = 0
if model_rank == 0:
    gmsh.clear()
    gmsh.option.setNumber("General.Terminal", 1)
    r_inner = 0.2
    r_outer = 0.25
    gap = 1E-5

    p1 = gmsh.model.occ.addPoint(0, 0, 0)
    p2 = gmsh.model.occ.addPoint(r_inner-gap, 0, 0)
    l1 = gmsh.model.occ.addLine(p1,p2,1)

    p3 = gmsh.model.occ.addPoint(r_inner+gap, 0, 0)
    p4 = gmsh.model.occ.addPoint(r_outer, 0, 0)
    l2 = gmsh.model.occ.addLine(p3,p4,2)
    gmsh.model.occ.synchronize()

    # Boundary Tags
    gmsh.model.addPhysicalGroup(0, [p1], 1) # Innest boundary
    gmsh.model.addPhysicalGroup(0, [p4], 2) # Outest boundary
    gmsh.model.addPhysicalGroup(0, [p2,p3], 0) # Outest boundary
    # Line Tags
    gmsh.model.addPhysicalGroup(1, [l1], tag=1) # Inner
    gmsh.model.addPhysicalGroup(1, [l2], tag=2) # Outer
    lc = 0.005 # 0.01
    gmsh.model.occ.synchronize()
    gmsh.option.setNumber("Mesh.MeshSizeMax", lc)
    gmsh.model.mesh.generate(1)
    # gmsh.fltk.run()

mesh, subdomains, facet_tags = gmshio.model_to_mesh(gmsh.model, mesh_comm, model_rank, gdim=gdim)
gmsh.clear()
gmsh.finalize()
MPI.COMM_WORLD.barrier()

says

Traceback (most recent call last):
  File "cg_mesh.py", line 40, in <module>
    mesh, subdomains, facet_tags = gmshio.model_to_mesh(gmsh.model, mesh_comm, model_rank, gdim=gdim)
  File "/home/ee331/Dev/Venvs/v060/lib/python3.8/dist-packages/dolfinx/io/gmshio.py", line 271, in model_to_mesh
    local_entities, local_values = _cpp.io.distribute_entity_data(
RuntimeError: Unrecognised cell type.

I think gmshio cannot read 1D physical tags of facets generated in gmsh.

6

That is why i pointed you to meshio for this, as I’ve only implemented it to read cell and facet tags.

7

Even if I use meshio I get the same error for this MWE;

import numpy as np
from mpi4py import MPI
from dolfinx.io import gmshio
import gmsh

gmsh.initialize()
gdim = 1
mesh_comm = MPI.COMM_WORLD
model_rank = 0
if model_rank == 0:
    gmsh.clear()
    gmsh.option.setNumber("General.Terminal", 1)
    r_inner = 0.2
    r_outer = 0.25
    gap = 1E-5

    p1 = gmsh.model.occ.addPoint(0, 0, 0)
    p2 = gmsh.model.occ.addPoint(r_inner-gap, 0, 0)
    l1 = gmsh.model.occ.addLine(p1,p2,1)

    p3 = gmsh.model.occ.addPoint(r_inner+gap, 0, 0)
    p4 = gmsh.model.occ.addPoint(r_outer, 0, 0)
    l2 = gmsh.model.occ.addLine(p3,p4,2)
    gmsh.model.occ.synchronize()

    # Boundary Tags
    gmsh.model.addPhysicalGroup(0, [p1], 1) # Innest boundary
    gmsh.model.addPhysicalGroup(0, [p4], 2) # Outest boundary
    gmsh.model.addPhysicalGroup(0, [p2,p3], 0) # Outest boundary
    # Line Tags
    gmsh.model.addPhysicalGroup(1, [l1], tag=1) # Inner
    gmsh.model.addPhysicalGroup(1, [l2], tag=2) # Outer
    lc = 0.005 # 0.01
    gmsh.model.occ.synchronize()
    gmsh.option.setNumber("Mesh.MeshSizeMax", lc)
    gmsh.model.mesh.generate(1)
    # gmsh.fltk.run()
    gmsh.write("mesh.msh")

# mesh, subdomains, facet_tags = gmshio.model_to_mesh(gmsh.model, mesh_comm, model_rank, gdim=gdim)
gmsh.clear()
gmsh.finalize()
MPI.COMM_WORLD.barrier()

import meshio
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] if prune_z else mesh.points
    out_mesh = meshio.Mesh(points=points, cells={cell_type: cells}, cell_data={"name_to_read":[cell_data]})
    return out_mesh

if MPI.COMM_WORLD.rank  == 0:
    # Read in mesh
    msh = meshio.read("mesh.msh")
   
    # Create and save one file for the mesh, and one file for the facets 
    line_mesh = create_mesh(msh, "line", prune_z=True)
    vertex_mesh = create_mesh(msh, "vertex", prune_z=True)
    meshio.write("mesh.xdmf", line_mesh)
    meshio.write("mt.xdmf", vertex_mesh)
    
from dolfinx.io import XDMFFile

with XDMFFile(MPI.COMM_WORLD, "mesh.xdmf", "r") as xdmf:
    mesh = xdmf.read_mesh(name="Grid")
    ct = xdmf.read_meshtags(mesh, name="Grid")
mesh.topology.create_connectivity(mesh.topology.dim, mesh.topology.dim-1)
with XDMFFile(MPI.COMM_WORLD, "mt.xdmf", "r") as xdmf:
    ft = xdmf.read_meshtags(mesh, name="Grid")

outputs;

  File "MWE_mesh.py", line 70, in <module>
    ft = xdmf.read_meshtags(mesh, name="Grid")
  File "/home/ee331/Dev/Venvs/v060/lib/python3.8/dist-packages/dolfinx/io/utils.py", line 177, in read_meshtags
    return super().read_meshtags(mesh, name, xpath)
RuntimeError: Unrecognised cell type.

It looks like I need to implement vertex tags and line tags using;

import numpy as np
def on_boundary(x):
    return np.isclose(x[0]), 0.2)
boundary_dofs = fem.locate_dofs_geometrical(V, on_boundary)

or something like that.

8

I’ve added a patch to DOLFINx which makes it possible to read vertex tags, see: Add point as cell type to be able to read in vertex tags by jorgensd · Pull Request #2549 · FEniCS/dolfinx · GitHub

1 Like
9

Does DOLFINx 0.6.0 support reading in vertex tags from gmsh? I have tried generating a 1D mesh using gmsh and pass it to gmshio.model_to_mesh function, but it still gives the error

in model_to_mesh
local_entities, local_values = _cpp.io.distribute_entity_data(
RuntimeError: Unrecognised cell type.
10

As I stated in february, this was fixed on the main branch: Add point as cell type to be able to read in vertex tags by jorgensd · Pull Request #2549 · FEniCS/dolfinx · GitHub
which did not make it into the 0.6.x release. This means that you have to use the main branch of DOLFINx, FFCx and Basix.