BiV ellipsoid

In this example we will create a simple bi-ventricular ellipsoid model.

from pathlib import Path
from mpi4py import MPI
import pyvista
import dolfinx
import numpy as np
import cardiac_geometries

geodir = Path("biv_ellipsoid")
if not geodir.exists():
    cardiac_geometries.mesh.biv_ellipsoid(outdir=geodir, create_fibers=True, fiber_space="P_1")

2025-12-12 19:49:44 [debug    ] Convert file biv_ellipsoid/biv_ellipsoid.msh to dolfin

Info    : Reading 'biv_ellipsoid/biv_ellipsoid.msh'...
Info    : 64 entities
Info    : 1472 nodes
Info    : 7421 elements
Info    : Done reading 'biv_ellipsoid/biv_ellipsoid.msh'

If the folder already exist, then we just load the geometry

geo = cardiac_geometries.geometry.Geometry.from_folder(
    comm=MPI.COMM_WORLD,
    folder=geodir,
)

Next we will use pyvista to plot the mesh

vtk_mesh = dolfinx.plot.vtk_mesh(geo.mesh, geo.mesh.topology.dim)
grid = pyvista.UnstructuredGrid(*vtk_mesh)
plotter = pyvista.Plotter()
plotter.add_mesh(grid, show_edges=True)
plotter.view_xy()
if not pyvista.OFF_SCREEN:
    plotter.show()
else:
    figure = plotter.screenshot("biv_mesh.png")

2025-12-12 19:49:52.857 (  11.177s) [    7F4E5C1B4140]vtkXOpenGLRenderWindow.:1458  WARN| bad X server connection. DISPLAY=

The facets of the mesh are also marked with a MeshTags object. We can access the facet tags with geo.ffun.values and the markers with geo.markers.

print(geo.markers)

{'LV_EPI_FW': [1, 2], 'RV_EPI_FW': [2, 2], 'BASE': [3, 2], 'LV_ENDO_FW': [4, 2], 'RV_ENDO_FW': [5, 2], 'RV_SEPTUM': [6, 2], 'LV_SEPTUM': [7, 2], 'Myocardium': [1, 3]}

Next lets plot the facet tags with pyvista.

assert geo.ffun is not None
vtk_bmesh = dolfinx.plot.vtk_mesh(geo.mesh, geo.ffun.dim, geo.ffun.indices)
bgrid = pyvista.UnstructuredGrid(*vtk_bmesh)
bgrid.cell_data["Facet tags"] = geo.ffun.values
bgrid.set_active_scalars("Facet tags")
p = pyvista.Plotter(window_size=[800, 800])
p.add_mesh(bgrid, show_edges=True)
if not pyvista.OFF_SCREEN:
    p.show()
else:
    figure = p.screenshot("facet_tags_biv.png")

# Now let us look at the fibers of the LV ellipsoid. The fibers are stored in the `f0` attribute of the `Geometry` object.
assert geo.f0 is not None
topology, cell_types, geometry = dolfinx.plot.vtk_mesh(geo.f0.function_space)
values = np.zeros((geometry.shape[0], 3), dtype=np.float64)
values[:, : len(geo.f0)] = geo.f0.x.array.real.reshape((geometry.shape[0], len(geo.f0)))
function_grid = pyvista.UnstructuredGrid(topology, cell_types, geometry)
function_grid["u"] = values
glyphs = function_grid.glyph(orient="u", factor=1.0)
grid = pyvista.UnstructuredGrid(*vtk_mesh)
plotter = pyvista.Plotter()
plotter.add_mesh(grid, style="wireframe", color="r")
plotter.add_mesh(glyphs)
plotter.view_xy()
if not pyvista.OFF_SCREEN:
    plotter.show()
else:
    fig_as_array = plotter.screenshot("fiber_biv.png")

Now we could try to redo the exercise but with create a thinner RV wall and move the RV a bit closer to the LV. We can also try with different fiber angles, such that the fibers are oriented more longitudinally at the endocardium and more circumferentially at the epicardium.

geodir2 = Path("biv_ellipsoid2")
if not geodir2.exists():
    cardiac_geometries.mesh.biv_ellipsoid(
        outdir=geodir2,
        create_fibers=True,
        fiber_space="P_1",
        fiber_angle_endo=80,
        fiber_angle_epi=-30,
        rv_wall_thickness=0.2,
        rv_offset_x=2.0,
    )

2025-12-12 19:49:57 [debug    ] Convert file biv_ellipsoid2/biv_ellipsoid.msh to dolfin

Info    : Reading 'biv_ellipsoid2/biv_ellipsoid.msh'...
Info    : 66 entities
Info    : 1487 nodes
Info    : 7509 elements
Info    : Done reading 'biv_ellipsoid2/biv_ellipsoid.msh'

If the folder already exist, then we just load the geometry

geo2 = cardiac_geometries.geometry.Geometry.from_folder(
    comm=MPI.COMM_WORLD,
    folder=geodir2,
)

Next we will use pyvista to plot the mesh

vtk_mesh = dolfinx.plot.vtk_mesh(geo2.mesh, geo2.mesh.topology.dim)
grid = pyvista.UnstructuredGrid(*vtk_mesh)
plotter = pyvista.Plotter()
plotter.add_mesh(grid, show_edges=True)
plotter.view_xy()
if not pyvista.OFF_SCREEN:
    plotter.show()
else:
    figure = plotter.screenshot("biv_mesh2.png")

assert geo2.ffun is not None
vtk_bmesh = dolfinx.plot.vtk_mesh(geo2.mesh, geo2.ffun.dim, geo2.ffun.indices)
bgrid = pyvista.UnstructuredGrid(*vtk_bmesh)
bgrid.cell_data["Facet tags"] = geo2.ffun.values
bgrid.set_active_scalars("Facet tags")
p = pyvista.Plotter(window_size=[800, 800])
p.add_mesh(bgrid, show_edges=True)
if not pyvista.OFF_SCREEN:
    p.show()
else:
    figure = p.screenshot("facet_tags2_biv.png")

assert geo2.f0 is not None
topology, cell_types, geometry = dolfinx.plot.vtk_mesh(geo2.f0.function_space)
values = np.zeros((geometry.shape[0], 3), dtype=np.float64)
values[:, : len(geo2.f0)] = geo2.f0.x.array.real.reshape((geometry.shape[0], len(geo2.f0)))
function_grid = pyvista.UnstructuredGrid(topology, cell_types, geometry)
function_grid["u"] = values
glyphs = function_grid.glyph(orient="u", factor=1.0)
grid = pyvista.UnstructuredGrid(*vtk_mesh)
plotter = pyvista.Plotter()
plotter.add_mesh(grid, style="wireframe", color="r")
plotter.add_mesh(glyphs)
plotter.view_xy()
if not pyvista.OFF_SCREEN:
    plotter.show()
else:
    fig_as_array = plotter.screenshot("fiber2_biv.png")