BiV Ellipsoid

BiV Ellipsoid#

In this example we will simulate an idealized bi-ventricular geometry. For this we will use cardiac-geometries to generate an idealized LV ellipsoid.

First lets do some imports

from pathlib import Path
from mpi4py import MPI
import numpy as np
import dolfinx
from dolfinx import log
import ldrb
import cardiac_geometries
import cardiac_geometries.geometry
import pulse

and lets turn on logging so that we can see more info from dolfinx

log.set_log_level(log.LogLevel.INFO)

Now we create the geometry using cardiac-geometries and save it to a folder called biv_ellipsoid. We will also create fiber orientations using the Laplace-Dirichlet Rule based (LDRB) algorithm, using the library fenicsx-ldrb package

outdir = Path("biv_ellipsoid")
outdir.mkdir(parents=True, exist_ok=True)
geodir = outdir / "geometry"
if not geodir.exists():
    geo = cardiac_geometries.mesh.biv_ellipsoid(outdir=geodir)
    markers = cardiac_geometries.mesh.transform_biv_markers(geo.markers)

    system = ldrb.dolfinx_ldrb(mesh=geo.mesh, ffun=geo.ffun, markers=markers, alpha_endo_lv=60, alpha_epi_lv=-60, beta_endo_lv=0, beta_epi_lv=0, fiber_space="P_2")
    cardiac_geometries.fibers.utils.save_microstructure(mesh=geo.mesh, functions=[system.f0, system.s0, system.n0], outdir=geodir)
