From d37b378cd73f26b486796fa280eadb7d6d6f8d1c Mon Sep 17 00:00:00 2001
From: Christoph Lehmann <christoph.lehmann@ufz.de>
Date: Sun, 21 Nov 2021 13:43:57 +0100
Subject: [PATCH] [PL] Unified local assembler factories for Taylor-Hood
elements
---
.../HydroMechanics/CreateLocalAssemblers.h | 14 +-
ProcessLib/StokesFlow/CreateLocalAssemblers.h | 18 +-
ProcessLib/StokesFlow/LocalDataInitializer.h | 152 ----------------
...calAssemblerFactoryTaylorHoodForOrderGeN.h | 171 ++++++++++++++++++
ProcessLib/Utils/LocalDataInitializerHM.h | 151 ----------------
5 files changed, 185 insertions(+), 321 deletions(-)
delete mode 100644 ProcessLib/StokesFlow/LocalDataInitializer.h
create mode 100644 ProcessLib/Utils/LocalAssemblerFactoryTaylorHoodForOrderGeN.h
delete mode 100644 ProcessLib/Utils/LocalDataInitializerHM.h
diff --git a/ProcessLib/HydroMechanics/CreateLocalAssemblers.h b/ProcessLib/HydroMechanics/CreateLocalAssemblers.h
index ddbb5937fbd..af29f10bdc9 100644
--- a/ProcessLib/HydroMechanics/CreateLocalAssemblers.h
+++ b/ProcessLib/HydroMechanics/CreateLocalAssemblers.h
@@ -13,7 +13,7 @@
#include "BaseLib/Logging.h"
#include "NumLib/DOF/LocalToGlobalIndexMap.h"
-#include "ProcessLib/Utils/LocalDataInitializerHM.h"
+#include "ProcessLib/Utils/LocalAssemblerFactoryTaylorHoodForOrderGeN.h"
namespace ProcessLib
{
@@ -31,21 +31,19 @@ void createLocalAssemblers(
std::vector<std::unique_ptr<LocalAssemblerInterface>>& local_assemblers,
ExtraCtorArgs&&... extra_ctor_args)
{
- // Shape matrices initializer
- using LocalDataInitializer =
- LocalDataInitializerHM<LocalAssemblerInterface,
- LocalAssemblerImplementation, GlobalDim,
- ExtraCtorArgs...>;
+ using LocAsmFactory = LocalAssemblerFactoryHM<LocalAssemblerInterface,
+ LocalAssemblerImplementation,
+ GlobalDim, ExtraCtorArgs...>;
DBUG("Create local assemblers.");
// Populate the vector of local assemblers.
local_assemblers.resize(mesh_elements.size());
- LocalDataInitializer initializer(dof_table);
+ LocAsmFactory factory(dof_table);
DBUG("Calling local assembler builder for all mesh elements.");
GlobalExecutor::transformDereferenced(
- initializer, mesh_elements, local_assemblers,
+ factory, mesh_elements, local_assemblers,
std::forward<ExtraCtorArgs>(extra_ctor_args)...);
}
diff --git a/ProcessLib/StokesFlow/CreateLocalAssemblers.h b/ProcessLib/StokesFlow/CreateLocalAssemblers.h
index 887b298963a..e4261afda20 100644
--- a/ProcessLib/StokesFlow/CreateLocalAssemblers.h
+++ b/ProcessLib/StokesFlow/CreateLocalAssemblers.h
@@ -12,8 +12,8 @@
#include <vector>
#include "BaseLib/Logging.h"
-#include "LocalDataInitializer.h"
#include "NumLib/DOF/LocalToGlobalIndexMap.h"
+#include "ProcessLib/Utils/LocalAssemblerFactoryTaylorHoodForOrderGeN.h"
namespace ProcessLib
{
@@ -31,21 +31,19 @@ void createLocalAssemblers(
std::vector<std::unique_ptr<LocalAssemblerInterface>>& local_assemblers,
ExtraCtorArgs&&... extra_ctor_args)
{
- // Shape matrices initializer
- using LocalDataInitializer =
- LocalDataInitializer<LocalAssemblerInterface,
- LocalAssemblerImplementation, GlobalDim,
- ExtraCtorArgs...>;
+ using LocAsmFactory =
+ LocalAssemblerFactoryStokes<LocalAssemblerInterface,
+ LocalAssemblerImplementation, GlobalDim,
+ ExtraCtorArgs...>;
DBUG("Create local assemblers.");
- // Populate the vector of local assemblers.
- local_assemblers.resize(mesh_elements.size());
- LocalDataInitializer initializer(dof_table);
+ LocAsmFactory factory(dof_table);
+ local_assemblers.resize(mesh_elements.size());
DBUG("Calling local assembler builder for all mesh elements.");
GlobalExecutor::transformDereferenced(
- initializer, mesh_elements, local_assemblers,
+ factory, mesh_elements, local_assemblers,
std::forward<ExtraCtorArgs>(extra_ctor_args)...);
}
diff --git a/ProcessLib/StokesFlow/LocalDataInitializer.h b/ProcessLib/StokesFlow/LocalDataInitializer.h
deleted file mode 100644
index e317b1698a6..00000000000
--- a/ProcessLib/StokesFlow/LocalDataInitializer.h
+++ /dev/null
@@ -1,152 +0,0 @@
-/**
- * \file
- * \copyright
- * Copyright (c) 2012-2021, OpenGeoSys Community (http://www.opengeosys.org)
- * Distributed under a Modified BSD License.
- * See accompanying file LICENSE.txt or
- * http://www.opengeosys.org/project/license
- *
- */
-
-#pragma once
-
-#include <functional>
-#include <memory>
-#include <type_traits>
-#include <typeindex>
-#include <typeinfo>
-#include <unordered_map>
-
-#include "MeshLib/Elements/Elements.h"
-#include "NumLib/DOF/LocalToGlobalIndexMap.h"
-#include "NumLib/Fem/Integration/GaussLegendreIntegrationPolicy.h"
-#include "ProcessLib/Utils/EnabledElements.h"
-
-namespace ProcessLib::StokesFlow
-{
-/// The LocalDataInitializer is a functor creating a local assembler data with
-/// corresponding to the mesh element type shape functions and calling
-/// initialization of the new local assembler data.
-/// For example for MeshLib::Quad a local assembler data with template argument
-/// NumLib::ShapeQuad4 is created.
-///
-/// \attention This is modified version of the ProcessLib::LocalDataInitializer
-/// class which does not include line or point elements. For the shape functions
-/// of order 2 (used for fluid velocity in StokesFlow) a shape function of order
-/// 1 will be used for the scalar variables.
-template <typename LocalAssemblerInterface,
- template <typename, typename, typename, int> class LocalAssemblerData,
- int GlobalDim, typename... ConstructorArgs>
-class LocalDataInitializer final
-{
- struct IsElementEnabled
- {
- template <typename ElementTraits>
- constexpr bool operator()(ElementTraits*) const
- {
- if constexpr (GlobalDim < ElementTraits::ShapeFunction::DIM)
- {
- return false;
- }
-
- return ElementTraits::Element::dimension >= 2 &&
- ElementTraits::ShapeFunction::ORDER == 2;
- }
- };
-
-public:
- using LADataIntfPtr = std::unique_ptr<LocalAssemblerInterface>;
-
- explicit LocalDataInitializer(
- NumLib::LocalToGlobalIndexMap const& dof_table)
- : _dof_table(dof_table)
- {
- using EnabledElementTraits =
- decltype(BaseLib::TMP::filter<EnabledElementTraitsLagrange>(
- std::declval<IsElementEnabled>()));
-
- BaseLib::TMP::foreach<EnabledElementTraits>(
- [this]<typename ET>(ET*)
- {
- using MeshElement = typename ET::Element;
- using ShapeFunction = typename ET::ShapeFunction;
- using LowerOrderShapeFunction =
- typename ET::LowerOrderShapeFunction;
-
- _builder[std::type_index(typeid(MeshElement))] =
- makeLocalAssemblerBuilder<ShapeFunction,
- LowerOrderShapeFunction>();
- });
- }
-
- /// Returns data pointer to the newly created local assembler data.
- ///
- /// \attention
- /// The index \c id is not necessarily the mesh item's id. Especially when
- /// having multiple meshes it will differ from the latter.
- LADataIntfPtr operator()(std::size_t const id,
- MeshLib::Element const& mesh_item,
- ConstructorArgs&&... args) const
- {
- auto const type_idx = std::type_index(typeid(mesh_item));
- auto const it = _builder.find(type_idx);
-
- if (it != _builder.end())
- {
- auto const num_local_dof = _dof_table.getNumberOfElementDOF(id);
- return it->second(mesh_item, num_local_dof,
- std::forward<ConstructorArgs>(args)...);
- }
- 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, or a mesh element order does not match shape "
- "function order given in the project file.",
- type_idx.name());
- }
-
-private:
- using LADataBuilder =
- std::function<LADataIntfPtr(MeshLib::Element const& e,
- std::size_t const local_matrix_size,
- ConstructorArgs&&...)>;
-
- template <typename ShapeFunction>
- using IntegrationMethod = typename NumLib::GaussLegendreIntegrationPolicy<
- typename ShapeFunction::MeshElement>::IntegrationMethod;
-
- template <typename ShapeFunctionDisplacement,
- typename ShapeFunctionPressure>
- using LAData =
- LocalAssemblerData<ShapeFunctionDisplacement, ShapeFunctionPressure,
- IntegrationMethod<ShapeFunctionDisplacement>,
- GlobalDim>;
-
- /// Generates a function that creates a new LocalAssembler of type
- /// LAData<ShapeFunctionDisplacement, ShapeFunctionPressure>. 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, typename LowerOrderShapeFunction>
- static LADataBuilder makeLocalAssemblerBuilder()
- {
- return [](MeshLib::Element const& e,
- std::size_t const local_matrix_size,
- ConstructorArgs&&... args)
- {
- return LADataIntfPtr{
- new LAData<ShapeFunction, LowerOrderShapeFunction>{
- e, local_matrix_size,
- std::forward<ConstructorArgs>(args)...}};
- };
- }
-
- /// Mapping of element types to local assembler constructors.
- std::unordered_map<std::type_index, LADataBuilder> _builder;
-
- NumLib::LocalToGlobalIndexMap const& _dof_table;
-};
-
-} // namespace ProcessLib::StokesFlow
diff --git a/ProcessLib/Utils/LocalAssemblerFactoryTaylorHoodForOrderGeN.h b/ProcessLib/Utils/LocalAssemblerFactoryTaylorHoodForOrderGeN.h
new file mode 100644
index 00000000000..6a819e69d4e
--- /dev/null
+++ b/ProcessLib/Utils/LocalAssemblerFactoryTaylorHoodForOrderGeN.h
@@ -0,0 +1,171 @@
+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2021, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#pragma once
+
+#include "EnabledElements.h"
+#include "GenericLocalDataInitializer.h"
+
+namespace ProcessLib
+{
+template <typename ShapeFunction,
+ typename LowerOrderShapeFunction,
+ typename LocalAssemblerInterface,
+ template <typename /* shp fct */,
+ typename /* lower order shp fct */,
+ typename /* int meth */,
+ int /* global dim */>
+ class LocalAssemblerImplementation,
+ int GlobalDim,
+ typename... ConstructorArgs>
+class LocalAssemblerBuilderFactoryTaylorHood
+{
+ using GLDI = GenericLocalDataInitializer<LocalAssemblerInterface,
+ ConstructorArgs...>;
+ using LocAsmIntfPtr = typename GLDI::LocAsmIntfPtr;
+ using LocAsmBuilder = typename GLDI::LocAsmBuilder;
+
+ using IntegrationMethod = typename NumLib::GaussLegendreIntegrationPolicy<
+ typename ShapeFunction::MeshElement>::IntegrationMethod;
+
+ using LocAsmImpl = LocalAssemblerImplementation<ShapeFunction,
+ LowerOrderShapeFunction,
+ IntegrationMethod,
+ GlobalDim>;
+
+ LocalAssemblerBuilderFactoryTaylorHood() = delete;
+
+public:
+ /// Generates a function that creates a new local assembler of type
+ /// \c LocAsmImpl.
+ static LocAsmBuilder create()
+ {
+ return [](MeshLib::Element const& e,
+ std::size_t const local_matrix_size,
+ ConstructorArgs&&... args)
+ {
+ return std::make_unique<LocAsmImpl>(
+ e, local_matrix_size, std::forward<ConstructorArgs>(args)...);
+ };
+ }
+};
+
+/**
+ * Local assembler factory for Taylor-Hood elements.
+ *
+ * Elements/shape functions must be of order greater than or equal to \c
+ * MinShapeFctOrder.
+ *
+ * If \c MinShapeFctOrder is 1, local assemblers are instantiated also for
+ * linear mesh elements. In this case we don't have Taylor-Hood elements for
+ * linear mesh elements. Instead, on linear mesh elements all shape functions
+ * will have the same order (namely 1).
+ */
+template <int MinShapeFctOrder,
+ typename LocalAssemblerInterface,
+ template <typename /* shp fct */,
+ typename /* lower order shp fct */,
+ typename /* int meth */,
+ int /* global dim */>
+ class LocalAssemblerImplementation,
+ int GlobalDim,
+ typename... ConstructorArgs>
+class LocalAssemblerFactoryTaylorHoodForOrderGeN final
+ : public ProcessLib::GenericLocalDataInitializer<LocalAssemblerInterface,
+ ConstructorArgs...>
+{
+ using Base =
+ ProcessLib::GenericLocalDataInitializer<LocalAssemblerInterface,
+ ConstructorArgs...>;
+
+ template <typename ShapeFunction, typename LowerOrderShapeFunction>
+ using LocAsmBuilderFactory =
+ LocalAssemblerBuilderFactoryTaylorHood<ShapeFunction,
+ LowerOrderShapeFunction,
+ LocalAssemblerInterface,
+ LocalAssemblerImplementation,
+ GlobalDim,
+ ConstructorArgs...>;
+
+ struct IsElementEnabled
+ {
+ template <typename ElementTraits>
+ constexpr bool operator()(ElementTraits*) const
+ {
+ if constexpr (GlobalDim < ElementTraits::ShapeFunction::DIM)
+ {
+ return false;
+ }
+
+ if constexpr (ElementTraits::Element::dimension < 2)
+ {
+ // in OGS Taylor-Hood elements are available in 2D and 3D only
+ return false;
+ }
+
+ return ElementTraits::ShapeFunction::ORDER >= MinShapeFctOrder;
+ }
+ };
+
+public:
+ explicit LocalAssemblerFactoryTaylorHoodForOrderGeN(
+ NumLib::LocalToGlobalIndexMap const& dof_table)
+ : Base(dof_table)
+ {
+ using EnabledElementTraits =
+ decltype(BaseLib::TMP::filter<EnabledElementTraitsLagrange>(
+ std::declval<IsElementEnabled>()));
+
+ BaseLib::TMP::foreach<EnabledElementTraits>(
+ [this]<typename ET>(ET*)
+ {
+ using Elt = typename ET::Element;
+ using Shp = typename ET::ShapeFunction;
+ using ShpLow = typename ET::LowerOrderShapeFunction;
+
+ Base::_builders[std::type_index(typeid(Elt))] =
+ LocalAssemblerBuilderFactoryTaylorHood<
+ Shp,
+ ShpLow,
+ LocalAssemblerInterface,
+ LocalAssemblerImplementation,
+ GlobalDim,
+ ConstructorArgs...>::create();
+ });
+ }
+};
+
+/// HM processes in OGS are defined for linear and higher order elements.
+template <typename LocalAssemblerInterface,
+ template <typename, typename, typename, int>
+ class LocalAssemblerImplementation,
+ int GlobalDim,
+ typename... ConstructorArgs>
+using LocalAssemblerFactoryHM =
+ LocalAssemblerFactoryTaylorHoodForOrderGeN<1,
+ LocalAssemblerInterface,
+ LocalAssemblerImplementation,
+ GlobalDim,
+ ConstructorArgs...>;
+
+/// Stokes flow in OGS is defined for higher order elements only.
+template <typename LocalAssemblerInterface,
+ template <typename, typename, typename, int>
+ class LocalAssemblerImplementation,
+ int GlobalDim,
+ typename... ConstructorArgs>
+using LocalAssemblerFactoryStokes =
+ LocalAssemblerFactoryTaylorHoodForOrderGeN<2,
+ LocalAssemblerInterface,
+ LocalAssemblerImplementation,
+ GlobalDim,
+ ConstructorArgs...>;
+
+} // namespace ProcessLib
diff --git a/ProcessLib/Utils/LocalDataInitializerHM.h b/ProcessLib/Utils/LocalDataInitializerHM.h
deleted file mode 100644
index ba4eb33b8f4..00000000000
--- a/ProcessLib/Utils/LocalDataInitializerHM.h
+++ /dev/null
@@ -1,151 +0,0 @@
-/**
- * \file
- * \copyright
- * Copyright (c) 2012-2021, OpenGeoSys Community (http://www.opengeosys.org)
- * Distributed under a Modified BSD License.
- * See accompanying file LICENSE.txt or
- * http://www.opengeosys.org/project/license
- *
- */
-
-#pragma once
-
-#include <functional>
-#include <memory>
-#include <type_traits>
-#include <typeindex>
-#include <typeinfo>
-#include <unordered_map>
-
-#include "EnabledElements.h"
-#include "MeshLib/Elements/Elements.h"
-#include "NumLib/DOF/LocalToGlobalIndexMap.h"
-#include "NumLib/Fem/FiniteElement/LowerDimShapeTable.h"
-#include "NumLib/Fem/Integration/GaussLegendreIntegrationPolicy.h"
-
-namespace ProcessLib
-{
-/// The LocalDataInitializerHM is a functor creating a local assembler data with
-/// corresponding to the mesh element type shape functions and calling
-/// initialization of the new local assembler data.
-/// For example for MeshLib::Quad a local assembler data with template argument
-/// NumLib::ShapeQuad4 is created.
-///
-/// \attention This is modified version of the ProcessLib::LocalDataInitializer
-/// class which does not include line or point elements. For the shape functions
-/// of order 2 (used for displacement) a shape function of order 1 will be used
-/// for the scalar variables.
-template <typename LocalAssemblerInterface,
- template <typename, typename, typename, int> class LocalAssemblerData,
- int GlobalDim, typename... ConstructorArgs>
-class LocalDataInitializerHM final
-{
- struct IsElementEnabled
- {
- template <typename ElementTraits>
- constexpr bool operator()(ElementTraits*) const
- {
- if constexpr (GlobalDim < ElementTraits::ShapeFunction::DIM)
- {
- return false;
- }
-
- // mechanics processes in OGS are defined in 2D and 3D only
- return ElementTraits::Element::dimension >= 2;
- }
- };
-
-public:
- using LADataIntfPtr = std::unique_ptr<LocalAssemblerInterface>;
-
- explicit LocalDataInitializerHM(
- NumLib::LocalToGlobalIndexMap const& dof_table)
- : _dof_table(dof_table)
- {
- using EnabledElementTraits =
- decltype(BaseLib::TMP::filter<EnabledElementTraitsLagrange>(
- std::declval<IsElementEnabled>()));
-
- BaseLib::TMP::foreach<EnabledElementTraits>(
- [this]<typename ET>(ET*)
- {
- using MeshElement = typename ET::Element;
- using Shp = typename ET::ShapeFunction;
- using ShpLow = typename ET::LowerOrderShapeFunction;
-
- _builder[std::type_index(typeid(MeshElement))] =
- makeLocalAssemblerBuilder<Shp, ShpLow>();
- });
- }
-
- /// Returns data pointer to the newly created local assembler data.
- ///
- /// \attention
- /// The index \c id is not necessarily the mesh item's id. Especially when
- /// having multiple meshes it will differ from the latter.
- LADataIntfPtr operator()(std::size_t const id,
- MeshLib::Element const& mesh_item,
- ConstructorArgs&&... args) const
- {
- auto const type_idx = std::type_index(typeid(mesh_item));
- auto const it = _builder.find(type_idx);
-
- if (it != _builder.end())
- {
- auto const num_local_dof = _dof_table.getNumberOfElementDOF(id);
- return it->second(mesh_item, num_local_dof,
- std::forward<ConstructorArgs>(args)...);
- }
- 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, or a mesh element order does not match shape "
- "function order given in the project file.",
- type_idx.name());
- }
-
-private:
- using LADataBuilder =
- std::function<LADataIntfPtr(MeshLib::Element const& e,
- std::size_t const local_matrix_size,
- ConstructorArgs&&...)>;
-
- template <typename ShapeFunction>
- using IntegrationMethod = typename NumLib::GaussLegendreIntegrationPolicy<
- typename ShapeFunction::MeshElement>::IntegrationMethod;
-
- template <typename ShapeFunctionDisplacement,
- typename ShapeFunctionPressure>
- using LAData =
- LocalAssemblerData<ShapeFunctionDisplacement, ShapeFunctionPressure,
- IntegrationMethod<ShapeFunctionDisplacement>,
- GlobalDim>;
-
- /// Generates a function that creates a new LocalAssembler of type
- /// LAData<ShapeFunctionDisplacement, ShapeFunctionPressure>. 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, typename LowerOrderShapeFunction>
- static LADataBuilder makeLocalAssemblerBuilder()
- {
- return [](MeshLib::Element const& e,
- std::size_t const local_matrix_size,
- ConstructorArgs&&... args)
- {
- return LADataIntfPtr{
- new LAData<ShapeFunction, LowerOrderShapeFunction>{
- e, local_matrix_size,
- std::forward<ConstructorArgs>(args)...}};
- };
- }
-
- /// Mapping of element types to local assembler constructors.
- std::unordered_map<std::type_index, LADataBuilder> _builder;
-
- NumLib::LocalToGlobalIndexMap const& _dof_table;
-};
-
-} // namespace ProcessLib
--
GitLab