From 97320b481363badde6e9964cfb842fec91d05a10 Mon Sep 17 00:00:00 2001 From: Patrick Jaap <patrick.jaap@tu-dresden.de> Date: Wed, 9 Feb 2022 10:26:47 +0100 Subject: [PATCH] Fix in mooneyrivlindensity.hh Square the determinant in the Ciarlet variant --- dune/elasticity/materials/mooneyrivlindensity.hh | 4 ++-- test/mooneyrivlintest.cc | 8 ++++---- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/dune/elasticity/materials/mooneyrivlindensity.hh b/dune/elasticity/materials/mooneyrivlindensity.hh index 7317242..2b4c35d 100644 --- a/dune/elasticity/materials/mooneyrivlindensity.hh +++ b/dune/elasticity/materials/mooneyrivlindensity.hh @@ -66,7 +66,7 @@ public: C = F^TF and B = FF^T have the same eigenvalues - we can use either one of them. There are three Mooney-Rivlin-Variants: - ciarlet: W(F) = mooneyrivlin_a*(normF)^2 + mooneyrivlin_b*(normFinv)^2*detF + mooneyrivlin_c*(detF)^2 - + ciarlet: W(F) = mooneyrivlin_a*(normF)^2 + mooneyrivlin_b*(normFinv)^2*detF^2 + mooneyrivlin_c*(detF)^2 - ((dim-1)*mooneyrivlin_a + mooneyrivlin_b + 2*mooneyrivlin_c)*ln(detF) log: W(F) = \sum_{i,j=0}^{i+j<=n} mooneyrivlin_ij * (I1 - 3)^i * (I2 - 3)^j + mooneyrivlin_k * ln(det(F))^2 square: W(F) = \sum_{i,j=0}^{i+j<=n} mooneyrivlin_ij * (I1 - 3)^i * (I2 - 3)^j + mooneyrivlin_k * 0.5 * (det(F) - 1)^2 @@ -99,7 +99,7 @@ public: gradientInverse.invert(); field_type frobeinusNormFInverseSquared = gradientInverse.frobenius_norm2(); using std::log; - return mooneyrivlin_a*frobeniusNormFsquared + mooneyrivlin_b*frobeinusNormFInverseSquared*detF + mooneyrivlin_c*detF*detF - ((dim-1)*mooneyrivlin_a + mooneyrivlin_b + 2*mooneyrivlin_c)*log(detF); + return mooneyrivlin_a*frobeniusNormFsquared + mooneyrivlin_b*frobeinusNormFInverseSquared*detF*detF + mooneyrivlin_c*detF*detF - ((dim-1)*mooneyrivlin_a + mooneyrivlin_b + 2*mooneyrivlin_c)*log(detF); } else { diff --git a/test/mooneyrivlintest.cc b/test/mooneyrivlintest.cc index 1134aa6..8af6cba 100644 --- a/test/mooneyrivlintest.cc +++ b/test/mooneyrivlintest.cc @@ -167,10 +167,10 @@ int main (int argc, char *argv[]) parameters["mooneyrivlin_a"] = "2.42e+6"; parameters["mooneyrivlin_b"] = "6.52e+6"; parameters["mooneyrivlin_c"] = "-7.34e+6"; - expectedEnergy = 83069427.2; - expectedGradientTwoNorm = 163210957; - expectedGradientInfinityNorm = 49284369.2; - expectedMatrixFrobeniusNorm = 5.92126562e+09; + expectedEnergy = 68753239.6; + expectedGradientTwoNorm = 31244643.1; + expectedGradientInfinityNorm = 9673572.39; + expectedMatrixFrobeniusNorm = 1.67660965e+09; int testCiarlet = assembleAndCompare(basis, parameters, x, expectedEnergy, expectedGradientTwoNorm, expectedGradientInfinityNorm, expectedMatrixFrobeniusNorm); return testCiarlet + testLog + testSquare; -- GitLab