[MaL] Extract fourth order rotation to KelvinVec.

Previously only used in the linear elastic orthotropic
model, but the transformation itself is more general
and is being used in other (not yet available) code
parts.

MaterialLib/SolidModels/LinearElasticOrthotropic.cpp

@@ -68,7 +68,14 @@ LinearElasticOrthotropic<DisplacementDim>::getElasticTensor(
     // clang-format on
 
     KelvinMatrixType<3> const C_ortho = S_ortho.inverse();
-    KelvinMatrixType<3> const Q = fourthOrderRotationMatrix(x);
+    auto const Q = [this, &x]() -> KelvinMatrixType<3> {
+        if (!_local_coordinate_system)
+        {
+            return MathLib::KelvinVector::KelvinMatrixType<3>::Identity();
+        }
+        return MathLib::KelvinVector::fourthOrderRotationMatrix(
+            _local_coordinate_system->transformation<3>(x));
+    }();
 
     // Rotate the forth-order tenser in Kelvin mapping with Q*C_ortho*Q^T and
     // return the top left corner block of size 4x4 for two-dimensional case or

MaterialLib/SolidModels/LinearElasticOrthotropic.h

@@ -172,47 +172,6 @@ public:
 
     MaterialProperties getMaterialProperties() const { return _mp; }
 
-private:
-    /// Rotation matrix of a forth order tensor in 3D.
-    MathLib::KelvinVector::KelvinMatrixType<3> fourthOrderRotationMatrix(
-        ParameterLib::SpatialPosition const& x) const
-    {
-        if (!_local_coordinate_system)
-        {
-            return MathLib::KelvinVector::KelvinMatrixType<3>::Identity();
-        }
-
-        // 1-based index access for convenience.
-        auto Q = [R = _local_coordinate_system->transformation<3>(x)](
-                     int const i, int const j) { return R(i - 1, j - 1); };
-
-        MathLib::KelvinVector::KelvinMatrixType<3> R;
-        R << Q(1, 1) * Q(1, 1), Q(1, 2) * Q(1, 2), Q(1, 3) * Q(1, 3),
-            std::sqrt(2) * Q(1, 1) * Q(1, 2), std::sqrt(2) * Q(1, 2) * Q(1, 3),
-            std::sqrt(2) * Q(1, 1) * Q(1, 3), Q(2, 1) * Q(2, 1),
-            Q(2, 2) * Q(2, 2), Q(2, 3) * Q(2, 3),
-            std::sqrt(2) * Q(2, 1) * Q(2, 2), std::sqrt(2) * Q(2, 2) * Q(2, 3),
-            std::sqrt(2) * Q(2, 1) * Q(2, 3), Q(3, 1) * Q(3, 1),
-            Q(3, 2) * Q(3, 2), Q(3, 3) * Q(3, 3),
-            std::sqrt(2) * Q(3, 1) * Q(3, 2), std::sqrt(2) * Q(3, 2) * Q(3, 3),
-            std::sqrt(2) * Q(3, 1) * Q(3, 3), std::sqrt(2) * Q(1, 1) * Q(2, 1),
-            std::sqrt(2) * Q(1, 2) * Q(2, 2), std::sqrt(2) * Q(1, 3) * Q(2, 3),
-            Q(1, 1) * Q(2, 2) + Q(1, 2) * Q(2, 1),
-            Q(1, 2) * Q(2, 3) + Q(1, 3) * Q(2, 2),
-            Q(1, 1) * Q(2, 3) + Q(1, 3) * Q(2, 1),
-            std::sqrt(2) * Q(2, 1) * Q(3, 1), std::sqrt(2) * Q(2, 2) * Q(3, 2),
-            std::sqrt(2) * Q(2, 3) * Q(3, 3),
-            Q(2, 1) * Q(3, 2) + Q(2, 2) * Q(3, 1),
-            Q(2, 2) * Q(3, 3) + Q(2, 3) * Q(3, 2),
-            Q(2, 1) * Q(3, 3) + Q(2, 3) * Q(3, 1),
-            std::sqrt(2) * Q(1, 1) * Q(3, 1), std::sqrt(2) * Q(1, 2) * Q(3, 2),
-            std::sqrt(2) * Q(1, 3) * Q(3, 3),
-            Q(1, 1) * Q(3, 2) + Q(1, 2) * Q(3, 1),
-            Q(1, 2) * Q(3, 3) + Q(1, 3) * Q(3, 2),
-            Q(1, 1) * Q(3, 3) + Q(1, 3) * Q(3, 1);
-        return R;
-    }
-
 protected:
     MaterialProperties _mp;
     boost::optional<ParameterLib::CoordinateSystem> const&

MathLib/KelvinVector.cpp

@@ -178,5 +178,37 @@ kelvinVectorToSymmetricTensor(Eigen::Matrix<double,
         "or 6, but a vector of size %d was given.",
         v.size());
 }
