Problem 3: Inflation and Active Contraction of a Ventricle

Problem 3: Inflation and Active Contraction of a Ventricle#

This example implements Problem 3 from the cardiac mechanics benchmark suite [Land et al. 2015].

Problem Description#

Geometry: The same truncated ellipsoid as in Problem 2.

Microstructure: A spatially varying fiber field.

  • \(\alpha_{endo} = +90^\circ\)

  • \(\alpha_{epi} = -90^\circ\)

  • The fiber angle varies linearly across the wall thickness.

Material: Transversely Isotropic Guccione material.

  • Constitutive parameters: \(C = 2.0\) kPa, \(b_f = 8.0\), \(b_t = 2.0\), \(b_{fs} = 4.0\).

  • Incompressible.

Active Contraction: An active stress \(T_a\) is applied along the fiber direction.

Loading Protocol:

  1. Inflation: Increase cavity pressure \(P\) to 15 kPa.

  2. Contraction: While maintaining pressure, increase active tension \(T_a\) to 60 kPa. (Note: The benchmark describes these as varying curves, here we usually perform a coupled ramp or sequential steps).

from pathlib import Path
from mpi4py import MPI
from dolfinx import log
import dolfinx
import numpy as np
import math
import cardiac_geometries
import pulse

log.set_log_level(log.LogLevel.INFO)

1. Geometry and Fibers#

We regenerate the mesh with the specific fiber angles required for Problem 3 (\(+90^\circ\) to \(-90^\circ\)).

geodir = Path("lv_ellipsoid-problem3")
comm = MPI.COMM_WORLD

if not geodir.exists():
    comm.barrier()
    cardiac_geometries.mesh.lv_ellipsoid(
        outdir=geodir,
        r_short_endo=7.0,
        r_short_epi=10.0,
        r_long_endo=17.0,
        r_long_epi=20.0,
        mu_apex_endo=-math.pi,
        mu_base_endo=-math.acos(5 / 17),
        mu_apex_epi=-math.pi,
        mu_base_epi=-math.acos(5 / 20),
        create_fibers=True,
        fiber_angle_endo=90,
        fiber_angle_epi=-90,
        fiber_space="Quadrature_6",  # High order quadrature space for accurate fiber representation
        comm=comm,
    )

2025-12-17 15:58:18 [debug    ] Convert file lv_ellipsoid-problem3/lv_ellipsoid.msh to dolfin

