10 #include <unordered_map> 14 #include <dune/common/fvector.hh> 15 #include <dune/common/hash.hh> 16 #include <dune/common/version.hh> 18 #include <dune/grid/common/geometry.hh> 19 #include <dune/grid/common/rangegenerators.hh> 21 #include <amdis/Output.hpp> 22 #include <amdis/common/ConcurrentCache.hpp> 23 #include <amdis/functions/FunctionFromCallable.hpp> 24 #include <amdis/operations/Assigner.hpp> 25 #include <amdis/typetree/Traversal.hpp> 33 template <
class LocalCoord>
34 std::size_t operator()(LocalCoord
const& coord)
const 37 for (std::size_t i = 0; i < coord.size(); ++i)
38 Dune::hash_combine(seed, coord[i]);
46 template <
class C,
class B>
47 DataTransfer<C,B>::DataTransfer(std::shared_ptr<B const> basis)
48 : basis_(
std::move(basis))
49 , mapper_(basis_->gridView(),
Dune::mcmgElementLayout())
54 template <
class C,
class B>
58 GridView gv = basis_->gridView();
59 LocalView lv = basis_->localView();
60 auto const& idSet = gv.grid().localIdSet();
62 Traversal::forEachLeafNode(lv.tree(), [&](
auto const& node,
auto const& tp) {
63 nodeDataTransfer_[tp].preAdaptInit(lv, coeff, node);
67 persistentContainer_.clear();
68 for (
const auto& e : elements(gv))
70 auto it = persistentContainer_.emplace(idSet.id(e), TypeTree::treeContainer<NodeElementData,true>(lv.tree()));
73 auto& treeContainer = it.first->second;
74 Traversal::forEachLeafNode(lv.tree(), [&](
auto const& node,
auto const& tp) {
75 nodeDataTransfer_[tp].cacheLocal(treeContainer[tp]);
83 auto maxLevel = gv.grid().maxLevel();
85 typename Grid::ctype
const checkInsideTolerance = sqrt(std::numeric_limits<typename Grid::ctype>::epsilon());
86 for (
auto const& e : elements(gv,
typename C::Traits::PartitionSet{}))
89 while (father.mightVanish() && father.hasFather())
91 father = father.father();
92 auto it = persistentContainer_.emplace(idSet.id(father), TypeTree::treeContainer<NodeElementData,true>(lv.tree()));
96 auto& treeContainer = it.first->second;
97 auto geo = father.geometry();
99 bool restrictLocalCompleted =
false;
100 auto hItEnd = father.hend(maxLevel);
101 for (
auto hIt = father.hbegin(maxLevel); hIt != hItEnd; ++hIt) {
105 auto const& child = *hIt;
106 auto search = persistentContainer_.find(idSet.id(child));
107 assert(search != persistentContainer_.end());
108 auto const& childContainer = search->second;
111 auto const& childGeo = child.geometry();
112 auto refTypeId = Dune::referenceElement(childGeo).type().id();
114 using BoolCoordPair = std::pair<bool, LocalCoordinate>;
115 using CacheImp = std::unordered_map<LocalCoordinate, BoolCoordPair, Impl::CoordHasher>;
120 auto xInChild = [&](LocalCoordinate
const& x) -> BoolCoordPair {
121 LocalCoordinate local = childGeo.local(geo.global(x));
124 bool isInside = Dune::Geo::Impl::checkInside(refTypeId, Geometry::mydimension, local, checkInsideTolerance);
125 return BoolCoordPair(isInside, std::move(local));
128 ChildCache childCache;
129 auto xInChildCached = [&](LocalCoordinate
const& x) -> BoolCoordPair {
130 return childCache.get(x, [&](LocalCoordinate
const& x) {
return xInChild(x); });
133 restrictLocalCompleted =
true;
134 Traversal::forEachLeafNode(lv.tree(), [&](
auto const& node,
auto const& tp) {
135 restrictLocalCompleted &=
136 nodeDataTransfer_[tp].restrictLocal(father, treeContainer[tp], xInChildCached,
137 childContainer[tp], init);
142 assert(restrictLocalCompleted);
149 template <
class C,
class B>
156 if (persistentContainer_.empty())
159 GridView gv = basis_->gridView();
160 LocalView lv = basis_->localView();
161 auto const& idSet = gv.grid().localIdSet();
162 Traversal::forEachLeafNode(lv.tree(), [&](
auto const& node,
auto const& tp) {
163 nodeDataTransfer_[tp].adaptInit(lv, coeff, node);
168 std::vector<bool> finished(mapper_.size(),
false);
169 for (
const auto& e : elements(gv,
typename C::Traits::PartitionSet{}))
171 auto index = mapper_.index(e);
175 auto it = persistentContainer_.find(idSet.id(e));
178 if (it != persistentContainer_.end()) {
180 auto const& treeContainer = it->second;
181 Traversal::forEachLeafNode(lv.tree(), [&](
auto const& node,
auto const& tp) {
182 nodeDataTransfer_[tp].copyLocal(treeContainer[tp]);
184 finished[index] =
true;
190 while (father.hasFather() && father.isNew())
191 father = father.father();
193 auto maxLevel = gv.grid().maxLevel();
194 auto fatherGeo = father.geometry();
197 auto father_it = persistentContainer_.find(idSet.id(father));
198 assert(father_it != persistentContainer_.end());
199 auto const& treeContainer = father_it->second;
201 auto hItEnd = father.hend(maxLevel);
202 for (
auto hIt = father.hbegin(maxLevel); hIt != hItEnd; ++hIt) {
206 auto const& child = *hIt;
210 auto xInFather = [&fatherGeo, childGeo = child.geometry()]
211 (LocalCoordinate
const& x) -> LocalCoordinate
213 return fatherGeo.local(childGeo.global(x));
216 Traversal::forEachLeafNode(lv.tree(), [&](
auto const& node,
auto const& tp) {
217 nodeDataTransfer_[tp].prolongLocal(father, treeContainer[tp], xInFather, init);
220 finished[mapper_.index(child)] =
true;
229 template <
class C,
class B>
232 persistentContainer_.clear();
239 template <
class Node,
class Container,
class Basis>
242 using T =
typename Container::value_type;
244 using Element =
typename Node::Element;
246 using LocalBasis =
typename Node::FiniteElement::Traits::LocalBasisType;
247 using LBRangeType =
typename LocalBasis::Traits::RangeType;
248 using LocalInterpolation =
typename Node::FiniteElement::Traits::LocalBasisType;
249 using LIDomainType =
typename LocalInterpolation::Traits::DomainType;
250 using LIRangeType =
typename LocalInterpolation::Traits::RangeType;
253 using NodeElementData = std::vector<T>;
259 void preAdaptInit(LocalView
const& lv, Container
const& coeff, Node
const& node)
263 fatherNode_ = std::make_unique<Node>(node);
264 constCoeff_ = &coeff;
275 constCoeff_->gather(*lv_, *node_, dofs);
293 template <
class Trafo>
294 bool restrictLocal(Element
const& father, NodeElementData& fatherDOFs, Trafo
const& trafo,
295 NodeElementData
const& childDOFs,
bool init);
299 void adaptInit(LocalView
const& lv, Container& coeff, Node
const& node)
303 fatherNode_ = std::make_unique<Node>(node);
326 template <
class Trafo>
327 void prolongLocal(Element
const& father, NodeElementData
const& fatherDOFs,
328 Trafo
const& trafo,
bool init);
331 LocalView
const* lv_ =
nullptr;
332 Node
const* node_ =
nullptr;
333 std::unique_ptr<Node> fatherNode_;
334 Container
const* constCoeff_ =
nullptr;
335 Container* coeff_ =
nullptr;
336 std::vector<bool> finishedDOFs_;
337 NodeElementData fatherDOFsTemp_;
341 template <
class N,
class C,
class B>
342 template <
class Trafo>
344 restrictLocal(Element
const& father, NodeElementData& fatherDOFs, Trafo
const& trafo,
345 NodeElementData
const& childDOFs,
bool init)
347 auto& fatherNode = *fatherNode_;
348 std::size_t currentDOF = 0;
352 bindTree(fatherNode, father);
354 auto const& childNode = *node_;
355 auto const& childFE = childNode.finiteElement();
356 auto const& fatherFE = fatherNode.finiteElement();
359 finishedDOFs_.assign(fatherFE.size(),
false);
360 fatherDOFsTemp_.assign(fatherFE.size(), 0);
363 auto evalLeaf = [&](LIDomainType
const& x) -> LIRangeType {
364 if (!finishedDOFs_[currentDOF])
366 auto const& insideLocal = trafo(x);
367 bool isInside = insideLocal.first;
370 auto const& local = insideLocal.second;
371 thread_local std::vector<LBRangeType> shapeValues;
372 childFE.localBasis().evaluateFunction(local, shapeValues);
374 assert(childDOFs.size() == shapeValues.size());
377 for (std::size_t i = 0; i < shapeValues.size(); ++i)
378 y += shapeValues[i] * childDOFs[i];
380 fatherDOFsTemp_[currentDOF] = T(y);
381 finishedDOFs_[currentDOF++] =
true;
385 return fatherDOFsTemp_[currentDOF++];
388 auto evalLeafFct = functionFromCallable<LIRangeType(LIDomainType)>(evalLeaf);
389 fatherFE.localInterpolation().interpolate(evalLeafFct, fatherDOFs);
392 return std::accumulate(finishedDOFs_.begin(), finishedDOFs_.end(),
true,
393 std::logical_and<bool>());
397 template <
class N,
class C,
class B>
398 template <
class Trafo>
401 Trafo
const& trafo,
bool init)
403 auto& fatherNode = *fatherNode_;
407 bindTree(fatherNode, father);
409 auto const& childNode = *node_;
412 auto evalFather = [&](LIDomainType
const& x) -> LIRangeType
414 thread_local std::vector<LBRangeType> shapeValues;
415 fatherNode.finiteElement().localBasis().evaluateFunction(trafo(x), shapeValues);
416 assert(shapeValues.size() == fatherDOFs.size());
419 for (std::size_t i = 0; i < shapeValues.size(); ++i)
420 y += shapeValues[i] * fatherDOFs[i];
425 auto const& childFE = childNode.finiteElement();
426 thread_local std::vector<T> childDOFs;
427 auto evalFatherFct = functionFromCallable<LIRangeType(LIDomainType)>(evalFather);
428 childFE.localInterpolation().interpolate(evalFatherFct, childDOFs);
void copyLocal(NodeElementData const &dofs) const
Copy already existing data to element bound to node_.
Definition: DataTransfer.inc.hpp:309
void adapt(Container &coeff) override
Definition: DataTransfer.inc.hpp:150
The class template ConcurrentCache describes an associative static container that allows the concurre...
Definition: ConcurrentCache.hpp:53
void prolongLocal(Element const &father, NodeElementData const &fatherDOFs, Trafo const &trafo, bool init)
Interpolate data from father onto the child element bound to node_ using the transformation trafo fro...
Definition: DataTransfer.inc.hpp:400
Definition: AdaptiveGrid.hpp:373
Definition: FieldMatVec.hpp:12
void postAdapt(Container &coeff) override
Definition: DataTransfer.inc.hpp:230
Definition: DataTransfer.hpp:74
Definition: AdaptBase.hpp:6
constexpr bool LocalView
A Dune::Functions::LocalView type.
Definition: Concepts.hpp:182
void adaptInit(LocalView const &lv, Container &coeff, Node const &node)
To be called once before copyLocal/prolongLocal are called within the adapt step. ...
Definition: DataTransfer.inc.hpp:299
void preAdapt(Container const &coeff, bool mightCoarsen) override
Definition: DataTransfer.inc.hpp:56
bool restrictLocal(Element const &father, NodeElementData &fatherDOFs, Trafo const &trafo, NodeElementData const &childDOFs, bool init)
Evaluate data on the child element bound to node_ and interpolate onto father entity using the coordi...
Definition: DataTransfer.inc.hpp:344
void cacheLocal(NodeElementData &dofs) const
Cache data on the element bound to node_.
Definition: DataTransfer.inc.hpp:273
void preAdaptInit(LocalView const &lv, Container const &coeff, Node const &node)
To be called once before cacheLocal/restrictLocal are called within the preAdapt step.
Definition: DataTransfer.inc.hpp:259
Definition: Assigner.hpp:7