OGS
ProcessLib::ConstraintDirichletBoundaryCondition Class Referencefinal

Detailed Description

The ConstraintDirichletBoundaryCondition class describes a Dirichlet-type boundary condition that is constant in space and time where the domain can shrink and grow within the simulation. The expected parameter in the passed configuration is "value" which, when not present defaults to zero.

Definition at line 25 of file ConstraintDirichletBoundaryCondition.h.

#include <ConstraintDirichletBoundaryCondition.h>

Inheritance diagram for ProcessLib::ConstraintDirichletBoundaryCondition:
[legend]
Collaboration diagram for ProcessLib::ConstraintDirichletBoundaryCondition:
[legend]

Public Member Functions

 ConstraintDirichletBoundaryCondition (ParameterLib::Parameter< double > const &parameter, NumLib::LocalToGlobalIndexMap const &dof_table_bulk, int const variable_id, int const component_id, MeshLib::Mesh const &bc_mesh, unsigned const integration_order, MeshLib::Mesh const &bulk_mesh, double const constraint_threshold, bool const lower, std::function< Eigen::Vector3d(std::size_t const, MathLib::Point3d const &, double const, std::vector< GlobalVector * > const &)> getFlux)
 
void preTimestep (double const t, std::vector< GlobalVector * > const &x, int const process_id) override
 
void getEssentialBCValues (const double t, const GlobalVector &x, NumLib::IndexValueVector< GlobalIndexType > &bc_values) const override
 Writes the values of essential BCs to bc_values.
 
- Public Member Functions inherited from ProcessLib::BoundaryCondition
virtual void applyNaturalBC (const double, std::vector< GlobalVector * > const &, int const, GlobalMatrix *, GlobalVector &, GlobalMatrix *)
 
virtual void postTimestep (const double, std::vector< GlobalVector * > const &, int const)
 
virtual ~BoundaryCondition ()=default
 

Private Attributes

ParameterLib::Parameter< double > const & _parameter
 
std::unique_ptr< NumLib::LocalToGlobalIndexMap_dof_table_boundary
 
int const _variable_id
 
int const _component_id
 
MeshLib::Mesh const & _bc_mesh
 
unsigned const _integration_order
 Integration order for integration over the lower-dimensional elements.
 
std::vector< std::pair< std::size_t, unsigned > > _bulk_ids
 
std::vector< double > _flux_values
 Stores the results of the flux computations per boundary element.
 
std::vector< std::unique_ptr< ConstraintDirichletBoundaryConditionLocalAssemblerInterface > > _local_assemblers
 Local assemblers for each boundary element.
 
double const _constraint_threshold
 
bool const _lower
 
MeshLib::Mesh const & _bulk_mesh
 
std::function< Eigen::Vector3d(std::size_t const, MathLib::Point3d const &, double const, std::vector< GlobalVector * > const &)> _getFlux
 The function _getFlux calculates the flux through the boundary element.
 

Constructor & Destructor Documentation

◆ ConstraintDirichletBoundaryCondition()

ProcessLib::ConstraintDirichletBoundaryCondition::ConstraintDirichletBoundaryCondition ( ParameterLib::Parameter< double > const & parameter,
NumLib::LocalToGlobalIndexMap const & dof_table_bulk,
int const variable_id,
int const component_id,
MeshLib::Mesh const & bc_mesh,
unsigned const integration_order,
MeshLib::Mesh const & bulk_mesh,
double const constraint_threshold,
bool const lower,
std::function< Eigen::Vector3d(std::size_t const, MathLib::Point3d const &, double const, std::vector< GlobalVector * > const &)> getFlux )
Parameters
parameterUsed for setting the values for the boundary condition.
dof_table_bulkThe bulk local to global index map is used to derive the local to global index map for the boundary.
variable_idThe variable id is needed to determine the global index.
component_idThe component id is needed to determine the global index.
bc_meshLower dimensional mesh the boundary condition is defined on. The bc_mesh must have the two PropertyVector objects 'bulk_element_ids' and 'bulk_node_ids' containing the corresponding information.
integration_orderOrder the order of integration used to compute the constraint.
bulk_meshThe FE mesh for the simulation.
constraint_thresholdThe threshold value used for the switch off/on decision.
lowerBoolean value used for the calculation of the constraint criterion, i.e., if lower is set to true the criterion 'calculated_value < constraint_threshold' is evaluated to switch on/off the boundary condition, else 'calculated_value > constraint_threshold' is evaluated.
getFluxThe function used for the flux calculation.
Note
The function has to be stored by value, else the process value is not captured properly.