Info    : Reading 'lv_ellipsoid-problem3/lv_ellipsoid.msh'...
Info    : 53 entities
Info    : 771 nodes
Info    : 4026 elements
Info    : Done reading 'lv_ellipsoid-problem3/lv_ellipsoid.msh'
[2025-12-17 15:58:18.204] [info] Extract basic topology: 11352->11352
[2025-12-17 15:58:18.204] [info] Build local dual graphs, re-order cells, and compute process boundary vertices.
[2025-12-17 15:58:18.204] [info] Build local part of mesh dual graph (mixed)
[2025-12-17 15:58:18.206] [info] GPS pseudo-diameter:(52) 22-2041
[2025-12-17 15:58:18.206] [info] Create topology (generalised)
[2025-12-17 15:58:18.207] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:18.207] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:18.207] [info] Compute ghost indices
[2025-12-17 15:58:18.207] [info] Computing communication graph edges (using PCX algorithm). Number of input edges: 0
[2025-12-17 15:58:18.207] [info] Finished graph edge discovery using PCX algorithm. Number of discovered edges 0
[2025-12-17 15:58:18.208] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:18.208] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:18.208] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:18.208] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:18.208] [info] Number of neighbourhood source ranks in distribute_to_postoffice: 0
[2025-12-17 15:58:18.208] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:18.208] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:18.208] [info] Neighbourhood destination ranks from post office in distribute_data (rank, num dests, num dests/mpi_size): 0, 0, 0
[2025-12-17 15:58:18.208] [info] Create Geometry (multiple)
[2025-12-17 15:58:18.208] [info] Got 1 dof layouts
[2025-12-17 15:58:18.208] [info] Checking required entities per dimension
[2025-12-17 15:58:18.208] [info] Cell type: 0 dofmap: 2838x4
[2025-12-17 15:58:18.208] [info] Global index computation
[2025-12-17 15:58:18.208] [info] Got 1 index_maps
[2025-12-17 15:58:18.208] [info] Get global indices
[2025-12-17 15:58:18.208] [info] Calling compute_local_to_global
[2025-12-17 15:58:18.208] [info] xdofs.size = 11352
[2025-12-17 15:58:18.208] [info] dofmap sizes = 11352 
[2025-12-17 15:58:18.208] [info] all_dofmaps.size = 11352
[2025-12-17 15:58:18.208] [info] nodes.size = 771
[2025-12-17 15:58:18.208] [info] Creating geometry with 1 dofmaps
[2025-12-17 15:58:18.208] [info] XDMF distribute entity data
[2025-12-17 15:58:18.209] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:18.209] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:18.209] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:18.209] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:18.209] [info] XDMF build map
[2025-12-17 15:58:18.210] [info] Requesting connectivity (3, 0) - (3, 0)
[2025-12-17 15:58:18.210] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-12-17 15:58:18.210] [info] Build list of mesh entity indices from the entity vertices.
[2025-12-17 15:58:18.211] [info] XDMF distribute entity data
[2025-12-17 15:58:18.211] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:18.211] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:18.211] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:18.211] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:18.211] [info] XDMF build map
[2025-12-17 15:58:18.212] [info] Computing mesh entities of dimension 2
[2025-12-17 15:58:18.213] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:18.213] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:18.214] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-12-17 15:58:18.214] [info] Computing mesh connectivity 2-3 from transpose.
[2025-12-17 15:58:18.214] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-12-17 15:58:18.214] [info] Build list of mesh entity indices from the entity vertices.
[2025-12-17 15:58:18.215] [info] XDMF distribute entity data
[2025-12-17 15:58:18.215] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:18.215] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:18.215] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:18.215] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:18.215] [info] XDMF build map
[2025-12-17 15:58:18.216] [info] Computing mesh entities of dimension 1
[2025-12-17 15:58:18.218] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:18.218] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:18.218] [info] Requesting connectivity (1, 0) - (3, 0)
[2025-12-17 15:58:18.218] [info] Computing mesh connectivity 1-3 from transpose.
[2025-12-17 15:58:18.218] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-12-17 15:58:18.218] [info] Build list of mesh entity indices from the entity vertices.
[2025-12-17 15:58:18.219] [info] XDMF distribute entity data
[2025-12-17 15:58:18.219] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:18.219] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:18.219] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:18.219] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:18.219] [info] XDMF build map
[2025-12-17 15:58:18.219] [info] Requesting connectivity (0, 0) - (3, 0)
[2025-12-17 15:58:18.219] [info] Computing mesh connectivity 0-3 from transpose.
[2025-12-17 15:58:18.219] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-12-17 15:58:18.219] [info] Build list of mesh entity indices from the entity vertices.
[2025-12-17 15:58:18.221] [info] Opened HDF5 file with id "72057594037927936"
[2025-12-17 15:58:18.221] [info] Adding mesh to node "/Xdmf/Domain"
[2025-12-17 15:58:18.221] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-12-17 15:58:18.221] [info] Adding geometry data to node "/Xdmf/Domain/Grid"
[2025-12-17 15:58:18.222] [info] XDMF: add meshtags (Cell tags)
[2025-12-17 15:58:18.222] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-12-17 15:58:18.222] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-12-17 15:58:18.222] [info] XDMF: add meshtags (Facet tags)
[2025-12-17 15:58:18.222] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-12-17 15:58:18.222] [info] Requesting connectivity (1, 0) - (3, 0)
[2025-12-17 15:58:18.222] [info] XDMF: add meshtags (Edge tags)
[2025-12-17 15:58:18.222] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-12-17 15:58:18.223] [info] Requesting connectivity (0, 0) - (3, 0)
[2025-12-17 15:58:18.223] [info] XDMF: add meshtags (Vertex tags)
[2025-12-17 15:58:18.223] [info] Adding topology data to node /Xdmf/Domain/Grid
[2025-12-17 15:58:18.225] [info] Checking required entities per dimension
[2025-12-17 15:58:18.225] [info] Cell type: 0 dofmap: 2838x4
[2025-12-17 15:58:18.225] [info] Global index computation
[2025-12-17 15:58:18.225] [info] Got 1 index_maps
[2025-12-17 15:58:18.225] [info] Get global indices

[2025-12-17 15:58:18.937] [info] Column ghost size increased from 0 to 0
[2025-12-17 15:58:18.944] [info] Checking required entities per dimension
[2025-12-17 15:58:18.944] [info] Cell type: 0 dofmap: 2838x24
[2025-12-17 15:58:18.944] [info] Global index computation
[2025-12-17 15:58:18.944] [info] Got 1 index_maps
[2025-12-17 15:58:18.945] [info] Get global indices

