19#include <dune/common/float_cmp.hh>
20#include <dune/common/hybridutilities.hh>
21#include <dune/common/indices.hh>
22#include <dune/functions/backends/istlvectorbackend.hh>
23#include <dune/functions/functionspacebases/boundarydofs.hh>
24#include <dune/functions/functionspacebases/flatmultiindex.hh>
25#include <dune/functions/functionspacebases/interpolate.hh>
26#include <dune/functions/functionspacebases/subspacebasis.hh>
30#include <autodiff/forward/real/real.hpp>
35template <
class K,
int N>
47template <Concepts::EigenVector T>
60 template <
typename Tree>
65 template <
typename Tree>
66 requires(Tree::isLeaf)
67 struct PreBasisInfo<Tree>
69 static constexpr std::size_t size = 0;
70 using NodalSolutionType = double;
73 template <
typename Tree>
74 requires(Tree::isPower)
75 struct PreBasisInfo<Tree>
77 static constexpr std::size_t size = Tree::degree();
78 using NodalSolutionType = Eigen::Vector<double, size>;
81 template <
typename Tree>
82 requires(Tree::isComposite)
83 struct PreBasisInfo<Tree>
85 using ChildTreeType = Tree::template Child<0>::Type;
86 static_assert(not ChildTreeType::isComposite,
87 "DirichletValues is not implemented to handle a composite basis within a composite basis.");
89 static constexpr std::size_t size = PreBasisInfo<ChildTreeType>::size;
90 using NodalSolutionType = PreBasisInfo<ChildTreeType>::NodalSolutionType;
106template <
typename B,
typename FC = std::vector<
bool>>
110 using Basis = std::remove_cvref_t<B>;
112 using BackendType =
decltype(Dune::Functions::istlVectorBackend(std::declval<FlagsType&>()));
121 dirichletFlagsBackend_{dirichletFlags_} {
122 dirichletFlagsBackend_.resize(basis_);
123 std::fill(dirichletFlags_.begin(), dirichletFlags_.end(),
false);
124 if constexpr (Tree::isComposite) {
125 Dune::Hybrid::forEach(
126 Dune::Hybrid::integralRange(Dune::index_constant<std::tuple_size_v<typename Basis::PreBasis::SubPreBases>>()),
128 static_assert(not Tree::template Child<i>::Type::isComposite,
129 "DirichletValues is not implemented to handle a composite basis within a composite basis.");
143 template <
typename F,
typename TreePath = Dune::TypeTree::Hybr
idTreePath<>>
145 using namespace Dune::Functions;
146 using SubSpaceLocalView =
147 typename std::remove_cvref_t<decltype(subspaceBasis(basis_, std::forward<TreePath>(treePath)))>
::LocalView;
150 auto lambda = [&](
auto&& indexGlobal) { f(dirichletFlagsBackend_, indexGlobal); };
151 Dune::Functions::forEachBoundaryDOF(subspaceBasis(basis_, std::forward<TreePath>(treePath)), lambda);
152 }
else if constexpr (Concepts::IsFunctorWithArgs<F, BackendType, int, SubSpaceLocalView>) {
153 auto lambda = [&](
auto&& localIndex,
auto&& localView) { f(dirichletFlagsBackend_, localIndex, localView); };
154 Dune::Functions::forEachBoundaryDOF(subspaceBasis(basis_, std::forward<TreePath>(treePath)), lambda);
155 }
else if constexpr (Concepts::IsFunctorWithArgs<F,
BackendType, int, SubSpaceLocalView,
156 typename Basis::GridView::Intersection>) {
157 auto lambda = [&](
auto&& localIndex,
auto&& localView,
auto&& intersection) {
158 f(dirichletFlagsBackend_, localIndex, localView, intersection);
160 Dune::Functions::forEachBoundaryDOF(subspaceBasis(basis_, std::forward<TreePath>(treePath)), lambda);
162 static_assert(Dune::AlwaysFalse<F>(),
"fixBoundaryDOFs: A function with this signature is not supported");
174 template <
typename F>
176 f(basis_, dirichletFlagsBackend_);
185 template <
typename MultiIndex>
186 requires(not std::integral<MultiIndex>)
188 dirichletFlagsBackend_[i] = flag;
199 requires(std::same_as<typename Basis::MultiIndex, Dune::Functions::FlatMultiIndex<size_t>>)
201 dirichletFlags_[i] = flag;
207 void reset() { std::fill(dirichletFlags_.begin(), dirichletFlags_.end(),
false); }
210 const auto&
basis()
const {
return basis_; }
213 template <
typename MultiIndex>
214 requires(not std::integral<MultiIndex>)
216 return dirichletFlagsBackend_[multiIndex];
221 requires(std::same_as<typename Basis::MultiIndex, Dune::Functions::FlatMultiIndex<size_t>>)
223 return dirichletFlags_[i];
227 auto fixedDOFsize()
const {
return std::ranges::count(dirichletFlags_,
true); }
230 auto size()
const {
return dirichletFlags_.size(); }
246 template <
typename F>
248 auto derivativeLambda = [&](
const auto& globalCoord,
const double& lambda) {
249 autodiff::real lambdaDual = lambda;
251 return derivative(f(globalCoord, lambdaDual));
253 dirichletFunctions_.push_back({f, derivativeLambda});
254 setInhomogeneousBoundaryConditionFlag(lambda);
263 for (Eigen::Index i = 0; i < xIh.size(); ++i)
278 Eigen::VectorXd inhomogeneousBoundaryVectorDummy;
279 inhomogeneousBoundaryVectorDummy.setZero(this->
size());
280 xIh.resizeLike(inhomogeneousBoundaryVectorDummy);
282 auto runInterpolateImpl = [&](
const auto&
basis) {
283 for (
auto& f : dirichletFunctions_) {
284 interpolate(
basis, inhomogeneousBoundaryVectorDummy,
285 [&](
const auto& globalCoord) {
return f.value(globalCoord, lambda); });
286 xIh += inhomogeneousBoundaryVectorDummy;
289 runInterpolate(runInterpolateImpl);
302 Eigen::VectorXd inhomogeneousBoundaryVectorDummy;
303 inhomogeneousBoundaryVectorDummy.setZero(this->
size());
304 xIh.resizeLike(inhomogeneousBoundaryVectorDummy);
306 auto runInterpolateImpl = [&](
const auto&
basis) {
307 for (
auto& f : dirichletFunctions_) {
308 interpolate(
basis, inhomogeneousBoundaryVectorDummy,
309 [&](
const auto& globalCoord) {
return f.derivative(globalCoord, lambda); });
310 xIh += inhomogeneousBoundaryVectorDummy;
313 runInterpolate(runInterpolateImpl);
320 struct DirichletFunctions
324 Signature derivative;
326 std::vector<DirichletFunctions> dirichletFunctions_;
328 void setInhomogeneousBoundaryConditionFlag(
double lambda) {
329 Eigen::VectorXd inhomogeneousBoundaryVectorDummy;
330 inhomogeneousBoundaryVectorDummy.setZero(this->
size());
332 for (
const std::size_t i : Dune::range(this->
size()))
333 if (Dune::FloatCmp::ne(inhomogeneousBoundaryVectorDummy[i], 0.0))
337 template <
typename F>
338 void runInterpolate(F&& f)
const {
339 if constexpr (Tree::isComposite)
340 f(subspaceBasis(basis_, Dune::Indices::_0));
Definition: assemblermanipulatorbuildingblocks.hh:22
Definition: utils/dirichletvalues.hh:34
Definition: utils/dirichletvalues.hh:36
A helper struct to derive the SizeType of the underlying container.
Definition: utils/dirichletvalues.hh:45
Eigen::Index SizeType
Definition: utils/dirichletvalues.hh:50
std::vector< bool >::size_type SizeType
Definition: utils/dirichletvalues.hh:56
Class for handling Dirichlet boundary conditions in Ikarus.
Definition: utils/dirichletvalues.hh:108
auto size() const
Definition: utils/dirichletvalues.hh:230
void setZeroAtConstrainedDofs(Eigen::VectorXd &xIh) const
Function to write zeros at constrained Dirichlet entries.
Definition: utils/dirichletvalues.hh:262
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:187
bool isConstrained(const MultiIndex &multiIndex) const
Definition: utils/dirichletvalues.hh:215
bool isConstrained(std::size_t i) const
Definition: utils/dirichletvalues.hh:220
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:198
std::remove_cvref_t< B > Basis
Definition: utils/dirichletvalues.hh:110
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:144
decltype(Dune::Functions::istlVectorBackend(std::declval< FlagsType & >())) BackendType
Definition: utils/dirichletvalues.hh:112
void evaluateInhomogeneousBoundaryConditionDerivative(Eigen::VectorXd &xIh, const double &lambda) const
Function to evaluate all stored inhomogeneous Dirichlet boundary derivative functions.
Definition: utils/dirichletvalues.hh:301
void storeInhomogeneousBoundaryCondition(F &&f, double lambda=1.0)
Function to insert a function for inhomogeneous Dirichlet boundary conditions.
Definition: utils/dirichletvalues.hh:247
LocalView::Tree Tree
Definition: utils/dirichletvalues.hh:116
static constexpr int worldDimension
Definition: utils/dirichletvalues.hh:114
void fixDOFs(F &&f)
Function to fix (set boolean values to true or false) degrees of freedom.
Definition: utils/dirichletvalues.hh:175
FC FlagsType
Definition: utils/dirichletvalues.hh:111
void evaluateInhomogeneousBoundaryCondition(Eigen::VectorXd &xIh, const double &lambda) const
Function to evaluate all stored inhomogeneous Dirichlet boundary functions.
Definition: utils/dirichletvalues.hh:277
Basis::LocalView LocalView
Definition: utils/dirichletvalues.hh:115
void reset()
Resets all degrees of freedom.
Definition: utils/dirichletvalues.hh:207
typename DeriveSizeType< FlagsType >::SizeType SizeType
Definition: utils/dirichletvalues.hh:113
const auto & basis() const
Definition: utils/dirichletvalues.hh:210
auto fixedDOFsize() const
Definition: utils/dirichletvalues.hh:227
DirichletValues(const B &basis)
Definition: utils/dirichletvalues.hh:119
static constexpr std::size_t numberOfChildrenAtNode
Definition: utils/dirichletvalues.hh:118
auto & container() const
Definition: utils/dirichletvalues.hh:233
typename Impl::PreBasisInfo< Tree >::NodalSolutionType NodalSolutionType
Definition: utils/dirichletvalues.hh:117
Concept defining the requirements for functors with arguments.
Definition: utils/concepts.hh:357