Cell Tags in the Wrong Order

I am working to create a mesh with subdomains (based on here), and I am visualizing it to make sure the subdomain assignments are correct.

The above example uses the model_to_mesh function from dolfinx.io.gmshio but this appears to have been removed (I get the error dolfinx.io is not a package and there isn’t reference to gmshio in the docs). Because of this, I switched to adapting some other examples where the gmsh geometry is directly used in FEniCSx.

When I specify the cell tags from gmsh in the pyvista code, the colors of the different subdomains are wrong. Looking at the cell tags, the number of cells associated with the different subdomains seems appropriate, so it looks like it might just be an ordering issue. Are there ways to ensure that the cell ordering is the same between the gmsh mesh and the ufl/FEniCSx mesh? I have the gmsh MWE below and what it generates, and what the coloring should look like (hand-colored).

import gmsh
import numpy as np
from mpi4py import MPI

import matplotlib.pyplot as plt

rank = MPI.COMM_WORLD.rank

gmsh.initialize()

gdim = 2  # Geometric dimension of the mesh
model_rank = 0
mesh_comm = MPI.COMM_WORLD

spacing = 0.01

body_dx = 1
body_dy = body_dx

floor_dx = 1.6
floor_dy = 0.1

indenter_dx = 0.5
indenter_dy = 0.5

if mesh_comm.rank == model_rank:
    floor = gmsh.model.occ.addRectangle(-floor_dx/2,0,0,floor_dx,floor_dy)
    
    body1 = gmsh.model.occ.addRectangle(-body_dx/2,floor_dy+spacing,0,body_dx,body_dy)
    
    indenter = gmsh.model.occ.addRectangle(-indenter_dx/2,body_dy+floor_dy+2*spacing,0,indenter_dx,indenter_dy)
    
    background = gmsh.model.occ.addRectangle(-floor_dx/2,0,0,floor_dx,indenter_dy+body_dy+floor_dy+2*spacing)
    
    gmsh.model.occ.synchronize()
    
    bodies = [(2,floor),(2,body1),(2,indenter)]
    gmsh.model.occ.synchronize()

    whole_domain = gmsh.model.occ.fragment([(2,background)],bodies)
    gmsh.model.occ.synchronize()
    
    background_surfaces = []
    other_surfaces = []
    for i,domain in enumerate(whole_domain[0]):
        com = gmsh.model.occ.getCenterOfMass(domain[0],domain[1])
        gmsh.model.addPhysicalGroup(domain[0], [domain[1]], tag=i)
        if np.isclose(com[1],1.6002/2):
            background_surfaces.append(domain)
        else:
            other_surfaces.append(domain)
    
    gmsh.model.mesh.generate(2)
    

from dolfinx import fem, mesh, plot, io, cpp, nls, la
import pyvista

if MPI.COMM_WORLD.rank == 0:
    # Get mesh geometry
    geometry_data = io.extract_gmsh_geometry(gmsh.model)
    # Get mesh topology for each element
    topology_data = io.extract_gmsh_topology_and_markers(gmsh.model)
    
if MPI.COMM_WORLD.rank == 0:
    # Extract the cell type and number of nodes per cell and broadcast
    # it to the other processors 
    gmsh_cell_type = list(topology_data.keys())[0]    
    properties = gmsh.model.mesh.getElementProperties(gmsh_cell_type)
    name, dim, order, num_nodes, local_coords, _ = properties
    cells = topology_data[gmsh_cell_type]["topology"]
    cell_id, num_nodes = MPI.COMM_WORLD.bcast([gmsh_cell_type, num_nodes], root=0)
else:        
    cell_id, num_nodes = MPI.COMM_WORLD.bcast([None, None], root=0)
    cells, geometry_data = np.empty([0, num_nodes]), np.empty([0, gdim])
    
ufl_domain = io.ufl_mesh_from_gmsh(cell_id, 2)
gmsh_cell_perm = io.cell_perm_gmsh(cpp.mesh.to_type(str(ufl_domain.ufl_cell())), num_nodes)
cells = cells[:, gmsh_cell_perm]

# Create distributed mesh
domain = mesh.create_mesh(MPI.COMM_WORLD, cells, geometry_data[:, :gdim], ufl_domain)

import pyvista
pyvista.set_jupyter_backend("pythreejs")
from dolfinx.plot import create_vtk_mesh

plotter = pyvista.Plotter()
topology, cell_types, geometry = create_vtk_mesh(domain, domain.topology.dim)
grid = pyvista.UnstructuredGrid(topology, cell_types, geometry )

num_local_cells = domain.topology.index_map(domain.topology.dim).size_local

grid.cell_data["Marker"] = topology_data[2]['cell_data']
grid.set_active_scalars("Marker")
actor = plotter.add_mesh(grid, show_edges=True)
plotter.view_xy()
if not pyvista.OFF_SCREEN:
    plotter.show()
else:
    pyvista.start_xvfb()
    cell_tag_fig = plotter.screenshot("cell_tags.png")

What the Code Generates:

What I am Looking to Generate:

What version of DOLFINx are you running? The dolfinx.Io.gmshio functions were added to main last week, and are therefore not available through conda or apt-get. You can obtain this version by using docker images, see: An overview of the FEniCS Project — FEniCSx tutorial or by installing from source.

Short answer is: no
Long answer: GMSH runs in serial, it generated cells numbered from 0-N-1.
DOLFINx is created to run in parallel, and therefore distributes the cells on each process, and reorders them for data locality.

2 Likes

I am running the JupyterLab Docker Image from the GitHub repo which seems to still not have gmshio and returns the same error I mentioned above. I will take a look at your Docker image. (EDIT: using your Docker image instead allowed for gmshio and I was able to directly port the example. Thanks!)

And thank you for the clarification about GMSH v. FEniCS!

Remember to pull the docker images docker pull dolfinx/lab to make sure that they are updated.

My image is based of dolfinx/lab, so it should be in there.
You can also print the git commit hash

import dolfinx
print(f"DOLFINx version: {dolfinx.__version__} based on GIT commit: {dolfinx.git_commit_hash} of https://github.com/FEniCS/dolfinx/")

to get which commit it is based on:
https://github.com/FEniCS/dolfinx/commit/insert-commit-hash-here