Merge pull request #27 from adtzlr/fix-space-piola-broadcasting

Fix `Space.piola()` for `ThirdOrderDeformation(C10=0.5)`

.readthedocs.yml

@@ -17,6 +17,7 @@ formats:
 
 python:
   install:
-    - requirements: docs/requirements.txt
     - method: pip
-      path: .
+      path: .
+      extra_requirements:
+        - docs

README.md

@@ -145,7 +145,7 @@ ax.set_title(rf"Mooney-Rivlin (C10={parameters['C10']}, C01={parameters['C01']})
 ![fig_lab-mr](https://github.com/adtzlr/hyperelastic/assets/5793153/1d4bb29b-885f-46d4-80dd-56e255b239eb)
 
 # License
-Hyperelastic - Constitutive hyperelastic material formulations for FElupe (C) 2023 Andreas Dutzler, Graz (Austria).
+Hyperelastic - Constitutive hyperelastic material formulations for FElupe (C) 2024 Andreas Dutzler, Graz (Austria).
 
 This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
 

docs/conf.py

@@ -10,7 +10,7 @@
 # import hyperelastic
 
 project = "Hyperelastic"
-copyright = "2023, Andreas Dutzler"
+copyright = "2024, Andreas Dutzler"
 author = "Andreas Dutzler"
 
 # -- General configuration ---------------------------------------------------

pyproject.toml

@@ -39,14 +39,19 @@ classifiers = [
   "Topic :: Utilities"
 ]
 dynamic = ["version"]
-requires-python = ">=3.7"
+requires-python = ">=3.8"
 dependencies = [
   "numpy",
   "scipy",
   "matplotlib",
 ]
 
 [project.optional-dependencies]
+docs = [
+    "pydata-sphinx-theme",
+    "sphinx-copybutton",
+    "myst-nb",
+]
 test = [
     "felupe",
     "sympy",

src/hyperelastic/__about__.py

@@ -1 +1 @@
-__version__ = "0.10.0"
+__version__ = "0.10.1"

src/hyperelastic/spaces/_space.py

@@ -58,6 +58,13 @@ def piola(self, F, S, detF=None, C4=None, invC=None):
 
                 # force --> deformed configuration
                 if self.force is None:
+                    # check if destination array stress must be broadcasted
+                    new_shape = np.broadcast_shapes(F.shape, stress.shape)
+
+                    if not np.allclose(new_shape, stress.shape):
+                        # broadcast stress and make a writable contiguous copy
+                        stress = np.broadcast_to(stress, new_shape).copy()
+
                     stress = self.einsum("iI...,IJ...->iJ...", F, stress, out=stress)
 
                 # area --> deformed configuration
@@ -83,6 +90,13 @@ def piola(self, F, S, detF=None, C4=None, invC=None):
 
                 # force --> deformed configuration
                 if self.force is None:
+                    # check if destination array elasticity must be broadcasted
+                    new_shape = np.broadcast_shapes((3, 3, *F.shape), elasticity.shape)
+
+                    if not np.allclose(new_shape, elasticity.shape):
+                        # broadcast elasticity and make a writable contiguous copy
+                        elasticity = np.broadcast_to(elasticity, new_shape).copy()
+
                     elasticity = self.einsum(
                         "iI...,kK...,IJKL...->iJkL...", F, F, elasticity, out=elasticity
                     )

tests/test_space.py

@@ -66,10 +66,19 @@ def test_deformation_invariants():
     fea(umat).evaluate(verbose=2)
 
 
+def test_deformation_invariants_broadcast():
+    model = hel.models.invariants.ThirdOrderDeformation(
+        C10=0.5
+    )
+    umat = hel.DeformationSpace(hel.InvariantsFramework(model))
+    fea(umat).evaluate(verbose=2)
+
+
 if __name__ == "__main__":
     test_distortional_stretches()
     test_distortional_invariants()
     test_dilatational_invariants()
     test_deformation_invariants()
+    test_deformation_invariants_broadcast()
     test_distortional_generalized_strain_invariants()
     test_distortional_generalized_deformation_invariants()