2025-11-03 11:20:06 [debug    ] Convert file biv_ellipsoid/geometry/biv_ellipsoid.msh to dolfin
Info    : Reading 'biv_ellipsoid/geometry/biv_ellipsoid.msh'...
Info    : 64 entities
Info    : 1472 nodes
Info    : 7421 elements
Info    : Done reading 'biv_ellipsoid/geometry/biv_ellipsoid.msh'
[2025-11-03 11:20:06.490] [info] Extract basic topology: 17924->17924
[2025-11-03 11:20:06.490] [info] Build local dual graph
[2025-11-03 11:20:06.490] [info] Build local part of mesh dual graph (mixed)
[2025-11-03 11:20:06.493] [info] GPS pseudo-diameter:(72) 2399-4155
[2025-11-03 11:20:06.494] [info] Create topology (single cell type)
[2025-11-03 11:20:06.494] [info] Create topology (generalised)
[2025-11-03 11:20:06.494] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:06.494] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:06.495] [info] Compute ghost indices
[2025-11-03 11:20:06.495] [info] Computing communication graph edges (using PCX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.495] [info] Finished graph edge discovery using PCX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.496] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.496] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.496] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.496] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.496] [info] Number of neighbourhood source ranks in distribute_to_postoffice: 0
[2025-11-03 11:20:06.496] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.496] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.496] [info] Neighbourhood destination ranks from post office in distribute_data (rank, num dests, num dests/mpi_size): 0, 0, 0
[2025-11-03 11:20:06.496] [info] Checking required entities per dimension
[2025-11-03 11:20:06.496] [info] Cell type: 0 dofmap: 4481x4
[2025-11-03 11:20:06.497] [info] Global index computation
[2025-11-03 11:20:06.497] [info] Got 1 index_maps
[2025-11-03 11:20:06.497] [info] Get global indices
[2025-11-03 11:20:06.497] [info] XDMF distribute entity data
[2025-11-03 11:20:06.497] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:06.497] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:06.497] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:06.497] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:06.498] [info] XDMF build map
[2025-11-03 11:20:06.500] [info] Requesting connectivity (3, 0) - (0, 0)
[2025-11-03 11:20:06.500] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-11-03 11:20:06.500] [info] Build list of mesh entity indices from the entity vertices.
[2025-11-03 11:20:06.502] [info] XDMF distribute entity data
[2025-11-03 11:20:06.502] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:06.502] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:06.502] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:06.502] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:06.502] [info] XDMF build map
[2025-11-03 11:20:06.504] [info] Computing mesh entities of dimension 2
[2025-11-03 11:20:06.506] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.506] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.506] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-11-03 11:20:06.506] [info] Computing mesh connectivity 2-3 from transpose.
[2025-11-03 11:20:06.506] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-11-03 11:20:06.506] [info] Build list of mesh entity indices from the entity vertices.
[2025-11-03 11:20:06.510] [info] Opened HDF5 file with id "72057594037927936"
[2025-11-03 11:20:06.510] [info] Adding mesh to node "/Xdmf/Domain"
[2025-11-03 11:20:06.510] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-11-03 11:20:06.510] [info] Adding geometry data to node "/Xdmf/Domain/Grid"
[2025-11-03 11:20:06.510] [info] XDMF: add meshtags (Cell tags)
[2025-11-03 11:20:06.510] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-11-03 11:20:06.511] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-11-03 11:20:06.511] [info] XDMF: add meshtags (Facet tags)
[2025-11-03 11:20:06.511] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-11-03 11:20:06.511] [info] Computing mesh entities of dimension 1
[2025-11-03 11:20:06.514] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.514] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.514] [info] Requesting connectivity (1, 0) - (3, 0)
[2025-11-03 11:20:06.514] [info] Computing mesh connectivity 1-3 from transpose.
[2025-11-03 11:20:06.514] [info] XDMF: add meshtags (Edge tags)
[2025-11-03 11:20:06.514] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-11-03 11:20:06.514] [info] Requesting connectivity (0, 0) - (3, 0)
[2025-11-03 11:20:06.514] [info] Computing mesh connectivity 0-3 from transpose.
[2025-11-03 11:20:06.514] [info] XDMF: add meshtags (Vertex tags)
[2025-11-03 11:20:06.514] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-11-03 11:20:06.515] [info] Opened HDF5 file with id "72057594037927937"
[2025-11-03 11:20:06.515] [info] Read topology data "Mesh" at /Xdmf/Domain
[2025-11-03 11:20:06.515] [info] HDF5 Read data rate: 6697.431106959365 MB/s
[2025-11-03 11:20:06.515] [info] IO permuting cells
[2025-11-03 11:20:06.516] [info] Read geometry data "Mesh" at /Xdmf/Domain
[2025-11-03 11:20:06.516] [info] HDF5 Read data rate: 3120.2967673555913 MB/s
[2025-11-03 11:20:06.519] [info] Using partitioner with 17924 cell data
[2025-11-03 11:20:06.519] [info] Compute partition of cells across ranks
[2025-11-03 11:20:06.519] [info] Building mesh dual graph
[2025-11-03 11:20:06.519] [info] Build local part of mesh dual graph (mixed)
[2025-11-03 11:20:06.521] [info] Build nonlocal part of mesh dual graph
[2025-11-03 11:20:06.521] [info] Graph edges (local: 14984, non-local: 0)
[2025-11-03 11:20:06.521] [info] Compute graph partition using PT-SCOTCH
[2025-11-03 11:20:06.522] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:06.522] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:06.522] [info] Extract basic topology: 17924->17924
[2025-11-03 11:20:06.522] [info] Build local dual graph
[2025-11-03 11:20:06.522] [info] Build local part of mesh dual graph (mixed)
[2025-11-03 11:20:06.524] [info] GPS pseudo-diameter:(72) 4478-1
[2025-11-03 11:20:06.525] [info] Create topology (single cell type)
[2025-11-03 11:20:06.525] [info] Create topology (generalised)
[2025-11-03 11:20:06.525] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:06.525] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:06.526] [info] Compute ghost indices
[2025-11-03 11:20:06.526] [info] Computing communication graph edges (using PCX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.526] [info] Finished graph edge discovery using PCX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.527] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.527] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.527] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.527] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.527] [info] Number of neighbourhood source ranks in distribute_to_postoffice: 0
[2025-11-03 11:20:06.527] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.527] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.527] [info] Neighbourhood destination ranks from post office in distribute_data (rank, num dests, num dests/mpi_size): 0, 0, 0
[2025-11-03 11:20:06.527] [info] Checking required entities per dimension
[2025-11-03 11:20:06.527] [info] Cell type: 0 dofmap: 4481x4
[2025-11-03 11:20:06.527] [info] Global index computation
[2025-11-03 11:20:06.527] [info] Got 1 index_maps
[2025-11-03 11:20:06.527] [info] Get global indices
[2025-11-03 11:20:06.528] [info] Requesting connectivity (3, 0) - (3, 0)
[2025-11-03 11:20:06.528] [info] XDMF read meshtags (Cell tags)
[2025-11-03 11:20:06.528] [info] Read topology data "Cell tags" at /Xdmf/Domain
[2025-11-03 11:20:06.528] [info] HDF5 Read data rate: 6681.826654240447 MB/s
[2025-11-03 11:20:06.528] [info] IO permuting cells
[2025-11-03 11:20:06.528] [info] HDF5 Read data rate: 1972.4881699130626 MB/s
[2025-11-03 11:20:06.528] [info] IO permuting cells
[2025-11-03 11:20:06.528] [info] XDMF distribute entity data
[2025-11-03 11:20:06.528] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:06.528] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:06.529] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:06.529] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:06.529] [info] XDMF build map
[2025-11-03 11:20:06.530] [info] XDMF create meshtags
[2025-11-03 11:20:06.530] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-11-03 11:20:06.530] [info] Build list of mesh entity indices from the entity vertices.
[2025-11-03 11:20:06.532] [info] Computing mesh entities of dimension 2
[2025-11-03 11:20:06.534] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.534] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.534] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-11-03 11:20:06.534] [info] Computing mesh connectivity 2-3 from transpose.
[2025-11-03 11:20:06.535] [info] XDMF read meshtags (Facet tags)
[2025-11-03 11:20:06.535] [info] Read topology data "Facet tags" at /Xdmf/Domain
[2025-11-03 11:20:06.535] [info] HDF5 Read data rate: 3630.1898441117455 MB/s
[2025-11-03 11:20:06.535] [info] IO permuting cells
[2025-11-03 11:20:06.535] [info] HDF5 Read data rate: 1171.5481171548117 MB/s
[2025-11-03 11:20:06.535] [info] IO permuting cells
[2025-11-03 11:20:06.535] [info] XDMF distribute entity data
[2025-11-03 11:20:06.535] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:06.535] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:06.535] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:06.535] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:06.535] [info] XDMF build map
[2025-11-03 11:20:06.536] [info] XDMF create meshtags
[2025-11-03 11:20:06.536] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-11-03 11:20:06.536] [info] Build list of mesh entity indices from the entity vertices.
[2025-11-03 11:20:06.539] [info] Computing mesh entities of dimension 1
[2025-11-03 11:20:06.541] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.541] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.542] [info] Requesting connectivity (1, 0) - (3, 0)
[2025-11-03 11:20:06.542] [info] Computing mesh connectivity 1-3 from transpose.
[2025-11-03 11:20:06.542] [info] XDMF read meshtags (Edge tags)
[2025-11-03 11:20:06.542] [info] Read topology data "Edge tags" at /Xdmf/Domain
[2025-11-03 11:20:06.542] [info] HDF5 Read data rate: 0 MB/s
[2025-11-03 11:20:06.542] [info] IO permuting cells
[2025-11-03 11:20:06.542] [info] HDF5 Read data rate: 0 MB/s
[2025-11-03 11:20:06.542] [info] IO permuting cells
[2025-11-03 11:20:06.542] [info] XDMF distribute entity data
[2025-11-03 11:20:06.542] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.542] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.542] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:06.542] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:06.542] [info] XDMF build map
[2025-11-03 11:20:06.543] [info] XDMF create meshtags
[2025-11-03 11:20:06.543] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-11-03 11:20:06.543] [info] Build list of mesh entity indices from the entity vertices.
[2025-11-03 11:20:06.544] [info] Requesting connectivity (0, 0) - (3, 0)
[2025-11-03 11:20:06.544] [info] Computing mesh connectivity 0-3 from transpose.
[2025-11-03 11:20:06.544] [info] XDMF read meshtags (Vertex tags)
[2025-11-03 11:20:06.544] [info] Read topology data "Vertex tags" at /Xdmf/Domain
[2025-11-03 11:20:06.544] [info] HDF5 Read data rate: 0 MB/s
[2025-11-03 11:20:06.544] [info] IO permuting cells
[2025-11-03 11:20:06.544] [info] HDF5 Read data rate: 0 MB/s
[2025-11-03 11:20:06.544] [info] IO permuting cells
[2025-11-03 11:20:06.544] [info] XDMF distribute entity data
[2025-11-03 11:20:06.544] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:06.544] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:06.544] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:06.544] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:06.545] [info] XDMF build map
[2025-11-03 11:20:06.545] [info] XDMF create meshtags
[2025-11-03 11:20:06.545] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-11-03 11:20:06.545] [info] Build list of mesh entity indices from the entity vertices.
[2025-11-03 11:20:06.547] [info] Checking required entities per dimension
[2025-11-03 11:20:06.547] [info] Cell type: 0 dofmap: 4481x4
[2025-11-03 11:20:06.548] [info] Global index computation
[2025-11-03 11:20:06.548] [info] Got 1 index_maps
[2025-11-03 11:20:06.548] [info] Get global indices
[2025-11-03 11:20:06.549] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-11-03 11:20:06.552] [info] Column ghost size increased from 0 to 0
[2025-11-03 11:20:06.564] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-11-03 11:20:06.565] [info] Requesting connectivity (2, 0) - (0, 0)
[2025-11-03 11:20:06.565] [info] Requesting connectivity (2, 0) - (0, 0)
[2025-11-03 11:20:06.565] [info] Requesting connectivity (0, 0) - (3, 0)
[2025-11-03 11:20:06.565] [info] Requesting connectivity (3, 0) - (0, 0)
[2025-11-03 11:20:06.569] [info] Column ghost size increased from 0 to 0
[2025-11-03 11:20:06.584] [info] Column ghost size increased from 0 to 0
[2025-11-03 11:20:06.594] [info] Column ghost size increased from 0 to 0
[2025-11-03 11:20:06.605] [info] Column ghost size increased from 0 to 0
[2025-11-03 11:20:06.614] [info] Column ghost size increased from 0 to 0
[2025-11-03 11:20:06.624] [info] Checking required entities per dimension
[2025-11-03 11:20:06.624] [info] Cell type: 0 dofmap: 4481x10
[2025-11-03 11:20:06.624] [info] Global index computation
[2025-11-03 11:20:06.625] [info] Got 2 index_maps
[2025-11-03 11:20:06.625] [info] Get global indices
[2025-11-03 11:20:06.626] [info] Checking required entities per dimension
[2025-11-03 11:20:06.626] [info] Cell type: 0 dofmap: 4481x10
[2025-11-03 11:20:06.627] [info] Global index computation
[2025-11-03 11:20:06.627] [info] Got 2 index_maps
[2025-11-03 11:20:06.627] [info] Get global indices
[2025-11-03 11:20:07.409] [info] Checking required entities per dimension
[2025-11-03 11:20:07.409] [info] Cell type: 0 dofmap: 4481x10
[2025-11-03 11:20:07.409] [info] Global index computation
[2025-11-03 11:20:07.409] [info] Got 2 index_maps
[2025-11-03 11:20:07.409] [info] Get global indices
[2025-11-03 11:20:07.411] [info] Checking required entities per dimension
[2025-11-03 11:20:07.411] [info] Cell type: 0 dofmap: 4481x10
[2025-11-03 11:20:07.412] [info] Global index computation
[2025-11-03 11:20:07.412] [info] Got 2 index_maps
[2025-11-03 11:20:07.412] [info] Get global indices
[2025-11-03 11:20:08.487] [info] Checking required entities per dimension
[2025-11-03 11:20:08.487] [info] Cell type: 0 dofmap: 4481x10
[2025-11-03 11:20:08.487] [info] Global index computation
[2025-11-03 11:20:08.487] [info] Got 2 index_maps
[2025-11-03 11:20:08.488] [info] Get global indices
[2025-11-03 11:20:08.489] [info] Checking required entities per dimension
[2025-11-03 11:20:08.489] [info] Cell type: 0 dofmap: 4481x10
[2025-11-03 11:20:08.490] [info] Global index computation
[2025-11-03 11:20:08.490] [info] Got 2 index_maps
[2025-11-03 11:20:08.490] [info] Get global indices
[2025-11-03 11:20:08.494] [info] Compute face permutations
[2025-11-03 11:20:08.494] [info] Computing permutations for face type 0
[2025-11-03 11:20:08.495] [info] Compute edge permutations

