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 dolfinx
from dolfinx import log
import fenicsx_pulse
import ldrb
import cardiac_geometries
import cardiac_geometries.geometry
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)
system = ldrb.dolfinx_ldrb(mesh=geo.mesh, ffun=geo.ffun, markers=geo.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-03-05 19:46:01 [debug ] Convert file biv_ellipsoid/geometry/biv_ellipsoid.msh to dolfin
Info : Reading 'biv_ellipsoid/geometry/biv_ellipsoid.msh'...
Info : 43 entities
Info : 664 nodes
Info : 3293 elements
Info : Done reading 'biv_ellipsoid/geometry/biv_ellipsoid.msh'
[2025-03-05 19:46:01.429] [info] Extract basic topology: 8668->8668
[2025-03-05 19:46:01.429] [info] Build local dual graph
[2025-03-05 19:46:01.429] [info] Build local part of mesh dual graph (mixed)
[2025-03-05 19:46:01.430] [info] GPS pseudo-diameter:(46) 1652-1941
[2025-03-05 19:46:01.430] [info] Create topology (single cell type)
[2025-03-05 19:46:01.430] [info] Create topology (generalised)
[2025-03-05 19:46:01.430] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:01.430] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:01.431] [info] Compute ghost indices
[2025-03-05 19:46:01.431] [info] Computing communication graph edges (using PCX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.431] [info] Finished graph edge discovery using PCX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.431] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.431] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.431] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.431] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.431] [info] Number of neighbourhood source ranks in distribute_to_postoffice: 0
[2025-03-05 19:46:01.431] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.431] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.431] [info] Neighbourhood destination ranks from post office in distribute_data (rank, num dests, num dests/mpi_size): 0, 0, 0
[2025-03-05 19:46:01.431] [info] Checking required entities per dimension
[2025-03-05 19:46:01.431] [info] Cell type: 0 dofmap: 2167x4
[2025-03-05 19:46:01.431] [info] Global index computation
[2025-03-05 19:46:01.431] [info] Got 1 index_maps
[2025-03-05 19:46:01.432] [info] Get global indices
[2025-03-05 19:46:01.432] [info] XDMF distribute entity data
[2025-03-05 19:46:01.432] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:01.432] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:01.432] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:01.432] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:01.432] [info] XDMF build map
[2025-03-05 19:46:01.433] [info] Requesting connectivity (3, 0) - (0, 0)
[2025-03-05 19:46:01.433] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-03-05 19:46:01.433] [info] Build list of mesh entity indices from the entity vertices.
[2025-03-05 19:46:01.434] [info] XDMF distribute entity data
[2025-03-05 19:46:01.434] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:01.434] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:01.434] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:01.434] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:01.434] [info] XDMF build map
[2025-03-05 19:46:01.435] [info] Computing mesh entities of dimension 2
[2025-03-05 19:46:01.436] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.436] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.436] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-03-05 19:46:01.436] [info] Computing mesh connectivity 2-3 from transpose.
[2025-03-05 19:46:01.436] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-03-05 19:46:01.436] [info] Build list of mesh entity indices from the entity vertices.
[2025-03-05 19:46:01.439] [info] Opened HDF5 file with id "72057594037927936"
[2025-03-05 19:46:01.439] [info] Adding mesh to node "/Xdmf/Domain"
[2025-03-05 19:46:01.439] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-03-05 19:46:01.439] [info] Adding geometry data to node "/Xdmf/Domain/Grid"
[2025-03-05 19:46:01.439] [info] XDMF: add meshtags (Cell tags)
[2025-03-05 19:46:01.439] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-03-05 19:46:01.439] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-03-05 19:46:01.439] [info] XDMF: add meshtags (Facet tags)
[2025-03-05 19:46:01.439] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-03-05 19:46:01.440] [info] Computing mesh entities of dimension 1
[2025-03-05 19:46:01.441] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.441] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.441] [info] Requesting connectivity (1, 0) - (3, 0)
[2025-03-05 19:46:01.441] [info] Computing mesh connectivity 1-3 from transpose.
[2025-03-05 19:46:01.441] [info] XDMF: add meshtags (Edge tags)
[2025-03-05 19:46:01.441] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-03-05 19:46:01.441] [info] Requesting connectivity (0, 0) - (3, 0)
[2025-03-05 19:46:01.441] [info] Computing mesh connectivity 0-3 from transpose.
[2025-03-05 19:46:01.441] [info] XDMF: add meshtags (Vertex tags)
[2025-03-05 19:46:01.441] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-03-05 19:46:01.442] [info] Opened HDF5 file with id "72057594037927937"
[2025-03-05 19:46:01.442] [info] Read topology data "Mesh" at /Xdmf/Domain
[2025-03-05 19:46:01.442] [info] HDF5 Read data rate: 4550.429818229543 MB/s
[2025-03-05 19:46:01.442] [info] IO permuting cells
[2025-03-05 19:46:01.442] [info] Read geometry data "Mesh" at /Xdmf/Domain
[2025-03-05 19:46:01.442] [info] HDF5 Read data rate: 1761.6626133097502 MB/s
[2025-03-05 19:46:01.444] [info] Using partitioner with 8668 cell data
[2025-03-05 19:46:01.444] [info] Compute partition of cells across ranks
[2025-03-05 19:46:01.444] [info] Building mesh dual graph
[2025-03-05 19:46:01.444] [info] Build local part of mesh dual graph (mixed)
[2025-03-05 19:46:01.445] [info] Build nonlocal part of mesh dual graph
[2025-03-05 19:46:01.445] [info] Graph edges (local: 7542, non-local: 0)
[2025-03-05 19:46:01.445] [info] Compute graph partition using PT-SCOTCH
[2025-03-05 19:46:01.446] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:01.446] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:01.446] [info] Extract basic topology: 8668->8668
[2025-03-05 19:46:01.446] [info] Build local dual graph
[2025-03-05 19:46:01.446] [info] Build local part of mesh dual graph (mixed)
[2025-03-05 19:46:01.447] [info] GPS pseudo-diameter:(46) 2166-0
[2025-03-05 19:46:01.447] [info] Create topology (single cell type)
[2025-03-05 19:46:01.447] [info] Create topology (generalised)
[2025-03-05 19:46:01.447] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:01.447] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:01.448] [info] Compute ghost indices
[2025-03-05 19:46:01.448] [info] Computing communication graph edges (using PCX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.448] [info] Finished graph edge discovery using PCX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.448] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.448] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.448] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.448] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.448] [info] Number of neighbourhood source ranks in distribute_to_postoffice: 0
[2025-03-05 19:46:01.448] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.448] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.448] [info] Neighbourhood destination ranks from post office in distribute_data (rank, num dests, num dests/mpi_size): 0, 0, 0
[2025-03-05 19:46:01.448] [info] Checking required entities per dimension
[2025-03-05 19:46:01.448] [info] Cell type: 0 dofmap: 2167x4
[2025-03-05 19:46:01.448] [info] Global index computation
[2025-03-05 19:46:01.448] [info] Got 1 index_maps
[2025-03-05 19:46:01.448] [info] Get global indices
[2025-03-05 19:46:01.449] [info] Requesting connectivity (3, 0) - (3, 0)
[2025-03-05 19:46:01.449] [info] XDMF read meshtags (Cell tags)
[2025-03-05 19:46:01.449] [info] Read topology data "Cell tags" at /Xdmf/Domain
[2025-03-05 19:46:01.449] [info] HDF5 Read data rate: 4488.575312317949 MB/s
[2025-03-05 19:46:01.449] [info] IO permuting cells
[2025-03-05 19:46:01.449] [info] HDF5 Read data rate: 1162.865575529917 MB/s
[2025-03-05 19:46:01.449] [info] IO permuting cells
[2025-03-05 19:46:01.449] [info] XDMF distribute entity data
[2025-03-05 19:46:01.449] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:01.449] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:01.449] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:01.449] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:01.449] [info] XDMF build map
[2025-03-05 19:46:01.450] [info] XDMF create meshtags
[2025-03-05 19:46:01.450] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-03-05 19:46:01.450] [info] Build list of mesh entity indices from the entity vertices.
[2025-03-05 19:46:01.451] [info] Computing mesh entities of dimension 2
[2025-03-05 19:46:01.452] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.452] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.452] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-03-05 19:46:01.452] [info] Computing mesh connectivity 2-3 from transpose.
[2025-03-05 19:46:01.452] [info] XDMF read meshtags (Facet tags)
[2025-03-05 19:46:01.452] [info] Read topology data "Facet tags" at /Xdmf/Domain
[2025-03-05 19:46:01.452] [info] HDF5 Read data rate: 2028.0675422138836 MB/s
[2025-03-05 19:46:01.452] [info] IO permuting cells
[2025-03-05 19:46:01.452] [info] HDF5 Read data rate: 635.8887477057743 MB/s
[2025-03-05 19:46:01.452] [info] IO permuting cells
[2025-03-05 19:46:01.452] [info] XDMF distribute entity data
[2025-03-05 19:46:01.452] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:01.452] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:01.452] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:01.452] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:01.453] [info] XDMF build map
[2025-03-05 19:46:01.453] [info] XDMF create meshtags
[2025-03-05 19:46:01.453] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-03-05 19:46:01.453] [info] Build list of mesh entity indices from the entity vertices.
[2025-03-05 19:46:01.454] [info] Computing mesh entities of dimension 1
[2025-03-05 19:46:01.455] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.455] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.455] [info] Requesting connectivity (1, 0) - (3, 0)
[2025-03-05 19:46:01.455] [info] Computing mesh connectivity 1-3 from transpose.
[2025-03-05 19:46:01.456] [info] XDMF read meshtags (Edge tags)
[2025-03-05 19:46:01.456] [info] Read topology data "Edge tags" at /Xdmf/Domain
[2025-03-05 19:46:01.456] [info] HDF5 Read data rate: 0 MB/s
[2025-03-05 19:46:01.456] [info] IO permuting cells
[2025-03-05 19:46:01.456] [info] HDF5 Read data rate: 0 MB/s
[2025-03-05 19:46:01.456] [info] IO permuting cells
[2025-03-05 19:46:01.456] [info] XDMF distribute entity data
[2025-03-05 19:46:01.456] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.456] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.456] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:01.456] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:01.456] [info] XDMF build map
[2025-03-05 19:46:01.456] [info] XDMF create meshtags
[2025-03-05 19:46:01.456] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-03-05 19:46:01.456] [info] Build list of mesh entity indices from the entity vertices.
[2025-03-05 19:46:01.457] [info] Requesting connectivity (0, 0) - (3, 0)
[2025-03-05 19:46:01.457] [info] Computing mesh connectivity 0-3 from transpose.
[2025-03-05 19:46:01.457] [info] XDMF read meshtags (Vertex tags)
[2025-03-05 19:46:01.457] [info] Read topology data "Vertex tags" at /Xdmf/Domain
[2025-03-05 19:46:01.457] [info] HDF5 Read data rate: 0 MB/s
[2025-03-05 19:46:01.457] [info] IO permuting cells
[2025-03-05 19:46:01.457] [info] HDF5 Read data rate: 0 MB/s
[2025-03-05 19:46:01.457] [info] IO permuting cells
[2025-03-05 19:46:01.457] [info] XDMF distribute entity data
[2025-03-05 19:46:01.457] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:01.457] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:01.457] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:01.457] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:01.457] [info] XDMF build map
[2025-03-05 19:46:01.458] [info] XDMF create meshtags
[2025-03-05 19:46:01.458] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-03-05 19:46:01.458] [info] Build list of mesh entity indices from the entity vertices.
[2025-03-05 19:46:01.460] [info] Checking required entities per dimension
[2025-03-05 19:46:01.460] [info] Cell type: 0 dofmap: 2167x4
[2025-03-05 19:46:01.461] [info] Global index computation
[2025-03-05 19:46:01.461] [info] Got 1 index_maps
[2025-03-05 19:46:01.461] [info] Get global indices
[2025-03-05 19:46:01.579] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-03-05 19:46:01.582] [info] Column ghost size increased from 0 to 0
[2025-03-05 19:46:01.587] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-03-05 19:46:01.587] [info] Requesting connectivity (2, 0) - (0, 0)
[2025-03-05 19:46:01.587] [info] Requesting connectivity (2, 0) - (0, 0)
[2025-03-05 19:46:01.587] [info] Requesting connectivity (0, 0) - (3, 0)
[2025-03-05 19:46:01.587] [info] Requesting connectivity (3, 0) - (0, 0)
[2025-03-05 19:46:01.590] [info] Column ghost size increased from 0 to 0
0.0 1.0
[2025-03-05 19:46:01.595] [info] Opened HDF5 file with id "72057594037927938"
[2025-03-05 19:46:01.595] [info] Adding mesh to node "/Xdmf/Domain"
[2025-03-05 19:46:01.595] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-03-05 19:46:01.596] [info] Adding geometry data to node "/Xdmf/Domain/Grid"
[2025-03-05 19:46:01.596] [info] Adding function to node "/Xdmf/Domain/Grid/Grid"
[2025-03-05 19:46:01.599] [info] Column ghost size increased from 0 to 0
[2025-03-05 19:46:01.605] [info] Column ghost size increased from 0 to 0
[2025-03-05 19:46:01.610] [info] Column ghost size increased from 0 to 0
[2025-03-05 19:46:01.615] [info] Column ghost size increased from 0 to 0
[2025-03-05 19:46:01.621] [info] Checking required entities per dimension
[2025-03-05 19:46:01.621] [info] Cell type: 0 dofmap: 2167x10
[2025-03-05 19:46:01.621] [info] Global index computation
[2025-03-05 19:46:01.621] [info] Got 2 index_maps
[2025-03-05 19:46:01.621] [info] Get global indices
[2025-03-05 19:46:01.623] [info] Checking required entities per dimension
[2025-03-05 19:46:01.623] [info] Cell type: 0 dofmap: 2167x10
[2025-03-05 19:46:01.623] [info] Global index computation
[2025-03-05 19:46:01.623] [info] Got 2 index_maps
[2025-03-05 19:46:01.623] [info] Get global indices
[2025-03-05 19:46:02.305] [info] Checking required entities per dimension
[2025-03-05 19:46:02.305] [info] Cell type: 0 dofmap: 2167x10
[2025-03-05 19:46:02.305] [info] Global index computation
[2025-03-05 19:46:02.305] [info] Got 2 index_maps
[2025-03-05 19:46:02.305] [info] Get global indices
[2025-03-05 19:46:02.307] [info] Checking required entities per dimension
[2025-03-05 19:46:02.307] [info] Cell type: 0 dofmap: 2167x10
[2025-03-05 19:46:02.307] [info] Global index computation
[2025-03-05 19:46:02.307] [info] Got 2 index_maps
[2025-03-05 19:46:02.307] [info] Get global indices
[2025-03-05 19:46:04.058] [info] Checking required entities per dimension
[2025-03-05 19:46:04.058] [info] Cell type: 0 dofmap: 2167x10
[2025-03-05 19:46:04.059] [info] Global index computation
[2025-03-05 19:46:04.059] [info] Got 2 index_maps
[2025-03-05 19:46:04.059] [info] Get global indices
[2025-03-05 19:46:04.061] [info] Checking required entities per dimension
[2025-03-05 19:46:04.061] [info] Cell type: 0 dofmap: 2167x10
[2025-03-05 19:46:04.061] [info] Global index computation
[2025-03-05 19:46:04.061] [info] Got 2 index_maps
[2025-03-05 19:46:04.061] [info] Get global indices
[2025-03-05 19:46:04.064] [info] Compute face permutations
[2025-03-05 19:46:04.064] [info] Computing permutations for face type 0
[2025-03-05 19:46:04.064] [info] Compute edge permutations
If the folder exist we just load it
geo = cardiac_geometries.geometry.Geometry.from_folder(
comm=MPI.COMM_WORLD,
folder=geodir,
)
[2025-03-05 19:46:04.076] [info] Opened HDF5 file with id "72057594037927939"
[2025-03-05 19:46:04.076] [info] Read topology data "Mesh" at /Xdmf/Domain
[2025-03-05 19:46:04.076] [info] HDF5 Read data rate: 2487.0525787246256 MB/s
[2025-03-05 19:46:04.076] [info] IO permuting cells
[2025-03-05 19:46:04.076] [info] Read geometry data "Mesh" at /Xdmf/Domain
[2025-03-05 19:46:04.076] [info] HDF5 Read data rate: 1014.3857415658815 MB/s
[2025-03-05 19:46:04.078] [info] Using partitioner with 8668 cell data
[2025-03-05 19:46:04.079] [info] Compute partition of cells across ranks
[2025-03-05 19:46:04.079] [info] Building mesh dual graph
[2025-03-05 19:46:04.079] [info] Build local part of mesh dual graph (mixed)
[2025-03-05 19:46:04.080] [info] Build nonlocal part of mesh dual graph
[2025-03-05 19:46:04.080] [info] Graph edges (local: 7542, non-local: 0)
[2025-03-05 19:46:04.080] [info] Compute graph partition using PT-SCOTCH
[2025-03-05 19:46:04.080] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:04.080] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:04.080] [info] Extract basic topology: 8668->8668
[2025-03-05 19:46:04.080] [info] Build local dual graph
[2025-03-05 19:46:04.080] [info] Build local part of mesh dual graph (mixed)
[2025-03-05 19:46:04.081] [info] GPS pseudo-diameter:(46) 2166-0
[2025-03-05 19:46:04.081] [info] Create topology (single cell type)
[2025-03-05 19:46:04.081] [info] Create topology (generalised)
[2025-03-05 19:46:04.081] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:04.081] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:04.082] [info] Compute ghost indices
[2025-03-05 19:46:04.082] [info] Computing communication graph edges (using PCX algorithm). Number of input edges: 0
[2025-03-05 19:46:04.082] [info] Finished graph edge discovery using PCX algorithm. Number of discovered edges 0
[2025-03-05 19:46:04.082] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:04.082] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:04.082] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:04.082] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:04.082] [info] Number of neighbourhood source ranks in distribute_to_postoffice: 0
[2025-03-05 19:46:04.082] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:04.082] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:04.082] [info] Neighbourhood destination ranks from post office in distribute_data (rank, num dests, num dests/mpi_size): 0, 0, 0
[2025-03-05 19:46:04.082] [info] Checking required entities per dimension
[2025-03-05 19:46:04.082] [info] Cell type: 0 dofmap: 2167x4
[2025-03-05 19:46:04.082] [info] Global index computation
[2025-03-05 19:46:04.082] [info] Got 1 index_maps
[2025-03-05 19:46:04.082] [info] Get global indices
[2025-03-05 19:46:04.084] [info] Requesting connectivity (3, 0) - (3, 0)
[2025-03-05 19:46:04.084] [info] XDMF read meshtags (Cell tags)
[2025-03-05 19:46:04.084] [info] Read topology data "Cell tags" at /Xdmf/Domain
[2025-03-05 19:46:04.084] [info] HDF5 Read data rate: 2934.078023186934 MB/s
[2025-03-05 19:46:04.084] [info] IO permuting cells
[2025-03-05 19:46:04.084] [info] HDF5 Read data rate: 807.8285181733457 MB/s
[2025-03-05 19:46:04.084] [info] IO permuting cells
[2025-03-05 19:46:04.084] [info] XDMF distribute entity data
[2025-03-05 19:46:04.084] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:04.084] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:04.084] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:04.084] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:04.084] [info] XDMF build map
[2025-03-05 19:46:04.085] [info] XDMF create meshtags
[2025-03-05 19:46:04.085] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-03-05 19:46:04.085] [info] Build list of mesh entity indices from the entity vertices.
[2025-03-05 19:46:04.086] [info] Computing mesh entities of dimension 2
[2025-03-05 19:46:04.087] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:04.087] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:04.087] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-03-05 19:46:04.087] [info] Computing mesh connectivity 2-3 from transpose.
[2025-03-05 19:46:04.087] [info] XDMF read meshtags (Facet tags)
[2025-03-05 19:46:04.087] [info] Read topology data "Facet tags" at /Xdmf/Domain
[2025-03-05 19:46:04.087] [info] HDF5 Read data rate: 1884.91316174932 MB/s
[2025-03-05 19:46:04.087] [info] IO permuting cells
[2025-03-05 19:46:04.087] [info] HDF5 Read data rate: 629.5778585406765 MB/s
[2025-03-05 19:46:04.087] [info] IO permuting cells
[2025-03-05 19:46:04.087] [info] XDMF distribute entity data
[2025-03-05 19:46:04.087] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:04.087] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:04.087] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:04.087] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:04.087] [info] XDMF build map
[2025-03-05 19:46:04.088] [info] XDMF create meshtags
[2025-03-05 19:46:04.088] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-03-05 19:46:04.088] [info] Build list of mesh entity indices from the entity vertices.
[2025-03-05 19:46:04.090] [info] Computing mesh entities of dimension 1
[2025-03-05 19:46:04.091] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:04.091] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:04.091] [info] Requesting connectivity (1, 0) - (3, 0)
[2025-03-05 19:46:04.091] [info] Computing mesh connectivity 1-3 from transpose.
[2025-03-05 19:46:04.091] [info] XDMF read meshtags (Edge tags)
[2025-03-05 19:46:04.091] [info] Read topology data "Edge tags" at /Xdmf/Domain
[2025-03-05 19:46:04.091] [info] HDF5 Read data rate: 0 MB/s
[2025-03-05 19:46:04.091] [info] IO permuting cells
[2025-03-05 19:46:04.091] [info] HDF5 Read data rate: 0 MB/s
[2025-03-05 19:46:04.091] [info] IO permuting cells
[2025-03-05 19:46:04.091] [info] XDMF distribute entity data
[2025-03-05 19:46:04.091] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:04.091] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:04.091] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:04.091] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:04.091] [info] XDMF build map
[2025-03-05 19:46:04.092] [info] XDMF create meshtags
[2025-03-05 19:46:04.092] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-03-05 19:46:04.092] [info] Build list of mesh entity indices from the entity vertices.
[2025-03-05 19:46:04.092] [info] Requesting connectivity (0, 0) - (3, 0)
[2025-03-05 19:46:04.092] [info] Computing mesh connectivity 0-3 from transpose.
[2025-03-05 19:46:04.093] [info] XDMF read meshtags (Vertex tags)
[2025-03-05 19:46:04.093] [info] Read topology data "Vertex tags" at /Xdmf/Domain
[2025-03-05 19:46:04.093] [info] HDF5 Read data rate: 0 MB/s
[2025-03-05 19:46:04.093] [info] IO permuting cells
[2025-03-05 19:46:04.093] [info] HDF5 Read data rate: 0 MB/s
[2025-03-05 19:46:04.093] [info] IO permuting cells
[2025-03-05 19:46:04.093] [info] XDMF distribute entity data
[2025-03-05 19:46:04.093] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-03-05 19:46:04.093] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-03-05 19:46:04.093] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-03-05 19:46:04.093] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-03-05 19:46:04.093] [info] XDMF build map
[2025-03-05 19:46:04.093] [info] XDMF create meshtags
[2025-03-05 19:46:04.093] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-03-05 19:46:04.093] [info] Build list of mesh entity indices from the entity vertices.
[2025-03-05 19:46:04.095] [info] Checking required entities per dimension
[2025-03-05 19:46:04.095] [info] Cell type: 0 dofmap: 2167x10
[2025-03-05 19:46:04.096] [info] Global index computation
[2025-03-05 19:46:04.096] [info] Got 2 index_maps
[2025-03-05 19:46:04.096] [info] Get global indices
[2025-03-05 19:46:04.097] [info] Compute face permutations
[2025-03-05 19:46:04.097] [info] Computing permutations for face type 0
[2025-03-05 19:46:04.097] [info] Compute edge permutations
[2025-03-05 19:46:04.104] [info] Checking required entities per dimension
[2025-03-05 19:46:04.104] [info] Cell type: 0 dofmap: 2167x10
[2025-03-05 19:46:04.104] [info] Global index computation
[2025-03-05 19:46:04.104] [info] Got 2 index_maps
[2025-03-05 19:46:04.104] [info] Get global indices
[2025-03-05 19:46:04.111] [info] Checking required entities per dimension
[2025-03-05 19:46:04.111] [info] Cell type: 0 dofmap: 2167x10
[2025-03-05 19:46:04.111] [info] Global index computation
[2025-03-05 19:46:04.111] [info] Got 2 index_maps
[2025-03-05 19:46:04.111] [info] Get global indices
In order to use the geometry with pulse we need to convert it to a fenicsx_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 = fenicsx_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 = fenicsx_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 = fenicsx_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: fenicsx_pulse.Compressibility = fenicsx_pulse.Incompressible()
else:
comp_model = fenicsx_pulse.Compressible()
Now we can assemble the CardiacModel
model = fenicsx_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 = fenicsx_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 = fenicsx_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 = fenicsx_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 = fenicsx_pulse.BoundaryConditions(neumann=(neumann_lv, neumann_rv), robin=(robin_per,))
create the problem
problem = fenicsx_pulse.StaticProblem(model=model, geometry=geometry, bcs=bcs, parameters={"base_bc": fenicsx_pulse.BaseBC.fixed})
[2025-03-05 19:46:04.165] [info] Checking required entities per dimension
[2025-03-05 19:46:04.165] [info] Cell type: 0 dofmap: 2167x10
[2025-03-05 19:46:04.165] [info] Global index computation
[2025-03-05 19:46:04.165] [info] Got 2 index_maps
[2025-03-05 19:46:04.165] [info] Get global indices
[2025-03-05 19:46:04.978] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-03-05 19:46:04.978] [info] Requesting connectivity (3, 0) - (2, 0)
[2025-03-05 19:46:04.978] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-03-05 19:46:04.978] [info] Requesting connectivity (3, 0) - (2, 0)
[2025-03-05 19:46:06.841] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-03-05 19:46:06.841] [info] Requesting connectivity (3, 0) - (2, 0)
[2025-03-05 19:46:06.842] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-03-05 19:46:06.842] [info] Requesting connectivity (3, 0) - (2, 0)
[2025-03-05 19:46:06.842] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-03-05 19:46:06.905] [info] Column ghost size increased from 0 to 0
and solve
problem.solve()
5
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:
pyvista.start_xvfb()
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-03-05 19:46:14.414] [info] Checking required entities per dimension
[2025-03-05 19:46:14.414] [info] Cell type: 0 dofmap: 2167x4
[2025-03-05 19:46:14.414] [info] Global index computation
[2025-03-05 19:46:14.414] [info] Got 1 index_maps
[2025-03-05 19:46:14.414] [info] Get global indices
error: XDG_RUNTIME_DIR is invalid or not set in the environment.