diff --git a/dune/elasticity/materials/mooneyrivlindensity.hh b/dune/elasticity/materials/mooneyrivlindensity.hh index 3cf4d5d591a863bcacc5c300688066d516c60cd1..d1dd3b15b82241b9c75a614efffdb7a7a9c482e0 100644 --- a/dune/elasticity/materials/mooneyrivlindensity.hh +++ b/dune/elasticity/materials/mooneyrivlindensity.hh @@ -25,7 +25,9 @@ public: mooneyrivlin_a = parameters.get<double>("mooneyrivlin_a"); mooneyrivlin_b = parameters.get<double>("mooneyrivlin_b"); mooneyrivlin_c = parameters.get<double>("mooneyrivlin_c"); - } else if (mooneyrivlin_energy == "log" or mooneyrivlin_energy == "square") { + } + else if (mooneyrivlin_energy == "log" or mooneyrivlin_energy == "square") + { mooneyrivlin_10 = parameters.get<double>("mooneyrivlin_10"); mooneyrivlin_01 = parameters.get<double>("mooneyrivlin_01"); mooneyrivlin_20 = parameters.get<double>("mooneyrivlin_20"); @@ -35,8 +37,10 @@ public: mooneyrivlin_03 = parameters.get<double>("mooneyrivlin_03"); mooneyrivlin_21 = parameters.get<double>("mooneyrivlin_21"); mooneyrivlin_12 = parameters.get<double>("mooneyrivlin_12"); - mooneyrivlin_k = parameters.get<double>("mooneyrivlin_k"); - } else { + mooneyrivlin_k = parameters.get<double>("mooneyrivlin_k"); + } + else + { DUNE_THROW(Exception, "Error: Selected mooneyrivlin implementation (" << mooneyrivlin_energy << ") not available!"); } } @@ -97,7 +101,10 @@ public: 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); - } else { // mooneyrivlin_energy is "log" or "square" + } + else + { + // mooneyrivlin_energy is "log" or "square" field_type a = pow(detF, 2.0/dim); field_type invariant1Minus3 = frobeniusNormFsquared/a - 3; field_type secondInvariantOfC = 0; @@ -106,19 +113,23 @@ public: secondInvariantOfC += C[i][i]*C[j][j] - C[i][j]*C[i][j]; field_type invariant2Minus3 = secondInvariantOfC/(a*a) - 3; strainEnergy = mooneyrivlin_10 * invariant1Minus3 + - mooneyrivlin_01 * invariant2Minus3 + - mooneyrivlin_20 * invariant1Minus3 * invariant1Minus3 + - mooneyrivlin_02 * invariant2Minus3 * invariant2Minus3 + - mooneyrivlin_11 * invariant1Minus3 * invariant2Minus3 + - mooneyrivlin_30 * invariant1Minus3 * invariant1Minus3 * invariant1Minus3 + - mooneyrivlin_21 * invariant1Minus3 * invariant1Minus3 * invariant2Minus3 + - mooneyrivlin_12 * invariant1Minus3 * invariant2Minus3 * invariant2Minus3 + - mooneyrivlin_03 * invariant2Minus3 * invariant2Minus3 * invariant2Minus3; - if (mooneyrivlin_energy == "log") { + mooneyrivlin_01 * invariant2Minus3 + + mooneyrivlin_20 * invariant1Minus3 * invariant1Minus3 + + mooneyrivlin_02 * invariant2Minus3 * invariant2Minus3 + + mooneyrivlin_11 * invariant1Minus3 * invariant2Minus3 + + mooneyrivlin_30 * invariant1Minus3 * invariant1Minus3 * invariant1Minus3 + + mooneyrivlin_21 * invariant1Minus3 * invariant1Minus3 * invariant2Minus3 + + mooneyrivlin_12 * invariant1Minus3 * invariant2Minus3 * invariant2Minus3 + + mooneyrivlin_03 * invariant2Minus3 * invariant2Minus3 * invariant2Minus3; + if (mooneyrivlin_energy == "log") + { using std::log; field_type logDetF = log(detF); return strainEnergy + 0.5 * mooneyrivlin_k* logDetF * logDetF; - } else { //mooneyrivlin_energy is "square" + } + else + { + //mooneyrivlin_energy is "square" field_type detFMinus1 = detF - 1; return strainEnergy + mooneyrivlin_k* detFMinus1 * detFMinus1; }