If the folder exist we just load it

[2025-11-03 11:20:08.511] [info] Opened HDF5 file with id "72057594037927938"
[2025-11-03 11:20:08.511] [info] Read topology data "Mesh" at /Xdmf/Domain
[2025-11-03 11:20:08.511] [info] HDF5 Read data rate: 3781.334880409272 MB/s
[2025-11-03 11:20:08.511] [info] IO permuting cells
[2025-11-03 11:20:08.511] [info] Read geometry data "Mesh" at /Xdmf/Domain
[2025-11-03 11:20:08.511] [info] HDF5 Read data rate: 2393.9825167717017 MB/s
[2025-11-03 11:20:08.515] [info] Using partitioner with 17924 cell data
[2025-11-03 11:20:08.515] [info] Compute partition of cells across ranks
[2025-11-03 11:20:08.515] [info] Building mesh dual graph
[2025-11-03 11:20:08.515] [info] Build local part of mesh dual graph (mixed)
[2025-11-03 11:20:08.517] [info] Build nonlocal part of mesh dual graph
[2025-11-03 11:20:08.517] [info] Graph edges (local: 14984, non-local: 0)
[2025-11-03 11:20:08.517] [info] Compute graph partition using PT-SCOTCH
[2025-11-03 11:20:08.517] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:08.517] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:08.518] [info] Extract basic topology: 17924->17924
[2025-11-03 11:20:08.518] [info] Build local dual graph
[2025-11-03 11:20:08.518] [info] Build local part of mesh dual graph (mixed)
[2025-11-03 11:20:08.520] [info] GPS pseudo-diameter:(72) 4478-1
[2025-11-03 11:20:08.521] [info] Create topology (single cell type)
[2025-11-03 11:20:08.521] [info] Create topology (generalised)
[2025-11-03 11:20:08.521] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:08.521] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:08.521] [info] Compute ghost indices
[2025-11-03 11:20:08.521] [info] Computing communication graph edges (using PCX algorithm). Number of input edges: 0
[2025-11-03 11:20:08.521] [info] Finished graph edge discovery using PCX algorithm. Number of discovered edges 0
[2025-11-03 11:20:08.522] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:08.522] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:08.523] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:08.523] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:08.523] [info] Number of neighbourhood source ranks in distribute_to_postoffice: 0
[2025-11-03 11:20:08.523] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:08.523] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:08.523] [info] Neighbourhood destination ranks from post office in distribute_data (rank, num dests, num dests/mpi_size): 0, 0, 0
[2025-11-03 11:20:08.523] [info] Checking required entities per dimension
[2025-11-03 11:20:08.523] [info] Cell type: 0 dofmap: 4481x4
[2025-11-03 11:20:08.523] [info] Global index computation
[2025-11-03 11:20:08.523] [info] Got 1 index_maps
[2025-11-03 11:20:08.523] [info] Get global indices
[2025-11-03 11:20:08.524] [info] Requesting connectivity (3, 0) - (3, 0)
[2025-11-03 11:20:08.524] [info] XDMF read meshtags (Cell tags)
[2025-11-03 11:20:08.524] [info] Read topology data "Cell tags" at /Xdmf/Domain
[2025-11-03 11:20:08.524] [info] HDF5 Read data rate: 4911.694183736385 MB/s
[2025-11-03 11:20:08.524] [info] IO permuting cells
[2025-11-03 11:20:08.524] [info] HDF5 Read data rate: 1250.1918113970844 MB/s
[2025-11-03 11:20:08.524] [info] IO permuting cells
[2025-11-03 11:20:08.524] [info] XDMF distribute entity data
[2025-11-03 11:20:08.525] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:08.525] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:08.525] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:08.525] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:08.525] [info] XDMF build map
[2025-11-03 11:20:08.527] [info] XDMF create meshtags
[2025-11-03 11:20:08.527] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-11-03 11:20:08.527] [info] Build list of mesh entity indices from the entity vertices.
[2025-11-03 11:20:08.529] [info] Computing mesh entities of dimension 2
[2025-11-03 11:20:08.531] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:08.531] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:08.532] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-11-03 11:20:08.532] [info] Computing mesh connectivity 2-3 from transpose.
[2025-11-03 11:20:08.532] [info] XDMF read meshtags (Facet tags)
[2025-11-03 11:20:08.532] [info] Read topology data "Facet tags" at /Xdmf/Domain
[2025-11-03 11:20:08.532] [info] HDF5 Read data rate: 4506.034868126955 MB/s
[2025-11-03 11:20:08.532] [info] IO permuting cells
[2025-11-03 11:20:08.532] [info] HDF5 Read data rate: 1399.1671624033315 MB/s
[2025-11-03 11:20:08.532] [info] IO permuting cells
[2025-11-03 11:20:08.532] [info] XDMF distribute entity data
[2025-11-03 11:20:08.532] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:08.532] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:08.532] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:08.532] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:08.532] [info] XDMF build map
[2025-11-03 11:20:08.534] [info] XDMF create meshtags
[2025-11-03 11:20:08.534] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-11-03 11:20:08.534] [info] Build list of mesh entity indices from the entity vertices.
[2025-11-03 11:20:08.536] [info] Computing mesh entities of dimension 1
[2025-11-03 11:20:08.538] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:08.538] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:08.539] [info] Requesting connectivity (1, 0) - (3, 0)
[2025-11-03 11:20:08.539] [info] Computing mesh connectivity 1-3 from transpose.
[2025-11-03 11:20:08.539] [info] XDMF read meshtags (Edge tags)
[2025-11-03 11:20:08.539] [info] Read topology data "Edge tags" at /Xdmf/Domain
[2025-11-03 11:20:08.539] [info] HDF5 Read data rate: 0 MB/s
[2025-11-03 11:20:08.539] [info] IO permuting cells
[2025-11-03 11:20:08.539] [info] HDF5 Read data rate: 0 MB/s
[2025-11-03 11:20:08.539] [info] IO permuting cells
[2025-11-03 11:20:08.539] [info] XDMF distribute entity data
[2025-11-03 11:20:08.539] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:08.539] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:08.539] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:08.539] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:08.539] [info] XDMF build map
[2025-11-03 11:20:08.540] [info] XDMF create meshtags
[2025-11-03 11:20:08.540] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-11-03 11:20:08.540] [info] Build list of mesh entity indices from the entity vertices.
[2025-11-03 11:20:08.541] [info] Requesting connectivity (0, 0) - (3, 0)
[2025-11-03 11:20:08.541] [info] Computing mesh connectivity 0-3 from transpose.
[2025-11-03 11:20:08.542] [info] XDMF read meshtags (Vertex tags)
[2025-11-03 11:20:08.542] [info] Read topology data "Vertex tags" at /Xdmf/Domain
[2025-11-03 11:20:08.542] [info] HDF5 Read data rate: 0 MB/s
[2025-11-03 11:20:08.542] [info] IO permuting cells
[2025-11-03 11:20:08.542] [info] HDF5 Read data rate: 0 MB/s
[2025-11-03 11:20:08.542] [info] IO permuting cells
[2025-11-03 11:20:08.542] [info] XDMF distribute entity data
[2025-11-03 11:20:08.542] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-11-03 11:20:08.542] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-11-03 11:20:08.542] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-11-03 11:20:08.542] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-11-03 11:20:08.542] [info] XDMF build map
[2025-11-03 11:20:08.543] [info] XDMF create meshtags
[2025-11-03 11:20:08.543] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-11-03 11:20:08.543] [info] Build list of mesh entity indices from the entity vertices.
[2025-11-03 11:20:08.545] [info] Checking required entities per dimension
[2025-11-03 11:20:08.545] [info] Cell type: 0 dofmap: 4481x10
[2025-11-03 11:20:08.545] [info] Global index computation
[2025-11-03 11:20:08.545] [info] Got 2 index_maps
[2025-11-03 11:20:08.545] [info] Get global indices
[2025-11-03 11:20:08.546] [info] Compute face permutations
[2025-11-03 11:20:08.546] [info] Computing permutations for face type 0
[2025-11-03 11:20:08.547] [info] Compute edge permutations
[2025-11-03 11:20:08.558] [info] Checking required entities per dimension
[2025-11-03 11:20:08.558] [info] Cell type: 0 dofmap: 4481x10
[2025-11-03 11:20:08.559] [info] Global index computation
[2025-11-03 11:20:08.559] [info] Got 2 index_maps
[2025-11-03 11:20:08.559] [info] Get global indices
[2025-11-03 11:20:08.570] [info] Checking required entities per dimension
[2025-11-03 11:20:08.570] [info] Cell type: 0 dofmap: 4481x10
[2025-11-03 11:20:08.570] [info] Global index computation
[2025-11-03 11:20:08.570] [info] Got 2 index_maps
[2025-11-03 11:20:08.570] [info] Get global indices

