Skip to content
Snippets Groups Projects
Commit f4e548e1 authored by Dmitri Naumov's avatar Dmitri Naumov Committed by GitHub
Browse files

Merge pull request #1520 from endJunction/ReduceNumberOfInstantiatedLocalAsms 

Reduce number of instantiated local asms
parents a1eefa79 00ac4eb4
No related branches found
No related tags found
No related merge requests found
Expand all Hide whitespace changes
Inline Side-by-side
ProcessLib/LIE/SmallDeformation/LocalAssembler/LocalDataInitializer.h
+ 126
− 96
View file @ f4e548e1
...@@ -14,13 +14,13 @@ ...@@ -14,13 +14,13 @@
#include <memory> #include <memory>
#include <typeindex> #include <typeindex>
#include <typeinfo> #include <typeinfo>
#include <type_traits>
#include <unordered_map> #include <unordered_map>
#include "MeshLib/Elements/Elements.h" #include "MeshLib/Elements/Elements.h"
#include "NumLib/DOF/LocalToGlobalIndexMap.h" #include "NumLib/DOF/LocalToGlobalIndexMap.h"
#include "NumLib/Fem/Integration/GaussIntegrationPolicy.h" #include "NumLib/Fem/Integration/GaussIntegrationPolicy.h"
#ifndef OGS_MAX_ELEMENT_DIM #ifndef OGS_MAX_ELEMENT_DIM
static_assert(false, "The macro OGS_MAX_ELEMENT_DIM is undefined."); static_assert(false, "The macro OGS_MAX_ELEMENT_DIM is undefined.");
#endif #endif
...@@ -47,7 +47,7 @@ static_assert(false, "The macro OGS_MAX_ELEMENT_ORDER is undefined."); ...@@ -47,7 +47,7 @@ static_assert(false, "The macro OGS_MAX_ELEMENT_ORDER is undefined.");
#ifdef OGS_ENABLE_ELEMENT_PRISM #ifdef OGS_ENABLE_ELEMENT_PRISM
#define ENABLED_ELEMENT_TYPE_PRISM 1u << 2 #define ENABLED_ELEMENT_TYPE_PRISM 1u << 2
#else #else
#define ENABLED_ELEMENT_TYPE_PRISM 0u #define ENABLED_ELEMENT_TYPE_PRISM 0u
#endif #endif
#ifdef OGS_ENABLE_ELEMENT_PYRAMID #ifdef OGS_ENABLE_ELEMENT_PYRAMID
...@@ -58,19 +58,18 @@ static_assert(false, "The macro OGS_MAX_ELEMENT_ORDER is undefined."); ...@@ -58,19 +58,18 @@ static_assert(false, "The macro OGS_MAX_ELEMENT_ORDER is undefined.");
// Dependent element types. // Dependent element types.
// Faces of tets, pyramids and prisms are triangles // Faces of tets, pyramids and prisms are triangles
#define ENABLED_ELEMENT_TYPE_TRI ((ENABLED_ELEMENT_TYPE_SIMPLEX) | (ENABLED_ELEMENT_TYPE_PYRAMID) \ #define ENABLED_ELEMENT_TYPE_TRI \
|(ENABLED_ELEMENT_TYPE_PRISM)) ((ENABLED_ELEMENT_TYPE_SIMPLEX) | (ENABLED_ELEMENT_TYPE_PYRAMID) | \
(ENABLED_ELEMENT_TYPE_PRISM))
// Faces of hexes, pyramids and prisms are quads // Faces of hexes, pyramids and prisms are quads
#define ENABLED_ELEMENT_TYPE_QUAD ((ENABLED_ELEMENT_TYPE_CUBOID) | (ENABLED_ELEMENT_TYPE_PYRAMID) \ #define ENABLED_ELEMENT_TYPE_QUAD \
|(ENABLED_ELEMENT_TYPE_PRISM)) ((ENABLED_ELEMENT_TYPE_CUBOID) | (ENABLED_ELEMENT_TYPE_PYRAMID) | \
(ENABLED_ELEMENT_TYPE_PRISM))
// All enabled element types // All enabled element types
#define OGS_ENABLED_ELEMENTS \ #define OGS_ENABLED_ELEMENTS \
( (ENABLED_ELEMENT_TYPE_SIMPLEX) \ ((ENABLED_ELEMENT_TYPE_SIMPLEX) | (ENABLED_ELEMENT_TYPE_CUBOID) | \
| (ENABLED_ELEMENT_TYPE_CUBOID ) \ (ENABLED_ELEMENT_TYPE_PYRAMID) | (ENABLED_ELEMENT_TYPE_PRISM))
| (ENABLED_ELEMENT_TYPE_PYRAMID) \
| (ENABLED_ELEMENT_TYPE_PRISM ) )
// Include only what is needed (Well, the conditions are not sharp). // Include only what is needed (Well, the conditions are not sharp).
#if OGS_ENABLED_ELEMENTS != 0 #if OGS_ENABLED_ELEMENTS != 0
...@@ -115,18 +114,19 @@ namespace LIE ...@@ -115,18 +114,19 @@ namespace LIE
{ {
namespace SmallDeformation namespace SmallDeformation
{ {
/// The LocalDataInitializer is a functor creating a local assembler data with /// The LocalDataInitializer is a functor creating a local assembler data with
/// corresponding to the mesh element type shape functions and calling /// corresponding to the mesh element type shape functions and calling
/// initialization of the new local assembler data. /// initialization of the new local assembler data.
/// For example for MeshLib::Quad a local assembler data with template argument /// For example for MeshLib::Quad a local assembler data with template argument
/// NumLib::ShapeQuad4 is created. /// NumLib::ShapeQuad4 is created.
template < template <typename LocalAssemblerInterface,
typename LocalAssemblerInterface, template <typename, typename, unsigned, int>
template <typename, typename, unsigned, int> class LocalAssemblerDataMatrix, class LocalAssemblerDataMatrix,
template <typename, typename, unsigned, int> class LocalAssemblerDataMatrixNearFracture, template <typename, typename, unsigned, int>
template <typename, typename, unsigned, int> class LocalAssemblerDataFracture, class LocalAssemblerDataMatrixNearFracture,
unsigned GlobalDim, int DisplacementDim, typename... ConstructorArgs> template <typename, typename, unsigned, int>
class LocalAssemblerDataFracture,
unsigned GlobalDim, int DisplacementDim, typename... ConstructorArgs>
class LocalDataInitializer final class LocalDataInitializer final
{ {
public: public:
...@@ -136,12 +136,12 @@ public: ...@@ -136,12 +136,12 @@ public:
NumLib::LocalToGlobalIndexMap const& dof_table) NumLib::LocalToGlobalIndexMap const& dof_table)
: _dof_table(dof_table) : _dof_table(dof_table)
{ {
//REMARKS: At the moment, only a 2D mesh with 1D elements are supported. // REMARKS: At the moment, only a 2D mesh with 1D elements are supported.
// /// Quads and Hexahedra /////////////////////////////////// // /// Quads and Hexahedra ///////////////////////////////////
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_QUAD) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_QUAD) != 0 && \
&& OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 1 OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 1
_builder[std::type_index(typeid(MeshLib::Quad))] = _builder[std::type_index(typeid(MeshLib::Quad))] =
makeLocalAssemblerBuilder<NumLib::ShapeQuad4>(); makeLocalAssemblerBuilder<NumLib::ShapeQuad4>();
#endif #endif
...@@ -154,8 +154,8 @@ public: ...@@ -154,8 +154,8 @@ public:
#endif #endif
*/ */
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_QUAD) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_QUAD) != 0 && \
&& OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 2 OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 2
_builder[std::type_index(typeid(MeshLib::Quad8))] = _builder[std::type_index(typeid(MeshLib::Quad8))] =
makeLocalAssemblerBuilder<NumLib::ShapeQuad8>(); makeLocalAssemblerBuilder<NumLib::ShapeQuad8>();
_builder[std::type_index(typeid(MeshLib::Quad9))] = _builder[std::type_index(typeid(MeshLib::Quad9))] =
...@@ -170,11 +170,10 @@ public: ...@@ -170,11 +170,10 @@ public:
#endif #endif
*/ */
// /// Simplices ////////////////////////////////////////////////
// /// Simplices //////////////////////////////////////////////// #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_TRI) != 0 && \
OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 1
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_TRI) != 0 \
&& OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 1
_builder[std::type_index(typeid(MeshLib::Tri))] = _builder[std::type_index(typeid(MeshLib::Tri))] =
makeLocalAssemblerBuilder<NumLib::ShapeTri3>(); makeLocalAssemblerBuilder<NumLib::ShapeTri3>();
#endif #endif
...@@ -187,8 +186,8 @@ public: ...@@ -187,8 +186,8 @@ public:
#endif #endif
*/ */
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_TRI) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_TRI) != 0 && \
&& OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 2 OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 2
_builder[std::type_index(typeid(MeshLib::Tri6))] = _builder[std::type_index(typeid(MeshLib::Tri6))] =
makeLocalAssemblerBuilder<NumLib::ShapeTri6>(); makeLocalAssemblerBuilder<NumLib::ShapeTri6>();
#endif #endif
...@@ -201,8 +200,7 @@ public: ...@@ -201,8 +200,7 @@ public:
#endif #endif
*/ */
// /// Prisms ////////////////////////////////////////////////////
// /// Prisms ////////////////////////////////////////////////////
/* /*
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PRISM) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PRISM) != 0 \
...@@ -218,7 +216,7 @@ public: ...@@ -218,7 +216,7 @@ public:
#endif #endif
*/ */
// /// Pyramids ////////////////////////////////////////////////// // /// Pyramids //////////////////////////////////////////////////
/* /*
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PYRAMID) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PYRAMID) != 0 \
...@@ -233,7 +231,7 @@ public: ...@@ -233,7 +231,7 @@ public:
makeLocalAssemblerBuilder<NumLib::ShapePyra13>(); makeLocalAssemblerBuilder<NumLib::ShapePyra13>();
#endif #endif
*/ */
// /// Lines /////////////////////////////////// // /// Lines ///////////////////////////////////
#if OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 1 #if OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 1
_builder[std::type_index(typeid(MeshLib::Line))] = _builder[std::type_index(typeid(MeshLib::Line))] =
...@@ -251,133 +249,165 @@ public: ...@@ -251,133 +249,165 @@ public:
/// \attention /// \attention
/// The index \c id is not necessarily the mesh item's id. Especially when /// The index \c id is not necessarily the mesh item's id. Especially when
/// having multiple meshes it will differ from the latter. /// having multiple meshes it will differ from the latter.
void operator()( void operator()(std::size_t const id,
std::size_t const id, MeshLib::Element const& mesh_item,
MeshLib::Element const& mesh_item, LADataIntfPtr& data_ptr,
LADataIntfPtr& data_ptr, ConstructorArgs&&... args) const
ConstructorArgs&&... args) const
{ {
auto const type_idx = std::type_index(typeid(mesh_item)); auto const type_idx = std::type_index(typeid(mesh_item));
auto const it = _builder.find(type_idx); auto const it = _builder.find(type_idx);
if (it != _builder.end()) { if (it != _builder.end())
{
auto const n_local_dof = _dof_table.getNumberOfElementDOF(id); auto const n_local_dof = _dof_table.getNumberOfElementDOF(id);
auto const n_global_components = _dof_table.getNumberOfElementComponents(id); auto const n_global_components =
const std::vector<std::size_t> varIDs(_dof_table.getElementVariableIDs(id)); _dof_table.getNumberOfElementComponents(id);
const std::vector<std::size_t> varIDs(
_dof_table.getElementVariableIDs(id));
std::vector<unsigned> dofIndex_to_localIndex; std::vector<unsigned> dofIndex_to_localIndex;
if (mesh_item.getDimension() < GlobalDim if (mesh_item.getDimension() < GlobalDim ||
|| n_global_components > GlobalDim) n_global_components > GlobalDim)
{ {
dofIndex_to_localIndex.resize(n_local_dof); dofIndex_to_localIndex.resize(n_local_dof);
unsigned dof_id = 0; unsigned dof_id = 0;
unsigned local_id = 0; unsigned local_id = 0;
for (auto i : varIDs) for (auto i : varIDs)
{ {
for (unsigned j=0; j<_dof_table.getNumberOfVariableComponents(i); j++) for (unsigned j = 0;
j < _dof_table.getNumberOfVariableComponents(i); j++)
{ {
auto& mss = _dof_table.getMeshSubsets(i, j); auto& mss = _dof_table.getMeshSubsets(i, j);
assert(mss.size() == 1); assert(mss.size() == 1);
auto mesh_id = mss.getMeshSubset(0).getMeshID(); auto mesh_id = mss.getMeshSubset(0).getMeshID();
for (unsigned k=0; k<mesh_item.getNumberOfNodes(); k++) for (unsigned k = 0; k < mesh_item.getNumberOfNodes();
k++)
{ {
MeshLib::Location l(mesh_id, MeshLib::Location l(mesh_id,
MeshLib::MeshItemType::Node, MeshLib::MeshItemType::Node,
mesh_item.getNodeIndex(k)); mesh_item.getNodeIndex(k));
auto global_index = _dof_table.getGlobalIndex(l, i, j); auto global_index =
_dof_table.getGlobalIndex(l, i, j);
if (global_index != NumLib::MeshComponentMap::nop) if (global_index != NumLib::MeshComponentMap::nop)
dofIndex_to_localIndex[dof_id++] = local_id; dofIndex_to_localIndex[dof_id++] = local_id;
local_id++; local_id++;
} }
} }
} }
} }
data_ptr = it->second(mesh_item, varIDs.size(), n_local_dof,
data_ptr = it->second( dofIndex_to_localIndex,
mesh_item, varIDs.size(), n_local_dof, dofIndex_to_localIndex, std::forward<ConstructorArgs>(args)...);
std::forward<ConstructorArgs>(args)...); }
} else { else
OGS_FATAL("You are trying to build a local assembler for an unknown mesh element type (%s)." {
" Maybe you have disabled this mesh element type in your build configuration.", OGS_FATAL(
"You are trying to build a local assembler for an unknown mesh "
"element type (%s)."
" Maybe you have disabled this mesh element type in your build "
"configuration.",
type_idx.name()); type_idx.name());
} }
} }
private: private:
using LADataBuilder = std::function<LADataIntfPtr( using LADataBuilder = std::function<LADataIntfPtr(
MeshLib::Element const& e, MeshLib::Element const& e,
std::size_t const n_variables, std::size_t const n_variables,
std::size_t const local_matrix_size, std::size_t const local_matrix_size,
std::vector<unsigned> const& dofIndex_to_localIndex, std::vector<unsigned> const& dofIndex_to_localIndex,
ConstructorArgs&&... ConstructorArgs&&...)>;
)>;
template <typename ShapeFunction> template <typename ShapeFunction>
using IntegrationMethod = typename NumLib::GaussIntegrationPolicy< using IntegrationMethod = typename NumLib::GaussIntegrationPolicy<
typename ShapeFunction::MeshElement>::IntegrationMethod; typename ShapeFunction::MeshElement>::IntegrationMethod;
template <typename ShapeFunction> template <typename ShapeFunction>
using LADataMatrix = using LADataMatrix =
LocalAssemblerDataMatrix<ShapeFunction, IntegrationMethod<ShapeFunction>, LocalAssemblerDataMatrix<ShapeFunction,
GlobalDim, DisplacementDim>; IntegrationMethod<ShapeFunction>, GlobalDim,
DisplacementDim>;
/// A helper forwarding to the correct version of makeLocalAssemblerBuilder
/// depending whether the global dimension is less than the shape function's
/// dimension or not.
template <typename ShapeFunction>
static LADataBuilder makeLocalAssemblerBuilder()
{
return makeLocalAssemblerBuilder<ShapeFunction>(
static_cast<std::integral_constant<
bool, (GlobalDim >= ShapeFunction::DIM)>*>(nullptr));
}
/// Mapping of element types to local assembler constructors.
std::unordered_map<std::type_index, LADataBuilder> _builder;
NumLib::LocalToGlobalIndexMap const& _dof_table;
template <typename ShapeFunction> template <typename ShapeFunction>
using LADataMatrixNearFracture = using LADataMatrixNearFracture =
LocalAssemblerDataMatrixNearFracture<ShapeFunction, IntegrationMethod<ShapeFunction>, LocalAssemblerDataMatrixNearFracture<ShapeFunction,
GlobalDim, DisplacementDim>; IntegrationMethod<ShapeFunction>,
GlobalDim, DisplacementDim>;
template <typename ShapeFunction> template <typename ShapeFunction>
using LAFractureData = using LAFractureData =
LocalAssemblerDataFracture<ShapeFunction, IntegrationMethod<ShapeFunction>, LocalAssemblerDataFracture<ShapeFunction,
GlobalDim, DisplacementDim>; IntegrationMethod<ShapeFunction>, GlobalDim,
DisplacementDim>;
/// Generates a function that creates a new LocalAssembler of type LAData<SHAPE_FCT> // local assembler builder implementations.
template<typename ShapeFct> private:
static LADataBuilder makeLocalAssemblerBuilder() /// Generates a function that creates a new LocalAssembler of type
/// LAData<ShapeFunction>. Only functions with shape function's dimension
/// less or equal to the global dimension are instantiated, e.g. following
/// combinations of shape functions and global dimensions: (Line2, 1),
/// (Line2, 2), (Line2, 3), (Hex20, 3) but not (Hex20, 2) or (Hex20, 1).
template <typename ShapeFunction>
static LADataBuilder makeLocalAssemblerBuilder(std::true_type*)
{ {
return [](MeshLib::Element const& e, return [](MeshLib::Element const& e,
std::size_t const n_variables, std::size_t const n_variables,
std::size_t const local_matrix_size, std::size_t const local_matrix_size,
std::vector<unsigned> const& dofIndex_to_localIndex, std::vector<unsigned> const& dofIndex_to_localIndex,
ConstructorArgs&&... args) ConstructorArgs&&... args) {
{
if (e.getDimension() == GlobalDim) if (e.getDimension() == GlobalDim)
{ {
if (dofIndex_to_localIndex.empty()) if (dofIndex_to_localIndex.empty())
{
return LADataIntfPtr{new LADataMatrix<ShapeFunction>{
e, local_matrix_size,
std::forward<ConstructorArgs>(args)...}};
}
else
{ {
return LADataIntfPtr{ return LADataIntfPtr{
new LADataMatrix<ShapeFct>{ new LADataMatrixNearFracture<ShapeFunction>{
e, local_matrix_size, e, n_variables, local_matrix_size,
std::forward<ConstructorArgs>(args)... dofIndex_to_localIndex,
}}; std::forward<ConstructorArgs>(args)...}};
} else {
return LADataIntfPtr{
new LADataMatrixNearFracture<ShapeFct>{
e, n_variables, local_matrix_size, dofIndex_to_localIndex,
std::forward<ConstructorArgs>(args)...
}};
} }
} }
return LADataIntfPtr{ return LADataIntfPtr{new LAFractureData<ShapeFunction>{
new LAFractureData<ShapeFct>{ e, local_matrix_size, dofIndex_to_localIndex,
e, local_matrix_size, dofIndex_to_localIndex, std::forward<ConstructorArgs>(args)...}};
std::forward<ConstructorArgs>(args)...
}};
}; };
} }
/// Mapping of element types to local assembler constructors. /// Returns nullptr for shape functions whose dimensions are less than the
std::unordered_map<std::type_index, LADataBuilder> _builder; /// global dimension.
template <typename ShapeFunction>
NumLib::LocalToGlobalIndexMap const& _dof_table; static LADataBuilder makeLocalAssemblerBuilder(std::false_type*)
{
return nullptr;
}
}; };
} // namespace SmallDeformation } // namespace SmallDeformationWithLIE
} // namespace LIE } // namespace ProcessLib
} // namespace ProcessLib } // namespace ProcessLib
#undef ENABLED_ELEMENT_TYPE_SIMPLEX #undef ENABLED_ELEMENT_TYPE_SIMPLEX
#undef ENABLED_ELEMENT_TYPE_CUBOID #undef ENABLED_ELEMENT_TYPE_CUBOID
......
ProcessLib/SmallDeformation/LocalDataInitializer.h
+ 98
− 79
View file @ f4e548e1
...@@ -14,13 +14,13 @@ ...@@ -14,13 +14,13 @@
#include <memory> #include <memory>
#include <typeindex> #include <typeindex>
#include <typeinfo> #include <typeinfo>
#include <type_traits>
#include <unordered_map> #include <unordered_map>
#include "MeshLib/Elements/Elements.h" #include "MeshLib/Elements/Elements.h"
#include "NumLib/DOF/LocalToGlobalIndexMap.h" #include "NumLib/DOF/LocalToGlobalIndexMap.h"
#include "NumLib/Fem/Integration/GaussIntegrationPolicy.h" #include "NumLib/Fem/Integration/GaussIntegrationPolicy.h"
#ifndef OGS_MAX_ELEMENT_DIM #ifndef OGS_MAX_ELEMENT_DIM
static_assert(false, "The macro OGS_MAX_ELEMENT_DIM is undefined."); static_assert(false, "The macro OGS_MAX_ELEMENT_DIM is undefined.");
#endif #endif
...@@ -47,7 +47,7 @@ static_assert(false, "The macro OGS_MAX_ELEMENT_ORDER is undefined."); ...@@ -47,7 +47,7 @@ static_assert(false, "The macro OGS_MAX_ELEMENT_ORDER is undefined.");
#ifdef OGS_ENABLE_ELEMENT_PRISM #ifdef OGS_ENABLE_ELEMENT_PRISM
#define ENABLED_ELEMENT_TYPE_PRISM 1u << 2 #define ENABLED_ELEMENT_TYPE_PRISM 1u << 2
#else #else
#define ENABLED_ELEMENT_TYPE_PRISM 0u #define ENABLED_ELEMENT_TYPE_PRISM 0u
#endif #endif
#ifdef OGS_ENABLE_ELEMENT_PYRAMID #ifdef OGS_ENABLE_ELEMENT_PYRAMID
...@@ -58,19 +58,18 @@ static_assert(false, "The macro OGS_MAX_ELEMENT_ORDER is undefined."); ...@@ -58,19 +58,18 @@ static_assert(false, "The macro OGS_MAX_ELEMENT_ORDER is undefined.");
// Dependent element types. // Dependent element types.
// Faces of tets, pyramids and prisms are triangles // Faces of tets, pyramids and prisms are triangles
#define ENABLED_ELEMENT_TYPE_TRI ((ENABLED_ELEMENT_TYPE_SIMPLEX) | (ENABLED_ELEMENT_TYPE_PYRAMID) \ #define ENABLED_ELEMENT_TYPE_TRI \
|(ENABLED_ELEMENT_TYPE_PRISM)) ((ENABLED_ELEMENT_TYPE_SIMPLEX) | (ENABLED_ELEMENT_TYPE_PYRAMID) | \
(ENABLED_ELEMENT_TYPE_PRISM))
// Faces of hexes, pyramids and prisms are quads // Faces of hexes, pyramids and prisms are quads
#define ENABLED_ELEMENT_TYPE_QUAD ((ENABLED_ELEMENT_TYPE_CUBOID) | (ENABLED_ELEMENT_TYPE_PYRAMID) \ #define ENABLED_ELEMENT_TYPE_QUAD \
|(ENABLED_ELEMENT_TYPE_PRISM)) ((ENABLED_ELEMENT_TYPE_CUBOID) | (ENABLED_ELEMENT_TYPE_PYRAMID) | \
(ENABLED_ELEMENT_TYPE_PRISM))
// All enabled element types // All enabled element types
#define OGS_ENABLED_ELEMENTS \ #define OGS_ENABLED_ELEMENTS \
( (ENABLED_ELEMENT_TYPE_SIMPLEX) \ ((ENABLED_ELEMENT_TYPE_SIMPLEX) | (ENABLED_ELEMENT_TYPE_CUBOID) | \
| (ENABLED_ELEMENT_TYPE_CUBOID ) \ (ENABLED_ELEMENT_TYPE_PYRAMID) | (ENABLED_ELEMENT_TYPE_PRISM))
| (ENABLED_ELEMENT_TYPE_PYRAMID) \
| (ENABLED_ELEMENT_TYPE_PRISM ) )
// Include only what is needed (Well, the conditions are not sharp). // Include only what is needed (Well, the conditions are not sharp).
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_SIMPLEX) != 0 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_SIMPLEX) != 0
...@@ -106,16 +105,14 @@ static_assert(false, "The macro OGS_MAX_ELEMENT_ORDER is undefined."); ...@@ -106,16 +105,14 @@ static_assert(false, "The macro OGS_MAX_ELEMENT_ORDER is undefined.");
namespace ProcessLib namespace ProcessLib
{ {
/// The LocalDataInitializer is a functor creating a local assembler data with /// The LocalDataInitializer is a functor creating a local assembler data with
/// corresponding to the mesh element type shape functions and calling /// corresponding to the mesh element type shape functions and calling
/// initialization of the new local assembler data. /// initialization of the new local assembler data.
/// For example for MeshLib::Quad a local assembler data with template argument /// For example for MeshLib::Quad a local assembler data with template argument
/// NumLib::ShapeQuad4 is created. /// NumLib::ShapeQuad4 is created.
template < template <typename LocalAssemblerInterface,
typename LocalAssemblerInterface, template <typename, typename, unsigned, int> class LocalAssemblerData,
template <typename, typename, unsigned, int> class LocalAssemblerData, unsigned GlobalDim, int DisplacementDim, typename... ConstructorArgs>
unsigned GlobalDim, int DisplacementDim, typename... ConstructorArgs>
class LocalDataInitializer final class LocalDataInitializer final
{ {
public: public:
...@@ -125,86 +122,84 @@ public: ...@@ -125,86 +122,84 @@ public:
NumLib::LocalToGlobalIndexMap const& dof_table) NumLib::LocalToGlobalIndexMap const& dof_table)
: _dof_table(dof_table) : _dof_table(dof_table)
{ {
// /// Quads and Hexahedra /////////////////////////////////// // /// Quads and Hexahedra ///////////////////////////////////
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_QUAD) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_QUAD) != 0 && \
&& OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 1 OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 1
_builder[std::type_index(typeid(MeshLib::Quad))] = _builder[std::type_index(typeid(MeshLib::Quad))] =
makeLocalAssemblerBuilder<NumLib::ShapeQuad4>(); makeLocalAssemblerBuilder<NumLib::ShapeQuad4>();
#endif #endif
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_CUBOID) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_CUBOID) != 0 && \
&& OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 1 OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 1
_builder[std::type_index(typeid(MeshLib::Hex))] = _builder[std::type_index(typeid(MeshLib::Hex))] =
makeLocalAssemblerBuilder<NumLib::ShapeHex8>(); makeLocalAssemblerBuilder<NumLib::ShapeHex8>();
#endif #endif
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_QUAD) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_QUAD) != 0 && \
&& OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 2 OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 2
_builder[std::type_index(typeid(MeshLib::Quad8))] = _builder[std::type_index(typeid(MeshLib::Quad8))] =
makeLocalAssemblerBuilder<NumLib::ShapeQuad8>(); makeLocalAssemblerBuilder<NumLib::ShapeQuad8>();
_builder[std::type_index(typeid(MeshLib::Quad9))] = _builder[std::type_index(typeid(MeshLib::Quad9))] =
makeLocalAssemblerBuilder<NumLib::ShapeQuad9>(); makeLocalAssemblerBuilder<NumLib::ShapeQuad9>();
#endif #endif
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_CUBOID) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_CUBOID) != 0 && \
&& OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 2 OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 2
_builder[std::type_index(typeid(MeshLib::Hex20))] = _builder[std::type_index(typeid(MeshLib::Hex20))] =
makeLocalAssemblerBuilder<NumLib::ShapeHex20>(); makeLocalAssemblerBuilder<NumLib::ShapeHex20>();
#endif #endif
// /// Simplices ////////////////////////////////////////////////
// /// Simplices //////////////////////////////////////////////// #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_TRI) != 0 && \
OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 1
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_TRI) != 0 \
&& OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 1
_builder[std::type_index(typeid(MeshLib::Tri))] = _builder[std::type_index(typeid(MeshLib::Tri))] =
makeLocalAssemblerBuilder<NumLib::ShapeTri3>(); makeLocalAssemblerBuilder<NumLib::ShapeTri3>();
#endif #endif
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_SIMPLEX) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_SIMPLEX) != 0 && \
&& OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 1 OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 1
_builder[std::type_index(typeid(MeshLib::Tet))] = _builder[std::type_index(typeid(MeshLib::Tet))] =
makeLocalAssemblerBuilder<NumLib::ShapeTet4>(); makeLocalAssemblerBuilder<NumLib::ShapeTet4>();
#endif #endif
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_TRI) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_TRI) != 0 && \
&& OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 2 OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 2
_builder[std::type_index(typeid(MeshLib::Tri6))] = _builder[std::type_index(typeid(MeshLib::Tri6))] =
makeLocalAssemblerBuilder<NumLib::ShapeTri6>(); makeLocalAssemblerBuilder<NumLib::ShapeTri6>();
#endif #endif
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_SIMPLEX) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_SIMPLEX) != 0 && \
&& OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 2 OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 2
_builder[std::type_index(typeid(MeshLib::Tet10))] = _builder[std::type_index(typeid(MeshLib::Tet10))] =
makeLocalAssemblerBuilder<NumLib::ShapeTet10>(); makeLocalAssemblerBuilder<NumLib::ShapeTet10>();
#endif #endif
// /// Prisms ////////////////////////////////////////////////////
// /// Prisms //////////////////////////////////////////////////// #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PRISM) != 0 && \
OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 1
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PRISM) != 0 \
&& OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 1
_builder[std::type_index(typeid(MeshLib::Prism))] = _builder[std::type_index(typeid(MeshLib::Prism))] =
makeLocalAssemblerBuilder<NumLib::ShapePrism6>(); makeLocalAssemblerBuilder<NumLib::ShapePrism6>();
#endif #endif
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PRISM) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PRISM) != 0 && \
&& OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 2 OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 2
_builder[std::type_index(typeid(MeshLib::Prism15))] = _builder[std::type_index(typeid(MeshLib::Prism15))] =
makeLocalAssemblerBuilder<NumLib::ShapePrism15>(); makeLocalAssemblerBuilder<NumLib::ShapePrism15>();
#endif #endif
// /// Pyramids ////////////////////////////////////////////////// // /// Pyramids //////////////////////////////////////////////////
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PYRAMID) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PYRAMID) != 0 && \
&& OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 1 OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 1
_builder[std::type_index(typeid(MeshLib::Pyramid))] = _builder[std::type_index(typeid(MeshLib::Pyramid))] =
makeLocalAssemblerBuilder<NumLib::ShapePyra5>(); makeLocalAssemblerBuilder<NumLib::ShapePyra5>();
#endif #endif
#if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PYRAMID) != 0 \ #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PYRAMID) != 0 && \
&& OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 2 OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 2
_builder[std::type_index(typeid(MeshLib::Pyramid13))] = _builder[std::type_index(typeid(MeshLib::Pyramid13))] =
makeLocalAssemblerBuilder<NumLib::ShapePyra13>(); makeLocalAssemblerBuilder<NumLib::ShapePyra13>();
#endif #endif
...@@ -215,67 +210,91 @@ public: ...@@ -215,67 +210,91 @@ public:
/// \attention /// \attention
/// The index \c id is not necessarily the mesh item's id. Especially when /// The index \c id is not necessarily the mesh item's id. Especially when
/// having multiple meshes it will differ from the latter. /// having multiple meshes it will differ from the latter.
void operator()( void operator()(std::size_t const id,
std::size_t const id, MeshLib::Element const& mesh_item,
MeshLib::Element const& mesh_item, LADataIntfPtr& data_ptr,
LADataIntfPtr& data_ptr, ConstructorArgs&&... args) const
ConstructorArgs&&... args) const
{ {
auto const type_idx = std::type_index(typeid(mesh_item)); auto const type_idx = std::type_index(typeid(mesh_item));
auto const it = _builder.find(type_idx); auto const it = _builder.find(type_idx);
if (it != _builder.end()) { if (it != _builder.end())
{
auto const num_local_dof = _dof_table.getNumberOfElementDOF(id); auto const num_local_dof = _dof_table.getNumberOfElementDOF(id);
data_ptr = it->second( data_ptr = it->second(mesh_item, num_local_dof,
mesh_item, num_local_dof, std::forward<ConstructorArgs>(args)...);
std::forward<ConstructorArgs>(args)...); }
} else { else
OGS_FATAL("You are trying to build a local assembler for an unknown mesh element type (%s)." {
" Maybe you have disabled this mesh element type in your build configuration.", OGS_FATAL(
"You are trying to build a local assembler for an unknown mesh "
"element type (%s)."
" Maybe you have disabled this mesh element type in your build "
"configuration.",
type_idx.name()); type_idx.name());
} }
} }
private: private:
using LADataBuilder = std::function<LADataIntfPtr( using LADataBuilder =
MeshLib::Element const& e, std::function<LADataIntfPtr(MeshLib::Element const& e,
std::size_t const local_matrix_size, std::size_t const local_matrix_size,
ConstructorArgs&&... ConstructorArgs&&...)>;
)>;
template <typename ShapeFunction> template <typename ShapeFunction>
using IntegrationMethod = typename NumLib::GaussIntegrationPolicy< using IntegrationMethod = typename NumLib::GaussIntegrationPolicy<
typename ShapeFunction::MeshElement>::IntegrationMethod; typename ShapeFunction::MeshElement>::IntegrationMethod;
template <typename ShapeFunction> template <typename ShapeFunction>
using LAData = using LAData =
LocalAssemblerData<ShapeFunction, IntegrationMethod<ShapeFunction>, LocalAssemblerData<ShapeFunction, IntegrationMethod<ShapeFunction>,
GlobalDim, DisplacementDim>; GlobalDim, DisplacementDim>;
/// Generates a function that creates a new LocalAssembler of type LAData<SHAPE_FCT> /// A helper forwarding to the correct version of makeLocalAssemblerBuilder
template<typename ShapeFct> /// depending whether the global dimension is less than the shape function's
/// dimension or not.
template <typename ShapeFunction>
static LADataBuilder makeLocalAssemblerBuilder() static LADataBuilder makeLocalAssemblerBuilder()
{ {
return [](MeshLib::Element const& e, return makeLocalAssemblerBuilder<ShapeFunction>(
std::size_t const local_matrix_size, static_cast<std::integral_constant<
ConstructorArgs&&... args) bool, (GlobalDim >= ShapeFunction::DIM)>*>(nullptr));
{
return LADataIntfPtr{
new LAData<ShapeFct>{
e, local_matrix_size,
std::forward<ConstructorArgs>(args)...
}};
};
} }
/// Mapping of element types to local assembler constructors. /// Mapping of element types to local assembler constructors.
std::unordered_map<std::type_index, LADataBuilder> _builder; std::unordered_map<std::type_index, LADataBuilder> _builder;
NumLib::LocalToGlobalIndexMap const& _dof_table; NumLib::LocalToGlobalIndexMap const& _dof_table;
};
} // namespace ProcessLib // local assembler builder implementations.
private:
/// Generates a function that creates a new LocalAssembler of type
/// LAData<ShapeFunction>. Only functions with shape function's dimension
/// less or equal to the global dimension are instantiated, e.g. following
/// combinations of shape functions and global dimensions: (Line2, 1),
/// (Line2, 2), (Line2, 3), (Hex20, 3) but not (Hex20, 2) or (Hex20, 1).
template <typename ShapeFunction>
static LADataBuilder makeLocalAssemblerBuilder(std::true_type*)
{
return [](MeshLib::Element const& e,
std::size_t const local_matrix_size,
ConstructorArgs&&... args) {
return LADataIntfPtr{new LAData<ShapeFunction>{
e, local_matrix_size, std::forward<ConstructorArgs>(args)...}};
};
}
/// Returns nullptr for shape functions whose dimensions are less than the
/// global dimension.
template <typename ShapeFunction>
static LADataBuilder makeLocalAssemblerBuilder(std::false_type*)
{
return nullptr;
}
};
} // namespace ProcessLib
#undef ENABLED_ELEMENT_TYPE_SIMPLEX #undef ENABLED_ELEMENT_TYPE_SIMPLEX
#undef ENABLED_ELEMENT_TYPE_CUBOID #undef ENABLED_ELEMENT_TYPE_CUBOID
......
ProcessLib/Utils/LocalDataInitializer.h
+ 171
− 148
View file @ f4e548e1
This diff is collapsed.
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment