Discrepancy in mesh rendering/viewing

I am trying to create a mesh comprising of three subdomains. On viewing the mesh with Gmsh it shows the correct rendering…but the PyVista engine is only tagging the surfaces and not the volumes. Not sure if it is a PyVista problem or the model generation.


import pyvista
from dolfinx.io import gmshio, XDMFFile
from dolfinx import mesh, fem, plot, default_scalar_type
from dolfinx.fem.petsc import LinearProblem
from mpi4py import MPI
import ufl,gmsh
import numpy as np

scale =1
L = 600.0*scale
W = 200.0*scale
H=50*scale
delta = W / L
airL=600*scale
airW=300*scale
airT=200*scale
airG=50*scale
model_rank = 0
mesh_comm = MPI.COMM_WORLD


def gmsh_cantilever (model:gmsh.model, name:str) ->gmsh.model:
    airbox=model.occ.addBox(0,0,0,airL,airW,airT,1)
    cantilever=model.occ.addBox(0,airW/2-W/2,airG+H,L,W,H,2)
    base=model.occ.addBox(0,airW/2-W/2,0,L,W,H,3)
    dom0,dom1=model.occ.fragment([(3,1)],[(3,i) for i in range(2,4)],removeTool=False,removeObject=True)
    model.occ.removeAllDuplicates()
    model.occ.synchronize()

    model.addPhysicalGroup(3,[dom0[-1][1]],tag=1001)
    model.setPhysicalName(3,1001,"Air")
    model.addPhysicalGroup(3,[dom0[0][1]],tag=1002)
    model.setPhysicalName(3,1002,"Cantilever")
    model.addPhysicalGroup(3,[dom0[1][1]],tag=int(1003))
    model.setPhysicalName(3,1003,"Base")
    lcar1=0.1*scale
    lcar2=0.05*scale
    model.mesh.generate(dim=3)
    gmsh.fltk.run()
    return model

gmsh.initialize()
model=gmsh.model()
model =gmsh_cantilever(model,"Structure")

domain, ct, facet_tags = gmshio.model_to_mesh(model, mesh_comm, model_rank, gdim=3)
gmsh.finalize()

material_tags = np.unique(ct.values)
airregion = ct.find(1001)
beamregion = ct.find(1002)
baseregion = ct.find(1003)

with XDMFFile(MPI.COMM_WORLD, "mt.xdmf", "w") as xdmf:
    xdmf.write_mesh(domain)
    xdmf.write_meshtags(ct, domain.geometry)

plotter = pyvista.Plotter()
tdim = domain.topology.dim
topology, cell_types, geometry = plot.vtk_mesh(domain)

grid = pyvista.UnstructuredGrid(topology, cell_types, geometry)
num_local_cells = domain.topology.index_map(tdim).size_local
grid.cell_data["Marker"] = ct.values[ct.indices]
grid.cell_data["Opacity"] = np.where(ct.values[ct.indices]==1001,0.5,1)

grid.set_active_scalars("Marker")
actor = plotter.add_mesh(grid, show_edges=True, opacity="Opacity", cmap="coolwarm")

plotter.view_xy()
plotter.show()

I’d guess that it may be due to the opacity not working as you’d expect.
You export to xdmf file, just check if the markers are correct there with paraview, where you can easily clip by the scalar marker.

Alternatively, consider trying to set up the same “clip by scalar” in pyvista.
We will not show you how to do that here, because this is a forum about FEniCS, not about pyvista. However, it should be pretty easy to find out what to do just by looking at their documentation.

Thank you,
You were absolutely correct. Somehow opacity is not working as advertised.

This is what I get when I clip it.

For reference this is the code I used I clip.

grid = pyvista.UnstructuredGrid(topology, cell_types, geometry)
num_local_cells = domain.topology.index_map(tdim).size_local
grid.cell_data["Marker"] = ct.values
cmodel=grid.clip('y', invert =False)
grid.set_active_scalars("Marker")
actor = plotter.add_mesh(grid, style='wireframe', cmap="coolwarm")
actor1 =plotter.add_mesh(cmodel, label='clipped')
plotter.view_xy()
plotter.show()