[2025-12-17 15:58:19.186] [info] Checking required entities per dimension
[2025-12-17 15:58:19.186] [info] Cell type: 0 dofmap: 2838x24
[2025-12-17 15:58:19.186] [info] Global index computation
[2025-12-17 15:58:19.186] [info] Got 1 index_maps
[2025-12-17 15:58:19.186] [info] Get global indices
[2025-12-17 15:58:19.199] [info] Compute face permutations
[2025-12-17 15:58:19.199] [info] Computing permutations for face type 0
[2025-12-17 15:58:19.200] [info] Compute edge permutations
[2025-12-17 15:58:19.219] [info] Opened HDF5 file with id "72057594037927938"
[2025-12-17 15:58:19.219] [info] Read topology data "Mesh" at /Xdmf/Domain
[2025-12-17 15:58:19.219] [info] HDF5 Read data rate: 1317.931154582922 MB/s
[2025-12-17 15:58:19.219] [info] IO permuting cells
[2025-12-17 15:58:19.219] [info] Read geometry data "Mesh" at /Xdmf/Domain
[2025-12-17 15:58:19.219] [info] HDF5 Read data rate: 1397.0554926387317 MB/s
[2025-12-17 15:58:19.222] [info] Using partitioner with cell data (1 cell types)
[2025-12-17 15:58:19.222] [info] Compute partition of cells across ranks
[2025-12-17 15:58:19.222] [info] Building mesh dual graph
[2025-12-17 15:58:19.222] [info] Build local part of mesh dual graph (mixed)
[2025-12-17 15:58:19.224] [info] Build nonlocal part of mesh dual graph
[2025-12-17 15:58:19.224] [info] Graph edges (local: 10252, non-local: 0)
[2025-12-17 15:58:19.224] [info] Compute graph partition using PT-SCOTCH
[2025-12-17 15:58:19.224] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.224] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.224] [info] Extract basic topology: 11352->11352
[2025-12-17 15:58:19.224] [info] Build local dual graphs, re-order cells, and compute process boundary vertices.
[2025-12-17 15:58:19.224] [info] Build local part of mesh dual graph (mixed)
[2025-12-17 15:58:19.226] [info] GPS pseudo-diameter:(52) 2837-0
[2025-12-17 15:58:19.226] [info] Create topology (generalised)
[2025-12-17 15:58:19.226] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.226] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.227] [info] Compute ghost indices
[2025-12-17 15:58:19.227] [info] Computing communication graph edges (using PCX algorithm). Number of input edges: 0
[2025-12-17 15:58:19.227] [info] Finished graph edge discovery using PCX algorithm. Number of discovered edges 0
[2025-12-17 15:58:19.227] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:19.227] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:19.228] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:19.228] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:19.228] [info] Number of neighbourhood source ranks in distribute_to_postoffice: 0
[2025-12-17 15:58:19.228] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:19.228] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:19.228] [info] Neighbourhood destination ranks from post office in distribute_data (rank, num dests, num dests/mpi_size): 0, 0, 0
[2025-12-17 15:58:19.228] [info] Create Geometry (multiple)
[2025-12-17 15:58:19.228] [info] Got 1 dof layouts
[2025-12-17 15:58:19.228] [info] Checking required entities per dimension
[2025-12-17 15:58:19.228] [info] Cell type: 0 dofmap: 2838x4
[2025-12-17 15:58:19.228] [info] Global index computation
[2025-12-17 15:58:19.228] [info] Got 1 index_maps
[2025-12-17 15:58:19.228] [info] Get global indices
[2025-12-17 15:58:19.228] [info] Calling compute_local_to_global
[2025-12-17 15:58:19.228] [info] xdofs.size = 11352
[2025-12-17 15:58:19.228] [info] dofmap sizes = 11352 
[2025-12-17 15:58:19.228] [info] all_dofmaps.size = 11352
[2025-12-17 15:58:19.228] [info] nodes.size = 771
[2025-12-17 15:58:19.228] [info] Creating geometry with 1 dofmaps
[2025-12-17 15:58:19.228] [info] Requesting connectivity (3, 0) - (3, 0)
[2025-12-17 15:58:19.228] [info] XDMF read meshtags (Cell tags)
[2025-12-17 15:58:19.228] [info] Read topology data "Cell tags" at /Xdmf/Domain
[2025-12-17 15:58:19.229] [info] HDF5 Read data rate: 3579.9432355723748 MB/s
[2025-12-17 15:58:19.229] [info] IO permuting cells
[2025-12-17 15:58:19.229] [info] HDF5 Read data rate: 1005.4025329908776 MB/s
[2025-12-17 15:58:19.229] [info] IO permuting cells
[2025-12-17 15:58:19.229] [info] XDMF distribute entity data
[2025-12-17 15:58:19.229] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.229] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.229] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.229] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.229] [info] XDMF build map
[2025-12-17 15:58:19.230] [info] XDMF create meshtags
[2025-12-17 15:58:19.230] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-12-17 15:58:19.230] [info] Build list of mesh entity indices from the entity vertices.
[2025-12-17 15:58:19.231] [info] Computing mesh entities of dimension 2
[2025-12-17 15:58:19.232] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:19.232] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:19.233] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-12-17 15:58:19.233] [info] Computing mesh connectivity 2-3 from transpose.
[2025-12-17 15:58:19.233] [info] XDMF read meshtags (Facet tags)
[2025-12-17 15:58:19.233] [info] Read topology data "Facet tags" at /Xdmf/Domain
[2025-12-17 15:58:19.233] [info] HDF5 Read data rate: 818.6046511627907 MB/s
[2025-12-17 15:58:19.233] [info] IO permuting cells
[2025-12-17 15:58:19.233] [info] HDF5 Read data rate: 607.4830871185973 MB/s
[2025-12-17 15:58:19.233] [info] IO permuting cells
[2025-12-17 15:58:19.233] [info] XDMF distribute entity data
[2025-12-17 15:58:19.233] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.233] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.233] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.233] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.233] [info] XDMF build map
[2025-12-17 15:58:19.234] [info] XDMF create meshtags
[2025-12-17 15:58:19.234] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-12-17 15:58:19.234] [info] Build list of mesh entity indices from the entity vertices.
[2025-12-17 15:58:19.235] [info] Computing mesh entities of dimension 1
[2025-12-17 15:58:19.237] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:19.237] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:19.237] [info] Requesting connectivity (1, 0) - (3, 0)
[2025-12-17 15:58:19.237] [info] Computing mesh connectivity 1-3 from transpose.
[2025-12-17 15:58:19.237] [info] XDMF read meshtags (Edge tags)
[2025-12-17 15:58:19.237] [info] Read topology data "Edge tags" at /Xdmf/Domain
[2025-12-17 15:58:19.237] [info] HDF5 Read data rate: 142.47256160695798 MB/s
[2025-12-17 15:58:19.237] [info] IO permuting cells
[2025-12-17 15:58:19.237] [info] HDF5 Read data rate: 58.9950265820614 MB/s
[2025-12-17 15:58:19.238] [info] IO permuting cells
[2025-12-17 15:58:19.238] [info] XDMF distribute entity data
[2025-12-17 15:58:19.238] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.238] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.238] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.238] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.238] [info] XDMF build map
[2025-12-17 15:58:19.238] [info] XDMF create meshtags
[2025-12-17 15:58:19.238] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-12-17 15:58:19.238] [info] Build list of mesh entity indices from the entity vertices.
[2025-12-17 15:58:19.239] [info] Requesting connectivity (0, 0) - (3, 0)
[2025-12-17 15:58:19.239] [info] Computing mesh connectivity 0-3 from transpose.
[2025-12-17 15:58:19.239] [info] XDMF read meshtags (Vertex tags)
[2025-12-17 15:58:19.239] [info] Read topology data "Vertex tags" at /Xdmf/Domain
[2025-12-17 15:58:19.239] [info] HDF5 Read data rate: 2.0499679692504804 MB/s
[2025-12-17 15:58:19.239] [info] IO permuting cells
[2025-12-17 15:58:19.239] [info] HDF5 Read data rate: 1.3766993632765445 MB/s
[2025-12-17 15:58:19.239] [info] IO permuting cells
[2025-12-17 15:58:19.239] [info] XDMF distribute entity data
[2025-12-17 15:58:19.239] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.239] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.239] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.239] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.239] [info] XDMF build map
[2025-12-17 15:58:19.240] [info] XDMF create meshtags
[2025-12-17 15:58:19.240] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-12-17 15:58:19.240] [info] Build list of mesh entity indices from the entity vertices.
[2025-12-17 15:58:19.242] [info] Checking required entities per dimension
[2025-12-17 15:58:19.242] [info] Cell type: 0 dofmap: 2838x24
[2025-12-17 15:58:19.242] [info] Global index computation
[2025-12-17 15:58:19.242] [info] Got 1 index_maps
[2025-12-17 15:58:19.242] [info] Get global indices
[2025-12-17 15:58:19.244] [info] Compute face permutations
[2025-12-17 15:58:19.244] [info] Computing permutations for face type 0
[2025-12-17 15:58:19.245] [info] Compute edge permutations
[2025-12-17 15:58:19.274] [info] Checking required entities per dimension
[2025-12-17 15:58:19.274] [info] Cell type: 0 dofmap: 2838x24
[2025-12-17 15:58:19.274] [info] Global index computation
[2025-12-17 15:58:19.274] [info] Got 1 index_maps
[2025-12-17 15:58:19.274] [info] Get global indices
[2025-12-17 15:58:19.305] [info] Checking required entities per dimension
[2025-12-17 15:58:19.305] [info] Cell type: 0 dofmap: 2838x24
[2025-12-17 15:58:19.305] [info] Global index computation
[2025-12-17 15:58:19.305] [info] Got 1 index_maps
[2025-12-17 15:58:19.305] [info] Get global indices

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