Now we need to redefine the markers to have so that facets on the endo- and epicardium combine both free wall and the septum.

markers = {"ENDO_LV": [1, 2], "ENDO_RV": [2, 2], "BASE": [3, 2], "EPI": [4, 2]}
marker_values = geo.ffun.values.copy()
marker_values[np.isin(geo.ffun.indices, geo.ffun.find(geo.markers["LV_ENDO_FW"][0]))] = markers["ENDO_LV"][0]
marker_values[np.isin(geo.ffun.indices, geo.ffun.find(geo.markers["LV_SEPTUM"][0]))] = markers["ENDO_LV"][0]
marker_values[np.isin(geo.ffun.indices, geo.ffun.find(geo.markers["RV_ENDO_FW"][0]))] = markers["ENDO_RV"][0]
marker_values[np.isin(geo.ffun.indices, geo.ffun.find(geo.markers["RV_SEPTUM"][0]))] = markers["ENDO_RV"][0]
marker_values[np.isin(geo.ffun.indices, geo.ffun.find(geo.markers["BASE"][0]))] = markers["BASE"][0]
marker_values[np.isin(geo.ffun.indices, geo.ffun.find(geo.markers["LV_EPI_FW"][0]))] = markers["EPI"][0]
marker_values[np.isin(geo.ffun.indices, geo.ffun.find(geo.markers["RV_EPI_FW"][0]))] = markers["EPI"][0]
geo.markers = markers
ffun = dolfinx.mesh.meshtags(
    geo.mesh,
    geo.ffun.dim,
    geo.ffun.indices,
    marker_values,
)
geo.ffun = ffun

