From 9d3a5516a86f798d01bd541af0e4d57a83a69393 Mon Sep 17 00:00:00 2001 From: Henrik Finsberg Date: Mon, 19 Jan 2026 08:20:59 +0100 Subject: [PATCH 1/2] Make sure we use deviatoric strains when using compressible models --- .cspell_dict.txt | 2 ++ src/pulse/cardiac_model.py | 71 +++++++++++++++++++++++++++++++++++-- tests/test_cardiac_model.py | 14 +++++--- 3 files changed, 79 insertions(+), 8 deletions(-) diff --git a/.cspell_dict.txt b/.cspell_dict.txt index 78f9a9c..795a5ff 100644 --- a/.cspell_dict.txt +++ b/.cspell_dict.txt @@ -19,6 +19,7 @@ Bucelli Bursi caplog cardiomyocytes +Cdev cellname Chiara CMAME @@ -41,6 +42,7 @@ Elife Elsevier endo Erlicher +Fdev Fedele fenics fenicsx diff --git a/src/pulse/cardiac_model.py b/src/pulse/cardiac_model.py index 212cf91..fa9a68a 100644 --- a/src/pulse/cardiac_model.py +++ b/src/pulse/cardiac_model.py @@ -66,13 +66,68 @@ def __post_init__(self): logger.debug(f" Compressibility: {type(self.compressibility).__name__}") logger.debug(f" Viscoelasticity: {type(self.viscoelasticity).__name__}") + def Cdev(self, C: ufl.core.expr.Expr) -> ufl.core.expr.Expr: + """Deviatoric part of the right Cauchy-Green deformation tensor. + + Parameters + ---------- + C : ufl.core.expr.Expr + Right Cauchy-Green deformation tensor + + Returns + ------- + ufl.core.expr.Expr + Deviatoric part of the right Cauchy-Green deformation tensor + """ + J = ufl.sqrt(ufl.det(C)) + dim = C.ufl_shape[0] + Cdev = J ** (-2.0 / dim) * C + return Cdev + + def Fdev(self, F: ufl.core.expr.Expr) -> ufl.core.expr.Expr: + """Deviatoric part of the deformation gradient. + + Parameters + ---------- + F : ufl.core.expr.Expr + Deformation gradient + + Returns + ------- + ufl.core.expr.Expr + Deviatoric part of the deformation gradient + """ + J = ufl.det(F) + dim = F.ufl_shape[0] + Fdev = J ** (-1.0 / dim) * F + return Fdev + def strain_energy( self, C: ufl.core.expr.Expr, C_dot: ufl.core.expr.Expr | None = None, ) -> ufl.core.expr.Expr: + """Total strain energy for the cardiac model. + + Parameters + ---------- + C : ufl.core.expr.Expr + Right Cauchy-Green deformation tensor + C_dot : ufl.core.expr.Expr | None, optional + Time derivative of the right Cauchy-Green deformation tensor, by default None + + Returns + ------- + ufl.core.expr.Expr + The total strain energy density + """ + if self.compressibility.is_compressible(): + Cdev = self.Cdev(C) + else: + Cdev = C + psi = ( - self.material.strain_energy(C) + self.material.strain_energy(Cdev) + self.active.strain_energy(C) + self.compressibility.strain_energy(C) ) @@ -87,8 +142,12 @@ def S( C_dot: ufl.core.expr.Expr | None = None, ) -> ufl.core.expr.Expr: """Cauchy stress for the cardiac model.""" + if self.compressibility.is_compressible(): + Cdev = self.Cdev(C) + else: + Cdev = C - S = self.material.S(C) + self.active.S(C) + self.compressibility.S(C) + S = self.material.S(Cdev) + self.active.S(C) + self.compressibility.S(C) if C_dot is not None: S += self.viscoelasticity.S(C_dot) return S @@ -99,7 +158,13 @@ def P( F_dot: ufl.core.expr.Expr | None = None, ) -> ufl.core.expr.Expr: """First Piola-Kirchhoff stress for the cardiac model.""" - P = self.material.P(F) + self.active.P(F) + self.compressibility.P(F) + + if self.compressibility.is_compressible(): + Fdev = self.Fdev(F) + else: + Fdev = F + + P = self.material.P(Fdev) + self.active.P(F) + self.compressibility.P(F) if F_dot is not None: P += self.viscoelasticity.P(F_dot) return P diff --git a/tests/test_cardiac_model.py b/tests/test_cardiac_model.py index decf4e0..854e4e5 100644 --- a/tests/test_cardiac_model.py +++ b/tests/test_cardiac_model.py @@ -52,7 +52,7 @@ def test_CardiacModel_HolzapfelOgden(comp_model_cls, isotropy, mesh, u): # ) # F = I + 0.1 I, C = 1.21 I, I4f = 1.21 - # J = det(F) = 1.1 ** 3 + # J = det(F) = 1.1 ** 3, J^{-2/3} = 1.1 ** -2 if isinstance(comp_model, pulse.compressibility.Incompressible): # psi = 0.5 * a * (I1 - 3) + 0.5 * a_f * (I4f - 1)**2 + p (J - 1) @@ -60,10 +60,14 @@ def test_CardiacModel_HolzapfelOgden(comp_model_cls, isotropy, mesh, u): # 1000.0 * (1.1 ** 3 - 1) = 668.05 assert math.isclose(value, 668.05) else: - # psi = 0.5 * a * (I1 - 3) + 0.5 * a_f * (I4f - 1)**2 + kappa * (J * ln(J) - J + 1) - # psi = 0.5 * 1000*1 * (3 * 1.21 - 3) + 0.5 * 1000.0 * 1.0 * (1.21 - 1)**2 + - # 1e6 * (1.1 ** 3 * math.log(1.1 ** 3) - 1.1 ** 3 + 1) = 49910.5979586692 - assert math.isclose(value, 49910.5979586692) + # J^{-2/3} = 1.1 ** -2 + # psi = 0.5 * a * (J^{-2/3} * I1 - 3) + 0.5 * a_f * (^{-2/3} * I4f - 1)**2 + # + kappa * (J * ln(J) - J + 1) + + # psi = 0.5 * 1000*1 * (1.1 ** -2 * 3 * 1.21 - 3) + # + 0.5 * 1000.0 * 1.0 * (1.1 ** -2 * 1.21 - 1)**2 + # +1e6 * (1.1 ** 3 * math.log(1.1 ** 3) - 1.1 ** 3 + 1) = 49573.547958669194 + assert math.isclose(value, 49573.547958669194) @pytest.mark.parametrize("isotropy", (pulse.active_stress.ActiveStressModels.transversely,)) From 93bbadcd4785dd348859ec1748dd076b64b74eea Mon Sep 17 00:00:00 2001 From: Henrik Finsberg Date: Mon, 19 Jan 2026 10:54:17 +0100 Subject: [PATCH 2/2] Remove deviatoric arguments from material models --- demo/geometries/biv_ellipsoid.py | 1 - src/pulse/cardiac_model.py | 58 +++------------------ src/pulse/kinematics.py | 38 ++++++++++++++ src/pulse/material_model.py | 23 ++++++-- src/pulse/material_models/guccione.py | 4 -- src/pulse/material_models/holzapfelogden.py | 20 ++++--- src/pulse/material_models/neo_hookean.py | 5 -- src/pulse/material_models/usyk.py | 4 -- tests/test_cardiac_model.py | 10 ++-- tests/test_material_models.py | 14 ++--- tests/test_static_problem_Holzapfel.py | 5 +- 11 files changed, 89 insertions(+), 93 deletions(-) diff --git a/demo/geometries/biv_ellipsoid.py b/demo/geometries/biv_ellipsoid.py index 4579de2..74ad8ce 100644 --- a/demo/geometries/biv_ellipsoid.py +++ b/demo/geometries/biv_ellipsoid.py @@ -157,7 +157,6 @@ def update_marker(old_name, new_id): # \Psi_{NH} = \frac{\mu}{2} (I_1 - 3) # $$ # -# `deviatoric=True` (default) means we use the isochoric invariant $\bar{I}_1 = J^{-2/3} I_1$. material = pulse.NeoHookean(mu=dolfinx.fem.Constant(geometry.mesh, dolfinx.default_scalar_type(15.0))) diff --git a/src/pulse/cardiac_model.py b/src/pulse/cardiac_model.py index fa9a68a..d641f06 100644 --- a/src/pulse/cardiac_model.py +++ b/src/pulse/cardiac_model.py @@ -11,6 +11,7 @@ import dolfinx import ufl +from . import kinematics from .viscoelasticity import NoneViscoElasticity logger = logging.getLogger(__name__) @@ -39,9 +40,9 @@ def register(self, p: dolfinx.fem.Function | None) -> None: ... class HyperElasticMaterial(Protocol): def strain_energy(self, C: ufl.core.expr.Expr) -> ufl.core.expr.Expr: ... - def P(self, F: ufl.core.expr.Expr) -> ufl.core.expr.Expr: ... + def P(self, F: ufl.core.expr.Expr, dev: bool) -> ufl.core.expr.Expr: ... - def S(self, C: ufl.core.expr.Expr) -> ufl.core.expr.Expr: ... + def S(self, C: ufl.core.expr.Expr, dev: bool) -> ufl.core.expr.Expr: ... class ViscoElasticity(Protocol): @@ -66,42 +67,6 @@ def __post_init__(self): logger.debug(f" Compressibility: {type(self.compressibility).__name__}") logger.debug(f" Viscoelasticity: {type(self.viscoelasticity).__name__}") - def Cdev(self, C: ufl.core.expr.Expr) -> ufl.core.expr.Expr: - """Deviatoric part of the right Cauchy-Green deformation tensor. - - Parameters - ---------- - C : ufl.core.expr.Expr - Right Cauchy-Green deformation tensor - - Returns - ------- - ufl.core.expr.Expr - Deviatoric part of the right Cauchy-Green deformation tensor - """ - J = ufl.sqrt(ufl.det(C)) - dim = C.ufl_shape[0] - Cdev = J ** (-2.0 / dim) * C - return Cdev - - def Fdev(self, F: ufl.core.expr.Expr) -> ufl.core.expr.Expr: - """Deviatoric part of the deformation gradient. - - Parameters - ---------- - F : ufl.core.expr.Expr - Deformation gradient - - Returns - ------- - ufl.core.expr.Expr - Deviatoric part of the deformation gradient - """ - J = ufl.det(F) - dim = F.ufl_shape[0] - Fdev = J ** (-1.0 / dim) * F - return Fdev - def strain_energy( self, C: ufl.core.expr.Expr, @@ -122,7 +87,7 @@ def strain_energy( The total strain energy density """ if self.compressibility.is_compressible(): - Cdev = self.Cdev(C) + Cdev = kinematics.Cdev(C) else: Cdev = C @@ -142,12 +107,9 @@ def S( C_dot: ufl.core.expr.Expr | None = None, ) -> ufl.core.expr.Expr: """Cauchy stress for the cardiac model.""" - if self.compressibility.is_compressible(): - Cdev = self.Cdev(C) - else: - Cdev = C + dev = self.compressibility.is_compressible() - S = self.material.S(Cdev) + self.active.S(C) + self.compressibility.S(C) + S = self.material.S(C, dev=dev) + self.active.S(C) + self.compressibility.S(C) if C_dot is not None: S += self.viscoelasticity.S(C_dot) return S @@ -158,13 +120,9 @@ def P( F_dot: ufl.core.expr.Expr | None = None, ) -> ufl.core.expr.Expr: """First Piola-Kirchhoff stress for the cardiac model.""" + dev = self.compressibility.is_compressible() - if self.compressibility.is_compressible(): - Fdev = self.Fdev(F) - else: - Fdev = F - - P = self.material.P(Fdev) + self.active.P(F) + self.compressibility.P(F) + P = self.material.P(F, dev=dev) + self.active.P(F) + self.compressibility.P(F) if F_dot is not None: P += self.viscoelasticity.P(F_dot) return P diff --git a/src/pulse/kinematics.py b/src/pulse/kinematics.py index 13f323f..a939a78 100644 --- a/src/pulse/kinematics.py +++ b/src/pulse/kinematics.py @@ -2,6 +2,44 @@ import ufl +def Fdev(F: ufl.core.expr.Expr) -> ufl.core.expr.Expr: + """Deviatoric part of the deformation gradient. + + Parameters + ---------- + F : ufl.core.expr.Expr + Deformation gradient + + Returns + ------- + ufl.core.expr.Expr + Deviatoric part of the deformation gradient + """ + J = ufl.det(F) + dim = F.ufl_shape[0] + Fdev = J ** (-1.0 / dim) * F + return Fdev + + +def Cdev(C: ufl.core.expr.Expr) -> ufl.core.expr.Expr: + """Deviatoric part of the right Cauchy-Green deformation tensor. + + Parameters + ---------- + C : ufl.core.expr.Expr + Right Cauchy-Green deformation tensor + + Returns + ------- + ufl.core.expr.Expr + Deviatoric part of the right Cauchy-Green deformation tensor + """ + J = ufl.sqrt(ufl.det(C)) + dim = C.ufl_shape[0] + Cdev = J ** (-2.0 / dim) * C + return Cdev + + # Second order identity tensor def SecondOrderIdentity(F: ufl.core.expr.Expr) -> ufl.core.expr.Expr: """Return identity with same dimension as input""" diff --git a/src/pulse/material_model.py b/src/pulse/material_model.py index f7be436..63e968a 100644 --- a/src/pulse/material_model.py +++ b/src/pulse/material_model.py @@ -72,13 +72,16 @@ def strain_energy(self, F: ufl.core.expr.Expr) -> ufl.core.expr.Expr: The strain energy density """ - def P(self, F: ufl.core.expr.Expr) -> ufl.core.expr.Expr: + def P(self, F: ufl.core.expr.Expr, dev: bool = True) -> ufl.core.expr.Expr: r"""First Piola-Kirchhoff stress tensor Parameters ---------- F : ufl.core.expr.Expr The deformation gradient + dev : bool + Whether to compute the stress for the deviatoric part only + This should be True for compressible materials Returns ------- @@ -95,22 +98,34 @@ def P(self, F: ufl.core.expr.Expr) -> ufl.core.expr.Expr: \mathbf{P} = \frac{\partial \Psi}{\partial \mathbf{F}} """ C = F.T * F - return ufl.diff(self.strain_energy(C), F) + if dev: + Cdev = kinematics.Cdev(C) + else: + Cdev = C + return ufl.diff(self.strain_energy(Cdev), F) def sigma(self, F: ufl.core.expr.Expr) -> ufl.core.expr.Expr: return kinematics.InversePiolaTransform(self.P(F), F) - def S(self, C: ufl.core.expr.Expr) -> ufl.core.expr.Expr: + def S(self, C: ufl.core.expr.Expr, dev: bool = True) -> ufl.core.expr.Expr: """Cauchy stress tensor for the hyperelastic material model. Parameters ---------- C : ufl.core.expr.Expr The right Cauchy-Green deformation tensor + dev : bool + Whether to compute the stress for the deviatoric part only + This should be True for compressible materials Returns ------- ufl.core.expr.Expr The Cauchy stress tensor + """ - return 2.0 * ufl.diff(self.strain_energy(C), C) + if dev: + Cdev = kinematics.Cdev(C) + else: + Cdev = C + return 2.0 * ufl.diff(self.strain_energy(Cdev), C) diff --git a/src/pulse/material_models/guccione.py b/src/pulse/material_models/guccione.py index 889993c..e146931 100644 --- a/src/pulse/material_models/guccione.py +++ b/src/pulse/material_models/guccione.py @@ -71,7 +71,6 @@ class Guccione(HyperElasticMaterial): bf: Variable = field(default_factory=lambda: Variable(8.0, "dimensionless")) bt: Variable = field(default_factory=lambda: Variable(2.0, "dimensionless")) bfs: Variable = field(default_factory=lambda: Variable(4.0, "dimensionless")) - deviatoric: bool = True def __post_init__(self): # Check that all values are positive @@ -130,9 +129,6 @@ def is_isotropic(self) -> bool: def _Q(self, C: ufl.core.expr.Expr) -> ufl.core.expr.Expr: dim = C.ufl_shape[0] - if self.deviatoric: - Jm23 = pow(ufl.det(C), -1.0 / dim) - C *= Jm23 # Make C deviatoric E = 0.5 * (C - ufl.Identity(dim)) bt = self.bt.to_base_units() diff --git a/src/pulse/material_models/holzapfelogden.py b/src/pulse/material_models/holzapfelogden.py index 15396b6..75d8d74 100644 --- a/src/pulse/material_models/holzapfelogden.py +++ b/src/pulse/material_models/holzapfelogden.py @@ -118,7 +118,6 @@ class HolzapfelOgden(HyperElasticMaterial): b_fs: Variable = field(default_factory=lambda: Variable(0.0, "dimensionless")) use_subplus: bool = field(default=True, repr=False) use_heaviside: bool = field(default=True, repr=False) - deviatoric: bool = field(default=True, repr=False) _W1_func: Invariant = field( init=False, @@ -331,15 +330,20 @@ def _W8fs(self, I8fs: ufl.core.expr.Expr) -> ufl.core.expr.Expr: return self._W8fs_func(I8fs) def strain_energy(self, C: ufl.core.expr.Expr) -> ufl.core.expr.Expr: - dim = C.ufl_shape[0] + r"""Strain energy density function - I1 = invariants.I1(C) + Parameters + ---------- + C : ufl.core.expr.Expr + Right Cauchy-Green deformation tensor + + Returns + ------- + ufl.core.expr.Expr + The strain energy density + """ - if self.deviatoric: - Jm23 = pow(ufl.det(C), -1.0 / dim) - C *= Jm23 # Make C deviatoric - Jm23 = pow(invariants.I3(C), -1 / dim) - I1 *= Jm23 # Convert to deviatoric I1 - see https://arxiv.org/pdf/2009.08754 + I1 = invariants.I1(C) I4f = self._I4f(C) I4s = self._I4s(C) diff --git a/src/pulse/material_models/neo_hookean.py b/src/pulse/material_models/neo_hookean.py index d536be0..5ed5880 100644 --- a/src/pulse/material_models/neo_hookean.py +++ b/src/pulse/material_models/neo_hookean.py @@ -34,7 +34,6 @@ class NeoHookean(HyperElasticMaterial): """ mu: Variable = field(default_factory=lambda: Variable(15.0, "kPa")) - deviatoric: bool = field(default=True) def __post_init__(self): if not isinstance(self.mu, Variable): @@ -61,10 +60,6 @@ def parameters(self) -> dict[str, Variable]: def strain_energy(self, C: ufl.core.expr.Expr) -> ufl.core.expr.Expr: I1 = invariants.I1(C) dim = C.ufl_shape[0] - if self.deviatoric: - # Deviatoric strain energy - Jm23 = pow(ufl.det(C), -1.0 / dim) - I1 *= Jm23 mu = self.mu.to_base_units() return 0.5 * mu * (I1 - dim) diff --git a/src/pulse/material_models/usyk.py b/src/pulse/material_models/usyk.py index c6422da..efff7b0 100644 --- a/src/pulse/material_models/usyk.py +++ b/src/pulse/material_models/usyk.py @@ -66,7 +66,6 @@ class Usyk(HyperElasticMaterial): bfs: Variable = field(default_factory=lambda: Variable(12.0, "dimensionless")) bfn: Variable = field(default_factory=lambda: Variable(3.0, "dimensionless")) bsn: Variable = field(default_factory=lambda: Variable(3.0, "dimensionless")) - deviatoric: bool = True def __post_init__(self): # Check that all values are positive @@ -141,9 +140,6 @@ def is_isotropic(self) -> bool: def _Q(self, C: ufl.core.expr.Expr) -> ufl.core.expr.Expr: dim = C.ufl_shape[0] - if self.deviatoric: - Jm23 = pow(ufl.det(C), -1.0 / dim) - C *= Jm23 # Make C deviatoric E = 0.5 * (C - ufl.Identity(dim)) bf = self.bf.to_base_units() diff --git a/tests/test_cardiac_model.py b/tests/test_cardiac_model.py index 854e4e5..d4ab96d 100644 --- a/tests/test_cardiac_model.py +++ b/tests/test_cardiac_model.py @@ -27,7 +27,6 @@ def test_CardiacModel_HolzapfelOgden(comp_model_cls, isotropy, mesh, u): b_f=0.0, a_fs=0.0, b_fs=0.0, - deviatoric=False, ) comp_model = comp_model_cls() active_model = pulse.ActiveStress(f0, isotropy=isotropy) @@ -66,7 +65,7 @@ def test_CardiacModel_HolzapfelOgden(comp_model_cls, isotropy, mesh, u): # psi = 0.5 * 1000*1 * (1.1 ** -2 * 3 * 1.21 - 3) # + 0.5 * 1000.0 * 1.0 * (1.1 ** -2 * 1.21 - 1)**2 - # +1e6 * (1.1 ** 3 * math.log(1.1 ** 3) - 1.1 ** 3 + 1) = 49573.547958669194 + # + 1e6 * (1.1 ** 3 * math.log(1.1 ** 3) - 1.1 ** 3 + 1) = 49573.547958669194 assert math.isclose(value, 49573.547958669194) @@ -78,7 +77,6 @@ def test_CardiacModel_HolzapfelOgden(comp_model_cls, isotropy, mesh, u): def test_CardiacModel_NeoHookean(comp_model_cls, isotropy, mesh, u): material = pulse.NeoHookean( mu=dolfinx.fem.Constant(mesh, dolfinx.default_scalar_type(15.0)), - deviatoric=False, ) f0 = dolfinx.fem.Constant(mesh, (1.0, 0.0, 0.0)) comp_model = comp_model_cls() @@ -99,7 +97,7 @@ def test_CardiacModel_NeoHookean(comp_model_cls, isotropy, mesh, u): if isinstance(comp_model, pulse.compressibility.Incompressible): assert math.isclose(value, 4725.331000000082) else: - assert math.isclose(value, 54298.54795867078) + assert math.isclose(value, 49573.54795867355) @pytest.mark.parametrize("isotropy", (pulse.active_stress.ActiveStressModels.transversely,)) @@ -112,7 +110,7 @@ def test_CardiacModel_Guccione(comp_model_cls, isotropy, mesh, u): s0 = dolfinx.fem.Constant(mesh, (0.0, 1.0, 0.0)) n0 = dolfinx.fem.Constant(mesh, (0.0, 0.0, 1.0)) material_params = pulse.Guccione.default_parameters() - material = pulse.Guccione(f0=f0, s0=s0, n0=n0, **material_params, deviatoric=False) + material = pulse.Guccione(f0=f0, s0=s0, n0=n0, **material_params) active_model = pulse.ActiveStress(f0, isotropy=isotropy) comp_model = comp_model_cls() model = pulse.CardiacModel( @@ -130,4 +128,4 @@ def test_CardiacModel_Guccione(comp_model_cls, isotropy, mesh, u): if isinstance(comp_model, pulse.compressibility.Incompressible): assert math.isclose(value, 141.78170311802802) else: - assert math.isclose(value, 49714.998661786354) + assert math.isclose(value, 49573.54795866912) diff --git a/tests/test_material_models.py b/tests/test_material_models.py index 14dd19d..62778d8 100644 --- a/tests/test_material_models.py +++ b/tests/test_material_models.py @@ -22,7 +22,7 @@ def test_holzapfel_ogden(params_func, expected_value, mesh, u) -> None: params = params_func() f0 = dolfinx.fem.Constant(mesh, (1.0, 0.0, 0.0)) s0 = dolfinx.fem.Constant(mesh, (0.0, 1.0, 0.0)) - model = pulse.HolzapfelOgden(f0=f0, s0=s0, **params, deviatoric=False) + model = pulse.HolzapfelOgden(f0=f0, s0=s0, **params) u.interpolate(lambda x: x / 10) F = pulse.kinematics.DeformationGradient(u) @@ -57,7 +57,7 @@ def test_holzapfel_ogden_raises_MissingModelAttribute(params, attr): def test_holzapfel_ogden_neohookean(u): - model = pulse.HolzapfelOgden(a=1.0, deviatoric=False) + model = pulse.HolzapfelOgden(a=1.0) u.interpolate(lambda x: x / 10) F = pulse.kinematics.DeformationGradient(u) C = F.T * F @@ -72,7 +72,7 @@ def test_holzapfel_ogden_neohookean(u): def test_holzapfel_ogden_pure_fiber(u, mesh): f0 = dolfinx.fem.Constant(mesh, (1.0, 0.0, 0.0)) - model = pulse.HolzapfelOgden(a_f=1.0, f0=f0, deviatoric=False) + model = pulse.HolzapfelOgden(a_f=1.0, f0=f0) u.interpolate(lambda x: x / 10) F = pulse.kinematics.DeformationGradient(u) C = F.T * F @@ -85,7 +85,7 @@ def test_holzapfel_ogden_pure_fiber(u, mesh): def test_holzapfel_ogden_pure_fiber_sheets(u, mesh): f0 = dolfinx.fem.Constant(mesh, (1.0, 0.0, 0.0)) - model = pulse.HolzapfelOgden(a_fs=1.0, f0=f0, s0=f0, deviatoric=False) + model = pulse.HolzapfelOgden(a_fs=1.0, f0=f0, s0=f0) u.interpolate(lambda x: x / 10) F = pulse.kinematics.DeformationGradient(u) C = F.T * F @@ -97,7 +97,7 @@ def test_holzapfel_ogden_pure_fiber_sheets(u, mesh): def test_neo_hookean(u, mesh): - model = pulse.NeoHookean(mu=1.0, deviatoric=False) + model = pulse.NeoHookean(mu=1.0) u.interpolate(lambda x: x / 10) F = pulse.kinematics.DeformationGradient(u) C = F.T * F @@ -128,7 +128,7 @@ def test_saint_venant_kirchhoff(u): def test_guccione_isotropic(u): C = 10.0 bf = bt = bfs = 1.0 - model = pulse.Guccione(C=C, bf=bf, bt=bt, bfs=bfs, deviatoric=False) + model = pulse.Guccione(C=C, bf=bf, bt=bt, bfs=bfs) assert model.is_isotropic() u.interpolate(lambda x: x / 10) @@ -151,7 +151,7 @@ def test_guccione_anisotropic(u, mesh): bf = 1.0 bt = 2.0 bfs = 3.0 - model = pulse.Guccione(C=C, bf=bf, bt=bt, bfs=bfs, f0=f0, s0=s0, n0=n0, deviatoric=False) + model = pulse.Guccione(C=C, bf=bf, bt=bt, bfs=bfs, f0=f0, s0=s0, n0=n0) assert not model.is_isotropic() u.interpolate(lambda x: x / 10) diff --git a/tests/test_static_problem_Holzapfel.py b/tests/test_static_problem_Holzapfel.py index a8e085d..6dd6584 100644 --- a/tests/test_static_problem_Holzapfel.py +++ b/tests/test_static_problem_Holzapfel.py @@ -112,10 +112,7 @@ def test_CompressibleProblem_and_boundary_conditions(mesh): ) f0 = dolfinx.fem.Constant(mesh, PETSc.ScalarType((1.0, 0.0, 0.0))) - material = pulse.NeoHookean( - mu=dolfinx.fem.Constant(mesh, PETSc.ScalarType(15.0)), - deviatoric=True, - ) + material = pulse.NeoHookean(mu=dolfinx.fem.Constant(mesh, PETSc.ScalarType(15.0))) Ta = dolfinx.fem.Constant(mesh, PETSc.ScalarType(0.0)) active_model = pulse.ActiveStress(f0, activation=Ta)