Definition at line 24 of file ConstraintDirichletBoundaryCondition.cpp.

34 : _parameter(parameter),
35 _variable_id(variable_id),
36 _component_id(component_id),
37 _bc_mesh(bc_mesh),
38 _integration_order(integration_order),
39 _constraint_threshold(constraint_threshold),
40 _lower(lower),
41 _bulk_mesh(bulk_mesh),
42 _getFlux(getFlux)
43{
44 if (variable_id >=
45 static_cast<int>(dof_table_bulk.getNumberOfVariables()) ||
46 component_id >=
47 dof_table_bulk.getNumberOfVariableComponents(variable_id))
48 {
50 "Variable id or component id too high. Actual values: ({:d}, "
51 "{:d}), maximum values: ({:d}, {:d}).",
52 variable_id, component_id, dof_table_bulk.getNumberOfVariables(),
53 dof_table_bulk.getNumberOfVariableComponents(variable_id));
54 }
55
56 std::vector<MeshLib::Node*> const& bc_nodes = _bc_mesh.getNodes();
57 DBUG(
58 "Found {:d} nodes for constraint Dirichlet BCs for the variable {:d} "
59 "and component {:d}",
60 bc_nodes.size(), variable_id, component_id);
61
62 MeshLib::MeshSubset bc_mesh_subset{_bc_mesh, bc_nodes};
63
64 // Create local DOF table from intersected mesh subsets for the given
65 // variable and component ids.
66 _dof_table_boundary = dof_table_bulk.deriveBoundaryConstrainedMap(
67 variable_id, {component_id}, std::move(bc_mesh_subset));
68
69 auto const& bc_elements(_bc_mesh.getElements());
70 _local_assemblers.resize(bc_elements.size());
71 _flux_values.resize(bc_elements.size());
72 // create _bulk_ids vector
73 auto const* bulk_element_ids = MeshLib::bulkElementIDs(_bc_mesh);
74 auto const* bulk_node_ids = MeshLib::bulkNodeIDs(_bc_mesh);
75 auto const& bulk_nodes = bulk_mesh.getNodes();
76
77 auto get_bulk_element_face_id =
78 [&](auto const bulk_element_id, MeshLib::Element const* bc_elem)
79 {
80 auto const* bulk_elem = _bulk_mesh.getElement(bulk_element_id);
81 std::array<MeshLib::Node const*, 3> nodes{
82 {bulk_nodes[(*bulk_node_ids)[bc_elem->getNode(0)->getID()]],
83 bulk_nodes[(*bulk_node_ids)[bc_elem->getNode(1)->getID()]],
84 bulk_nodes[(*bulk_node_ids)[bc_elem->getNode(2)->getID()]]}};
85 return bulk_elem->identifyFace(nodes.data());
86 };
87
88 _bulk_ids.reserve(bc_elements.size());
89 std::transform(
90 begin(bc_elements), end(bc_elements), std::back_inserter(_bulk_ids),
91 [&](auto const* bc_element)
92 {
93 auto const bulk_element_id =
94 (*bulk_element_ids)[bc_element->getID()];
95 return std::make_pair(
96 bulk_element_id,
97 get_bulk_element_face_id(bulk_element_id, bc_element));
98 });
99
100 const int shape_function_order = 1;
101
103 ConstraintDirichletBoundaryConditionLocalAssembler>(
105 shape_function_order, _local_assemblers,
106 NumLib::IntegrationOrder{_integration_order},
108}
#define OGS_FATAL(...)
Definition Error.h:26
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:30
virtual unsigned identifyFace(Node const *nodes[3]) const =0
Returns the ID of a face given an array of nodes.
A subset of nodes on a single mesh.
Definition MeshSubset.h:26
bool isAxiallySymmetric() const
Definition Mesh.h:137
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition Mesh.h:106
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition Mesh.h:109
unsigned getDimension() const
Returns the dimension of the mesh (determined by the maximum dimension over all elements).
Definition Mesh.h:88
const Element * getElement(std::size_t idx) const
Get the element with the given index.
Definition Mesh.h:94
std::function< Eigen::Vector3d(std::size_t const, MathLib::Point3d const &, double const, std::vector< GlobalVector * > const &)> _getFlux
The function _getFlux calculates the flux through the boundary element.
std::vector< std::unique_ptr< ConstraintDirichletBoundaryConditionLocalAssemblerInterface > > _local_assemblers
Local assemblers for each boundary element.
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _dof_table_boundary
std::vector< std::pair< std::size_t, unsigned > > _bulk_ids
std::vector< double > _flux_values
Stores the results of the flux computations per boundary element.
unsigned const _integration_order
Integration order for integration over the lower-dimensional elements.
PropertyVector< std::size_t > const * bulkElementIDs(Mesh const &mesh)
Definition Mesh.cpp:300
PropertyVector< std::size_t > const * bulkNodeIDs(Mesh const &mesh)
Definition Mesh.cpp:292
void createLocalAssemblers(const unsigned dimension, std::vector< MeshLib::Element * > const &mesh_elements, NumLib::LocalToGlobalIndexMap const &dof_table, const unsigned shapefunction_order, std::vector< std::unique_ptr< LocalAssemblerInterface > > &local_assemblers, NumLib::IntegrationOrder const integration_order, ExtraCtorArgs &&... extra_ctor_args)