Scale the geometry to be in meters

geo.mesh.geometry.x[:] *= 1.4e-2

In order to use the geometry with pulse we need to convert it to a pulse.Geometry object. We can do this by using the from_cardiac_geometries method. We also specify that we want to use a quadrature degree of 4

geometry = pulse.Geometry.from_cardiac_geometries(geo, metadata={"quadrature_degree": 4})

Next we create the material object, and we will use the transversely isotropic version of the Neo Hookean model

material = pulse.NeoHookean(mu=dolfinx.fem.Constant(geometry.mesh, dolfinx.default_scalar_type(15.0)))
# and use an active stress approach
Setting mu to c_6 kPa
Ta = dolfinx.fem.Constant(geometry.mesh, dolfinx.default_scalar_type(0.0))
active_model = pulse.ActiveStress(geo.f0, activation=Ta)
Activation is not a Variable, defaulting to kPa

Now we will also implement two different versions, one where we use a compressible model and one where we use an incompressible model. To do this we will introduce a flag

incompressible = False

And in both cases we will use different compressible models, mechanics problems and different ways to get the displacement.

if incompressible:
    comp_model: pulse.Compressibility = pulse.Incompressible()
else:
    comp_model = pulse.Compressible()

Now we can assemble the CardiacModel

model = pulse.CardiacModel(
    material=material,
    active=active_model,
    compressibility=comp_model,
)