[2025-12-17 15:58:19.340] [info] Opened HDF5 file with id "72057594037927939"
[2025-12-17 15:58:19.340] [info] Read topology data "Mesh" at /Xdmf/Domain
[2025-12-17 15:58:19.340] [info] HDF5 Read data rate: 2363.645827911093 MB/s
[2025-12-17 15:58:19.340] [info] IO permuting cells
[2025-12-17 15:58:19.340] [info] Read geometry data "Mesh" at /Xdmf/Domain
[2025-12-17 15:58:19.340] [info] HDF5 Read data rate: 1345.1584763012504 MB/s
[2025-12-17 15:58:19.343] [info] Using partitioner with cell data (1 cell types)
[2025-12-17 15:58:19.343] [info] Compute partition of cells across ranks
[2025-12-17 15:58:19.343] [info] Building mesh dual graph
[2025-12-17 15:58:19.343] [info] Build local part of mesh dual graph (mixed)
[2025-12-17 15:58:19.344] [info] Build nonlocal part of mesh dual graph
[2025-12-17 15:58:19.344] [info] Graph edges (local: 10252, non-local: 0)
[2025-12-17 15:58:19.344] [info] Compute graph partition using PT-SCOTCH
[2025-12-17 15:58:19.344] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.344] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.345] [info] Extract basic topology: 11352->11352
[2025-12-17 15:58:19.345] [info] Build local dual graphs, re-order cells, and compute process boundary vertices.
[2025-12-17 15:58:19.345] [info] Build local part of mesh dual graph (mixed)
[2025-12-17 15:58:19.346] [info] GPS pseudo-diameter:(52) 2837-0
[2025-12-17 15:58:19.347] [info] Create topology (generalised)
[2025-12-17 15:58:19.347] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.347] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.347] [info] Compute ghost indices
[2025-12-17 15:58:19.347] [info] Computing communication graph edges (using PCX algorithm). Number of input edges: 0
[2025-12-17 15:58:19.347] [info] Finished graph edge discovery using PCX algorithm. Number of discovered edges 0
[2025-12-17 15:58:19.348] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:19.348] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:19.348] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:19.348] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:19.348] [info] Number of neighbourhood source ranks in distribute_to_postoffice: 0
[2025-12-17 15:58:19.348] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:19.348] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:19.348] [info] Neighbourhood destination ranks from post office in distribute_data (rank, num dests, num dests/mpi_size): 0, 0, 0
[2025-12-17 15:58:19.348] [info] Create Geometry (multiple)
[2025-12-17 15:58:19.348] [info] Got 1 dof layouts
[2025-12-17 15:58:19.348] [info] Checking required entities per dimension
[2025-12-17 15:58:19.348] [info] Cell type: 0 dofmap: 2838x4
[2025-12-17 15:58:19.348] [info] Global index computation
[2025-12-17 15:58:19.348] [info] Got 1 index_maps
[2025-12-17 15:58:19.348] [info] Get global indices
[2025-12-17 15:58:19.348] [info] Calling compute_local_to_global
[2025-12-17 15:58:19.348] [info] xdofs.size = 11352
[2025-12-17 15:58:19.348] [info] dofmap sizes = 11352 
[2025-12-17 15:58:19.348] [info] all_dofmaps.size = 11352
[2025-12-17 15:58:19.348] [info] nodes.size = 771
[2025-12-17 15:58:19.348] [info] Creating geometry with 1 dofmaps
[2025-12-17 15:58:19.349] [info] Requesting connectivity (3, 0) - (3, 0)
[2025-12-17 15:58:19.349] [info] XDMF read meshtags (Cell tags)
[2025-12-17 15:58:19.349] [info] Read topology data "Cell tags" at /Xdmf/Domain
[2025-12-17 15:58:19.349] [info] HDF5 Read data rate: 3690.8071202145816 MB/s
[2025-12-17 15:58:19.349] [info] IO permuting cells
[2025-12-17 15:58:19.349] [info] HDF5 Read data rate: 1033.879781420765 MB/s
[2025-12-17 15:58:19.349] [info] IO permuting cells
[2025-12-17 15:58:19.349] [info] XDMF distribute entity data
[2025-12-17 15:58:19.349] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.349] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.349] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.349] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.349] [info] XDMF build map
[2025-12-17 15:58:19.350] [info] XDMF create meshtags
[2025-12-17 15:58:19.350] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-12-17 15:58:19.350] [info] Build list of mesh entity indices from the entity vertices.
[2025-12-17 15:58:19.351] [info] Computing mesh entities of dimension 2
[2025-12-17 15:58:19.353] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:19.353] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:19.353] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-12-17 15:58:19.353] [info] Computing mesh connectivity 2-3 from transpose.
[2025-12-17 15:58:19.353] [info] XDMF read meshtags (Facet tags)
[2025-12-17 15:58:19.353] [info] Read topology data "Facet tags" at /Xdmf/Domain
[2025-12-17 15:58:19.353] [info] HDF5 Read data rate: 1702.2374105358178 MB/s
[2025-12-17 15:58:19.353] [info] IO permuting cells
[2025-12-17 15:58:19.353] [info] HDF5 Read data rate: 580.0922874093606 MB/s
[2025-12-17 15:58:19.353] [info] IO permuting cells
[2025-12-17 15:58:19.353] [info] XDMF distribute entity data
[2025-12-17 15:58:19.353] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.353] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.353] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.353] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.353] [info] XDMF build map
[2025-12-17 15:58:19.354] [info] XDMF create meshtags
[2025-12-17 15:58:19.354] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-12-17 15:58:19.354] [info] Build list of mesh entity indices from the entity vertices.
[2025-12-17 15:58:19.355] [info] Computing mesh entities of dimension 1
[2025-12-17 15:58:19.357] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 0
[2025-12-17 15:58:19.357] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 0
[2025-12-17 15:58:19.357] [info] Requesting connectivity (1, 0) - (3, 0)
[2025-12-17 15:58:19.357] [info] Computing mesh connectivity 1-3 from transpose.
[2025-12-17 15:58:19.357] [info] XDMF read meshtags (Edge tags)
[2025-12-17 15:58:19.357] [info] Read topology data "Edge tags" at /Xdmf/Domain
[2025-12-17 15:58:19.357] [info] HDF5 Read data rate: 149.94006756020485 MB/s
[2025-12-17 15:58:19.357] [info] IO permuting cells
[2025-12-17 15:58:19.357] [info] HDF5 Read data rate: 57.61179031987942 MB/s
[2025-12-17 15:58:19.357] [info] IO permuting cells
[2025-12-17 15:58:19.357] [info] XDMF distribute entity data
[2025-12-17 15:58:19.357] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.357] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.357] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.357] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.357] [info] XDMF build map
[2025-12-17 15:58:19.358] [info] XDMF create meshtags
[2025-12-17 15:58:19.358] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-12-17 15:58:19.358] [info] Build list of mesh entity indices from the entity vertices.
[2025-12-17 15:58:19.359] [info] Requesting connectivity (0, 0) - (3, 0)
[2025-12-17 15:58:19.359] [info] Computing mesh connectivity 0-3 from transpose.
[2025-12-17 15:58:19.359] [info] XDMF read meshtags (Vertex tags)
[2025-12-17 15:58:19.359] [info] Read topology data "Vertex tags" at /Xdmf/Domain
[2025-12-17 15:58:19.359] [info] HDF5 Read data rate: 1.9863438857852267 MB/s
[2025-12-17 15:58:19.359] [info] IO permuting cells
[2025-12-17 15:58:19.359] [info] HDF5 Read data rate: 1.4257708073427198 MB/s
[2025-12-17 15:58:19.359] [info] IO permuting cells
[2025-12-17 15:58:19.359] [info] XDMF distribute entity data
[2025-12-17 15:58:19.359] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.359] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.359] [info] Computing communication graph edges (using NBX algorithm). Number of input edges: 1
[2025-12-17 15:58:19.359] [info] Finished graph edge discovery using NBX algorithm. Number of discovered edges 1
[2025-12-17 15:58:19.359] [info] XDMF build map
[2025-12-17 15:58:19.359] [info] XDMF create meshtags
[2025-12-17 15:58:19.359] [info] Building MeshTags object from tagged entities (defined by vertices).
[2025-12-17 15:58:19.359] [info] Build list of mesh entity indices from the entity vertices.
[2025-12-17 15:58:19.362] [info] Checking required entities per dimension
[2025-12-17 15:58:19.362] [info] Cell type: 0 dofmap: 2838x24
[2025-12-17 15:58:19.362] [info] Global index computation
[2025-12-17 15:58:19.363] [info] Got 1 index_maps
[2025-12-17 15:58:19.363] [info] Get global indices
[2025-12-17 15:58:19.365] [info] Compute face permutations
[2025-12-17 15:58:19.365] [info] Computing permutations for face type 0
[2025-12-17 15:58:19.365] [info] Compute edge permutations

