From bd11829a0c4855cde088d048fa3eff21a2e2a789 Mon Sep 17 00:00:00 2001
From: Andreas Dutzler <a.dutzler@gmail.com>
Date: Sat, 6 Jan 2024 08:29:02 +0100
Subject: [PATCH] Update test_generalized.py

enhance the test
---
 tests/test_generalized.py | 43 +++++++++++++++++++++++----------------
 1 file changed, 26 insertions(+), 17 deletions(-)

diff --git a/tests/test_generalized.py b/tests/test_generalized.py
index 22b7fa8..313dd80 100644
--- a/tests/test_generalized.py
+++ b/tests/test_generalized.py
@@ -7,30 +7,39 @@
 def test_generalized():
     ϵ = np.sqrt(np.finfo(np.float64).eps)
     λ = np.append(1, np.linspace(1, 2, 21) + ϵ)
+    ux = ("uniaxial", 3, 1 / np.sqrt(λ), 1 / np.sqrt(λ))
+    bx = ("biaxial", 6, λ, 1 / λ**2)
     z = np.zeros_like(λ)
 
-    # uniaxial incompressible deformation
-    defgrad = np.array([[λ, z, z], [z, 1 / np.sqrt(λ), z], [z, z, 1 / np.sqrt(λ)]])
+    for title, initial_modulus, λ2, λ3 in [ux, bx]:
+        # create a new figure
+        plt.figure()
 
-    # loop over strain exponents
-    for k in [-0.5, 0.001, 2.0, 4.2]:
-        tod = hyperelastic.models.invariants.ThirdOrderDeformation(
-            strain=False, C10=0.3, C01=0.2, C20=-0.1, C30=0.05, C11=0.1
-        )
-        fun = hyperelastic.models.generalized.deformation
-        fwg = hyperelastic.GeneralizedInvariantsFramework(tod, fun=fun, exponent=k)
+        # uniaxial incompressible deformation
+        defgrad = np.array([[λ, z, z], [z, λ2, z], [z, z, λ3]])
 
-        umat = hyperelastic.DeformationSpace(fwg)
-        P = umat.gradient([defgrad, None])[0]
+        # loop over strain exponents
+        for k in [-0.5, 0.001, 2.0, 4.2]:
+            tod = hyperelastic.models.invariants.ThirdOrderDeformation(
+                strain=False, C10=0.3, C01=0.2, C20=-0.1, C30=0.05, C11=0.1
+            )
+            fun = hyperelastic.models.generalized.deformation
+            fwg = hyperelastic.GeneralizedInvariantsFramework(tod, fun=fun, exponent=k)
 
-        force = P[0, 0] - P[2, 2] * 1 / (λ * np.sqrt(λ))
-        dforce = np.diff(force)[0] / ϵ
+            umat = hyperelastic.DeformationSpace(fwg)
+            P = umat.gradient([defgrad, None])[0]
 
-        plt.plot(λ, force, label=f"k={k}")
-        plt.ylim(0, 5)
-        plt.legend()
+            force = P[0, 0] - P[2, 2] * λ3 / λ
+            dforce = np.diff(force)[0] / ϵ
 
-        assert np.isclose(dforce, 3)
+            plt.plot(λ, force, label=f"k={k}")
+            plt.xlabel("stretch $\lambda_1$ $\longrightarrow$")
+            plt.ylabel("force per undeformed area $N_1/A$ $\longrightarrow$")
+            plt.title(title)
+            plt.ylim(0, 5)
+            plt.legend()
+
+            assert np.isclose(dforce, initial_modulus, atol=1e-4)
 
 
 if __name__ == "__main__":