diff --git a/src/fesa/elements/hex8_high_low.cpp b/src/fesa/elements/hex8_high_low.cpp index 6d5842b..93c7acd 100644 --- a/src/fesa/elements/hex8_high_low.cpp +++ b/src/fesa/elements/hex8_high_low.cpp @@ -204,6 +204,18 @@ Result Hex8Element::EvaluateHighLow( for (std::size_t a = 0; a < 8; ++a) for (std::size_t i = 0; i < 3; ++i) integrated_fg[a][i] += fg[a][i] * HighLow{point.weight}; + // Reuse identical first products within this trial point only. The second + // multiplication and the original left association remain unchanged. + std::array deviatoric_fg{}, shear_fg{}; + Tensor shear_fft{}; + for (std::size_t a = 0; a < 8; ++a) + for (std::size_t i = 0; i < 3; ++i) { + deviatoric_fg[a][i] = deviatoric_coefficient * fg[a][i]; + shear_fg[a][i] = HighLow{mu_} * fg[a][i]; + } + for (std::size_t i = 0; i < 3; ++i) + for (std::size_t k = 0; k < 3; ++k) + shear_fft[i][k] = HighLow{mu_} * fft[i][k]; for (std::size_t a = 0; a < 8; ++a) for (std::size_t b = 0; b < 8; ++b) { const auto dot = Dot(gradients[a], gradients[b]); @@ -216,22 +228,27 @@ Result Hex8Element::EvaluateHighLow( for (std::size_t k = 0; k < 3; ++k) { const auto row = 3U * a + i; const auto column = 3U * b + k; - const auto mat = deviatoric_coefficient * fg[a][i] * fg[b][k] + - HighLow{mu_} * fg[b][i] * fg[a][k] + - HighLow{mu_} * fft[i][k] * dot; + const auto mat = deviatoric_fg[a][i] * fg[b][k] + + shear_fg[b][i] * fg[a][k] + + shear_fft[i][k] * dot; material[row][column] += mat * HighLow{point.weight}; if (i == k) geometric[row][column] += geo * HighLow{point.weight}; } } } - if (nonlinear || !paired_linear_tangent_ready_) + if (nonlinear || !paired_linear_tangent_ready_) { + std::array bulk_fg{}; + for (std::size_t a = 0; a < 8; ++a) + for (std::size_t i = 0; i < 3; ++i) + bulk_fg[a][i] = bulk_coefficient * integrated_fg[a][i]; for (std::size_t a = 0; a < 8; ++a) for (std::size_t b = 0; b < 8; ++b) for (std::size_t i = 0; i < 3; ++i) for (std::size_t k = 0; k < 3; ++k) material[3 * a + i][3 * b + k] += - bulk_coefficient * integrated_fg[a][i] * integrated_fg[b][k]; + bulk_fg[a][i] * integrated_fg[b][k]; + } if (!IsNormalized(energy)) return Result::Failure( Failure("Hex8 high+low energy is unsupported or nonfinite."));