[2025-12-17 15:58:19.394] [info] Checking required entities per dimension
[2025-12-17 15:58:19.394] [info] Cell type: 0 dofmap: 2838x24
[2025-12-17 15:58:19.394] [info] Global index computation
[2025-12-17 15:58:19.394] [info] Got 1 index_maps
[2025-12-17 15:58:19.395] [info] Get global indices
[2025-12-17 15:58:19.425] [info] Checking required entities per dimension
[2025-12-17 15:58:19.425] [info] Cell type: 0 dofmap: 2838x24
[2025-12-17 15:58:19.425] [info] Global index computation
[2025-12-17 15:58:19.425] [info] Got 1 index_maps
[2025-12-17 15:58:19.425] [info] Get global indices

geometry = pulse.HeartGeometry.from_cardiac_geometries(geo, metadata={"quadrature_degree": 6})

[2025-12-17 15:58:19.459] [info] Requesting connectivity (2, 0) - (3, 0)

2. Constitutive Model#

Material: Transversely Isotropic Guccione. Active: Active Stress model \(S_{active} = T_a (\mathbf{f}_0 \otimes \mathbf{f}_0)\).

material_params = {
    "C": dolfinx.fem.Constant(geometry.mesh, dolfinx.default_scalar_type(2.0)),
    "bf": dolfinx.fem.Constant(geometry.mesh, dolfinx.default_scalar_type(8.0)),
    "bt": dolfinx.fem.Constant(geometry.mesh, dolfinx.default_scalar_type(2.0)),
    "bfs": dolfinx.fem.Constant(geometry.mesh, dolfinx.default_scalar_type(4.0)),
}

