diff --git a/MathLib/MathTools.h b/MathLib/MathTools.h
index 40d09f3ad688feaae4f4e996992a808e844e0f4c..be8b7016db7b829002f9f74ed0a4b49f0fe35241 100644
--- a/MathLib/MathTools.h
+++ b/MathLib/MathTools.h
@@ -166,6 +166,16 @@ T fastpow (T base, std::size_t exp)
return result;
}
+/**
+ * Template metaprogramming, compile-time version of pow() for integral
+ * exponents.
+ */
+template <typename T>
+inline constexpr T pow(T const x, unsigned const y)
+{
+ return (y == 0) ? 1 : x * pow(x, y - 1);
+};
+
} // namespace
#endif /* MATHTOOLS_H_ */
diff --git a/MeshLib/Elements/Cell.cpp b/MeshLib/Elements/Cell.cpp
index e9a549e94d061ce3930bffe01e8c6d86240c0546..357d89784e7130879808e24f37103fe6c61d7dbe 100644
--- a/MeshLib/Elements/Cell.cpp
+++ b/MeshLib/Elements/Cell.cpp
@@ -19,8 +19,7 @@
namespace MeshLib {
-const unsigned Cell::dimension = 3u;
-
+constexpr unsigned Cell::dimension;
/*
Cell::Cell(Node** nodes, MeshElemType type, unsigned value)
: Element(nodes, type, value)
diff --git a/MeshLib/Elements/Cell.h b/MeshLib/Elements/Cell.h
index dfe4b9b9946df75c0f6c4626927103dc697b1084..150a46b0ea540be4afae3c640d830e5b6891cc2c 100644
--- a/MeshLib/Elements/Cell.h
+++ b/MeshLib/Elements/Cell.h
@@ -27,7 +27,7 @@ class Cell : public Element
{
public:
/// Constant: Dimension of this mesh element
- static const unsigned dimension;
+ static constexpr unsigned dimension = 3;
/// Returns the length, area or volume of a 1D, 2D or 3D element
double getContent() const { return _volume; };
diff --git a/MeshLib/Elements/Face.cpp b/MeshLib/Elements/Face.cpp
index f741b6e4baaa80568b89efbe4d543ccf9eeb5524..dfd145af46b5f57517925f5ad21e790bf7de8d6a 100644
--- a/MeshLib/Elements/Face.cpp
+++ b/MeshLib/Elements/Face.cpp
@@ -21,7 +21,7 @@
namespace MeshLib {
-const unsigned Face::dimension = 2u;
+constexpr unsigned Face::dimension;
/*
Face::Face(Node** nodes, MeshElemType type, unsigned value)
diff --git a/MeshLib/Elements/Face.h b/MeshLib/Elements/Face.h
index b338d7a2694f279cdb27b2ea5163fdde91d52d7f..1bc0b8c452a27919470ffc2b0b655221e6006320 100644
--- a/MeshLib/Elements/Face.h
+++ b/MeshLib/Elements/Face.h
@@ -32,7 +32,7 @@ class Face : public Element
{
public:
/// Constant: Dimension of this mesh element
- static const unsigned dimension;
+ static constexpr unsigned dimension = 2;
/// Get the area of this 2d element.
virtual double getArea() const { return _area; };
diff --git a/NumLib/Fem/CoordinatesMapping/ShapeMatrices.h b/NumLib/Fem/CoordinatesMapping/ShapeMatrices.h
index 51ce2877d444171e8edb72314048f23ab7d4bd52..19788db1829fa333e5315e6f08eeb6f86dd4e674 100644
--- a/NumLib/Fem/CoordinatesMapping/ShapeMatrices.h
+++ b/NumLib/Fem/CoordinatesMapping/ShapeMatrices.h
@@ -53,6 +53,12 @@ struct ShapeMatrices
JacobianType invJ; ///< Inverse matrix of the Jacobian
DShapeType dNdx; ///< Matrix of gradient of shape functions in physical coordinates, dN(r)/dx
+ /** The default constructor is used by fixed-size (at compile-time)
+ * matrix/vector types where no resizing of the matrices/vectors is required
+ * in this constructor.
+ */
+ ShapeMatrices() = default;
+
/**
* Initialize matrices and vectors
*
diff --git a/NumLib/Fem/FiniteElement/C0IsoparametricElements.h b/NumLib/Fem/FiniteElement/C0IsoparametricElements.h
index f05d76d5509e828b4a0e79b66794a623069959c5..2c94faf678b67c84c422f77c4c7654d5b0adbe0d 100644
--- a/NumLib/Fem/FiniteElement/C0IsoparametricElements.h
+++ b/NumLib/Fem/FiniteElement/C0IsoparametricElements.h
@@ -22,36 +22,34 @@
#include "NumLib/Fem/ShapeFunction/ShapeTri3.h"
#include "NumLib/Fem/ShapeFunction/ShapeQuad4.h"
#include "NumLib/Fem/ShapeFunction/ShapeHex8.h"
-#include "NumLib/Fem/Integration/IntegrationGaussRegular.h"
-#include "NumLib/Fem/Integration/IntegrationGaussTri.h"
#include "TemplateIsoparametric.h"
namespace NumLib
{
-template <class T_SHAPE_VECTOR, class T_DSHAPE_MATRIX, class T_JACOBIAN_MATRIX>
+template <template <typename> class T_SHAPE_MATRIX_POLICY>
struct FeLINE2
{
- typedef TemplateIsoparametric<MeshLib::Line, ShapeLine2, IntegrationGaussRegular<1>, T_SHAPE_VECTOR, T_DSHAPE_MATRIX, T_JACOBIAN_MATRIX> type;
+ typedef TemplateIsoparametric<ShapeLine2, T_SHAPE_MATRIX_POLICY<ShapeLine2>> type;
};
-template <class T_SHAPE_VECTOR, class T_DSHAPE_MATRIX, class T_JACOBIAN_MATRIX>
+template <template <typename> class T_SHAPE_MATRIX_POLICY>
struct FeTRI3
{
- typedef TemplateIsoparametric<MeshLib::Tri, ShapeTri3, IntegrationGaussTri, T_SHAPE_VECTOR, T_DSHAPE_MATRIX, T_JACOBIAN_MATRIX> type;
+ typedef TemplateIsoparametric<ShapeTri3, T_SHAPE_MATRIX_POLICY<ShapeTri3>> type;
};
-template <class T_SHAPE_VECTOR, class T_DSHAPE_MATRIX, class T_JACOBIAN_MATRIX>
+template <template <typename> class T_SHAPE_MATRIX_POLICY>
struct FeQUAD4
{
- typedef TemplateIsoparametric<MeshLib::Quad, ShapeQuad4, IntegrationGaussRegular<2>, T_SHAPE_VECTOR, T_DSHAPE_MATRIX, T_JACOBIAN_MATRIX> type;
+ typedef TemplateIsoparametric<ShapeQuad4, T_SHAPE_MATRIX_POLICY<ShapeQuad4>> type;
};
-template <class T_SHAPE_VECTOR, class T_DSHAPE_MATRIX, class T_JACOBIAN_MATRIX>
+template <template <typename> class T_SHAPE_MATRIX_POLICY>
struct FeHEX8
{
- typedef TemplateIsoparametric<MeshLib::Hex, ShapeHex8, IntegrationGaussRegular<3>, T_SHAPE_VECTOR, T_DSHAPE_MATRIX, T_JACOBIAN_MATRIX> type;
+ typedef TemplateIsoparametric<ShapeHex8, T_SHAPE_MATRIX_POLICY<ShapeHex8>> type;
};
} // NumLib
diff --git a/NumLib/Fem/FiniteElement/TemplateIsoparametric.h b/NumLib/Fem/FiniteElement/TemplateIsoparametric.h
index 3659dc5c85ff7d3206451e619ccd073eedf7b260..7d3166f0e2a40f33cdcbe63613a938152d0e1639 100644
--- a/NumLib/Fem/FiniteElement/TemplateIsoparametric.h
+++ b/NumLib/Fem/FiniteElement/TemplateIsoparametric.h
@@ -26,32 +26,28 @@ namespace NumLib
/**
* \brief Template class for isoparametric elements
*
- * \tparam T_ELEMENT Mesh element class
- * \tparam T_SHAPE Shape function
- * \tparam T_INTEGRAL Integration method
- * \tparam T_NODAL_VECTOR Nodal vector class
- * \tparam T_DIM_NODAL_MATRIX Matrix class for a size of dim * nnodes
- * \tparam T_DIM_MATRIX Matrix class for a size of dim * dim
+ * \tparam ShapeFunctionType_ The shape function type.
+ * \tparam ShapeMatrixTypes_ An aggregate of shape matrix types.
*/
template <
- class T_MESH_ELEMENT,
- class T_SHAPE,
- class T_INTEGRAL,
- class T_NODAL_VECTOR,
- class T_DIM_NODAL_MATRIX,
- class T_DIM_MATRIX
+ class ShapeFunctionType_,
+ class ShapeMatrixTypes_
>
class TemplateIsoparametric
{
public:
- typedef T_MESH_ELEMENT MeshElementType;
- typedef T_SHAPE ShapeFunctionType;
- typedef T_INTEGRAL IntegrationMethod;
- typedef T_NODAL_VECTOR NodalVectorType;
- typedef T_DIM_NODAL_MATRIX DimNodalMatrixType;
- typedef T_DIM_MATRIX DimMatrixType;
- typedef ShapeMatrices<NodalVectorType, DimNodalMatrixType, DimMatrixType> ShapeMatricesType;
- typedef NaturalCoordinatesMapping<MeshElementType, ShapeFunctionType, ShapeMatricesType> NaturalCoordsMappingType;
+ using ShapeFunctionType = ShapeFunctionType_;
+
+ /// Coordinate mapping matrices type.
+ using ShapeMatrices = typename ShapeMatrixTypes_::ShapeMatrices;
+
+ /// Type of the underlying geometrical element.
+ using MeshElementType = typename ShapeFunctionType_::MeshElement;
+
+ /// Natural coordinates mapping tools specialization for specific
+ /// MeshElement and ShapeFunction types.
+ using NaturalCoordsMappingType = NaturalCoordinatesMapping<
+ MeshElementType, ShapeFunctionType, ShapeMatrices>;
/**
* Constructor without specifying a mesh element. setMeshElement() must be called afterwards.
@@ -88,7 +84,7 @@ public:
* @param natural_pt position in natural coordinates
* @param shape evaluated shape function matrices
*/
- void computeShapeFunctions(const double *natural_pt, ShapeMatricesType &shape) const
+ void computeShapeFunctions(const double *natural_pt, ShapeMatrices &shape) const
{
NaturalCoordsMappingType::computeShapeMatrices(*_ele, natural_pt, shape);
}
@@ -101,7 +97,7 @@ public:
* @param shape evaluated shape function matrices
*/
template <ShapeMatrixType T_SHAPE_MATRIX_TYPE>
- void computeShapeFunctions(const double *natural_pt, ShapeMatricesType &shape) const
+ void computeShapeFunctions(const double *natural_pt, ShapeMatrices &shape) const
{
NaturalCoordsMappingType::template computeShapeMatrices<T_SHAPE_MATRIX_TYPE>(*_ele, natural_pt, shape);
}
diff --git a/NumLib/Fem/Integration/GaussIntegrationPolicy.h b/NumLib/Fem/Integration/GaussIntegrationPolicy.h
new file mode 100644
index 0000000000000000000000000000000000000000..5f8707fae27ae28e4bc2df801cd1eee505cc88bc
--- /dev/null
+++ b/NumLib/Fem/Integration/GaussIntegrationPolicy.h
@@ -0,0 +1,44 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2014, 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 NUMLIB_GAUSSINTEGRATIONPOLICY_H_
+#define NUMLIB_GAUSSINTEGRATIONPOLICY_H_
+
+#include "NumLib/Fem/Integration/IntegrationGaussRegular.h"
+#include "NumLib/Fem/Integration/IntegrationGaussTri.h"
+
+namespace NumLib
+{
+
+/// An integration policy providing an integration method suitable for the given
+/// mesh element.
+/// Gauss-Legendre integration method is used. The default implementation is
+/// choosing the regular-element integration, for other elements a
+/// specialization must be provided, as for example for triangles
+/// (\see GaussIntegrationPolicy<MeshLib::Tri>).
+/// The integration method depends on the dimension of the element and correctly
+/// chosen number and placement of the integration points within the element.
+template <typename MeshElement_>
+struct GaussIntegrationPolicy
+{
+ using MeshElement = MeshElement_;
+ using IntegrationMethod =
+ NumLib::IntegrationGaussRegular<MeshElement::dimension>;
+};
+
+template<>
+struct GaussIntegrationPolicy<MeshLib::Tri>
+{
+ using MeshElement = MeshLib::Tri;
+ using IntegrationMethod = NumLib::IntegrationGaussTri;
+};
+
+} // namespace NumLib
+
+#endif // NUMLIB_GAUSSINTEGRATIONPOLICY_H_
diff --git a/NumLib/Fem/ShapeFunction/ShapeHex8.h b/NumLib/Fem/ShapeFunction/ShapeHex8.h
index 4d07b1d3e796ec901b3e7fc090fea36a4b65603d..a04325c26e4fde2312c157d7c91d9c8f6ea5f26b 100644
--- a/NumLib/Fem/ShapeFunction/ShapeHex8.h
+++ b/NumLib/Fem/ShapeFunction/ShapeHex8.h
@@ -10,6 +10,8 @@
#ifndef SHAPEHEX8_H_
#define SHAPEHEX8_H_
+#include "MeshLib/Elements/Hex.h"
+
namespace NumLib
{
@@ -55,6 +57,10 @@ public:
*/
template <class T_X, class T_N>
static void computeGradShapeFunction(const T_X &r, T_N &dN);
+
+ using MeshElement = MeshLib::Hex;
+ static constexpr std::size_t DIM = MeshElement::dimension;
+ static constexpr std::size_t NPOINTS = MeshElement::n_all_nodes;
};
}
diff --git a/NumLib/Fem/ShapeFunction/ShapeLine2.h b/NumLib/Fem/ShapeFunction/ShapeLine2.h
index e372e6abb7dccd45a363ed4186f3963d8da66be7..3f7974f2010b07cd1f0640b4d5e34b899079678c 100644
--- a/NumLib/Fem/ShapeFunction/ShapeLine2.h
+++ b/NumLib/Fem/ShapeFunction/ShapeLine2.h
@@ -10,6 +10,8 @@
#ifndef SHAPELINE2_H_
#define SHAPELINE2_H_
+#include "MeshLib/Elements/Line.h"
+
namespace NumLib
{
@@ -41,6 +43,10 @@ public:
*/
template <class T_X, class T_N>
static void computeGradShapeFunction(const T_X &r, T_N &dN);
+
+ using MeshElement = MeshLib::Line;
+ static constexpr std::size_t DIM = MeshElement::dimension;
+ static constexpr std::size_t NPOINTS = MeshElement::n_all_nodes;
};
}
diff --git a/NumLib/Fem/ShapeFunction/ShapeQuad4.h b/NumLib/Fem/ShapeFunction/ShapeQuad4.h
index bacc6ac41ef1a4a50c310eab7c091dedfa41efb9..5da243edddb1628678dd3c146ab7d0211e8add81 100644
--- a/NumLib/Fem/ShapeFunction/ShapeQuad4.h
+++ b/NumLib/Fem/ShapeFunction/ShapeQuad4.h
@@ -13,6 +13,8 @@
#ifndef SHAPEQUAD4_H_
#define SHAPEQUAD4_H_
+#include "MeshLib/Elements/Quad.h"
+
namespace NumLib
{
@@ -49,6 +51,10 @@ public:
*/
template <class T_X, class T_N>
static void computeGradShapeFunction(const T_X &r, T_N &dN);
+
+ using MeshElement = MeshLib::Quad;
+ static constexpr std::size_t DIM = MeshElement::dimension;
+ static constexpr std::size_t NPOINTS = MeshElement::n_all_nodes;
};
}
diff --git a/NumLib/Fem/ShapeFunction/ShapeTri3.h b/NumLib/Fem/ShapeFunction/ShapeTri3.h
index eaca0f915eea2c5dd289a55575156288a7696b43..d5148e3fdd7b0b8c19be305091140a04384ad5a5 100644
--- a/NumLib/Fem/ShapeFunction/ShapeTri3.h
+++ b/NumLib/Fem/ShapeFunction/ShapeTri3.h
@@ -15,6 +15,8 @@
#ifndef SHAPETRI3_H_
#define SHAPETRI3_H_
+#include "MeshLib/Elements/Tri.h"
+
namespace NumLib
{
@@ -49,6 +51,10 @@ public:
*/
template <class T_X, class T_N>
static void computeGradShapeFunction(const T_X &/*r*/, T_N &dN6);
+
+ using MeshElement = MeshLib::Tri;
+ static constexpr std::size_t DIM = MeshElement::dimension;
+ static constexpr std::size_t NPOINTS = MeshElement::n_all_nodes;
};
}
diff --git a/NumLib/Fem/ShapeMatrixPolicy.h b/NumLib/Fem/ShapeMatrixPolicy.h
new file mode 100644
index 0000000000000000000000000000000000000000..36ab839dc153d62032b1044d37964a516b159977
--- /dev/null
+++ b/NumLib/Fem/ShapeMatrixPolicy.h
@@ -0,0 +1,74 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2014, 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 SHAPE_MATRIX_POLICY_H_
+#define SHAPE_MATRIX_POLICY_H_
+
+#include "NumLib/Fem/CoordinatesMapping/ShapeMatrices.h"
+
+#ifdef OGS_USE_EIGEN
+#include <Eigen/Eigen>
+
+/// An implementation of ShapeMatrixPolicy using fixed size (compile-time) eigen
+/// matrices and vectors.
+template <typename ShapeFunction>
+struct EigenFixedShapeMatrixPolicy
+{
+ template <int N, int M>
+ using _MatrixType = Eigen::Matrix<double, N, M, Eigen::RowMajor>;
+ template <int N>
+ using _VectorType = Eigen::Matrix<double, N, 1>;
+
+ using NodalMatrixType = _MatrixType<ShapeFunction::NPOINTS, ShapeFunction::NPOINTS>;
+ using NodalVectorType = _VectorType<ShapeFunction::NPOINTS>;
+ using DimNodalMatrixType = _MatrixType<ShapeFunction::DIM, ShapeFunction::NPOINTS>;
+ using DimMatrixType = _MatrixType<ShapeFunction::DIM, ShapeFunction::DIM>;
+
+ using ShapeMatrices =
+ NumLib::ShapeMatrices<
+ NodalVectorType,
+ DimNodalMatrixType,
+ DimMatrixType>;
+};
+
+/// An implementation of ShapeMatrixPolicy using dynamic size eigen matrices and
+/// vectors.
+template <typename ShapeFunction>
+struct EigenDynamicShapeMatrixPolicy
+{
+ // Dynamic size local matrices are much slower in allocation than their
+ // fixed counterparts.
+
+ using _MatrixType =
+ Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>;
+ using _VectorType =
+ Eigen::Matrix<double, Eigen::Dynamic, 1>;
+
+ using NodalMatrixType = _MatrixType;
+ using NodalVectorType = _VectorType;
+ using DimNodalMatrixType = _MatrixType;
+ using DimMatrixType = _MatrixType;
+
+ using ShapeMatrices =
+ NumLib::ShapeMatrices<
+ NodalVectorType,
+ DimNodalMatrixType,
+ DimMatrixType>;
+};
+
+/// Default choice of the ShapeMatrixPolicy.
+template <typename ShapeFunction>
+using ShapeMatrixPolicyType = EigenFixedShapeMatrixPolicy<ShapeFunction>;
+
+#endif // OGS_USE_EIGEN
+
+//static_assert(std::is_class<ShapeMatrixPolicyType<>::value,
+ //"ShapeMatrixPolicyType was not defined.");
+
+#endif // SHAPE_MATRIX_POLICY_H_
diff --git a/Tests/NumLib/FeTestData/TestFeHEX8.h b/Tests/NumLib/FeTestData/TestFeHEX8.h
index c2c12df6ac4bb9487507e8f27e9a9ed34038c404..f738e26adbb57212263f88444d1e0ddba62a0179 100644
--- a/Tests/NumLib/FeTestData/TestFeHEX8.h
+++ b/Tests/NumLib/FeTestData/TestFeHEX8.h
@@ -22,12 +22,12 @@ namespace FeTestData
class TestFeHEX8
{
public:
+ using ShapeFunction = NumLib::ShapeHex8;
+
// Fe type information
- template <class T_MATRIX_TYPES>
- struct FeType
- {
- typedef NumLib::FeHEX8<typename T_MATRIX_TYPES::NodalVectorType, typename T_MATRIX_TYPES::DimNodalMatrixType, typename T_MATRIX_TYPES::DimMatrixType> type;
- };
+ template <template <typename> class ShapeMatrixPolicy_>
+ using FeType = NumLib::FeHEX8<ShapeMatrixPolicy_>;
+
typedef MeshLib::Hex MeshElementType;
static const unsigned dim = 3; //MeshElementType::dimension;
static const unsigned e_nnodes = MeshElementType::n_all_nodes;
diff --git a/Tests/NumLib/FeTestData/TestFeLINE2.h b/Tests/NumLib/FeTestData/TestFeLINE2.h
index f3d674b39fad45460409580087e3e1d56d9f49bf..80b6cff559ad5ed456a189459d1f7218f68a1aaa 100644
--- a/Tests/NumLib/FeTestData/TestFeLINE2.h
+++ b/Tests/NumLib/FeTestData/TestFeLINE2.h
@@ -22,12 +22,12 @@ namespace FeTestData
class TestFeLINE2
{
public:
+ using ShapeFunction = NumLib::ShapeLine2;
+
// Fe type information
- template <class T_MATRIX_TYPES>
- struct FeType
- {
- typedef NumLib::FeLINE2<typename T_MATRIX_TYPES::NodalVectorType, typename T_MATRIX_TYPES::DimNodalMatrixType, typename T_MATRIX_TYPES::DimMatrixType> type;
- };
+ template <template <typename> class ShapeMatrixPolicy_>
+ using FeType = NumLib::FeLINE2<ShapeMatrixPolicy_>;
+
typedef MeshLib::Line MeshElementType;
static const unsigned dim = MeshElementType::dimension;
static const unsigned e_nnodes = MeshElementType::n_all_nodes;
diff --git a/Tests/NumLib/FeTestData/TestFeQUAD4.h b/Tests/NumLib/FeTestData/TestFeQUAD4.h
index afa0c64bb7aece24ef030e2607cd790a3bb6f0b2..bc9517a90b4a3960f85ac21c634369289b1d4cc4 100644
--- a/Tests/NumLib/FeTestData/TestFeQUAD4.h
+++ b/Tests/NumLib/FeTestData/TestFeQUAD4.h
@@ -22,12 +22,12 @@ namespace FeTestData
class TestFeQUAD4
{
public:
+ using ShapeFunction = NumLib::ShapeQuad4;
+
// Fe type information
- template <class T_MATRIX_TYPES>
- struct FeType
- {
- typedef NumLib::FeQUAD4<typename T_MATRIX_TYPES::NodalVectorType, typename T_MATRIX_TYPES::DimNodalMatrixType, typename T_MATRIX_TYPES::DimMatrixType> type;
- };
+ template <template <typename> class ShapeMatrixPolicy_>
+ using FeType = NumLib::FeQUAD4<ShapeMatrixPolicy_>;
+
typedef MeshLib::Quad MeshElementType;
static const unsigned dim = 2; //MeshElementType::dimension;
static const unsigned e_nnodes = MeshElementType::n_all_nodes;
diff --git a/Tests/NumLib/FeTestData/TestFeTRI3.h b/Tests/NumLib/FeTestData/TestFeTRI3.h
index f800e5093ed8c264311c45781c485211adba3dbb..2d40ae6232f2824f8a6bbeef9321e9bcb89a95d1 100644
--- a/Tests/NumLib/FeTestData/TestFeTRI3.h
+++ b/Tests/NumLib/FeTestData/TestFeTRI3.h
@@ -22,12 +22,12 @@ namespace FeTestData
class TestFeTRI3
{
public:
+ using ShapeFunction = NumLib::ShapeTri3;
+
// Fe type information
- template <class T_MATRIX_TYPES>
- struct FeType
- {
- typedef NumLib::FeTRI3<typename T_MATRIX_TYPES::NodalVectorType, typename T_MATRIX_TYPES::DimNodalMatrixType, typename T_MATRIX_TYPES::DimMatrixType> type;
- };
+ template <template <typename> class ShapeMatrixPolicy_>
+ using FeType = NumLib::FeTRI3<ShapeMatrixPolicy_>;
+
typedef MeshLib::Tri MeshElementType;
static const unsigned dim = 2; //MeshElementType::dimension;
static const unsigned e_nnodes = MeshElementType::n_all_nodes;
diff --git a/Tests/NumLib/TestFe.cpp b/Tests/NumLib/TestFe.cpp
index aee03d1c5a499ec7d52b7f83bab9ffa501decf63..81d0e72d6a82acbb89fe29b4317835fa09a2b734 100644
--- a/Tests/NumLib/TestFe.cpp
+++ b/Tests/NumLib/TestFe.cpp
@@ -17,6 +17,8 @@
#include "NumLib/Fem/CoordinatesMapping/ShapeMatrices.h"
#include "NumLib/Fem/FiniteElement/C0IsoparametricElements.h"
+#include "NumLib/Fem/Integration/GaussIntegrationPolicy.h"
+#include "NumLib/Fem/ShapeMatrixPolicy.h"
#include "FeTestData/TestFeLINE2.h"
#include "FeTestData/TestFeTRI3.h"
@@ -30,71 +32,60 @@ using namespace FeTestData;
namespace
{
-#ifdef OGS_USE_EIGEN
-template <typename T_FE>
-struct EigenFixedMatrixTypes
-{
- typedef Eigen::Matrix<double, T_FE::e_nnodes, T_FE::e_nnodes, Eigen::RowMajor> NodalMatrixType;
- typedef Eigen::Matrix<double, T_FE::e_nnodes, 1> NodalVectorType;
- typedef Eigen::Matrix<double, T_FE::dim, T_FE::e_nnodes, Eigen::RowMajor> DimNodalMatrixType;
- typedef Eigen::Matrix<double, T_FE::dim, T_FE::dim, Eigen::RowMajor> DimMatrixType;
-};
-
-struct EigenDynamicMatrixTypes
-{
- typedef Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> NodalMatrixType;
- typedef NodalMatrixType DimNodalMatrixType;
- typedef NodalMatrixType DimMatrixType;
- typedef Eigen::VectorXd NodalVectorType;
-};
-
-typedef EigenDynamicMatrixTypes DefaultMatrixType;
-#endif // OGS_USE_EIGEN
// test cases
-template <class T_FE_TYPE, class T_MAT_TYPES=DefaultMatrixType>
+template <class TestFeType_, template <typename> class ShapeMatrixPolicy_>
struct TestCase
{
- typedef T_FE_TYPE T_FE;
- typedef T_MAT_TYPES T_MATRIX_TYPES;
+ typedef TestFeType_ TestFeType;
+ typedef ShapeMatrixPolicy_<typename TestFeType::ShapeFunction> ShapeMatrixTypes;
+ template <typename X>
+ using ShapeMatrixPolicy = ShapeMatrixPolicy_<X>;
};
typedef ::testing::Types<
- TestCase<TestFeLINE2>,
- TestCase<TestFeTRI3>,
- TestCase<TestFeQUAD4>,
- TestCase<TestFeHEX8>
#ifdef OGS_USE_EIGEN
- ,TestCase<TestFeLINE2, EigenFixedMatrixTypes<TestFeLINE2> >
- ,TestCase<TestFeTRI3, EigenFixedMatrixTypes<TestFeTRI3> >
- ,TestCase<TestFeQUAD4, EigenFixedMatrixTypes<TestFeQUAD4> >
- ,TestCase<TestFeHEX8, EigenFixedMatrixTypes<TestFeHEX8> >
+ TestCase<TestFeLINE2, EigenDynamicShapeMatrixPolicy>,
+ TestCase<TestFeTRI3, EigenDynamicShapeMatrixPolicy>,
+ TestCase<TestFeQUAD4, EigenDynamicShapeMatrixPolicy>,
+ TestCase<TestFeHEX8, EigenDynamicShapeMatrixPolicy>,
+ TestCase<TestFeLINE2, EigenFixedShapeMatrixPolicy>,
+ TestCase<TestFeTRI3, EigenFixedShapeMatrixPolicy>,
+ TestCase<TestFeQUAD4, EigenFixedShapeMatrixPolicy>,
+ TestCase<TestFeHEX8, EigenFixedShapeMatrixPolicy>
#endif
> TestTypes;
}
template <class T>
-class NumLibFemIsoTest : public ::testing::Test, public T::T_FE
+class NumLibFemIsoTest : public ::testing::Test, public T::TestFeType
{
public:
- typedef typename T::T_MATRIX_TYPES T_MATRIX_TYPES;
- typedef typename T::T_FE T_FE;
+ typedef typename T::ShapeMatrixTypes ShapeMatrixTypes;
+ typedef typename T::TestFeType TestFeType;
// Matrix types
- typedef T_MATRIX_TYPES MatrixType;
- typedef typename T_MATRIX_TYPES::NodalMatrixType NodalMatrix;
- typedef typename T_MATRIX_TYPES::NodalVectorType NodalVector;
- typedef typename T_MATRIX_TYPES::DimNodalMatrixType DimNodalMatrix;
- typedef typename T_MATRIX_TYPES::DimMatrixType DimMatrix;
+ typedef typename ShapeMatrixTypes::NodalMatrixType NodalMatrix;
+ typedef typename ShapeMatrixTypes::NodalVectorType NodalVector;
+ typedef typename ShapeMatrixTypes::DimNodalMatrixType DimNodalMatrix;
+ typedef typename ShapeMatrixTypes::DimMatrixType DimMatrix;
+
// Finite element type
- typedef typename T_FE::template FeType<T_MATRIX_TYPES>::type::type FeType;
+ template <typename X>
+ using ShapeMatrixPolicy = typename T::template ShapeMatrixPolicy<X>;
+ typedef typename TestFeType::template FeType<ShapeMatrixPolicy>::type FeType;
+
// Shape matrix data type
- typedef typename FeType::ShapeMatricesType ShapeMatricesType;
- typedef typename T_FE::MeshElementType MeshElementType;
+ typedef typename ShapeMatrixTypes::ShapeMatrices ShapeMatricesType;
+ typedef typename TestFeType::MeshElementType MeshElementType;
+
+ static const unsigned dim = TestFeType::dim;
+ static const unsigned e_nnodes = TestFeType::e_nnodes;
+ static const unsigned n_sample_pt_order2 = TestFeType::n_sample_pt_order2;
+ static const unsigned n_sample_pt_order3 = TestFeType::n_sample_pt_order3;
+
+ using IntegrationMethod =
+ typename NumLib::GaussIntegrationPolicy<MeshElementType>::IntegrationMethod;
- static const unsigned dim = T_FE::dim;
- static const unsigned e_nnodes = T_FE::e_nnodes;
- static const unsigned n_sample_pt_order2 = T_FE::n_sample_pt_order2;
- static const unsigned n_sample_pt_order3 = T_FE::n_sample_pt_order3;
public:
NumLibFemIsoTest() :
@@ -135,7 +126,7 @@ class NumLibFemIsoTest : public ::testing::Test, public T::T_FE
DimMatrix D;
NodalMatrix expectedM;
NodalMatrix expectedK;
- typename FeType::IntegrationMethod integration_method;
+ IntegrationMethod integration_method;
std::vector<const MeshLib::Node*> vec_nodes;
std::vector<const MeshElementType*> vec_eles;