version 0.4.4
utils/dirichletvalues.hh
Go to the documentation of this file.
1// SPDX-FileCopyrightText: 2021-2025 The Ikarus Developers ikarus@ibb.uni-stuttgart.de
2// SPDX-License-Identifier: LGPL-3.0-or-later
3
10#pragma once
11#include <algorithm>
12#include <concepts>
13#include <cstddef>
14#include <functional>
15#include <type_traits>
16#include <utility>
17#include <vector>
18
19#include <dune/common/float_cmp.hh>
20#include <dune/functions/backends/istlvectorbackend.hh>
21#include <dune/functions/functionspacebases/boundarydofs.hh>
22#include <dune/functions/functionspacebases/flatmultiindex.hh>
23#include <dune/functions/functionspacebases/interpolate.hh>
24#include <dune/functions/functionspacebases/subspacebasis.hh>
25
26#include <Eigen/Core>
27
28#include <autodiff/forward/real/real.hpp>
29
31
32namespace Dune {
33template <class K, int N>
35}
36
37namespace Ikarus {
38
42template <class>
44
45template <Concepts::EigenVector T>
47{
48 using SizeType = Eigen::Index;
49};
50
51template <>
52struct DeriveSizeType<std::vector<bool>>
53{
54 using SizeType = std::vector<bool>::size_type;
55};
56
67template <typename B, typename FC = std::vector<bool>>
69{
70public:
71 using Basis = std::remove_cvref_t<B>;
72 using FlagsType = FC;
73 static constexpr int worldDimension = Basis::GridView::dimensionworld;
74 using BackendType = decltype(Dune::Functions::istlVectorBackend(std::declval<FlagsType&>()));
76 explicit DirichletValues(const B& basis)
77 : basis_{basis},
78 dirichletFlagsBackend_{dirichletFlags_} {
79 dirichletFlagsBackend_.resize(basis_);
80 std::fill(dirichletFlags_.begin(), dirichletFlags_.end(), false);
81 }
82
92 template <typename F, typename TreePath = Dune::TypeTree::HybridTreePath<>>
93 void fixBoundaryDOFs(F&& f, TreePath&& treePath = {}) {
94 using namespace Dune::Functions;
95 using SubSpaceLocalView =
96 typename std::remove_cvref_t<decltype(subspaceBasis(basis_, std::forward<TreePath>(treePath)))>::LocalView;
97
99 auto lambda = [&](auto&& indexGlobal) { f(dirichletFlagsBackend_, indexGlobal); };
100 Dune::Functions::forEachBoundaryDOF(subspaceBasis(basis_, std::forward<TreePath>(treePath)), lambda);
101 } else if constexpr (Concepts::IsFunctorWithArgs<F, BackendType, int, SubSpaceLocalView>) {
102 auto lambda = [&](auto&& localIndex, auto&& localView) { f(dirichletFlagsBackend_, localIndex, localView); };
103 Dune::Functions::forEachBoundaryDOF(subspaceBasis(basis_, std::forward<TreePath>(treePath)), lambda);
104 } else if constexpr (Concepts::IsFunctorWithArgs<F, BackendType, int, SubSpaceLocalView,
105 typename Basis::GridView::Intersection>) {
106 auto lambda = [&](auto&& localIndex, auto&& localView, auto&& intersection) {
107 f(dirichletFlagsBackend_, localIndex, localView, intersection);
108 };
109 Dune::Functions::forEachBoundaryDOF(subspaceBasis(basis_, std::forward<TreePath>(treePath)), lambda);
110 } else {
111 static_assert(Dune::AlwaysFalse<F>(), "fixBoundaryDOFs: A function with this signature is not supported");
112 }
113 }
114
123 template <typename F>
124 void fixDOFs(F&& f) {
125 f(basis_, dirichletFlagsBackend_);
126 }
127
134 template <typename MultiIndex>
135 requires(not std::integral<MultiIndex>)
136 void setSingleDOF(const MultiIndex i, bool flag) {
137 dirichletFlagsBackend_[i] = flag;
138 }
139
147 void setSingleDOF(std::size_t i, bool flag)
148 requires(std::same_as<typename Basis::MultiIndex, Dune::Functions::FlatMultiIndex<size_t>>)
149 {
150 dirichletFlags_[i] = flag;
151 }
152
156 void reset() { std::fill(dirichletFlags_.begin(), dirichletFlags_.end(), false); }
157
158 /* \brief Returns the local basis object */
159 const auto& basis() const { return basis_; }
160
161 /* \brief Returns a boolean values, if the give multiIndex is constrained */
162 template <typename MultiIndex>
163 requires(not std::integral<MultiIndex>)
164 [[nodiscard]] bool isConstrained(const MultiIndex& multiIndex) const {
165 return dirichletFlagsBackend_[multiIndex];
166 }
167
168 /* \brief Returns a boolean values, if the i-th degree of freedom is constrained */
169 [[nodiscard]] bool isConstrained(std::size_t i) const
170 requires(std::same_as<typename Basis::MultiIndex, Dune::Functions::FlatMultiIndex<size_t>>)
171 {
172 return dirichletFlags_[i];
173 }
174
175 /* \brief Returns how many degrees of freedoms are fixed */
176 auto fixedDOFsize() const { return std::ranges::count(dirichletFlags_, true); }
177
178 /* \brief Returns how many degrees of freedom there are */
179 auto size() const { return dirichletFlags_.size(); }
180
181 /* \brief Returns the underlying dof flag container */
182 auto& container() const { return dirichletFlags_; }
183
195 template <typename F>
196 void storeInhomogeneousBoundaryCondition(F&& f, double lambda = 1.0) {
197 auto derivativeLambda = [&](const auto& globalCoord, const double& lambda) {
198 autodiff::real lambdaDual = lambda;
199 lambdaDual[1] = 1; // Setting the derivative in lambda direction to 1
200 return derivative(f(globalCoord, lambdaDual));
201 };
202 dirichletFunctions_.push_back({f, derivativeLambda});
203 setInhomogeneousBoundaryConditionFlag(lambda);
204 }
205
211 void setZeroAtConstrainedDofs(Eigen::VectorXd& xIh) const {
212 for (Eigen::Index i = 0; i < xIh.size(); ++i)
213 if (this->isConstrained(i))
214 xIh[i] = 0.0;
215 }
216
226 void evaluateInhomogeneousBoundaryCondition(Eigen::VectorXd& xIh, const double& lambda) const {
227 Eigen::VectorXd inhomogeneousBoundaryVectorDummy;
228 inhomogeneousBoundaryVectorDummy.setZero(this->size());
229 xIh.resizeLike(inhomogeneousBoundaryVectorDummy);
230 xIh.setZero();
231 for (auto& f : dirichletFunctions_) {
232 interpolate(basis_, inhomogeneousBoundaryVectorDummy,
233 [&](const auto& globalCoord) { return f.value(globalCoord, lambda); });
234 xIh += inhomogeneousBoundaryVectorDummy;
235 }
236 }
237
247 void evaluateInhomogeneousBoundaryConditionDerivative(Eigen::VectorXd& xIh, const double& lambda) const {
248 Eigen::VectorXd inhomogeneousBoundaryVectorDummy;
249 inhomogeneousBoundaryVectorDummy.setZero(this->size());
250 xIh.resizeLike(inhomogeneousBoundaryVectorDummy);
251 xIh.setZero();
252 for (auto& f : dirichletFunctions_) {
253 interpolate(basis_, inhomogeneousBoundaryVectorDummy,
254 [&](const auto& globalCoord) { return f.derivative(globalCoord, lambda); });
255 xIh += inhomogeneousBoundaryVectorDummy;
256 }
257 }
258
259private:
260 Basis basis_;
261 FlagsType dirichletFlags_;
262 BackendType dirichletFlagsBackend_;
263 struct DirichletFunctions
264 {
265 using Signature = std::function<Eigen::Vector<double, worldDimension>(
266 const Dune::FieldVector<double, worldDimension>&, const double&)>;
267 Signature value;
268 Signature derivative;
269 };
270 std::vector<DirichletFunctions> dirichletFunctions_;
271
272 void setInhomogeneousBoundaryConditionFlag(double lambda) {
273 Eigen::VectorXd inhomogeneousBoundaryVectorDummy;
274 inhomogeneousBoundaryVectorDummy.setZero(this->size());
275 this->evaluateInhomogeneousBoundaryCondition(inhomogeneousBoundaryVectorDummy, lambda);
276 for (const std::size_t i : Dune::range(this->size()))
277 if (Dune::FloatCmp::ne(inhomogeneousBoundaryVectorDummy[i], 0.0))
278 this->setSingleDOF(i, true);
279 }
280};
281
282} // namespace Ikarus
Definition: assemblermanipulatorbuildingblocks.hh:22
Definition: utils/dirichletvalues.hh:32
Definition: utils/dirichletvalues.hh:34
A helper struct to derive the SizeType of the underlying container.
Definition: utils/dirichletvalues.hh:43
Eigen::Index SizeType
Definition: utils/dirichletvalues.hh:48
std::vector< bool >::size_type SizeType
Definition: utils/dirichletvalues.hh:54
Class for handling Dirichlet boundary conditions in Ikarus.
Definition: utils/dirichletvalues.hh:69
auto size() const
Definition: utils/dirichletvalues.hh:179
void setZeroAtConstrainedDofs(Eigen::VectorXd &xIh) const
Function to write zeros at constrained Dirichlet entries.
Definition: utils/dirichletvalues.hh:211
void setSingleDOF(const MultiIndex i, bool flag)
Fixes and unfixes (set boolean value to true or false) a specific degree of freedom.
Definition: utils/dirichletvalues.hh:136
bool isConstrained(const MultiIndex &multiIndex) const
Definition: utils/dirichletvalues.hh:164
bool isConstrained(std::size_t i) const
Definition: utils/dirichletvalues.hh:169
void setSingleDOF(std::size_t i, bool flag)
Fixes or unfixes (set boolean value to true or false) a specific degree of freedom.
Definition: utils/dirichletvalues.hh:147
std::remove_cvref_t< B > Basis
Definition: utils/dirichletvalues.hh:71
void fixBoundaryDOFs(F &&f, TreePath &&treePath={})
Function to fix (set boolean values to true or false) degrees of freedom on the boundary.
Definition: utils/dirichletvalues.hh:93
decltype(Dune::Functions::istlVectorBackend(std::declval< FlagsType & >())) BackendType
Definition: utils/dirichletvalues.hh:74
void evaluateInhomogeneousBoundaryConditionDerivative(Eigen::VectorXd &xIh, const double &lambda) const
Function to evaluate all stored inhomogeneous Dirichlet boundary derivative functions.
Definition: utils/dirichletvalues.hh:247
void storeInhomogeneousBoundaryCondition(F &&f, double lambda=1.0)
Function to insert a function for inhomogeneous Dirichlet boundary conditions.
Definition: utils/dirichletvalues.hh:196
static constexpr int worldDimension
Definition: utils/dirichletvalues.hh:73
void fixDOFs(F &&f)
Function to fix (set boolean values to true or false) degrees of freedom.
Definition: utils/dirichletvalues.hh:124
FC FlagsType
Definition: utils/dirichletvalues.hh:72
void evaluateInhomogeneousBoundaryCondition(Eigen::VectorXd &xIh, const double &lambda) const
Function to evaluate all stored inhomogeneous Dirichlet boundary functions.
Definition: utils/dirichletvalues.hh:226
void reset()
Resets all degrees of freedom.
Definition: utils/dirichletvalues.hh:156
typename DeriveSizeType< FlagsType >::SizeType SizeType
Definition: utils/dirichletvalues.hh:75
const auto & basis() const
Definition: utils/dirichletvalues.hh:159
auto fixedDOFsize() const
Definition: utils/dirichletvalues.hh:176
DirichletValues(const B &basis)
Definition: utils/dirichletvalues.hh:76
auto & container() const
Definition: utils/dirichletvalues.hh:182
Concept defining the requirements for functors with arguments.
Definition: utils/concepts.hh:357
Several concepts.