# Note: f0, s0, n0 come from the geometry object where they were loaded from files
material = pulse.Guccione(f0=geo.f0, s0=geo.s0, n0=geo.n0, **material_params)

Setting C to c_3 kPa

Ta = dolfinx.fem.Constant(geometry.mesh, dolfinx.default_scalar_type(0.0))
active_model = pulse.ActiveStress(geo.f0, activation=pulse.Variable(Ta, "kPa"))

comp_model = pulse.Incompressible()

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

3. Boundary Conditions#

  • Base: Fixed.

  • Pressure: Endocardial load.

pressure = dolfinx.fem.Constant(geometry.mesh, dolfinx.default_scalar_type(0.0))
neumann = pulse.NeumannBC(
    traction=pulse.Variable(pressure, "kPa"),
    marker=geo.markers["ENDO"][0],
)

bcs = pulse.BoundaryConditions(neumann=(neumann,))

4. Solution Process#

We simultaneously ramp up Pressure (\(0 \to 15\) kPa) and Active Tension (\(0 \to 60\) kPa). This represents a path in the P-Ta space.

problem = pulse.StaticProblem(
    model=model, geometry=geometry, bcs=bcs, parameters={"base_bc": pulse.BaseBC.fixed},
)

[2025-12-17 15:58:19.496] [info] Checking required entities per dimension
[2025-12-17 15:58:19.496] [info] Cell type: 0 dofmap: 2838x10
[2025-12-17 15:58:19.497] [info] Global index computation
[2025-12-17 15:58:19.497] [info] Got 2 index_maps
[2025-12-17 15:58:19.497] [info] Get global indices
[2025-12-17 15:58:19.498] [info] Checking required entities per dimension
[2025-12-17 15:58:19.498] [info] Cell type: 0 dofmap: 2838x4
[2025-12-17 15:58:19.498] [info] Global index computation
[2025-12-17 15:58:19.498] [info] Got 1 index_maps
[2025-12-17 15:58:19.498] [info] Get global indices

