3 #ifndef DUNE_FUNCTIONS_FUNCTIONSPACEBASES_INTERPOLATE_HH 4 #define DUNE_FUNCTIONS_FUNCTIONSPACEBASES_INTERPOLATE_HH 9 #include <dune/common/exceptions.hh> 10 #include <dune/common/bitsetvector.hh> 12 #include <dune/typetree/childextraction.hh> 14 #include <dune/localfunctions/common/virtualinterface.hh> 25 #include <dune/typetree/traversal.hh> 26 #include <dune/typetree/visitor.hh> 33 struct AllTrueBitSetVector
37 bool test(
int i)
const {
return true; }
46 const AllTrueBitSetVector& operator[](
const I&)
const 52 void resize(
const SP&)
const 59 template <
class B,
class T,
class NTRE,
class HV,
class LF,
class HBV>
60 class LocalInterpolateVisitor
61 :
public TypeTree::TreeVisitor
62 ,
public TypeTree::DynamicTraversal
68 using LocalIndexSet =
typename B::LocalIndexSet;
69 using MultiIndex =
typename LocalIndexSet::MultiIndex;
71 using LocalFunction = LF;
75 using HierarchicVector = HV;
76 using HierarchicBitVector = HBV;
78 using NodeToRangeEntry = NTRE;
80 using GridView =
typename Basis::GridView;
81 using Element =
typename GridView::template Codim<0>::Entity;
83 using LocalDomain =
typename Element::Geometry::LocalCoordinate;
85 using GlobalDomain =
typename Element::Geometry::GlobalCoordinate;
87 using CoefficientBlock =
typename std::decay<decltype(std::declval<HierarchicVector>()[std::declval<MultiIndex>()])>::type;
88 using BitVectorBlock =
typename std::decay<decltype(std::declval<HierarchicBitVector>()[std::declval<MultiIndex>()])>::type;
90 LocalInterpolateVisitor(
const B& basis, HV& coeff,
const HBV& bitVector,
const LF& localF,
const LocalIndexSet& localIndexSet,
const NodeToRangeEntry& nodeToRangeEntry) :
93 bitVector_(bitVector),
94 localIndexSet_(localIndexSet),
95 nodeToRangeEntry_(nodeToRangeEntry)
97 static_assert(Dune::Functions::Concept::isCallable<LocalFunction, LocalDomain>(),
"Function passed to LocalInterpolateVisitor does not model the Callable<LocalCoordinate> concept");
100 template<
typename Node,
typename TreePath>
101 void pre(Node& node, TreePath treePath)
104 template<
typename Node,
typename TreePath>
105 void post(Node& node, TreePath treePath)
108 template<
typename Node,
typename TreePath>
109 void leaf(Node& node, TreePath treePath)
111 using FiniteElement =
typename Node::FiniteElement;
112 using FiniteElementRange =
typename FiniteElement::Traits::LocalBasisType::Traits::RangeType;
113 using FiniteElementRangeField =
typename FiniteElement::Traits::LocalBasisType::Traits::RangeFieldType;
114 using FunctionBaseClass =
typename Dune::LocalFiniteElementFunctionBase<FiniteElement>::type;
121 auto localFj = [&](
const LocalDomain& x){
122 const auto& y = localF_(x);
123 const auto& y_node = nodeToRangeEntry_(node, y);
124 using FunctionRange =
typename std::decay<decltype(y_node)>::type;
130 auto interpolationCoefficients = std::vector<FiniteElementRangeField>();
132 auto&& fe = node.finiteElement();
139 for(j=0; j<blockSize; ++j)
143 fe.localInterpolation().interpolate(FunctionFromCallable(localFj), interpolationCoefficients);
144 for (
size_t i=0; i<fe.localBasis().size(); ++i)
146 auto multiIndex = localIndexSet_.index(node.localIndex(i));
147 const auto& bitVectorBlock = bitVector_[multiIndex];
152 auto&& vectorBlock = vector_[multiIndex];
162 HierarchicVector& vector_;
163 const LocalFunction& localF_;
164 const HierarchicBitVector& bitVector_;
165 const LocalIndexSet& localIndexSet_;
166 const NodeToRangeEntry& nodeToRangeEntry_;
192 template <
class B,
class... TreeIndices,
class NTRE,
class C,
class F,
class BV>
193 void interpolateTreeSubset(
const B& basis,
const TypeTree::HybridTreePath<TreeIndices...>& treePath, C&& coeff,
const F& f,
const NTRE& nodeToRangeEntry,
const BV& bv)
195 using GridView =
typename B::GridView;
196 using Element =
typename GridView::template Codim<0>::Entity;
198 using Tree =
typename std::decay<decltype(TypeTree::child(basis.localView().tree(),treePath))>::type;
200 using GlobalDomain =
typename Element::Geometry::GlobalCoordinate;
202 static_assert(Dune::Functions::Concept::isCallable<F, GlobalDomain>(),
"Function passed to interpolate does not model the Callable<GlobalCoordinate> concept");
204 auto&& gridView = basis.gridView();
206 auto&& bitVector = makeHierarchicVectorForMultiIndex<typename B::MultiIndex>(bv);
207 auto&& vector = makeHierarchicVectorForMultiIndex<typename B::MultiIndex>(coeff);
218 auto localView = basis.localView();
219 auto localIndexSet = basis.localIndexSet();
221 for (
const auto& e : elements(gridView))
224 localIndexSet.bind(localView);
227 auto&& subTree = TypeTree::child(localView.tree(),treePath);
229 Imp::LocalInterpolateVisitor<B, Tree, NTRE, decltype(vector), decltype(localF), decltype(bitVector)> localInterpolateVisitor(basis, vector, bitVector, localF, localIndexSet, nodeToRangeEntry);
230 TypeTree::applyToTree(subTree,localInterpolateVisitor);
235 template <
class B,
class... TreeIndices,
class C,
class F,
class BV>
236 void interpolateTreeSubset(
const B& basis,
const TypeTree::HybridTreePath<TreeIndices...>& treePath, C&& coeff,
const F& f,
const BV& bitVector)
242 template <
class B,
class... TreeIndices,
class NTRE,
class C,
class F>
243 void interpolateTree(
const B& basis,
const TypeTree::HybridTreePath<TreeIndices...>& treePath, C&& coeff,
const F& f,
const NTRE& nodeToRangeEntry)
249 template <
class B,
class... TreeIndices,
class C,
class F>
250 void interpolateTree(
const B& basis,
const TypeTree::HybridTreePath<TreeIndices...>& treePath, C&& coeff,
const F& f)
273 template <
class B,
class... TreeIndices,
class C,
class F,
class BV>
274 void interpolate(
const B& basis,
const TypeTree::HybridTreePath<TreeIndices...>& treePath, C&& coeff,
const F& f,
const BV& bitVector)
282 constexpr std::true_type isHybridTreePath(
const TypeTree::HybridTreePath<T...>&)
286 constexpr std::false_type isHybridTreePath(
const T&)
290 constexpr
auto isHybridTreePath() -> decltype(isHybridTreePath(std::declval<std::decay_t<T>>()))
311 template <
class B,
class C,
class F,
class BV,
312 std::enable_if_t<not Imp::isHybridTreePath<C>(),
int> = 0>
313 void interpolate(
const B& basis, C&& coeff,
const F& f,
const BV& bitVector)
315 auto root = Dune::TypeTree::hybridTreePath();
335 template <
class B,
class C,
class F>
338 interpolate (basis, Dune::TypeTree::hybridTreePath(), coeff, f, Imp::AllTrueBitSetVector());
356 template <
class B,
class... TreeIndices,
class C,
class F>
357 void interpolate(
const B& basis,
const TypeTree::HybridTreePath<TreeIndices...>& treePath, C&& coeff,
const F& f)
359 interpolate (basis, treePath, coeff, f, Imp::AllTrueBitSetVector());
365 #endif // DUNE_FUNCTIONS_FUNCTIONSPACEBASES_INTERPOLATE_HH void localFunction(DiscreteGlobalBasisFunction< TT... > &&t)=delete
Construction of local functions from a temporary DiscreteGlobalBasisFunction (forbidden) ...
static auto size(VV &&v) -> decltype(v.size())
Definition: flatvectorbackend.hh:41
std::decay< F >::type makeGridViewFunction(F &&f, const GridView &gridView)
Construct a function modeling GridViewFunction from function and grid view.
Definition: gridviewfunction.hh:68
void interpolateTreeSubset(const B &basis, const TypeTree::HybridTreePath< TreeIndices... > &treePath, C &&coeff, const F &f, const NTRE &nodeToRangeEntry, const BV &bv)
Interpolate given function in discrete function space.
Definition: interpolate.hh:193
static auto getEntry(VV &&v, const Index &i) -> decltype(v[i])
Definition: flatvectorbackend.hh:25
void interpolate(const B &basis, const TypeTree::HybridTreePath< TreeIndices... > &treePath, C &&coeff, const F &f, const BV &bitVector)
Interpolate given function in discrete function space.
Definition: interpolate.hh:274
void interpolateTree(const B &basis, const TypeTree::HybridTreePath< TreeIndices... > &treePath, C &&coeff, const F &f, const NTRE &nodeToRangeEntry)
Definition: interpolate.hh:243
DefaultNodeToRangeMap< Tree > makeDefaultNodeToRangeMap(const Tree &tree)
Definition: defaultnodetorangemap.hh:106
Definition: polynomial.hh:7
SizeInfo< Basis > sizeInfo(const Basis &basis)
Definition: sizeinfo.hh:69
Definition: functionfromcallable.hh:20