AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
DiscreteFunction.inc.hpp
1 #pragma once
2 
3 #include <type_traits>
4 
5 #include <dune/grid/utility/hierarchicsearch.hh>
6 
7 #include <amdis/functions/Interpolate.hpp>
8 #include <amdis/gridfunctions/GridFunction.hpp>
9 #include <amdis/linearalgebra/VectorFacade.hpp>
10 
11 namespace AMDiS {
12 
13 // Evaluate DiscreteFunction in global coordinates
14 template <class C, class GB, class TP, class R>
15 typename DiscreteFunction<C const,GB,TP,R>::Range DiscreteFunction<C const,GB,TP,R>::
16  operator()(Domain const& x) const
17 {
18  using Grid = typename GlobalBasis::GridView::Grid;
19  using IS = typename GlobalBasis::GridView::IndexSet;
20 
21  auto const& gv = this->basis().gridView();
22  Dune::HierarchicSearch<Grid,IS> hsearch{gv.grid(), gv.indexSet()};
23 
24  auto element = hsearch.findEntity(x);
25  auto geometry = element.geometry();
26  auto localFct = localFunction(*this);
27  localFct.bind(element);
28  return localFct(geometry.local(x));
29 }
30 
31 
32 // Interpolation of GridFunction to DOFVector
33 template <class C, class GB, class TP, class R>
34  template <class Expr, class Tag>
36  interpolate_noalias(Expr&& expr, Tag strategy)
37 {
38  auto const& basis = this->basis();
39  auto const& treePath = this->treePath();
40 
41  auto&& gf = makeGridFunction(FWD(expr), basis.gridView());
42 
43  if (std::is_same_v<Tag, tag::average>) {
44  VectorType_t<short,Coefficients> counter(basis);
45  AMDiS::interpolate(basis, coefficients(), gf, treePath, counter);
46 
47  coefficients().forEach([&counter](auto dof, auto& coeff)
48  {
49  coeff /= std::max(double(counter.at(dof)), 1.0);
50  });
51  } else {
52  AMDiS::interpolate(basis, coefficients(), gf, treePath);
53  }
54 }
55 
56 
57 // Interpolation of GridFunction to DOFVector
58 template <class C, class GB, class TP, class R>
59  template <class Expr, class Tag>
61  interpolate(Expr&& expr, Tag strategy)
62 {
63  // create temporary copy of data
64  Coefficients tmp(coefficients());
65 
66  Self tmpView{tmp, this->basis(), this->treePath()};
67  tmpView.interpolate_noalias(FWD(expr), strategy);
68 
69  // move data from temporary vector into stored DOFVector
70  coefficients() = std::move(tmp);
71 }
72 
73 } // end namespace AMDiS
void interpolate(Basis const &basis, Vec &vec, GF const &gf, TP const &tp_, C &&c_, BV &&bv_, Assign &&a_)
Interpolate given function in discrete function space.
Definition: Interpolate.hpp:120
void interpolate(Expr &&expr, Tag strategy={})
Interpolation of GridFunction to DOFVector.
Definition: DiscreteFunction.inc.hpp:61
decltype(auto) makeGridFunction(PreGridFct const &preGridFct, GridView const &gridView)
Generator for Gridfunctions from Expressions (PreGridfunctions)
Definition: GridFunction.hpp:154
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
A mutable view on the subspace of a DOFVector,.
Definition: DiscreteFunction.hpp:36
typename Impl::VectorTypeImpl< T, Facade >::type VectorType_t
Type transformation changing the value type of the vector.
Definition: VectorFacade.hpp:369
void interpolate_noalias(Expr &&expr, Tag strategy={})
Interpolation of GridFunction to DOFVector, assuming that there is no reference to this DOFVector in ...
Definition: DiscreteFunction.inc.hpp:36