10 #include <unordered_map> 14 #include <dune/common/fvector.hh> 15 #include <dune/common/hash.hh> 17 #include <dune/grid/common/geometry.hh> 18 #include <dune/grid/common/rangegenerators.hh> 20 #include <amdis/Output.hpp> 21 #include <amdis/common/ConcurrentCache.hpp> 22 #include <amdis/functions/FunctionFromCallable.hpp> 23 #include <amdis/operations/Assigner.hpp> 24 #include <amdis/typetree/Traversal.hpp> 32 template <
class LocalCoord>
33 std::size_t operator()(LocalCoord
const& coord)
const 36 for (std::size_t i = 0; i < coord.size(); ++i)
37 Dune::hash_combine(seed, coord[i]);
45 template <
class C,
class B>
46 DataTransfer<C,B>::DataTransfer(std::shared_ptr<B const> basis)
47 : basis_(
std::move(basis))
48 , mapper_(basis_->gridView().grid(),
Dune::mcmgElementLayout())
53 template <
class C,
class B>
57 GridView gv = basis_->gridView();
58 LocalView lv = basis_->localView();
59 auto const& idSet = gv.grid().localIdSet();
61 Traversal::forEachLeafNode(lv.tree(), [&](
auto const& node,
auto const& tp) {
62 nodeDataTransfer_[tp].preAdaptInit(lv, coeff, node);
66 persistentContainer_.clear();
67 for (
const auto& e : elements(gv))
69 auto it = persistentContainer_.emplace(idSet.id(e), TypeTree::treeContainer<NodeElementData,true>(lv.tree()));
72 auto& treeContainer = it.first->second;
73 Traversal::forEachLeafNode(lv.tree(), [&](
auto const& node,
auto const& tp) {
74 nodeDataTransfer_[tp].cacheLocal(treeContainer[tp]);
82 auto maxLevel = gv.grid().maxLevel();
84 typename Grid::ctype
const checkInsideTolerance = sqrt(std::numeric_limits<typename Grid::ctype>::epsilon());
85 for (
auto const& e : elements(gv,
typename C::Traits::PartitionSet{}))
88 while (father.mightVanish() && father.hasFather())
90 father = father.father();
91 auto it = persistentContainer_.emplace(idSet.id(father), TypeTree::treeContainer<NodeElementData,true>(lv.tree()));
95 auto& treeContainer = it.first->second;
96 auto geo = father.geometry();
98 bool restrictLocalCompleted =
false;
99 auto hItEnd = father.hend(maxLevel);
100 for (
auto hIt = father.hbegin(maxLevel); hIt != hItEnd; ++hIt) {
104 auto const& child = *hIt;
105 auto search = persistentContainer_.find(idSet.id(child));
106 assert(search != persistentContainer_.end());
107 auto const& childContainer = search->second;
110 auto const& childGeo = child.geometry();
111 auto refTypeId = Dune::referenceElement(childGeo).type().id();
113 using BoolCoordPair = std::pair<bool, LocalCoordinate>;
114 using CacheImp = std::unordered_map<LocalCoordinate, BoolCoordPair, Impl::CoordHasher>;
119 auto xInChild = [&](LocalCoordinate
const& x) -> BoolCoordPair {
120 LocalCoordinate local = childGeo.local(geo.global(x));
123 bool isInside = Dune::Geo::Impl::checkInside(refTypeId, Geometry::coorddimension, local, checkInsideTolerance);
124 return BoolCoordPair(isInside, std::move(local));
127 ChildCache childCache;
128 auto xInChildCached = [&](LocalCoordinate
const& x) -> BoolCoordPair {
129 return childCache.get(x, [&](LocalCoordinate
const& x) {
return xInChild(x); });
132 restrictLocalCompleted =
true;
133 Traversal::forEachLeafNode(lv.tree(), [&](
auto const& node,
auto const& tp) {
134 restrictLocalCompleted &=
135 nodeDataTransfer_[tp].restrictLocal(father, treeContainer[tp], xInChildCached,
136 childContainer[tp], init);
141 assert(restrictLocalCompleted);
148 template <
class C,
class B>
155 if (persistentContainer_.empty())
158 GridView gv = basis_->gridView();
159 LocalView lv = basis_->localView();
160 auto const& idSet = gv.grid().localIdSet();
161 Traversal::forEachLeafNode(lv.tree(), [&](
auto const& node,
auto const& tp) {
162 nodeDataTransfer_[tp].adaptInit(lv, coeff, node);
166 std::vector<bool> finished(mapper_.size(),
false);
167 for (
const auto& e : elements(gv,
typename C::Traits::PartitionSet{}))
169 auto index = mapper_.index(e);
173 auto it = persistentContainer_.find(idSet.id(e));
176 if (it != persistentContainer_.end()) {
178 auto const& treeContainer = it->second;
179 Traversal::forEachLeafNode(lv.tree(), [&](
auto const& node,
auto const& tp) {
180 nodeDataTransfer_[tp].copyLocal(treeContainer[tp]);
182 finished[index] =
true;
188 while (father.hasFather() && father.isNew())
189 father = father.father();
191 auto maxLevel = gv.grid().maxLevel();
192 auto fatherGeo = father.geometry();
195 auto father_it = persistentContainer_.find(idSet.id(father));
196 assert(father_it != persistentContainer_.end());
197 auto const& treeContainer = father_it->second;
199 auto hItEnd = father.hend(maxLevel);
200 for (
auto hIt = father.hbegin(maxLevel); hIt != hItEnd; ++hIt) {
204 auto const& child = *hIt;
208 auto xInFather = [&fatherGeo, childGeo = child.geometry()]
209 (LocalCoordinate
const& x) -> LocalCoordinate
211 return fatherGeo.local(childGeo.global(x));
214 Traversal::forEachLeafNode(lv.tree(), [&](
auto const& node,
auto const& tp) {
215 nodeDataTransfer_[tp].prolongLocal(father, treeContainer[tp], xInFather, init);
218 finished[mapper_.index(child)] =
true;
227 template <
class C,
class B>
230 persistentContainer_.clear();
237 template <
class Node,
class Container,
class Basis>
240 using T =
typename Container::value_type;
242 using Element =
typename Node::Element;
244 using LocalBasis =
typename Node::FiniteElement::Traits::LocalBasisType;
245 using LBRangeType =
typename LocalBasis::Traits::RangeType;
246 using LocalInterpolation =
typename Node::FiniteElement::Traits::LocalBasisType;
247 using LIDomainType =
typename LocalInterpolation::Traits::DomainType;
248 using LIRangeType =
typename LocalInterpolation::Traits::RangeType;
251 using NodeElementData = std::vector<T>;
257 void preAdaptInit(LocalView
const& lv, Container
const& coeff, Node
const& node)
261 fatherNode_ = std::make_unique<Node>(node);
262 constCoeff_ = &coeff;
273 constCoeff_->gather(*lv_, *node_, dofs);
291 template <
class Trafo>
292 bool restrictLocal(Element
const& father, NodeElementData& fatherDOFs, Trafo
const& trafo,
293 NodeElementData
const& childDOFs,
bool init);
297 void adaptInit(LocalView
const& lv, Container& coeff, Node
const& node)
301 fatherNode_ = std::make_unique<Node>(node);
324 template <
class Trafo>
325 void prolongLocal(Element
const& father, NodeElementData
const& fatherDOFs,
326 Trafo
const& trafo,
bool init);
329 LocalView
const* lv_ =
nullptr;
330 Node
const* node_ =
nullptr;
331 std::unique_ptr<Node> fatherNode_;
332 Container
const* constCoeff_ =
nullptr;
333 Container* coeff_ =
nullptr;
334 std::vector<bool> finishedDOFs_;
335 NodeElementData fatherDOFsTemp_;
339 template <
class N,
class C,
class B>
340 template <
class Trafo>
342 restrictLocal(Element
const& father, NodeElementData& fatherDOFs, Trafo
const& trafo,
343 NodeElementData
const& childDOFs,
bool init)
345 auto& fatherNode = *fatherNode_;
346 std::size_t currentDOF = 0;
350 bindTree(fatherNode, father);
352 auto const& childNode = *node_;
353 auto const& childFE = childNode.finiteElement();
354 auto const& fatherFE = fatherNode.finiteElement();
357 finishedDOFs_.assign(fatherFE.size(),
false);
358 fatherDOFsTemp_.assign(fatherFE.size(), 0);
361 auto evalLeaf = [&](LIDomainType
const& x) -> LIRangeType {
362 if (!finishedDOFs_[currentDOF])
364 auto const& insideLocal = trafo(x);
365 bool isInside = insideLocal.first;
368 auto const& local = insideLocal.second;
369 thread_local std::vector<LBRangeType> shapeValues;
370 childFE.localBasis().evaluateFunction(local, shapeValues);
372 assert(childDOFs.size() == shapeValues.size());
375 for (std::size_t i = 0; i < shapeValues.size(); ++i)
376 y += shapeValues[i] * childDOFs[i];
378 fatherDOFsTemp_[currentDOF] = T(y);
379 finishedDOFs_[currentDOF++] =
true;
383 return fatherDOFsTemp_[currentDOF++];
386 auto evalLeafFct = functionFromCallable<LIRangeType(LIDomainType)>(evalLeaf);
387 fatherFE.localInterpolation().interpolate(evalLeafFct, fatherDOFs);
390 return std::accumulate(finishedDOFs_.begin(), finishedDOFs_.end(),
true,
391 std::logical_and<bool>());
395 template <
class N,
class C,
class B>
396 template <
class Trafo>
399 Trafo
const& trafo,
bool init)
401 auto& fatherNode = *fatherNode_;
405 bindTree(fatherNode, father);
407 auto const& childNode = *node_;
410 auto evalFather = [&](LIDomainType
const& x) -> LIRangeType
412 thread_local std::vector<LBRangeType> shapeValues;
413 fatherNode.finiteElement().localBasis().evaluateFunction(trafo(x), shapeValues);
414 assert(shapeValues.size() == fatherDOFs.size());
417 for (std::size_t i = 0; i < shapeValues.size(); ++i)
418 y += shapeValues[i] * fatherDOFs[i];
423 auto const& childFE = childNode.finiteElement();
424 thread_local std::vector<T> childDOFs;
425 auto evalFatherFct = functionFromCallable<LIRangeType(LIDomainType)>(evalFather);
426 childFE.localInterpolation().interpolate(evalFatherFct, childDOFs);
void copyLocal(NodeElementData const &dofs) const
Copy already existing data to element bound to node_.
Definition: DataTransfer.inc.hpp:307
void adapt(Container &coeff) override
Definition: DataTransfer.inc.hpp:149
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:398
Definition: AdaptiveGrid.hpp:373
Definition: FieldMatVec.hpp:12
void postAdapt(Container &coeff) override
Definition: DataTransfer.inc.hpp:228
Definition: DataTransfer.hpp:74
Contains all classes needed for solving linear and non linear equation systems.
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:297
void preAdapt(Container const &coeff, bool mightCoarsen) override
Definition: DataTransfer.inc.hpp:55
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:342
void cacheLocal(NodeElementData &dofs) const
Cache data on the element bound to node_.
Definition: DataTransfer.inc.hpp:271
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:257
void init(int &argc, char **&argv, std::string const &initFileName="")
Initialized the Environment for MPI.
Definition: AMDiS.hpp:29
Definition: Assigner.hpp:7