We will add a pressure on the LV endocarium

lvp = dolfinx.fem.Constant(geometry.mesh, dolfinx.default_scalar_type(0.0))
neumann_lv = pulse.NeumannBC(traction=lvp, marker=geometry.markers["ENDO_LV"][0])
Traction is not a Variable, defaulting to kPa

and on the RV endocardium

rvp = dolfinx.fem.Constant(geometry.mesh, dolfinx.default_scalar_type(0.0))
neumann_rv = pulse.NeumannBC(traction=lvp, marker=geometry.markers["ENDO_RV"][0])
Traction is not a Variable, defaulting to kPa

We will also add a Robin type spring on the epicardial surface to mimic the pericardium.

pericardium = dolfinx.fem.Constant(geometry.mesh, dolfinx.default_scalar_type(1.0))
robin_per = pulse.RobinBC(value=pericardium, marker=geometry.markers["EPI"][0])
Value is not a Variable, defaulting to Pa / m

We collect all the boundary conditions

bcs = pulse.BoundaryConditions(neumann=(neumann_lv, neumann_rv), robin=(robin_per,))

create the problem

problem = pulse.StaticProblem(model=model, geometry=geometry, bcs=bcs, parameters={"base_bc": pulse.BaseBC.fixed})
[2025-11-03 11:20:08.639] [info] Checking required entities per dimension
[2025-11-03 11:20:08.639] [info] Cell type: 0 dofmap: 4481x10
[2025-11-03 11:20:08.640] [info] Global index computation
[2025-11-03 11:20:08.640] [info] Got 2 index_maps
[2025-11-03 11:20:08.640] [info] Get global indices
[2025-11-03 11:20:09.500] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-11-03 11:20:09.500] [info] Requesting connectivity (3, 0) - (2, 0)
[2025-11-03 11:20:09.500] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-11-03 11:20:09.500] [info] Requesting connectivity (3, 0) - (2, 0)
[2025-11-03 11:20:13.001] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-11-03 11:20:13.001] [info] Requesting connectivity (3, 0) - (2, 0)
[2025-11-03 11:20:13.001] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-11-03 11:20:13.001] [info] Requesting connectivity (3, 0) - (2, 0)
[2025-11-03 11:20:13.001] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-11-03 11:20:13.147] [info] Column ghost size increased from 0 to 0