[2025-12-17 15:58:20.350] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-12-17 15:58:20.350] [info] Requesting connectivity (3, 0) - (2, 0)
[2025-12-17 15:58:20.350] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-12-17 15:58:20.350] [info] Requesting connectivity (3, 0) - (2, 0)

[2025-12-17 15:58:27.479] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-12-17 15:58:27.479] [info] Requesting connectivity (3, 0) - (2, 0)
[2025-12-17 15:58:27.479] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-12-17 15:58:27.479] [info] Requesting connectivity (3, 0) - (2, 0)

[2025-12-17 15:58:27.965] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-12-17 15:58:27.965] [info] Requesting connectivity (3, 0) - (2, 0)

[2025-12-17 15:58:28.494] [info] Requesting connectivity (2, 0) - (3, 0)
[2025-12-17 15:58:28.583] [info] Column ghost size increased from 0 to 0

Initial solve

problem.solve()

1

Ramping

target_pressure = 15.0
target_Ta = 60.0
steps = 20

pressures = np.linspace(0, target_pressure, steps)
tensions = np.linspace(0, target_Ta, steps)

vtx = dolfinx.io.VTXWriter(geometry.mesh.comm, "problem3.bp", [problem.u], engine="BP4")
vtx.write(0.0)

for i, (p, ta) in enumerate(zip(pressures, tensions)):
    if i == 0:
        continue
    print(f"Step {i}/{steps}: Pressure={p:.2f} kPa, Ta={ta:.2f} kPa")

    pressure.value = p
    Ta.value = ta

    problem.solve()
    vtx.write(float(i))
