Klotz curve#
Inlate a geometry to a pressure using different material models, and compare with the Klotz curve.
import math
import matplotlib.pyplot as plt
import numpy as np
import pulse
try:
from dolfin_adjoint import Constant, DirichletBC
except ImportError:
from dolfin import Constant, DirichletBC
ED_pressure = 1.6 # kPa
def setup_material(material_model):
"""
Choose parameters based on
Hadjicharalambous, Myrianthi, et al. "Analysis of passive
cardiac constitutive laws for parameter estimation using 3D
tagged MRI." Biomechanics and modeling in mechanobiology 14.4
(2015): 807-828.
These parameters did not really match the Klotz curve here.
Perhaps they did some more tuning?
"""
if material_model == "guccione":
matparams = pulse.Guccione.default_parameters()
matparams["C"] = 0.18 # kPa
matparams["bf"] = 27.75
matparams["bt"] = 5.37
matparams["bfs"] = 2.445
material = pulse.Guccione(
parameters=matparams,
f0=geometry.f0,
s0=geometry.s0,
n0=geometry.n0,
)
elif material_model == "neo_hookean":
matparams = pulse.NeoHookean.default_parameters()
matparams["mu"] = 10.0 # kPa
material = pulse.NeoHookean(parameters=matparams)
elif material_model == "holzapfel_ogden":
matparams = pulse.HolzapfelOgden.default_parameters()
matparams["a"] = 4.0 # kPa
matparams["a_f"] = 10.0 # kPa
matparams["b"] = 5.0
matparams["b_f"] = 5.0
material = pulse.HolzapfelOgden(
parameters=matparams,
f0=geometry.f0,
s0=geometry.s0,
n0=geometry.n0,
)
return material
def klotz_curve(ED_pressure):
"""
EDPVR based on Klotz curve
Klotz, Stefan, et al. "Single-beat estimation of end-diastolic
pressure-volume relationship: a novel method with potential for
noninvasive application." American Journal of Physiology-Heart and
Circulatory Physiology 291.1 (2006): H403-H412.
"""
# Some point at the EDPVR line
Vm = 148.663
Pm = ED_pressure
# Some constants
An = 27.8
Bn = 2.76
# kpa to mmhg
Pm = Pm * 760 / 101.325
V0 = Vm * (0.6 - 0.006 * Pm)
V30 = V0 + (Vm - V0) / (Pm / An) ** (1.0 / Bn)
beta = math.log(Pm / 30.0) / math.log(Vm / V30)
alpha = 30.0 / V30**beta
# Unloaded volume (not used here)
# P_V0 = alpha * V0 ** beta
vs = [V0]
ps = [0.0]
for p in np.linspace(1.0, 12.0):
vi = (p / alpha) ** (1.0 / beta)
vs.append(vi)
ps.append(p * 101.325 / 760) # Convert from mmhg to kPa
return vs, ps
def fix_basal_plane(W):
V = W if W.sub(0).num_sub_spaces() == 0 else W.sub(0)
bc = DirichletBC(
V,
Constant((0.0, 0.0, 0.0)),
geometry.ffun,
geometry.markers["BASE"][0],
)
return bc
geometry = pulse.HeartGeometry.from_file(pulse.mesh_paths["simple_ellipsoid"])
geometry.mesh.coordinates()[:] *= 4.7
# geometry = pulse.geometries.prolate_ellipsoid_geometry(mesh_size_factor=3.0)
# geometry.mesh.coordinates()[:] *= 0.37
2024-09-05 11:51:29,937 [760] INFO pulse.geometry_utils:
Load mesh from h5
dirichlet_bc = [fix_basal_plane]
lvp = Constant(0.0)
lv_marker = geometry.markers["ENDO"][0]
lv_pressure = pulse.NeumannBC(traction=lvp, marker=lv_marker, name="lv")
neumann_bc = [lv_pressure]
bcs = pulse.BoundaryConditions(dirichlet=dirichlet_bc, neumann=neumann_bc)
fig, ax = plt.subplots()
for material_model in ["neo_hookean", "guccione", "holzapfel_ogden"]:
material = setup_material(material_model)
problem = pulse.MechanicsProblem(geometry, material, bcs)
pressures = [0.0]
volumes = [geometry.cavity_volume()]
for p in np.linspace(0, ED_pressure, 10)[1:]:
pulse.iterate.iterate(problem, lvp, p)
pressures.append(p)
volumes.append(geometry.cavity_volume(u=problem.state.split()[0]))
ax.plot(volumes, pressures, label=" ".join(material_model.split("_")))
# Reset pressure
lvp.assign(Constant(0.0))
2024-09-05 11:51:30,072 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:51:30,073 [760] INFO pulse.iterate: Current control: f_20 = 0.000
2024-09-05 11:51:30,073 [760] INFO pulse.iterate: Target: 0.178
2024-09-05 11:51:30,177 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:51:30,178 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:51:30,178 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:51:30,179 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:51:30,798 [760] DEBUG UFL_LEGACY: Blocks of each mode:
99 full
2024-09-05 11:51:31,012 [760] DEBUG UFL_LEGACY: Blocks of each mode:
18 full
2024-09-05 11:52:09,650 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:52:09,651 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:52:09,652 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:52:09,652 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:52:09,791 [760] DEBUG UFL_LEGACY: Blocks of each mode:
10 full
2024-09-05 11:52:09,983 [760] DEBUG UFL_LEGACY: Blocks of each mode:
3 full
2024-09-05 11:52:12,619 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:52:12,620 [760] INFO pulse.iterate: Current control: f_20 = 0.178
2024-09-05 11:52:12,621 [760] INFO pulse.iterate: Target: 0.356
2024-09-05 11:52:14,482 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:52:14,483 [760] INFO pulse.iterate: Current control: f_20 = 0.356
2024-09-05 11:52:14,483 [760] INFO pulse.iterate: Target: 0.533
2024-09-05 11:52:16,337 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:52:16,338 [760] INFO pulse.iterate: Current control: f_20 = 0.533
2024-09-05 11:52:16,338 [760] INFO pulse.iterate: Target: 0.711
2024-09-05 11:52:18,200 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:52:18,200 [760] INFO pulse.iterate: Current control: f_20 = 0.711
2024-09-05 11:52:18,201 [760] INFO pulse.iterate: Target: 0.889
2024-09-05 11:52:20,062 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:52:20,063 [760] INFO pulse.iterate: Current control: f_20 = 0.889
2024-09-05 11:52:20,064 [760] INFO pulse.iterate: Target: 1.067
2024-09-05 11:52:21,923 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:52:21,923 [760] INFO pulse.iterate: Current control: f_20 = 1.067
2024-09-05 11:52:21,924 [760] INFO pulse.iterate: Target: 1.244
2024-09-05 11:52:23,789 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:52:23,790 [760] INFO pulse.iterate: Current control: f_20 = 1.244
2024-09-05 11:52:23,790 [760] INFO pulse.iterate: Target: 1.422
2024-09-05 11:52:25,652 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:52:25,653 [760] INFO pulse.iterate: Current control: f_20 = 1.422
2024-09-05 11:52:25,654 [760] INFO pulse.iterate: Target: 1.600
2024-09-05 11:52:27,579 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:52:27,579 [760] INFO pulse.iterate: Current control: f_20 = 0.000
2024-09-05 11:52:27,580 [760] INFO pulse.iterate: Target: 0.178
2024-09-05 11:52:27,685 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:52:27,686 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:52:27,687 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:52:27,687 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:52:27,688 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:52:27,689 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:52:27,689 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:52:28,439 [760] DEBUG UFL_LEGACY: Blocks of each mode:
99 full
2024-09-05 11:52:28,648 [760] DEBUG UFL_LEGACY: Blocks of each mode:
18 full
2024-09-05 11:53:41,803 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:53:41,804 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:53:41,805 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:53:41,805 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:53:41,806 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:53:41,806 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:53:41,807 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:53:42,000 [760] DEBUG UFL_LEGACY: Blocks of each mode:
10 full
2024-09-05 11:53:42,192 [760] DEBUG UFL_LEGACY: Blocks of each mode:
3 full
2024-09-05 11:53:49,016 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:53:49,017 [760] INFO pulse.iterate: Current control: f_20 = 0.178
2024-09-05 11:53:49,017 [760] INFO pulse.iterate: Target: 0.356
2024-09-05 11:53:54,114 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:53:54,115 [760] INFO pulse.iterate: Current control: f_20 = 0.356
2024-09-05 11:53:54,116 [760] INFO pulse.iterate: Target: 0.533
2024-09-05 11:53:59,227 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:53:59,228 [760] INFO pulse.iterate: Current control: f_20 = 0.533
2024-09-05 11:53:59,228 [760] INFO pulse.iterate: Target: 0.711
2024-09-05 11:54:03,842 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:54:03,843 [760] INFO pulse.iterate: Current control: f_20 = 0.711
2024-09-05 11:54:03,843 [760] INFO pulse.iterate: Target: 0.889
2024-09-05 11:54:08,500 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:54:08,501 [760] INFO pulse.iterate: Current control: f_20 = 0.889
2024-09-05 11:54:08,501 [760] INFO pulse.iterate: Target: 1.067
2024-09-05 11:54:12,612 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:54:12,613 [760] INFO pulse.iterate: Current control: f_20 = 1.067
2024-09-05 11:54:12,613 [760] INFO pulse.iterate: Target: 1.244
2024-09-05 11:54:16,730 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:54:16,731 [760] INFO pulse.iterate: Current control: f_20 = 1.244
2024-09-05 11:54:16,732 [760] INFO pulse.iterate: Target: 1.422
2024-09-05 11:54:20,823 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:54:20,823 [760] INFO pulse.iterate: Current control: f_20 = 1.422
2024-09-05 11:54:20,824 [760] INFO pulse.iterate: Target: 1.600
2024-09-05 11:54:25,006 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:54:25,007 [760] INFO pulse.iterate: Current control: f_20 = 0.000
2024-09-05 11:54:25,007 [760] INFO pulse.iterate: Target: 0.178
2024-09-05 11:54:25,099 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:54:25,100 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:54:25,101 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:54:25,101 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:54:25,102 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:54:25,102 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:54:25,103 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:54:25,104 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:54:25,104 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:54:25,105 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:54:25,106 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:54:25,847 [760] DEBUG UFL_LEGACY: Blocks of each mode:
99 full
2024-09-05 11:54:26,056 [760] DEBUG UFL_LEGACY: Blocks of each mode:
18 full
2024-09-05 11:55:29,450 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:55:29,451 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:55:29,452 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:55:29,452 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:55:29,453 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:55:29,453 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:55:29,453 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:55:29,454 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:55:29,455 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:55:29,455 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:55:29,455 [760] INFO UFL_LEGACY: Adjusting missing element cell to tetrahedron.
2024-09-05 11:55:29,627 [760] DEBUG UFL_LEGACY: Blocks of each mode:
10 full
2024-09-05 11:55:29,819 [760] DEBUG UFL_LEGACY: Blocks of each mode:
3 full
2024-09-05 11:55:34,809 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:55:34,810 [760] INFO pulse.iterate: Current control: f_20 = 0.178
2024-09-05 11:55:34,810 [760] INFO pulse.iterate: Target: 0.356
2024-09-05 11:55:39,224 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:55:39,225 [760] INFO pulse.iterate: Current control: f_20 = 0.356
2024-09-05 11:55:39,225 [760] INFO pulse.iterate: Target: 0.533
2024-09-05 11:55:43,172 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:55:43,172 [760] INFO pulse.iterate: Current control: f_20 = 0.533
2024-09-05 11:55:43,173 [760] INFO pulse.iterate: Target: 0.711
2024-09-05 11:55:47,110 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:55:47,111 [760] INFO pulse.iterate: Current control: f_20 = 0.711
2024-09-05 11:55:47,111 [760] INFO pulse.iterate: Target: 0.889
2024-09-05 11:55:51,529 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:55:51,530 [760] INFO pulse.iterate: Current control: f_20 = 0.889
2024-09-05 11:55:51,531 [760] INFO pulse.iterate: Target: 1.067
2024-09-05 11:55:55,959 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:55:55,959 [760] INFO pulse.iterate: Current control: f_20 = 1.067
2024-09-05 11:55:55,960 [760] INFO pulse.iterate: Target: 1.244
2024-09-05 11:55:59,423 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:55:59,423 [760] INFO pulse.iterate: Current control: f_20 = 1.244
2024-09-05 11:55:59,424 [760] INFO pulse.iterate: Target: 1.422
2024-09-05 11:56:03,866 [760] INFO pulse.iterate: Iterating to target control (f_20)...
2024-09-05 11:56:03,867 [760] INFO pulse.iterate: Current control: f_20 = 1.422
2024-09-05 11:56:03,867 [760] INFO pulse.iterate: Target: 1.600
vs, ps = klotz_curve(ED_pressure)
ax.plot(vs, ps, linestyle="--", label="Klotz curve")
ax.legend(loc="best")
ax.set_xlabel("Volume (ml)")
ax.set_ylabel("Pressure (kPa)")
plt.show()
# plt.savefig("klotz_curve.png")