OGS 6.2.1-97-g73d1aeda3
LocalDataInitializer.h
Go to the documentation of this file.
1 
10 #pragma once
11 
12 #include <functional>
13 #include <memory>
14 #include <type_traits>
15 #include <typeindex>
16 #include <typeinfo>
17 #include <unordered_map>
18 
22 
23 #ifndef OGS_MAX_ELEMENT_DIM
24 static_assert(false, "The macro OGS_MAX_ELEMENT_DIM is undefined.");
25 #endif
26 
27 #ifndef OGS_MAX_ELEMENT_ORDER
28 static_assert(false, "The macro OGS_MAX_ELEMENT_ORDER is undefined.");
29 #endif
30 
31 // The following macros decide which element types will be compiled, i.e.
32 // which element types will be available for use in simulations.
33 
34 #ifdef OGS_ENABLE_ELEMENT_SIMPLEX
35 #define ENABLED_ELEMENT_TYPE_SIMPLEX 1u
36 #else
37 #define ENABLED_ELEMENT_TYPE_SIMPLEX 0u
38 #endif
39 
40 #ifdef OGS_ENABLE_ELEMENT_CUBOID
41 #define ENABLED_ELEMENT_TYPE_CUBOID 1u << 1
42 #else
43 #define ENABLED_ELEMENT_TYPE_CUBOID 0u
44 #endif
45 
46 #ifdef OGS_ENABLE_ELEMENT_PRISM
47 #define ENABLED_ELEMENT_TYPE_PRISM 1u << 2
48 #else
49 #define ENABLED_ELEMENT_TYPE_PRISM 0u
50 #endif
51 
52 #ifdef OGS_ENABLE_ELEMENT_PYRAMID
53 #define ENABLED_ELEMENT_TYPE_PYRAMID 1u << 3
54 #else
55 #define ENABLED_ELEMENT_TYPE_PYRAMID 0u
56 #endif
57 
58 // Dependent element types.
59 // Faces of tets, pyramids and prisms are triangles
60 #define ENABLED_ELEMENT_TYPE_TRI \
61  ((ENABLED_ELEMENT_TYPE_SIMPLEX) | (ENABLED_ELEMENT_TYPE_PYRAMID) | \
62  (ENABLED_ELEMENT_TYPE_PRISM))
63 // Faces of hexes, pyramids and prisms are quads
64 #define ENABLED_ELEMENT_TYPE_QUAD \
65  ((ENABLED_ELEMENT_TYPE_CUBOID) | (ENABLED_ELEMENT_TYPE_PYRAMID) | \
66  (ENABLED_ELEMENT_TYPE_PRISM))
67 
68 // All enabled element types
69 #define OGS_ENABLED_ELEMENTS \
70  ((ENABLED_ELEMENT_TYPE_SIMPLEX) | (ENABLED_ELEMENT_TYPE_CUBOID) | \
71  (ENABLED_ELEMENT_TYPE_PYRAMID) | (ENABLED_ELEMENT_TYPE_PRISM))
72 
73 // Include only what is needed (Well, the conditions are not sharp).
74 #if OGS_ENABLED_ELEMENTS != 0
77 #endif
78 
79 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_SIMPLEX) != 0
82 #endif
83 
84 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_TRI) != 0
87 #endif
88 
89 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_CUBOID) != 0
92 #endif
93 
94 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_QUAD) != 0
98 #endif
99 
100 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PRISM) != 0
103 #endif
104 
105 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PYRAMID) != 0
108 #endif
109 
110 namespace ProcessLib
111 {
112 namespace LIE
113 {
114 namespace SmallDeformation
115 {
121 template <typename LocalAssemblerInterface,
122  template <typename, typename, int> class LocalAssemblerDataMatrix,
123  template <typename, typename, int>
124  class LocalAssemblerDataMatrixNearFracture,
125  template <typename, typename, int> class LocalAssemblerDataFracture,
126  int GlobalDim, typename... ConstructorArgs>
128 {
129 public:
130  using LADataIntfPtr = std::unique_ptr<LocalAssemblerInterface>;
131 
133  NumLib::LocalToGlobalIndexMap const& dof_table)
134  : _dof_table(dof_table)
135  {
136  // REMARKS: At the moment, only a 2D mesh with 1D elements are
137  // supported.
138 
139  // /// Quads and Hexahedra ///////////////////////////////////
140 
141 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_QUAD) != 0 && \
142  OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 1
143  _builder[std::type_index(typeid(MeshLib::Quad))] =
144  makeLocalAssemblerBuilder<NumLib::ShapeQuad4>();
145 #endif
146 
147 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_CUBOID) != 0 && \
148  OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 1
149  _builder[std::type_index(typeid(MeshLib::Hex))] =
150  makeLocalAssemblerBuilder<NumLib::ShapeHex8>();
151 #endif
152 
153 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_QUAD) != 0 && \
154  OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 2
155  _builder[std::type_index(typeid(MeshLib::Quad8))] =
156  makeLocalAssemblerBuilder<NumLib::ShapeQuad8>();
157  _builder[std::type_index(typeid(MeshLib::Quad9))] =
158  makeLocalAssemblerBuilder<NumLib::ShapeQuad9>();
159 #endif
160 
161 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_CUBOID) != 0 && \
162  OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 2
163  _builder[std::type_index(typeid(MeshLib::Hex20))] =
164  makeLocalAssemblerBuilder<NumLib::ShapeHex20>();
165 #endif
166 
167  // /// Simplices ////////////////////////////////////////////////
168 
169 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_TRI) != 0 && \
170  OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 1
171  _builder[std::type_index(typeid(MeshLib::Tri))] =
172  makeLocalAssemblerBuilder<NumLib::ShapeTri3>();
173 #endif
174 
175 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_SIMPLEX) != 0 && \
176  OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 1
177  _builder[std::type_index(typeid(MeshLib::Tet))] =
178  makeLocalAssemblerBuilder<NumLib::ShapeTet4>();
179 #endif
180 
181 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_TRI) != 0 && \
182  OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 2
183  _builder[std::type_index(typeid(MeshLib::Tri6))] =
184  makeLocalAssemblerBuilder<NumLib::ShapeTri6>();
185 #endif
186 
187 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_SIMPLEX) != 0 && \
188  OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 2
189  _builder[std::type_index(typeid(MeshLib::Tet10))] =
190  makeLocalAssemblerBuilder<NumLib::ShapeTet10>();
191 #endif
192 
193  // /// Prisms ////////////////////////////////////////////////////
194 
195 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PRISM) != 0 && \
196  OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 1
197  _builder[std::type_index(typeid(MeshLib::Prism))] =
198  makeLocalAssemblerBuilder<NumLib::ShapePrism6>();
199 #endif
200 
201 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PRISM) != 0 && \
202  OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 2
203  _builder[std::type_index(typeid(MeshLib::Prism15))] =
204  makeLocalAssemblerBuilder<NumLib::ShapePrism15>();
205 #endif
206 
207  // /// Pyramids //////////////////////////////////////////////////
208 
209 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PYRAMID) != 0 && \
210  OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 1
211  _builder[std::type_index(typeid(MeshLib::Pyramid))] =
212  makeLocalAssemblerBuilder<NumLib::ShapePyra5>();
213 #endif
214 
215 #if (OGS_ENABLED_ELEMENTS & ENABLED_ELEMENT_TYPE_PYRAMID) != 0 && \
216  OGS_MAX_ELEMENT_DIM >= 3 && OGS_MAX_ELEMENT_ORDER >= 2
217  _builder[std::type_index(typeid(MeshLib::Pyramid13))] =
218  makeLocalAssemblerBuilder<NumLib::ShapePyra13>();
219 #endif
220  // /// Lines ///////////////////////////////////
221 
222 #if OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 1
223  _builder[std::type_index(typeid(MeshLib::Line))] =
224  makeLocalAssemblerBuilder<NumLib::ShapeLine2>();
225 #endif
226 
227 #if OGS_MAX_ELEMENT_DIM >= 2 && OGS_MAX_ELEMENT_ORDER >= 2
228  _builder[std::type_index(typeid(MeshLib::Line3))] =
229  makeLocalAssemblerBuilder<NumLib::ShapeLine3>();
230 #endif
231  }
232 
238  void operator()(std::size_t const id,
239  MeshLib::Element const& mesh_item,
240  LADataIntfPtr& data_ptr,
241  ConstructorArgs&&... args) const
242  {
243  auto const type_idx = std::type_index(typeid(mesh_item));
244  auto const it = _builder.find(type_idx);
245 
246  if (it == _builder.end())
247  {
248  OGS_FATAL(
249  "You are trying to build a local assembler for an unknown mesh "
250  "element type (%s)."
251  " Maybe you have disabled this mesh element type in your build "
252  "configuration or this process requires higher order elements.",
253  type_idx.name());
254  }
255 
256  auto const n_local_dof = _dof_table.getNumberOfElementDOF(id);
257  auto const n_global_components =
259  auto const varIDs = _dof_table.getElementVariableIDs(id);
260 
261  std::vector<unsigned> dofIndex_to_localIndex;
262  if (mesh_item.getDimension() < GlobalDim ||
263  n_global_components > GlobalDim)
264  {
265  dofIndex_to_localIndex.resize(n_local_dof);
266  unsigned dof_id = 0;
267  unsigned local_id = 0;
268  for (auto i : varIDs)
269  {
270  for (int j = 0; j < _dof_table.getNumberOfVariableComponents(i);
271  j++)
272  {
273  auto const& ms = _dof_table.getMeshSubset(i, j);
274  auto const mesh_id = ms.getMeshID();
275  for (unsigned k = 0; k < mesh_item.getNumberOfNodes(); k++)
276  {
277  MeshLib::Location l(mesh_id,
279  mesh_item.getNodeIndex(k));
280  auto global_index = _dof_table.getGlobalIndex(l, i, j);
281  if (global_index != NumLib::MeshComponentMap::nop)
282  {
283  dofIndex_to_localIndex[dof_id++] = local_id;
284  }
285  local_id++;
286  }
287  }
288  }
289  }
290 
291  data_ptr = it->second(mesh_item, varIDs.size(), n_local_dof,
292  dofIndex_to_localIndex,
293  std::forward<ConstructorArgs>(args)...);
294  }
295 
296 private:
297  using LADataBuilder = std::function<LADataIntfPtr(
298  MeshLib::Element const& e,
299  std::size_t const n_variables,
300  std::size_t const local_matrix_size,
301  std::vector<unsigned> const& dofIndex_to_localIndex,
302  ConstructorArgs&&...)>;
303 
304  template <typename ShapeFunction>
306  typename ShapeFunction::MeshElement>::IntegrationMethod;
307 
308  template <typename ShapeFunction>
309  using LADataMatrix =
310  LocalAssemblerDataMatrix<ShapeFunction,
312 
316  template <typename ShapeFunction>
318  {
319  return makeLocalAssemblerBuilder<ShapeFunction>(
320  static_cast<std::integral_constant<
321  bool, (GlobalDim >= ShapeFunction::DIM)>*>(nullptr));
322  }
323 
325  std::unordered_map<std::type_index, LADataBuilder> _builder;
326 
328 
329  template <typename ShapeFunction>
330  using LADataMatrixNearFracture = LocalAssemblerDataMatrixNearFracture<
331  ShapeFunction, IntegrationMethod<ShapeFunction>, GlobalDim>;
332 
333  template <typename ShapeFunction>
334  using LAFractureData =
335  LocalAssemblerDataFracture<ShapeFunction,
336  IntegrationMethod<ShapeFunction>, GlobalDim>;
337 
338  // local assembler builder implementations.
339 private:
345  template <typename ShapeFunction>
346  static LADataBuilder makeLocalAssemblerBuilder(std::true_type* /*unused*/)
347  {
348  return [](MeshLib::Element const& e,
349  std::size_t const n_variables,
350  std::size_t const local_matrix_size,
351  std::vector<unsigned> const& dofIndex_to_localIndex,
352  ConstructorArgs&&... args) {
353  if (e.getDimension() == GlobalDim)
354  {
355  if (dofIndex_to_localIndex.empty())
356  {
358  e, local_matrix_size,
359  std::forward<ConstructorArgs>(args)...}};
360  }
361 
362  return LADataIntfPtr{
364  e, n_variables, local_matrix_size,
365  dofIndex_to_localIndex,
366  std::forward<ConstructorArgs>(args)...}};
367  }
369  e, n_variables, local_matrix_size, dofIndex_to_localIndex,
370  std::forward<ConstructorArgs>(args)...}};
371  };
372  }
373 
376  template <typename ShapeFunction>
377  static LADataBuilder makeLocalAssemblerBuilder(std::false_type* /*unused*/)
378  {
379  return nullptr;
380  }
381 };
382 
383 } // namespace SmallDeformation
384 } // namespace LIE
385 } // namespace ProcessLib
386 
387 #undef ENABLED_ELEMENT_TYPE_SIMPLEX
388 #undef ENABLED_ELEMENT_TYPE_CUBOID
389 #undef ENABLED_ELEMENT_TYPE_PYRAMID
390 #undef ENABLED_ELEMENT_TYPE_PRISM
391 #undef ENABLED_ELEMENT_TYPE_TRI
392 #undef ENABLED_ELEMENT_TYPE_QUAD
393 #undef OGS_ENABLED_ELEMENTS
static LADataBuilder makeLocalAssemblerBuilder(std::false_type *)
std::vector< int > getElementVariableIDs(std::size_t const mesh_item_id) const
static LADataBuilder makeLocalAssemblerBuilder(std::true_type *)
std::size_t getMeshID() const
return this mesh ID
Definition: MeshSubset.h:71
std::unique_ptr< LocalAssemblerInterface > LADataIntfPtr
std::size_t getNumberOfElementDOF(std::size_t const mesh_item_id) const
int getNumberOfVariableComponents(int variable_id) const
virtual unsigned getDimension() const =0
Get dimension of the mesh element.
std::function< LADataIntfPtr(MeshLib::Element const &e, std::size_t const n_variables, std::size_t const local_matrix_size, std::vector< unsigned > const &dofIndex_to_localIndex, ConstructorArgs &&...)> LADataBuilder
LocalAssemblerDataMatrix< ShapeFunction, IntegrationMethod< ShapeFunction >, GlobalDim > LADataMatrix
LocalDataInitializer(NumLib::LocalToGlobalIndexMap const &dof_table)
std::size_t getNodeIndex(unsigned i) const
Definition: Element.cpp:182
GlobalIndexType getGlobalIndex(MeshLib::Location const &l, int const variable_id, int const component_id) const
void operator()(std::size_t const id, MeshLib::Element const &mesh_item, LADataIntfPtr &data_ptr, ConstructorArgs &&... args) const
MeshLib::MeshSubset const & getMeshSubset(int const variable_id, int const component_id) const
std::size_t getNumberOfElementComponents(std::size_t const mesh_item_id) const
std::unordered_map< std::type_index, LADataBuilder > _builder
Mapping of element types to local assembler constructors.
#define OGS_FATAL(fmt,...)
Definition: Error.h:63
virtual unsigned getNumberOfNodes() const =0
Returns the number of all nodes including both linear and nonlinear nodes.
typename NumLib::GaussLegendreIntegrationPolicy< typename ShapeFunction::MeshElement >::IntegrationMethod IntegrationMethod
LocalAssemblerDataMatrixNearFracture< ShapeFunction, IntegrationMethod< ShapeFunction >, GlobalDim > LADataMatrixNearFracture
static NUMLIB_EXPORT GlobalIndexType const nop
LocalAssemblerDataFracture< ShapeFunction, IntegrationMethod< ShapeFunction >, GlobalDim > LAFractureData