AMDiS 2.10
The Adaptive Multi-Dimensional Simulation Toolbox
SimpleInterpolator.hpp
1#pragma once
2
3#include <vector>
4
5#include <amdis/common/FakeContainer.hpp>
6#include <amdis/common/FieldMatVec.hpp>
7#include <amdis/functions/EntitySet.hpp>
8#include <amdis/functions/HierarchicNodeToRangeMap.hpp>
9#include <amdis/functions/NodeIndices.hpp>
10#include <amdis/linearalgebra/Traits.hpp>
11#include <amdis/operations/Assigner.hpp>
12#include <amdis/typetree/Traversal.hpp>
13
14namespace AMDiS
15{
16 namespace tag
17 {
18 struct assign {};
19 }
20
21 template <class Basis,
22 class TreePath = Dune::TypeTree::HybridTreePath<>>
24 {
25 template <class Coeff, class GridFct, class BitVector>
26 void operator()(Coeff& coeff, GridFct const& gf, BitVector const& bitVec) const
27 {
28 // Obtain a local view of the gridFunction
29 auto lf = localFunction(gf);
30 auto localView = basis_.localView();
31
32 std::vector<typename Coeff::value_type> localCoeff;
33
34 for (const auto& e : entitySet(basis_))
35 {
36 localView.bind(e);
37 lf.bind(e);
38
39 auto&& subTree = Dune::TypeTree::child(localView.tree(),treePath_);
40 Traversal::forEachLeafNode(subTree, [&](auto const& node, auto const& tp)
41 {
42 auto bitVecRange = mappedRangeView(Dune::range(node.size()), [&](std::size_t i) -> bool {
43 return bitVec[localView.index(node.localIndex(i))];
44 });
45
46 // check whether there is a local contribution
47 std::vector<bool> mask(bitVecRange.begin(), bitVecRange.end());
48 if (std::all_of(mask.begin(), mask.end(), [](bool m) { return !m; }))
49 return;
50
51 // compute local interpolation coefficients
52 node.finiteElement().localInterpolation().interpolate(HierarchicNodeWrapper{tp,lf}, localCoeff);
53
54 // assign local coefficients to global vector
55 coeff.scatter(localView, node, localCoeff, mask, Assigner::assign{});
56 });
57
58 lf.unbind();
59 }
60 coeff.finish();
61 }
62
63 template <class Coeff, class GridFct>
64 void operator()(Coeff& coeff, GridFct const& gf) const
65 {
66 (*this)(coeff, gf, FakeContainer<bool,true>{});
67 }
68
69 SimpleInterpolator(Basis const& basis, TreePath treePath = {})
70 : basis_{basis}
71 , treePath_{treePath}
72 {}
73
74 Basis const& basis_;
75 TreePath treePath_ = {};
76 };
77
78
79 template <class Tag>
81
82 template <>
83 struct InterpolatorFactory<tag::assign>
84 {
85 template <class Basis,
86 class TreePath = Dune::TypeTree::HybridTreePath<>>
87 static auto create(Basis const& basis, TreePath treePath = {})
88 {
89 return SimpleInterpolator<Basis,TreePath>{basis, treePath};
90 }
91 };
92
93} // end namespace AMDiS
A container-like data-structure not storing anything and with empty implementations in many container...
Definition: FakeContainer.hpp:36
Definition: Assigner.hpp:8
Definition: HierarchicNodeToRangeMap.hpp:42
Definition: SimpleInterpolator.hpp:80
Definition: SimpleInterpolator.hpp:24
Definition: SimpleInterpolator.hpp:18