From 96c96de7eec42a830b6be7f3b9aa54c2c35d45da Mon Sep 17 00:00:00 2001
From: Dmitri Naumov <github@naumov.de>
Date: Wed, 1 Jun 2016 00:35:35 +0000
Subject: [PATCH] [PL] Implement BMatrices and Kelvin vec/mat types.
---
ProcessLib/Deformation/BMatrixPolicy.h | 94 ++++++++++++++++++++++++++
ProcessLib/Deformation/LinearBMatrix.h | 58 ++++++++++++++++
2 files changed, 152 insertions(+)
create mode 100644 ProcessLib/Deformation/BMatrixPolicy.h
create mode 100644 ProcessLib/Deformation/LinearBMatrix.h
diff --git a/ProcessLib/Deformation/BMatrixPolicy.h b/ProcessLib/Deformation/BMatrixPolicy.h
new file mode 100644
index 00000000000..00b33e4f442
--- /dev/null
+++ b/ProcessLib/Deformation/BMatrixPolicy.h
@@ -0,0 +1,94 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef PROCESSLIB_FEM_BMATRIXPOLICY_H_
+#define PROCESSLIB_FEM_BMATRIXPOLICY_H_
+
+#include "NumLib/Fem/ShapeMatrixPolicy.h"
+
+namespace ProcessLib
+{
+/// Kelvin vector dimensions for given displacement dimension.
+template <int DisplacementDim>
+struct KelvinVectorDimensions;
+
+template <>
+struct KelvinVectorDimensions<2>
+{
+ static int const value = 4;
+};
+
+template <>
+struct KelvinVectorDimensions<3>
+{
+ static int const value = 6;
+};
+
+//
+// Kelvin vector and matrix templates for given displacement dimension.
+//
+
+/// \TODO Maybe better to move the KelvinVector/MatrixType into
+/// MaterialLib/SolidModels/KelvinVector.h.
+
+/// Kelvin vector type for given displacement dimension.
+/// \note The Eigen vector is always a fixed size vector in contrast to the
+/// BMatrixPolicyType::StressVectorType.
+template <int DisplacementDim>
+using KelvinVectorType =
+ Eigen::Matrix<double, KelvinVectorDimensions<DisplacementDim>::value, 1,
+ Eigen::ColMajor>;
+
+/// Kelvin matrix type for given displacement dimension.
+template <int DisplacementDim>
+using KelvinMatrixType =
+ Eigen::Matrix<double, KelvinVectorDimensions<DisplacementDim>::value,
+ KelvinVectorDimensions<DisplacementDim>::value,
+ Eigen::RowMajor>;
+
+/// An implementation of B-Matrix policy using same matrix and vector types
+/// (fixed size or dynamic) as in the ShapeMatrixPolicyType.
+template <typename ShapeFunction, unsigned DisplacementDim>
+class BMatrixPolicyType
+{
+private:
+ /// Reusing the ShapeMatrixPolicy vector type.
+ template <int N>
+ using VectorType =
+ typename ShapeMatrixPolicyType<ShapeFunction,
+ DisplacementDim>::template VectorType<N>;
+
+ /// Reusing the ShapeMatrixPolicy matrix type.
+ template <int N, int M>
+ using MatrixType = typename ShapeMatrixPolicyType<
+ ShapeFunction, DisplacementDim>::template MatrixType<N, M>;
+
+ // Dimensions of specific b-matrix for n-points and displacement dimension.
+ static int const _number_of_dof = ShapeFunction::NPOINTS * DisplacementDim;
+ static int const _kelvin_vector_size =
+ KelvinVectorDimensions<DisplacementDim>::value;
+
+public:
+ using StiffnessMatrixType = MatrixType<_number_of_dof, _number_of_dof>;
+
+ /// Rhs residual
+ using NodalForceVectorType = VectorType<_number_of_dof>;
+
+ /// Sigma
+ using StressVectorType = VectorType<_kelvin_vector_size>;
+
+ /// C
+ using ModulusMatrixType =
+ MatrixType<_kelvin_vector_size, _kelvin_vector_size>;
+
+ using BMatrixType = MatrixType<_kelvin_vector_size, _number_of_dof>;
+};
+} // namespace ProcessLib
+
+#endif // PROCESSLIB_FEM_BMATRIXPOLICY_H_
diff --git a/ProcessLib/Deformation/LinearBMatrix.h b/ProcessLib/Deformation/LinearBMatrix.h
new file mode 100644
index 00000000000..093c222a1ad
--- /dev/null
+++ b/ProcessLib/Deformation/LinearBMatrix.h
@@ -0,0 +1,58 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef PROCESSLIB_FEM_DEFORMATION_LINEARBMATRIX_H_
+#define PROCESSLIB_FEM_DEFORMATION_LINEARBMATRIX_H_
+
+#include <cmath>
+
+namespace ProcessLib
+{
+namespace LinearBMatrix
+{
+/// Fills a B-matrix based on given shape function dN/dx values.
+template <int DisplacementDim, int NPOINTS, typename DNDX_Type,
+ typename BMatrixType>
+void computeBMatrix(DNDX_Type const& dNdx, BMatrixType& b_matrix)
+{
+ static_assert(0 < DisplacementDim && DisplacementDim <= 3,
+ "LinearBMatrix::computeBMatrix: DisplacementDim must be in "
+ "range [1,3].");
+
+ b_matrix.setZero();
+
+ switch (DisplacementDim)
+ {
+ case 3:
+ for (int i = 0; i < NPOINTS; ++i)
+ {
+ b_matrix(2, 2 * NPOINTS + i) = dNdx(2, i);
+ b_matrix(4, NPOINTS + i) = dNdx(2, i) / std::sqrt(2);
+ b_matrix(4, 2 * NPOINTS + i) = dNdx(1, i) / std::sqrt(2);
+ b_matrix(5, i) = dNdx(2, i) / std::sqrt(2);
+ b_matrix(5, 2 * NPOINTS + i) = dNdx(0, i) / std::sqrt(2);
+ }
+ // no break for fallthrough.
+ case 2:
+ for (int i = 0; i < NPOINTS; ++i)
+ {
+ b_matrix(1, NPOINTS + i) = dNdx(1, i);
+ b_matrix(3, i) = dNdx(1, i) / std::sqrt(2);
+ b_matrix(3, NPOINTS + i) = dNdx(0, i) / std::sqrt(2);
+ b_matrix(0, i) = dNdx(0, i);
+ }
+ default:
+ break;
+ }
+}
+
+} // namespace LinearBMatrix
+} // namespace ProcessLib
+
+#endif // PROCESSLIB_FEM_DEFORMATION_LINEARBMATRIX_H_
--
GitLab