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>
28#include <autodiff/forward/real/real.hpp>
33template <
class K,
int N>
45template <Concepts::EigenVector T>
67template <
typename B,
typename FC = std::vector<
bool>>
71 using Basis = std::remove_cvref_t<B>;
74 using BackendType =
decltype(Dune::Functions::istlVectorBackend(std::declval<FlagsType&>()));
78 dirichletFlagsBackend_{dirichletFlags_} {
79 dirichletFlagsBackend_.resize(basis_);
80 std::fill(dirichletFlags_.begin(), dirichletFlags_.end(),
false);
92 template <
typename F,
typename TreePath = Dune::TypeTree::Hybr
idTreePath<>>
94 using namespace Dune::Functions;
95 using SubSpaceLocalView =
96 typename std::remove_cvref_t<decltype(subspaceBasis(basis_, std::forward<TreePath>(treePath)))>::LocalView;
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);
109 Dune::Functions::forEachBoundaryDOF(subspaceBasis(basis_, std::forward<TreePath>(treePath)), lambda);
111 static_assert(Dune::AlwaysFalse<F>(),
"fixBoundaryDOFs: A function with this signature is not supported");
123 template <
typename F>
125 f(basis_, dirichletFlagsBackend_);
134 template <
typename MultiIndex>
135 requires(not std::integral<MultiIndex>)
137 dirichletFlagsBackend_[i] = flag;
148 requires(std::same_as<typename Basis::MultiIndex, Dune::Functions::FlatMultiIndex<size_t>>)
150 dirichletFlags_[i] = flag;
156 void reset() { std::fill(dirichletFlags_.begin(), dirichletFlags_.end(),
false); }
159 const auto&
basis()
const {
return basis_; }
162 template <
typename MultiIndex>
163 requires(not std::integral<MultiIndex>)
165 return dirichletFlagsBackend_[multiIndex];
170 requires(std::same_as<typename Basis::MultiIndex, Dune::Functions::FlatMultiIndex<size_t>>)
172 return dirichletFlags_[i];
176 auto fixedDOFsize()
const {
return std::ranges::count(dirichletFlags_,
true); }
179 auto size()
const {
return dirichletFlags_.size(); }
195 template <
typename F>
197 auto derivativeLambda = [&](
const auto& globalCoord,
const double& lambda) {
198 autodiff::real lambdaDual = lambda;
200 return derivative(f(globalCoord, lambdaDual));
202 dirichletFunctions_.push_back({f, derivativeLambda});
203 setInhomogeneousBoundaryConditionFlag(lambda);
212 for (Eigen::Index i = 0; i < xIh.size(); ++i)
227 Eigen::VectorXd inhomogeneousBoundaryVectorDummy;
228 inhomogeneousBoundaryVectorDummy.setZero(this->
size());
229 xIh.resizeLike(inhomogeneousBoundaryVectorDummy);
231 for (
auto& f : dirichletFunctions_) {
232 interpolate(basis_, inhomogeneousBoundaryVectorDummy,
233 [&](
const auto& globalCoord) {
return f.value(globalCoord, lambda); });
234 xIh += inhomogeneousBoundaryVectorDummy;
248 Eigen::VectorXd inhomogeneousBoundaryVectorDummy;
249 inhomogeneousBoundaryVectorDummy.setZero(this->
size());
250 xIh.resizeLike(inhomogeneousBoundaryVectorDummy);
252 for (
auto& f : dirichletFunctions_) {
253 interpolate(basis_, inhomogeneousBoundaryVectorDummy,
254 [&](
const auto& globalCoord) {
return f.derivative(globalCoord, lambda); });
255 xIh += inhomogeneousBoundaryVectorDummy;
263 struct DirichletFunctions
265 using Signature = std::function<Eigen::Vector<double, worldDimension>(
268 Signature derivative;
270 std::vector<DirichletFunctions> dirichletFunctions_;
272 void setInhomogeneousBoundaryConditionFlag(
double lambda) {
273 Eigen::VectorXd inhomogeneousBoundaryVectorDummy;
274 inhomogeneousBoundaryVectorDummy.setZero(this->
size());
276 for (
const std::size_t i : Dune::range(this->
size()))
277 if (Dune::FloatCmp::ne(inhomogeneousBoundaryVectorDummy[i], 0.0))
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