References _bc_mesh, _dof_table_boundary, _flux_values, _local_assemblers, MeshLib::bulkElementIDs(), MeshLib::bulkNodeIDs(), DBUG(), NumLib::LocalToGlobalIndexMap::deriveBoundaryConstrainedMap(), MeshLib::Mesh::getElements(), MeshLib::Mesh::getNodes(), NumLib::LocalToGlobalIndexMap::getNumberOfVariableComponents(), and NumLib::LocalToGlobalIndexMap::getNumberOfVariables().

Member Function Documentation

◆ getEssentialBCValues()

void ProcessLib::ConstraintDirichletBoundaryCondition::getEssentialBCValues ( const double ,
const GlobalVector & ,
NumLib::IndexValueVector< GlobalIndexType > &  ) const
overridevirtual

Writes the values of essential BCs to bc_values.

Reimplemented from ProcessLib::BoundaryCondition.

Definition at line 127 of file ConstraintDirichletBoundaryCondition.cpp.

130{
132
133 bc_values.ids.clear();
134 bc_values.values.clear();
135
136 std::vector<std::pair<GlobalIndexType, double>> tmp_bc_values;
137
138 auto isFlux = [&](const std::size_t element_id)
139 {
140 return _lower ? _flux_values[element_id] < _constraint_threshold
141 : _flux_values[element_id] > _constraint_threshold;
142 };
143
144 for (auto const* boundary_element : _bc_mesh.getElements())
145 {
146 // check if the boundary element is active
147 if (isFlux(boundary_element->getID()))
148 {
149 continue;
150 }
151
152 // loop over the boundary element nodes and set the dirichlet value
153 unsigned const number_nodes = boundary_element->getNumberOfNodes();
154 for (unsigned i = 0; i < number_nodes; ++i)
155 {
156 auto const id = boundary_element->getNode(i)->getID();
157 pos.setAll(id, boundary_element->getID(), {}, {});
158
160 id);
161 // TODO: that might be slow, but only done once
162 const auto g_idx = _dof_table_boundary->getGlobalIndex(
165 {
166 continue;
167 }
168 // For the DDC approach (e.g. with PETSc option), the negative
169 // index of g_idx means that the entry by that index is a ghost one,
170 // which should be dropped. Especially for PETSc routines
171 // MatZeroRows and MatZeroRowsColumns, which are called to apply the
172 // Dirichlet BC, the negative index is not accepted like other
173 // matrix or vector PETSc routines. Therefore, the following
174 // if-condition is applied.
175 if (g_idx >= 0)
176 {
177 tmp_bc_values.emplace_back(g_idx, _parameter(t, pos).front());
178 }
179 }
180 }
181
182 if (tmp_bc_values.empty())
183 {
184 DBUG("The domain on which the boundary is defined is empty");
185 return;
186 }
187
188 // store unique pairs of node id and value with respect to the node id in
189 // the bc_values vector. The values related to the same node id are
190 // averaged.
191 // first: sort the (node id, value) pairs according to the node id
192 std::sort(tmp_bc_values.begin(), tmp_bc_values.end(),
193 [](std::pair<GlobalIndexType, double> const& a,
194 std::pair<GlobalIndexType, double> const& b)
195 { return a.first < b.first; });
196 // second: average the values over equal node id ranges
197 unsigned cnt = 1;
198 GlobalIndexType current_id = tmp_bc_values.begin()->first;
199 double sum = tmp_bc_values.begin()->second;
200 for (auto const& tmp_bc_value : tmp_bc_values)
201 {
202 if (tmp_bc_value.first == current_id)
203 {
204 cnt++;
205 sum += tmp_bc_value.second;
206 }
207 else
208 {
209 bc_values.ids.emplace_back(current_id);
210 bc_values.values.emplace_back(sum / cnt);
211 cnt = 1;
212 sum = tmp_bc_value.second;
213 current_id = tmp_bc_value.first;
214 }
215 }
216 bc_values.ids.emplace_back(current_id);
217 bc_values.values.emplace_back(sum / cnt);
218
219 DBUG("Found {:d} constraint dirichlet boundary condition values.",
220 bc_values.ids.size());
221}
GlobalMatrix::IndexType GlobalIndexType
std::size_t getID() const
Get id of the mesh.
Definition Mesh.h:121
static constexpr NUMLIB_EXPORT GlobalIndexType const nop
void setAll(std::optional< std::size_t > const &node_id, std::optional< std::size_t > const &element_id, std::optional< unsigned > const &integration_point, std::optional< MathLib::Point3d > const &coordinates)