and solve

problem.solve()
1

Now let us inflate the two ventricles and save the displacement

vtx = dolfinx.io.VTXWriter(geometry.mesh.comm, outdir / "biv_displacement.bp", [problem.u], engine="BP4")
vtx.write(0.0)
i = 1
for plv in [0.1]: #, 0.5, 1.0, 2.0]:
    print(f"plv: {plv}")
    lvp.value = plv
    rvp.value = plv * 0.2
    problem.solve()
    vtx.write(float(i))
    i += 1
plv: 0.1

and then apply an active tension

for ta in [0.1] : #, 1.0, 5.0, 10.0]:
    print(f"ta: {ta}")
    Ta.value = ta
    problem.solve()
    vtx.write(float(i))
    i += 1
ta: 0.1
vtx.close()
try:
    import pyvista
except ImportError:
    print("Pyvista is not installed")
else:
    V = dolfinx.fem.functionspace(geometry.mesh, ("Lagrange", 1, (geometry.mesh.geometry.dim,)))
    uh = dolfinx.fem.Function(V)
    uh.interpolate(problem.u)
    # Create plotter and pyvista grid
    p = pyvista.Plotter()
    topology, cell_types, geometry = dolfinx.plot.vtk_mesh(V)
    grid = pyvista.UnstructuredGrid(topology, cell_types, geometry)

    # Attach vector values to grid and warp grid by vector
    grid["u"] = uh.x.array.reshape((geometry.shape[0], 3))
    actor_0 = p.add_mesh(grid, style="wireframe", color="k")
    warped = grid.warp_by_vector("u", factor=1.5)
    actor_1 = p.add_mesh(warped, show_edges=True)
    p.show_axes()
    if not pyvista.OFF_SCREEN:
        p.show()
    else:
        figure_as_array = p.screenshot("biv_ellipsoid_pressure.png")
[2025-11-03 11:20:19.108] [info] Checking required entities per dimension
[2025-11-03 11:20:19.108] [info] Cell type: 0 dofmap: 4481x4
[2025-11-03 11:20:19.108] [info] Global index computation
[2025-11-03 11:20:19.108] [info] Got 1 index_maps
[2025-11-03 11:20:19.108] [info] Get global indices
2025-11-03 11:20:19.700 (   0.824s) [    7F077014E140]vtkXOpenGLRenderWindow.:1458  WARN| bad X server connection. DISPLAY=:99.0