+
+KelvinMatrixType<3> fourthOrderRotationMatrix(
+    Eigen::Matrix3d const& transformation)
+{
+    // 1-based index access for convenience.
+    auto Q = [&](int const i, int const j) {
+        return transformation(i - 1, j - 1);
+    };
+
+    MathLib::KelvinVector::KelvinMatrixType<3> R;
+    R << Q(1, 1) * Q(1, 1), Q(1, 2) * Q(1, 2), Q(1, 3) * Q(1, 3),
+        std::sqrt(2) * Q(1, 1) * Q(1, 2), std::sqrt(2) * Q(1, 2) * Q(1, 3),
+        std::sqrt(2) * Q(1, 1) * Q(1, 3), Q(2, 1) * Q(2, 1), Q(2, 2) * Q(2, 2),
+        Q(2, 3) * Q(2, 3), std::sqrt(2) * Q(2, 1) * Q(2, 2),
+        std::sqrt(2) * Q(2, 2) * Q(2, 3), std::sqrt(2) * Q(2, 1) * Q(2, 3),
+        Q(3, 1) * Q(3, 1), Q(3, 2) * Q(3, 2), Q(3, 3) * Q(3, 3),
+        std::sqrt(2) * Q(3, 1) * Q(3, 2), std::sqrt(2) * Q(3, 2) * Q(3, 3),
+        std::sqrt(2) * Q(3, 1) * Q(3, 3), std::sqrt(2) * Q(1, 1) * Q(2, 1),
+        std::sqrt(2) * Q(1, 2) * Q(2, 2), std::sqrt(2) * Q(1, 3) * Q(2, 3),
+        Q(1, 1) * Q(2, 2) + Q(1, 2) * Q(2, 1),
+        Q(1, 2) * Q(2, 3) + Q(1, 3) * Q(2, 2),
+        Q(1, 1) * Q(2, 3) + Q(1, 3) * Q(2, 1), std::sqrt(2) * Q(2, 1) * Q(3, 1),
+        std::sqrt(2) * Q(2, 2) * Q(3, 2), std::sqrt(2) * Q(2, 3) * Q(3, 3),
+        Q(2, 1) * Q(3, 2) + Q(2, 2) * Q(3, 1),
+        Q(2, 2) * Q(3, 3) + Q(2, 3) * Q(3, 2),
+        Q(2, 1) * Q(3, 3) + Q(2, 3) * Q(3, 1), std::sqrt(2) * Q(1, 1) * Q(3, 1),
+        std::sqrt(2) * Q(1, 2) * Q(3, 2), std::sqrt(2) * Q(1, 3) * Q(3, 3),
+        Q(1, 1) * Q(3, 2) + Q(1, 2) * Q(3, 1),
+        Q(1, 2) * Q(3, 3) + Q(1, 3) * Q(3, 2),
+        Q(1, 1) * Q(3, 3) + Q(1, 3) * Q(3, 1);
+    return R;
+}
 }  // namespace KelvinVector
 }  // namespace MathLib

MathLib/KelvinVector.h

@@ -205,6 +205,13 @@ symmetricTensorToKelvinVector(Eigen::MatrixBase<Derived> const& v)
     }
     return result;
 }
+
+/// Rotation tensor for Kelvin mapped vectors and tensors. It is meant to be
+/// used for rotation of stress/strain tensors epsilon:Q and tangent stiffness
+/// tensors Q*C*Q^t.
+KelvinMatrixType<3> fourthOrderRotationMatrix(
+    Eigen::Matrix3d const& transformation);
+
 }  // namespace KelvinVector
 }  // namespace MathLib