References DBUG(), NumLib::IndexValueVector< typename >::ids, MeshLib::Node, NumLib::MeshComponentMap::nop, ParameterLib::SpatialPosition::setAll(), and NumLib::IndexValueVector< typename >::values.

◆ preTimestep()

void ProcessLib::ConstraintDirichletBoundaryCondition::preTimestep ( double const t,
std::vector< GlobalVector * > const & x,
int const process_id )
overridevirtual

Reimplemented from ProcessLib::BoundaryCondition.

Definition at line 110 of file ConstraintDirichletBoundaryCondition.cpp.

113{
114 DBUG(
115 "ConstraintDirichletBoundaryCondition::preTimestep: computing flux "
116 "constraints");
117 for (auto const id : _bc_mesh.getElements() | MeshLib::views::ids)
118 {
119 _flux_values[id] = _local_assemblers[id]->integrate(
120 x, t,
121 [this](std::size_t const element_id, MathLib::Point3d const& pnt,
122 double const t, std::vector<GlobalVector*> const& x)
123 { return _getFlux(element_id, pnt, t, x); });
124 }
125}
constexpr ranges::views::view_closure ids
For an element of a range view return its id.
Definition Mesh.h:225

References DBUG(), and MeshLib::views::ids.

Member Data Documentation

◆ _bc_mesh

MeshLib::Mesh const& ProcessLib::ConstraintDirichletBoundaryCondition::_bc_mesh
private

