7#include <amdis/DataTransfer.hpp>
8#include <amdis/functions/EntitySet.hpp>
9#include <amdis/gridfunctions/CoarsenedGridFunction.hpp>
10#include <amdis/gridfunctions/DiscreteFunction.hpp>
12#include <dune/grid/common/capabilities.hh>
32template <
class Basis,
class Container>
35 using CoefficientVector = std::vector<typename Container::value_type>;
36 using EntityId =
typename Basis::GridView::Grid::LocalIdSet::IdType;
37 using PersistentStorage = std::map<EntityId, CoefficientVector>;
38 using Domain =
typename Basis::GridView::template Codim<0>::Entity::Geometry::LocalCoordinate;
41 PersistentStorage persistentStorage_;
46 static_assert(Dune::Capabilities::hasSingleGeometryType<typename Basis::GridView::Grid>::v);
50 void preAdapt(Basis
const& basis, Container
const& coeff,
bool mightCoarsen)
52 CoefficientVector interpolationCoeff;
53 auto localView = basis.localView();
54 auto interpol = [&](
auto const& lf,
auto& localCoeff) {
55 localCoeff.resize(localView.size());
56 Traversal::forEachLeafNode(localView.tree(), [&](
auto const& node,
auto const& tp)
58 auto lf_at_tp = HierarchicNodeWrapper{tp,lf};
59 node.finiteElement().localInterpolation().interpolate(lf_at_tp, interpolationCoeff);
60 assert(node.size() == interpolationCoeff.size());
61 for (std::size_t i = 0; i < node.size(); ++i)
62 localCoeff[node.localIndex(i)] = interpolationCoeff[i];
67 CoarsenedGridFunction c_gf{basis.gridView(), coeffFct};
68 auto c_lf = localFunction(c_gf);
69 auto const& idSet = basis.gridView().grid().localIdSet();
70 for (
auto const& e : entitySet(basis))
75 coeff.gather(localView, persistentStorage_[idSet.id(e)]);
79 if (e.mightVanish()) {
81 while (father.hasFather()) {
82 father = father.father();
83 auto [it,inserted] = persistentStorage_.try_emplace(idSet.id(father));
86 interpol(c_lf, it->second);
89 if (!inserted || !father.mightVanish())
98 void adapt(Basis
const& basis, Container& coeff)
100 coeff.init(basis,
false);
101 auto localView = basis.localView();
103 CoefficientVector interpolationCoeff;
104 auto interpol = [&](
auto const& local,
auto const& localCoeff) {
105 Traversal::forEachLeafNode(localView.tree(), [&](
auto const& node,
auto const& tp)
107 auto const& localFE = node.finiteElement();
108 auto const& localBasis = localFE.localBasis();
109 auto const& localInterpolation = localFE.localInterpolation();
111 using LocalBasis = TYPEOF(localBasis);
112 using RangeType = typename LocalBasis::Traits::RangeType;
113 std::vector<RangeType> shapeValues;
114 localInterpolation.interpolate([&](auto const& x) {
115 localBasis.evaluateFunction(local(x), shapeValues);
116 auto range = localCoeff[node.localIndex(0)] * shapeValues[0];
117 assert(node.size() == shapeValues.size());
118 for (std::size_t i = 1; i < shapeValues.size(); ++i)
119 range.axpy(localCoeff[node.localIndex(i)], shapeValues[i]);
121 }, interpolationCoeff);
126 auto const& idSet = basis.gridView().grid().localIdSet();
127 for (
const auto& e : entitySet(basis))
133 std::function<Domain(Domain
const&)> local = [](Domain
const& x) {
return x; };
134 [[maybe_unused]]
bool found =
false;
135 while (father.hasFather()) {
137 local = [local, geo=father.geometryInFather()](Domain
const& x) {
138 return geo.global(local(x));
140 father = father.father();
143 auto it = persistentStorage_.find(idSet.id(father));
144 if (it != persistentStorage_.end()) {
145 interpol(local, it->second);
153 auto it = persistentStorage_.find(idSet.id(e));
154 test_exit(it != persistentStorage_.end(),
155 "Element should be in persistent storage but isn't.");
156 coeff.scatter(localView, it->second, Assigner::assign{});
166 persistentStorage_.clear();
174 template <
class Basis,
class T,
template <
class>
class Impl,
175 REQUIRES(Concepts::GlobalBasis<Basis>)>
The base class for data transfer classes.
Definition: DataTransfer.hpp:68
A mutable view on the subspace of a DOFVector,.
Definition: DiscreteFunction.hpp:45
Definition: SimpleDataTransfer.hpp:34
void preAdapt(Basis const &basis, Container const &coeff, bool mightCoarsen)
Create a persistent storage of the coefficients.
Definition: SimpleDataTransfer.hpp:50
void postAdapt(Container &coeff)
Is called after the preAdapt step of the grid.
Definition: SimpleDataTransfer.hpp:164
void adapt(Basis const &basis, Container &coeff)
Transfer the content of the coefficients to the new basis.
Definition: SimpleDataTransfer.hpp:98
The basic container that stores a base vector and a corresponding basis.
Definition: VectorFacade.hpp:39
Definition: Assigner.hpp:8
Definition: DataTransfer.hpp:133
Base tag type for all data transfer tags.
Definition: DataTransfer.hpp:121
Tag that referrs to the SimpleDataTransfer.
Definition: SimpleDataTransfer.hpp:27