vtx.close()

Step 1/20: Pressure=0.79 kPa, Ta=3.16 kPa

Step 2/20: Pressure=1.58 kPa, Ta=6.32 kPa

Step 3/20: Pressure=2.37 kPa, Ta=9.47 kPa

Step 4/20: Pressure=3.16 kPa, Ta=12.63 kPa

Step 5/20: Pressure=3.95 kPa, Ta=15.79 kPa

Step 6/20: Pressure=4.74 kPa, Ta=18.95 kPa

Step 7/20: Pressure=5.53 kPa, Ta=22.11 kPa

Step 8/20: Pressure=6.32 kPa, Ta=25.26 kPa

Step 9/20: Pressure=7.11 kPa, Ta=28.42 kPa

Step 10/20: Pressure=7.89 kPa, Ta=31.58 kPa

Step 11/20: Pressure=8.68 kPa, Ta=34.74 kPa

Step 12/20: Pressure=9.47 kPa, Ta=37.89 kPa

Step 13/20: Pressure=10.26 kPa, Ta=41.05 kPa

Step 14/20: Pressure=11.05 kPa, Ta=44.21 kPa

Step 15/20: Pressure=11.84 kPa, Ta=47.37 kPa

Step 16/20: Pressure=12.63 kPa, Ta=50.53 kPa

Step 17/20: Pressure=13.42 kPa, Ta=53.68 kPa

Step 18/20: Pressure=14.21 kPa, Ta=56.84 kPa

Step 19/20: Pressure=15.00 kPa, Ta=60.00 kPa

5. Post-processing#

Compute apex position

U = dolfinx.fem.Function(problem.u.function_space)
U.interpolate(lambda x: (x[0], x[1], x[2]))

Use utility to evaluate at the specific vertex tag for the apex (ENDOPT/EPIPT generated by cardiac-geometries)

endo_apex_coord = pulse.utils.evaluate_at_vertex_tag(U, geo.vfun, geo.markers["ENDOPT"][0])
u_endo_apex = pulse.utils.evaluate_at_vertex_tag(problem.u, geo.vfun, geo.markers["ENDOPT"][0])
endo_apex_pos = pulse.utils.gather_broadcast_array(geo.mesh.comm, endo_apex_coord + u_endo_apex)
print(f"\nGet longitudinal position of endocardial apex: {endo_apex_pos[0, 0]:4f} mm")

Get longitudinal position of endocardial apex: -11.992595 mm

epi_apex_coord = pulse.utils.evaluate_at_vertex_tag(U, geo.vfun, geo.markers["EPIPT"][0])
u_epi_apex = pulse.utils.evaluate_at_vertex_tag(problem.u, geo.vfun, geo.markers["EPIPT"][0])
epi_apex_pos = pulse.utils.gather_broadcast_array(geo.mesh.comm, epi_apex_coord + u_epi_apex)
print(f"\nGet longitudinal position of epicardial apex: {epi_apex_pos[0, 0]:4f} mm")

Get longitudinal position of epicardial apex: -15.228496 mm

Visualization

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)

    p = pyvista.Plotter()
    topology, cell_types, geometry = dolfinx.plot.vtk_mesh(V)
    grid = pyvista.UnstructuredGrid(topology, cell_types, geometry)
    grid["u"] = uh.x.array.reshape((geometry.shape[0], 3))

    p.add_mesh(grid, style="wireframe", color="k", opacity=0.3, label="Reference")
    warped = grid.warp_by_vector("u", factor=1.0)
    p.add_mesh(warped, show_edges=True, color="blue", label="Contracted")

    p.show_axes()
    if not pyvista.OFF_SCREEN:
        p.show()
    else:
        p.screenshot("problem3.png")

[2025-12-17 15:59:47.687] [info] Checking required entities per dimension
[2025-12-17 15:59:47.687] [info] Cell type: 0 dofmap: 2838x4
[2025-12-17 15:59:47.687] [info] Global index computation
[2025-12-17 15:59:47.687] [info] Got 1 index_maps
[2025-12-17 15:59:47.687] [info] Get global indices

2025-12-17 15:59:48.299 (   0.846s) [    7F4635B2D140]vtkXOpenGLRenderWindow.:1458  WARN| bad X server connection. DISPLAY=:99.0
Contents