3 #include <amdis/Output.hpp> 4 #include <amdis/common/DerivativeTraits.hpp> 5 #include <amdis/common/FieldMatVec.hpp> 6 #include <amdis/common/Math.hpp> 7 #include <amdis/common/RecursiveForEach.hpp> 8 #include <amdis/functions/HierarchicNodeToRangeMap.hpp> 9 #include <amdis/functions/NodeIndices.hpp> 10 #include <amdis/functions/Order.hpp> 11 #include <amdis/typetree/FiniteElementType.hpp> 12 #include <amdis/utility/LocalToGlobalAdapter.hpp> 14 #include <dune/common/ftraits.hh> 20 template <
class Coeff,
class LocalView,
class LocalCoeff,
21 class = decltype(std::declval<Coeff>().gather(std::declval<LocalView>(), std::declval<LocalCoeff&>()))>
22 void gather(Coeff
const& coeff, LocalView
const& localView, LocalCoeff& localCoeff, Dune::PriorityTag<2>)
24 coeff.gather(localView, localCoeff);
28 template <
class Coeff,
class LocalView,
class LocalCoeff,
29 class = decltype(std::declval<Coeff>()[std::declval<typename LocalView::MultiIndex>()])>
30 void gather(Coeff
const& coeff, LocalView
const& localView, LocalCoeff& localCoeff, Dune::PriorityTag<1>)
32 localCoeff.resize(localView.size());
33 auto it = localCoeff.begin();
39 template <
class Coeff,
class LocalView,
class LocalCoeff,
40 class = decltype(std::declval<Coeff>()[std::integral_constant<std::size_t,0>{}])>
41 void gather(Coeff
const& coeff, LocalView
const& localView, LocalCoeff& localCoeff, Dune::PriorityTag<0>)
43 auto backend = Dune::Functions::istlVectorBackend(coeff);
44 localCoeff.resize(localView.size());
45 auto it = localCoeff.begin();
53 template <
class C,
class GB,
class TP,
class R>
57 using DomainType =
typename DiscreteFunction::Domain;
58 using RangeType =
typename DiscreteFunction::Range;
64 using Domain =
typename EntitySet::LocalCoordinate;
65 using Range = DerivativeRange<Type>;
67 enum { hasDerivative = std::is_same<Type, tag::value>::value };
71 using Element =
typename EntitySet::Element;
72 using Geometry =
typename Element::Geometry;
78 LocalFunction(std::shared_ptr<LV> localView, TP
const& treePath, C
const& coefficients, Type type)
79 : localView_(
std::move(localView))
81 , coefficients_(coefficients)
86 void bind(Element
const& element)
88 localView_->bind(element);
89 geometry_.emplace(element.geometry());
90 Impl::gather(coefficients_, *localView_, localCoefficients_, Dune::PriorityTag<4>{});
103 Range operator() (
const Domain& local)
const 106 if constexpr (std::is_same_v<Type, tag::value>)
107 return evaluateFunction(local);
108 else if constexpr (std::is_same_v<Type, tag::gradient>)
109 return evaluateJacobian(local);
110 else if constexpr (std::is_same_v<Type, tag::partial>)
111 return evaluatePartialDerivative(local);
112 else if constexpr (std::is_same_v<Type, tag::divergence>) {
113 static constexpr
bool isVectorNode
114 = SubTree::isPower && SubTree::degree() == GridView::dimensionworld;
115 static_assert(isVectorNode);
117 return evaluateDivergence(local);
124 template <
class DerType>
127 static_assert(std::is_same_v<Type, tag::value>);
128 return {localView_, treePath_, coefficients_, type};
143 return localView_->element();
148 template <
class LocalCoordinate>
149 auto evaluateFunction (
const LocalCoordinate& local)
const;
151 template <
class LocalCoordinate>
152 auto evaluateJacobian (
const LocalCoordinate& local)
const;
154 template <
class LocalCoordinate>
155 auto evaluateDivergence (
const LocalCoordinate& local)
const;
157 template <
class LocalCoordinate>
158 auto evaluatePartialDerivative (
const LocalCoordinate& local)
const;
161 Geometry
const& geometry()
const 168 SubTree
const& subTree()
const 170 return Dune::TypeTree::child(localView_->tree(), treePath_);
176 class = decltype(std::declval<LV>().treeCache())>
177 auto const& subTreeCache(Dune::PriorityTag<1>)
const 179 return Dune::TypeTree::child(localView_->treeCache(), treePath_);
184 auto subTreeCache(Dune::PriorityTag<0>)
const 189 decltype(
auto) subTreeCache()
const 191 return subTreeCache<LocalView>(Dune::PriorityTag<2>());
195 std::shared_ptr<LocalView> localView_;
197 Coefficients
const& coefficients_;
201 std::optional<Geometry> geometry_;
202 std::vector<Coefficient> localCoefficients_;
207 template <
class C,
class GB,
class TP,
class R>
208 template <
class Type>
209 template <
class LocalCoordinate>
214 Recursive::forEach(y, [](
auto& v) { v = 0; });
216 Traversal::forEachLeafNode(this->subTreeCache(), [&](
auto const& node,
auto const& tp)
218 auto const& shapeFunctionValues = node.localBasisValuesAt(local);
219 std::size_t size = node.finiteElement().size();
222 auto re = Dune::Functions::flatVectorView(nodeToRangeMap_(node, tp, y));
224 for (std::size_t i = 0; i < size; ++i) {
226 auto c = Dune::Functions::flatVectorView(localCoefficients_[node.localIndex(i)]);
229 auto v = Dune::Functions::flatVectorView(shapeFunctionValues[i]);
231 std::size_t dimC = c.size();
232 std::size_t dimV = v.size();
233 assert(dimC*dimV == std::size_t(re.size()));
234 for(std::size_t j = 0; j < dimC; ++j) {
236 for(std::size_t k = 0; k < dimV; ++k)
237 re[j*dimV + k] += c_j*v[k];
246 template <
class C,
class GB,
class TP,
class R>
247 template <
class Type>
248 template <
class LocalCoordinate>
253 Recursive::forEach(dy, [](
auto& v) { v = 0; });
255 Traversal::forEachLeafNode(this->subTreeCache(), [&](
auto const& node,
auto const& tp)
258 auto const& gradients = localBasis.
gradientsAt(local);
261 auto re = Dune::Functions::flatVectorView(nodeToRangeMap_(node, tp, dy));
263 for (std::size_t i = 0; i < localBasis.
size(); ++i) {
265 auto c = Dune::Functions::flatVectorView(localCoefficients_[node.localIndex(i)]);
268 auto grad = Dune::Functions::flatVectorView(gradients[i]);
270 std::size_t dimC = c.size();
271 std::size_t dimV = grad.size();
272 assert(dimC*dimV == std::size_t(re.size()));
273 for(std::size_t j = 0; j < dimC; ++j) {
275 for(std::size_t k = 0; k < dimV; ++k)
276 re[j*dimV + k] += c_j*grad[k];
285 template <
class C,
class GB,
class TP,
class R>
286 template <
class Type>
287 template <
class LocalCoordinate>
291 static_assert(static_size_v<Range> == 1,
"Divergence implemented for scalar coefficients only.");
294 Recursive::forEach(dy, [](
auto& v) { v = 0; });
296 auto&& node = this->subTreeCache();
299 auto const& gradients = localBasis.
gradientsAt(local);
301 auto re = Dune::Functions::flatVectorView(dy);
302 assert(
int(re.size()) == 1);
303 for (std::size_t i = 0; i < localBasis.size(); ++i) {
304 auto grad = Dune::Functions::flatVectorView(gradients[i]);
306 assert(
int(grad.size()) == GridView::dimensionworld);
307 for (std::size_t j = 0; j < GridView::dimensionworld; ++j)
308 re[0] += localCoefficients_[node.child(j).localIndex(i)] * grad[j];
315 template <
class C,
class GB,
class TP,
class R>
316 template <
class Type>
317 template <
class LocalCoordinate>
321 std::size_t comp = this->type_.comp;
324 Recursive::forEach(dy, [](
auto& v) { v = 0; });
326 Traversal::forEachLeafNode(this->subTreeCache(), [&](
auto const& node,
auto const& tp)
329 auto const& partial = localBasis.
partialsAt(local, comp);
332 auto re = Dune::Functions::flatVectorView(nodeToRangeMap_(node, tp, dy));
334 for (std::size_t i = 0; i < localBasis.
size(); ++i) {
336 auto c = Dune::Functions::flatVectorView(localCoefficients_[node.localIndex(i)]);
339 auto d_comp = Dune::Functions::flatVectorView(partial[i]);
341 std::size_t dimC = c.size();
342 std::size_t dimV = d_comp.size();
343 assert(dimC*dimV == std::size_t(re.size()));
344 for(std::size_t j = 0; j < dimC; ++j) {
346 for(std::size_t k = 0; k < dimV; ++k)
347 re[j*dimV + k] += c_j*d_comp[k];
LocalFunction(std::shared_ptr< LV > localView, TP const &treePath, C const &coefficients, Type type)
Constructor. Stores a copy of the DiscreteFunction.
Definition: DiscreteLocalFunction.inc.hpp:78
Convert a simple (scalar) local basis into a global basis.
Definition: LocalToGlobalAdapter.hpp:64
auto makeNodeCache(Node const &node)
Construct a new local-basis cache from a basis-node.
Definition: NodeCache.hpp:35
Definition: FieldMatVec.hpp:12
AMDiS::LocalView< Self > LocalView
Type of the local view on the restriction of the basis to a single element.
Definition: GlobalBasis.hpp:70
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
A mutable view on the subspace of a DOFVector,.
Definition: DiscreteFunction.hpp:36
std::size_t size() const
Return the number of local basis functions.
Definition: LocalToGlobalAdapter.hpp:96
A simple node to range map using the nested tree indices.
Definition: HierarchicNodeToRangeMap.hpp:23
auto order(F const &f) -> decltype(&F::operator(), f.order())
polynomial order of functions
Definition: Order.hpp:11
auto const & partialsAt(typename Traits::DomainLocal const &x, std::size_t comp) const
Definition: LocalToGlobalAdapter.hpp:176
LocalFunction< DerType > makeDerivative(DerType type) const
Create a LocalFunction representing the derivative.
Definition: DiscreteLocalFunction.inc.hpp:125
void bind(Element const &element)
Bind the LocalView to the element.
Definition: DiscreteLocalFunction.inc.hpp:86
Definition: DerivativeTraits.hpp:29
void unbind()
Unbind the LocalView from the element.
Definition: DiscreteLocalFunction.inc.hpp:95
auto order() const -> decltype(AMDiS::order(std::declval< SubTree >()))
The polynomial degree of the LocalFunctions basis-functions.
Definition: DiscreteLocalFunction.inc.hpp:132
Element const & localContext() const
Return the bound element.
Definition: DiscreteLocalFunction.inc.hpp:140
auto nodeIndices(LocalView const &localView, Node const &node)
Returns a range over (flat) DOF indices on a node, given by the localView.
Definition: NodeIndices.hpp:13
Definition: DiscreteLocalFunction.inc.hpp:55
auto const & gradientsAt(typename Traits::DomainLocal const &x) const
Definition: LocalToGlobalAdapter.hpp:147