dune-pdelab  2.5-dev
interpolate.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 
4 #ifndef DUNE_PDELAB_GRIDFUNCTIONSPACE_INTERPOLATE_HH
5 #define DUNE_PDELAB_GRIDFUNCTIONSPACE_INTERPOLATE_HH
6 
7 #include <vector>
8 #include <utility>
9 
10 #include <dune/common/exceptions.hh>
11 
12 #include <dune/localfunctions/common/interfaceswitch.hh>
13 #include <dune/localfunctions/common/virtualinterface.hh>
14 #include <dune/functions/common/functionfromcallable.hh>
15 
16 #include <dune/typetree/typetree.hh>
17 #include <dune/typetree/pairtraversal.hh>
18 
23 
24 namespace Dune {
25  namespace PDELab {
26 
30 
31  // Backend for standard local interpolation
33  {
34  template<typename FE, typename ElemFunction, typename XL>
35  void interpolate(const FE &fe, const ElemFunction &elemFunction,
36  XL &xl) const
37  {
38  FiniteElementInterfaceSwitch<FE>::interpolation(fe).
39  interpolate(elemFunction,xl);
40  }
41  };
42 
43  namespace {
44 
45  template<typename IB, typename LF, typename XG>
46  struct InterpolateLeafFromScalarVisitor
47  : public TypeTree::TreeVisitor
48  , public TypeTree::DynamicTraversal
49  {
50 
51  template<typename LFS, typename TreePath>
52  void leaf(const LFS& lfs, TreePath treePath) const
53  {
54  std::vector<typename XG::ElementType> xl(lfs.size());
55 
56  // call interpolate for the basis
57  ib.interpolate(lfs.finiteElement(), lf, xl);
58 
59  // write coefficients into local vector
60  xg.write_sub_container(lfs,xl);
61  }
62 
63  InterpolateLeafFromScalarVisitor(const IB& ib_, const LF& lf_, XG& xg_)
64  : ib(ib_)
65  , lf(lf_)
66  , xg(xg_)
67  {}
68 
69  const IB& ib;
70  const LF& lf;
71  XG& xg;
72 
73  };
74 
75 
76  template<typename IB, typename LF, typename XG>
77  struct InterpolateLeafFromVectorVisitor
78  : public TypeTree::TreeVisitor
79  , public TypeTree::DynamicTraversal
80  {
81  using Domain = typename Functions::SignatureTraits<LF>::Domain;
82  using Range = typename Functions::SignatureTraits<LF>::Range;
83  using RangeField = typename FieldTraits<Range>::field_type;
84 
85  template<typename LFS, typename TreePath>
86  void leaf(const LFS& lfs, TreePath treePath) const
87  {
88  std::vector<typename XG::ElementType> xl(lfs.size());
89 
90  using LFSRange = typename LFS::Traits::FiniteElement::Traits::LocalBasisType::Traits::RangeType;
91  static_assert(std::is_convertible<LFSRange, typename FieldTraits< LFSRange >::field_type>::value,
92  "only interpolation into scalar leaf function spaces is implemented");
93 
94  // call interpolate for the basis
95  auto f = [&](const Domain& x) -> RangeField { return lf(x)[index]; };
96 
97  using LocalFunction = typename Dune::Functions::FunctionFromCallable<RangeField(Domain), decltype(f), Dune::Function<Domain,RangeField> >;
98  LocalFunction fnkt(f);
99 
100  ib.interpolate(lfs.finiteElement(), fnkt, xl);
101 
102  // write coefficients into local vector
103  xg.write_sub_container(lfs,xl);
104 
105  // increment index
106  index++;
107  }
108 
109  InterpolateLeafFromVectorVisitor(const IB& ib_, const LF& lf_, XG& xg_)
110  : ib(ib_)
111  , lf(lf_)
112  , index(0)
113  , xg(xg_)
114  {}
115 
116  const IB& ib;
117  const LF& lf;
118  mutable std::size_t index;
119  XG& xg;
120 
121  };
122 
123 
124  template<typename IB, typename E, typename XG>
125  struct InterpolateVisitor
126  : public TypeTree::TreePairVisitor
127  , public TypeTree::DynamicTraversal
128  {
129 
130  template<typename F, typename LFS, typename TreePath>
131  typename std::enable_if<F::isLeaf && LFS::isLeaf>::type
132  leaf(const F& f, const LFS& lfs, TreePath treePath) const
133  {
134  std::vector<typename XG::ElementType> xl(lfs.size());
135  // call interpolate for the basis
136  using Domain = typename Functions::SignatureTraits<F>::Domain;
137  using Range = typename Functions::SignatureTraits<F>::Range;
138  using LocalFunction = typename Dune::Functions::FunctionFromCallable<Range(Domain), F, Dune::Function<Domain,Range> >;
139  LocalFunction lf(f);
140  ib.interpolate(lfs.finiteElement(), lf, xl);
141  // write coefficients into local vector
142  xg.write_sub_container(lfs,xl);
143  }
144 
145  // interpolate PowerLFS from scalar-valued function
146  template<typename F, typename LFS, typename TreePath,
147  typename Range = typename Functions::SignatureTraits<F>::Range>
148  typename std::enable_if<F::isLeaf &&
149  std::is_convertible<Range, typename FieldTraits< Range >::field_type>::value &&
150  (!LFS::isLeaf)>::type
151  leaf(const F& f, const LFS& lfs, TreePath treePath) const
152  {
153  // call interpolate for the basis
154  using Domain = typename Functions::SignatureTraits<F>::Domain;
155  using LocalFunction = typename Dune::Functions::FunctionFromCallable<Range(Domain), F, Dune::Function<Domain,Range> >;
156  LocalFunction lf(f);
157  TypeTree::applyToTree(lfs,InterpolateLeafFromScalarVisitor<IB,LocalFunction,XG>(ib,lf,xg));
158 
159  }
160 
161  // interpolate PowerLFS from vector-valued function
162  template<typename F, typename LFS, typename TreePath,
163  typename Range = typename Functions::SignatureTraits<F>::Range>
164  typename std::enable_if<F::isLeaf &&
165  (!std::is_convertible<Range, typename FieldTraits< Range >::field_type>::value) &&
166  (!LFS::isLeaf)>::type
167  leaf(const F& f, const LFS& lfs, TreePath treePath) const
168  {
169  static_assert(TypeTree::TreeInfo<LFS>::leafCount == Range::dimension,
170  "Number of leaves and dimension of range type " \
171  "must match for automatic interpolation of " \
172  "vector-valued function");
173 
174  TypeTree::applyToTree(lfs,InterpolateLeafFromVectorVisitor<IB,F,XG>(ib,f,xg));
175  }
176 
177  InterpolateVisitor(IB ib_, const E& e_, XG& xg_)
178  : ib(ib_)
179  , e(e_)
180  , xg(xg_)
181  {}
182 
183  private:
184  IB ib;
185  const E& e;
186  XG& xg;
187  };
188 
189  } // anonymous namespace
190 
192 
203  template<typename F, typename GFS, typename XG>
204  void interpolate (const F& f, const GFS& gfs, XG& xg)
205  {
206  // this is the leaf version now
207 
208  // get some types
209  using EntitySet = typename GFS::Traits::EntitySet;
210  using Element = typename EntitySet::Element;
211 
212  auto entity_set = gfs.entitySet();
213 
214  // make local function
215  auto lf = makeLocalFunctionTree(f, gfs.gridView());
216 
217  // make local function space
218  typedef LocalFunctionSpace<GFS> LFS;
219  LFS lfs(gfs);
220  typedef LFSIndexCache<LFS> LFSCache;
221  LFSCache lfs_cache(lfs);
222  typedef typename XG::template LocalView<LFSCache> XView;
223 
224  XView x_view(xg);
225 
226  // loop once over the grid
227  for (const auto& element : elements(entity_set))
228  {
229  // bind local function space to element
230  lf.bind(element);
231  lfs.bind(element);
232  lfs_cache.update();
233  x_view.bind(lfs_cache);
234 
235  // call interpolate
236  TypeTree::applyToTreePair(lf,lfs,InterpolateVisitor<InterpolateBackendStandard,Element,XView>(InterpolateBackendStandard(),element,x_view));
237 
238  x_view.unbind();
239  lf.unbind();
240  }
241 
242  x_view.detach();
243  }
244 
246  } // namespace PDELab
247 } // namespace Dune
248 
249 #endif // DUNE_PDELAB_GRIDFUNCTIONSPACE_INTERPOLATE_HH
lfsindexcache.hh
gridfunctionspace.hh
index
std::size_t index
Definition: interpolate.hh:118
lf
const LF & lf
Definition: interpolate.hh:70
Dune
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::PDELab::InterpolateBackendStandard
Definition: interpolate.hh:32
ib
const IB & ib
Definition: interpolate.hh:69
xg
XG & xg
Definition: interpolate.hh:71
localfunction.hh
e
const Entity & e
Definition: localfunctionspace.hh:120
function.hh
Dune::PDELab::LocalFunctionSpace< GFS >
value
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
Dune::PDELab::InterpolateBackendStandard::interpolate
void interpolate(const FE &fe, const ElemFunction &elemFunction, XL &xl) const
Definition: interpolate.hh:35
Dune::PDELab::interpolate
void interpolate(const F &f, const GFS &gfs, XG &xg)
interpolation from a given grid function
Definition: interpolate.hh:204
Dune::PDELab::makeLocalFunctionTree
auto makeLocalFunctionTree(const F &f, const GV &gv) -> Imp::LocalFunctionLeafNodeWrapper< decltype(localFunction(f)) >
turn a given function tree into a local function tree
Definition: localfunction.hh:125
Dune::PDELab::LFSIndexCache< LFS >