Vector of (lower-dimensional) boundary elements on which the boundary condition is defined.

Definition at line 82 of file ConstraintDirichletBoundaryCondition.h.

Referenced by ConstraintDirichletBoundaryCondition().

◆ _bulk_ids

std::vector<std::pair<std::size_t, unsigned> > ProcessLib::ConstraintDirichletBoundaryCondition::_bulk_ids
private

The first item of the pair is the element id in the bulk mesh, the second item is the face id of the bulk element that is part of the boundary

Definition at line 90 of file ConstraintDirichletBoundaryCondition.h.

◆ _bulk_mesh

MeshLib::Mesh const& ProcessLib::ConstraintDirichletBoundaryCondition::_bulk_mesh
private

The mesh _bulk_mesh is the discretized domain the process(es) are defined on. It is needed to get values for the constraint calculation.

Definition at line 112 of file ConstraintDirichletBoundaryCondition.h.

◆ _component_id

int const ProcessLib::ConstraintDirichletBoundaryCondition::_component_id
private

Definition at line 78 of file ConstraintDirichletBoundaryCondition.h.

◆ _constraint_threshold

double const ProcessLib::ConstraintDirichletBoundaryCondition::_constraint_threshold
private

The threshold value used to the switch off/on the Dirichlet-type boundary condition.

Definition at line 102 of file ConstraintDirichletBoundaryCondition.h.

◆ _dof_table_boundary

std::unique_ptr<NumLib::LocalToGlobalIndexMap> ProcessLib::ConstraintDirichletBoundaryCondition::_dof_table_boundary
private

Local dof table, a subset of the global one restricted to the participating number of elements of the boundary condition.

Definition at line 75 of file ConstraintDirichletBoundaryCondition.h.

Referenced by ConstraintDirichletBoundaryCondition().

◆ _flux_values

std::vector<double> ProcessLib::ConstraintDirichletBoundaryCondition::_flux_values
private

Stores the results of the flux computations per boundary element.

Definition at line 93 of file ConstraintDirichletBoundaryCondition.h.

Referenced by ConstraintDirichletBoundaryCondition().

◆ _getFlux

std::function<Eigen::Vector3d(std::size_t const, MathLib::Point3d const&, double const, std::vector<GlobalVector*> const&)> ProcessLib::ConstraintDirichletBoundaryCondition::_getFlux
private

The function _getFlux calculates the flux through the boundary element.

Definition at line 118 of file ConstraintDirichletBoundaryCondition.h.

◆ _integration_order

unsigned const ProcessLib::ConstraintDirichletBoundaryCondition::_integration_order
private

Integration order for integration over the lower-dimensional elements.

Definition at line 85 of file ConstraintDirichletBoundaryCondition.h.

◆ _local_assemblers

std::vector<std::unique_ptr< ConstraintDirichletBoundaryConditionLocalAssemblerInterface> > ProcessLib::ConstraintDirichletBoundaryCondition::_local_assemblers
private

Local assemblers for each boundary element.

Definition at line 98 of file ConstraintDirichletBoundaryCondition.h.

Referenced by ConstraintDirichletBoundaryCondition().

◆ _lower

bool const ProcessLib::ConstraintDirichletBoundaryCondition::_lower
private

The boolean value lower is used for the calculation of the constraint criterion, i.e., if lower is set to true the criterion 'calculated_value < constraint_threshold' is evaluated to switch on/off the boundary condition, else 'calculated_value > constraint_threshold' is evaluated.

Definition at line 108 of file ConstraintDirichletBoundaryCondition.h.

◆ _parameter

ParameterLib::Parameter<double> const& ProcessLib::ConstraintDirichletBoundaryCondition::_parameter
private

Definition at line 71 of file ConstraintDirichletBoundaryCondition.h.

◆ _variable_id

int const ProcessLib::ConstraintDirichletBoundaryCondition::_variable_id
private

Definition at line 77 of file ConstraintDirichletBoundaryCondition.h.


The documentation for this class was generated from the following files: