dune-pdelab  2.5-dev
gridfunctionspaceutilities.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 #ifndef DUNE_PDELAB_GRIDFUNCTIONSPACE_GRIDFUNCTIONSPACEUTILITIES_HH
4 #define DUNE_PDELAB_GRIDFUNCTIONSPACE_GRIDFUNCTIONSPACEUTILITIES_HH
5 
6 #include <cstdlib>
7 #include<vector>
8 
9 #include<dune/common/exceptions.hh>
10 #include <dune/common/fvector.hh>
11 
12 #include <dune/localfunctions/common/interfaceswitch.hh>
13 
14 #include"../common/function.hh"
16 #include"gridfunctionspace.hh"
19 
20 namespace Dune {
21  namespace PDELab {
22  template <typename ES>
24  : public TypeTree::TreeVisitor
25  , public TypeTree::DynamicTraversal
26  {
27 
28  template<typename LeafGFS, typename TreePath>
29  void leaf(const LeafGFS& leaf_gfs, TreePath tp)
30  {
31  leaf_gfs._es = es;
32  }
33 
34  explicit set_entity_set_visitor(const ES& es_) : es(es_)
35  {
36  }
37 
38  const ES& es;
39  };
40 
41 
45 
46  //===============================================================
47  // output: convert grid function space to discrete grid function
48  //===============================================================
49 
50 
71  template<typename T, typename X>
73  : public TypeTree::LeafNode
75  GridFunctionTraits<
76  typename T::Traits::GridViewType,
77  typename BasisInterfaceSwitch<
78  typename FiniteElementInterfaceSwitch<
79  typename T::Traits::FiniteElementType
80  >::Basis
81  >::RangeField,
82  BasisInterfaceSwitch<
83  typename FiniteElementInterfaceSwitch<
84  typename T::Traits::FiniteElementType
85  >::Basis
86  >::dimRange,
87  typename BasisInterfaceSwitch<
88  typename FiniteElementInterfaceSwitch<
89  typename T::Traits::FiniteElementType
90  >::Basis
91  >::Range
92  >,
93  DiscreteGridFunction<T,X>
94  >
95  {
96  typedef T GFS;
97 
98  typedef typename Dune::BasisInterfaceSwitch<
99  typename FiniteElementInterfaceSwitch<
100  typename T::Traits::FiniteElementType
101  >::Basis
102  > BasisSwitch;
103  typedef GridFunctionInterface<
105  typename T::Traits::GridViewType,
106  typename BasisSwitch::RangeField,
107  BasisSwitch::dimRange,
108  typename BasisSwitch::Range
109  >,
111  > BaseT;
112 
113  public:
114  typedef typename BaseT::Traits Traits;
115 
121  DiscreteGridFunction (const GFS& gfs, const X& x_)
122  : pgfs(stackobject_to_shared_ptr(gfs))
123  , lfs(gfs)
124  , lfs_cache(lfs)
125  , x_view(x_)
126  , xl(gfs.maxLocalSize())
127  , yb(gfs.maxLocalSize())
128  {
129  }
130 
136  DiscreteGridFunction (std::shared_ptr<const GFS> gfs, std::shared_ptr<const X> x_)
137  : pgfs(gfs)
138  , lfs(*gfs)
139  , lfs_cache(lfs)
140  , x_view(*x_)
141  , xl(gfs->maxLocalSize())
142  , yb(gfs->maxLocalSize())
143  , px(x_) // FIXME: The LocalView should handle a shared_ptr correctly!
144  {
145  }
146 
147  // Evaluate
148  inline void evaluate (const typename Traits::ElementType& e,
149  const typename Traits::DomainType& x,
150  typename Traits::RangeType& y) const
151  {
152  typedef FiniteElementInterfaceSwitch<
154  > FESwitch;
155  lfs.bind(e);
156  lfs_cache.update();
157  x_view.bind(lfs_cache);
158  x_view.read(xl);
159  x_view.unbind();
160  FESwitch::basis(lfs.finiteElement()).evaluateFunction(x,yb);
161  y = 0;
162  for (unsigned int i=0; i<yb.size(); i++)
163  {
164  y.axpy(xl[i],yb[i]);
165  }
166  }
167 
169  inline const typename Traits::GridViewType& getGridView () const
170  {
171  return pgfs->gridView();
172  }
173 
174  private:
175 
178  typedef typename X::template ConstLocalView<LFSCache> XView;
179 
180  std::shared_ptr<GFS const> pgfs;
181  mutable LFS lfs;
182  mutable LFSCache lfs_cache;
183  mutable XView x_view;
184  mutable std::vector<typename Traits::RangeFieldType> xl;
185  mutable std::vector<typename Traits::RangeType> yb;
186  std::shared_ptr<const X> px; // FIXME: dummy pointer to make sure we take ownership of X
187  };
188 
198  template<typename T, typename X>
200  public GridFunctionInterface<
201  GridFunctionTraits<
202  typename T::Traits::GridViewType,
203  typename JacobianToCurl<typename T::Traits::FiniteElementType::
204  Traits::LocalBasisType::Traits::JacobianType>::CurlField,
205  JacobianToCurl<typename T::Traits::FiniteElementType::Traits::LocalBasisType::
206  Traits::JacobianType>::dimCurl,
207  typename JacobianToCurl<typename T::Traits::FiniteElementType::
208  Traits::LocalBasisType::Traits::JacobianType>::Curl
209  >,
210  DiscreteGridFunctionCurl<T,X>
211  >
212  {
213  typedef T GFS;
214  typedef typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
215  JacobianType Jacobian;
217 
218  public:
219  typedef GridFunctionTraits<
220  typename T::Traits::GridViewType,
221  typename J2C::CurlField, J2C::dimCurl, typename J2C::Curl
223 
229  DiscreteGridFunctionCurl(const GFS& gfs, const X& x_)
230  : pgfs(stackobject_to_shared_ptr(gfs))
231  , lfs(gfs)
232  , lfs_cache(lfs)
233  , x_view(x_)
234  , xl(gfs.maxLocalSize())
235  , jacobian(gfs.maxLocalSize())
236  , yb(gfs.maxLocalSize())
237  , px(stackobject_to_shared_ptr(x_))
238  {}
239 
240  // Evaluate
241  void evaluate (const typename Traits::ElementType& e,
242  const typename Traits::DomainType& x,
243  typename Traits::RangeType& y) const
244  {
245  static const J2C& j2C = J2C();
246 
247  lfs.bind();
248  lfs_cache.update();
249  x_view.bind(lfs_cache);
250  x_view.read(xl);
251  x_view.unbind();
252 
253  lfs.finiteElement().basis().evaluateJacobian(x,jacobian);
254 
255  y = 0;
256  for (std::size_t i=0; i < lfs.size(); i++) {
257  j2C(jacobian[i], yb);
258  y.axpy(xl[i], yb);
259  }
260  }
261 
263  const typename Traits::GridViewType& getGridView() const
264  { return pgfs->gridView(); }
265 
266  private:
268  BaseT;
271  typedef typename X::template ConstLocalView<LFSCache> XView;
272 
273  std::shared_ptr<GFS const> pgfs;
274  mutable LFS lfs;
275  mutable LFSCache lfs_cache;
276  mutable XView x_view;
277  mutable std::vector<typename Traits::RangeFieldType> xl;
278  mutable std::vector<Jacobian> jacobian;
279  mutable std::vector<typename Traits::RangeType> yb;
280  std::shared_ptr<const X> px; // FIXME: dummy pointer to make sure we take ownership of X
281  };
282 
284 
296  template<typename GV, typename RangeFieldType, int dimRangeOfBasis>
298  static_assert(AlwaysFalse<GV>::value,
299  "DiscreteGridFunctionCurl (and friends) work in 2D "
300  "and 3D only");
301  };
303 
309  template<typename GV, typename RangeFieldType>
311  : public GridFunctionTraits<GV,
312  RangeFieldType, 2,
313  FieldVector<RangeFieldType, 2> >
314  {
315  static_assert(GV::dimensionworld == 2,
316  "World dimension of grid must be 2 for the curl of a "
317  "scalar (1D) quantity");
318  };
320 
325  template<typename GV, typename RangeFieldType>
327  : public GridFunctionTraits<GV,
328  RangeFieldType, 1,
329  FieldVector<RangeFieldType, 1> >
330  {
331  static_assert(GV::dimensionworld == 2,
332  "World dimension of grid must be 2 for the curl of a"
333  "2D quantity");
334  };
336 
341  template<typename GV, typename RangeFieldType>
343  : public GridFunctionTraits<GV,
344  RangeFieldType, 3,
345  FieldVector<RangeFieldType, 3> >
346  {
347  static_assert(GV::dimensionworld == 3,
348  "World dimension of grid must be 3 for the curl of a"
349  "3D quantity");
350  };
351 
354 
368  template<typename T, typename X>
370  : public GridFunctionInterface<
371  DiscreteGridFunctionCurlTraits<
372  typename T::Traits::GridViewType,
373  typename T::Traits::FiniteElementType::Traits::
374  LocalBasisType::Traits::RangeFieldType,
375  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
376  dimRange>,
377  DiscreteGridFunctionGlobalCurl<T,X> >
378  {
379  public:
381  typename T::Traits::GridViewType,
382  typename T::Traits::FiniteElementType::Traits::
383  LocalBasisType::Traits::RangeFieldType,
384  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
385  dimRange> Traits;
386 
387  private:
388  typedef T GFS;
389  typedef GridFunctionInterface<
390  Traits,
392  typedef typename T::Traits::FiniteElementType::Traits::
393  LocalBasisType::Traits LBTraits;
394 
395  public:
401  DiscreteGridFunctionGlobalCurl (const GFS& gfs, const X& x_)
402  : pgfs(stackobject_to_shared_ptr(gfs))
403  , lfs(gfs)
404  , lfs_cache(lfs)
405  , x_view(x_)
406  , xl(gfs.maxLocalSize())
407  , J(gfs.maxLocalSize())
408  , px(stackobject_to_shared_ptr(x_))
409  {}
410 
411 
412  // Evaluate
413  inline void evaluate (const typename Traits::ElementType& e,
414  const typename Traits::DomainType& x,
415  typename Traits::RangeType& y) const
416  {
417  lfs.bind(e);
418  lfs_cache.update();
419  x_view.bind(lfs_cache);
420  x_view.read(xl);
421  x_view.unbind();
422 
423  lfs.finiteElement().localBasis().
424  evaluateJacobianGlobal(x,J,e.geometry());
425  y = 0;
426  for (unsigned int i=0; i<J.size(); i++)
427  // avoid a "case label value exceeds maximum value for type"
428  // warning: since dimRange is an anonymous enum, its type may
429  // contain only the values 0 and 1, resulting in a warning.
430  switch(unsigned(Traits::dimRange)) {
431  case 1:
432  y[0] += xl[i] * J[i][0][1];
433  y[1] += xl[i] * -J[i][0][0];
434  break;
435  case 2:
436  y[0] += xl[i]*(J[i][1][0] - J[i][0][1]);
437  break;
438  case 3:
439  y[0] += xl[i]*(J[i][2][1] - J[i][1][2]);
440  y[1] += xl[i]*(J[i][0][2] - J[i][2][0]);
441  y[2] += xl[i]*(J[i][1][0] - J[i][0][1]);
442  break;
443  default:
444  //how did that pass all the static asserts?
445  std::abort();
446  }
447  }
448 
450  inline const typename Traits::GridViewType& getGridView () const
451  {
452  return pgfs->gridView();
453  }
454 
455  private:
458  typedef typename X::template ConstLocalView<LFSCache> XView;
459 
460  std::shared_ptr<GFS const> pgfs;
461  mutable LFS lfs;
462  mutable LFSCache lfs_cache;
463  mutable XView x_view;
464  mutable std::vector<typename Traits::RangeFieldType> xl;
465  mutable std::vector<typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::JacobianType> J;
466  std::shared_ptr<const X> px; // FIXME: dummy pointer to make sure we take ownership of X
467  };
468 
471 
479  template<typename T, typename X>
481  : public GridFunctionInterface<
482  GridFunctionTraits<
483  typename T::Traits::GridViewType,
484  typename T::Traits::FiniteElementType::Traits::LocalBasisType
485  ::Traits::RangeFieldType,
486  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
487  ::dimDomain,
488  FieldVector<
489  typename T::Traits::FiniteElementType::Traits
490  ::LocalBasisType::Traits::RangeFieldType,
491  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
492  ::dimDomain> >,
493  DiscreteGridFunctionGradient<T,X> >
494  {
495  typedef T GFS;
496  typedef typename GFS::Traits::FiniteElementType::Traits::
497  LocalBasisType::Traits LBTraits;
498 
499  public:
500  typedef GridFunctionTraits<
501  typename GFS::Traits::GridViewType,
502  typename LBTraits::RangeFieldType,
503  LBTraits::dimDomain,
504  FieldVector<
505  typename LBTraits::RangeFieldType,
506  LBTraits::dimDomain> > Traits;
507 
508  private:
509  typedef GridFunctionInterface<
510  Traits,
512 
513  public:
519  DiscreteGridFunctionGradient (const GFS& gfs, const X& x_)
520  : pgfs(stackobject_to_shared_ptr(gfs))
521  , lfs(gfs)
522  , lfs_cache(lfs)
523  , x_view(x_)
524  , xl(lfs.size())
525  { }
526 
527  // Evaluate
528  inline void evaluate (const typename Traits::ElementType& e,
529  const typename Traits::DomainType& x,
530  typename Traits::RangeType& y) const
531  {
532  // get and bind local functions space
533  lfs.bind(e);
534  lfs_cache.update();
535  x_view.bind(lfs_cache);
536 
537  // get local coefficients
538  xl.resize(lfs.size());
539  x_view.read(xl);
540  x_view.unbind();
541 
542  // get Jacobian of geometry
543  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
544  JgeoIT = e.geometry().jacobianInverseTransposed(x);
545 
546  // get local Jacobians/gradients of the shape functions
547  std::vector<typename LBTraits::JacobianType> J(lfs.size());
548  lfs.finiteElement().localBasis().evaluateJacobian(x,J);
549 
550  typename Traits::RangeType gradphi;
551  y = 0;
552  for(unsigned int i = 0; i < lfs.size(); ++i) {
553  // compute global gradient of shape function i
554  gradphi = 0;
555  JgeoIT.umv(J[i][0], gradphi);
556 
557  // sum up global gradients, weighting them with the appropriate coeff
558  y.axpy(xl[i], gradphi);
559  }
560 
561  }
562 
564  inline const typename Traits::GridViewType& getGridView () const
565  {
566  return pgfs->gridView();
567  }
568 
569  private:
572  typedef typename X::template ConstLocalView<LFSCache> XView;
573 
574  std::shared_ptr<GFS const> pgfs;
575  mutable LFS lfs;
576  mutable LFSCache lfs_cache;
577  mutable XView x_view;
578  mutable std::vector<typename Traits::RangeFieldType> xl;
579  };
580 
585  template<typename T, typename X>
587  : public GridFunctionInterface<
588  GridFunctionTraits<
589  typename T::Traits::GridViewType,
590  typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
591  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
592  typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType
593  >,
594  DiscreteGridFunctionPiola<T,X>
595  >
596  {
597  typedef T GFS;
598 
599  typedef GridFunctionInterface<
601  typename T::Traits::GridViewType,
602  typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
603  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
604  typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType
605  >,
607  > BaseT;
608 
609  public:
610  typedef typename BaseT::Traits Traits;
611 
616  DiscreteGridFunctionPiola (const GFS& gfs, const X& x_)
617  : pgfs(stackobject_to_shared_ptr(gfs))
618  , lfs(gfs)
619  , lfs_cache(lfs)
620  , x_view(x_)
621  , xl(pgfs->maxLocalSize())
622  , yb(pgfs->maxLocalSize())
623  {
624  }
625 
626  inline void evaluate (const typename Traits::ElementType& e,
627  const typename Traits::DomainType& x,
628  typename Traits::RangeType& y) const
629  {
630  // evaluate shape function on the reference element as before
631  lfs.bind(e);
632  lfs_cache.update();
633  x_view.bind(lfs_cache);
634  x_view.read(xl);
635  x_view.unbind();
636 
637  lfs.finiteElement().localBasis().evaluateFunction(x,yb);
638  typename Traits::RangeType yhat;
639  yhat = 0;
640  for (unsigned int i=0; i<yb.size(); i++)
641  yhat.axpy(xl[i],yb[i]);
642 
643  // apply Piola transformation
644  typename Traits::ElementType::Geometry::JacobianInverseTransposed
645  J = e.geometry().jacobianInverseTransposed(x);
646  J.invert();
647  y = 0;
648  J.umtv(yhat,y);
649  y /= J.determinant();
650  }
651 
653  inline const typename Traits::GridViewType& getGridView () const
654  {
655  return pgfs->gridView();
656  }
657 
658  private:
659 
662  typedef typename X::template ConstLocalView<LFSCache> XView;
663 
664  std::shared_ptr<GFS const> pgfs;
665  mutable LFS lfs;
666  mutable LFSCache lfs_cache;
667  mutable XView x_view;
668  mutable std::vector<typename Traits::RangeFieldType> xl;
669  mutable std::vector<typename Traits::RangeType> yb;
670 
671  };
672 
684 #ifdef __clang__
685  // clang is too stupid to correctly apply the constexpr qualifier of staticDegree in this context
687 #else
689 #endif
691  : public GridFunctionInterface<
692  GridFunctionTraits<
693  typename T::Traits::GridViewType,
694  typename T::template Child<0>::Type::Traits::FiniteElementType
695  ::Traits::LocalBasisType::Traits::RangeFieldType,
696  dimR,
697  Dune::FieldVector<
698  typename T::template Child<0>::Type::Traits::FiniteElementType
699  ::Traits::LocalBasisType::Traits::RangeFieldType,
700  dimR
701  >
702  >,
703  VectorDiscreteGridFunction<T,X>
704  >,
705  public TypeTree::LeafNode
706  {
707  typedef T GFS;
708 
709  typedef GridFunctionInterface<
711  typename T::Traits::GridViewType,
712  typename T::template Child<0>::Type::Traits::FiniteElementType
713  ::Traits::LocalBasisType::Traits::RangeFieldType,
714  dimR,
715  Dune::FieldVector<
716  typename T::template Child<0>::Type::Traits::FiniteElementType
717  ::Traits::LocalBasisType::Traits::RangeFieldType,
718  dimR
719  >
720  >,
722  > BaseT;
723 
724  public:
725  typedef typename BaseT::Traits Traits;
726  typedef typename T::template Child<0>::Type ChildType;
727  typedef typename ChildType::Traits::FiniteElementType
728  ::Traits::LocalBasisType::Traits::RangeFieldType RF;
729  typedef typename ChildType::Traits::FiniteElementType
730  ::Traits::LocalBasisType::Traits::RangeType RT;
731 
733 
738  VectorDiscreteGridFunction(const GFS& gfs, const X& x_,
739  std::size_t start = 0)
740  : pgfs(stackobject_to_shared_ptr(gfs))
741  , lfs(gfs)
742  , lfs_cache(lfs)
743  , x_view(x_)
744  , xl(gfs.maxLocalSize())
745  , yb(gfs.maxLocalSize())
746  {
747  for(std::size_t i = 0; i < dimR; ++i)
748  remap[i] = i + start;
749  }
750 
752 
763  template<class Remap>
764  VectorDiscreteGridFunction(const GFS& gfs, const X& x_,
765  const Remap &remap_)
766  : pgfs(stackobject_to_shared_ptr(gfs))
767  , lfs(gfs)
768  , lfs_cache(lfs)
769  , x_view(x_)
770  , xl(gfs.maxLocalSize())
771  , yb(gfs.maxLocalSize())
772  , px(stackobject_to_shared_ptr(x_))
773  {
774  for(std::size_t i = 0; i < dimR; ++i)
775  remap[i] = remap_[i];
776  }
777 
778  inline void evaluate (const typename Traits::ElementType& e,
779  const typename Traits::DomainType& x,
780  typename Traits::RangeType& y) const
781  {
782  lfs.bind(e);
783  lfs_cache.update();
784  x_view.bind(lfs_cache);
785  x_view.read(xl);
786  x_view.unbind();
787  for (unsigned int k=0; k < dimR; k++)
788  {
789  lfs.child(remap[k]).finiteElement().localBasis().
790  evaluateFunction(x,yb);
791  y[k] = 0.0;
792  for (unsigned int i=0; i<yb.size(); i++)
793  y[k] += xl[lfs.child(remap[k]).localIndex(i)]*yb[i];
794  }
795  }
796 
798  inline const typename Traits::GridViewType& getGridView () const
799  {
800  return pgfs->gridView();
801  }
802 
803  private:
806  typedef typename X::template ConstLocalView<LFSCache> XView;
807 
808  std::shared_ptr<GFS const> pgfs;
809  std::size_t remap[dimR];
810  mutable LFS lfs;
811  mutable LFSCache lfs_cache;
812  mutable XView x_view;
813  mutable std::vector<RF> xl;
814  mutable std::vector<RT> yb;
815  std::shared_ptr<const X> px; // FIXME: dummy pointer to make sure we take ownership of X
816  };
817 
823  template<typename T, typename X>
825  : public GridFunctionInterface<
826  GridFunctionTraits<
827  typename T::Traits::GridViewType,
828  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
829  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain,
830  TypeTree::StaticDegree<T>::value,
831  Dune::FieldMatrix<
832  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
833  TypeTree::StaticDegree<T>::value,
834  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain
835  >
836  >,
837  VectorDiscreteGridFunctionGradient<T,X>
838  >,
839  public TypeTree::LeafNode
840  {
841  typedef T GFS;
842 
843  typedef GridFunctionInterface<
845  typename T::Traits::GridViewType,
846  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
847  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain,
849  Dune::FieldMatrix<
850  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
851  TypeTree::StaticDegree<T>::value,
852  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain>
853  >,
855  > BaseT;
856 
857  public:
858  typedef typename BaseT::Traits Traits;
859  typedef typename T::template Child<0>::Type ChildType;
860  typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits;
861 
862  typedef typename LBTraits::RangeFieldType RF;
863  typedef typename LBTraits::JacobianType JT;
864 
865  VectorDiscreteGridFunctionGradient (const GFS& gfs, const X& x_)
866  : pgfs(stackobject_to_shared_ptr(gfs))
867  , lfs(gfs)
868  , lfs_cache(lfs)
869  , x_view(x_)
870  , xl(gfs.maxLocalSize())
871  , J(gfs.maxLocalSize())
872  {
873  }
874 
875  inline void evaluate(const typename Traits::ElementType& e,
876  const typename Traits::DomainType& x,
877  typename Traits::RangeType& y) const
878  {
879  // get and bind local functions space
880  lfs.bind(e);
881  lfs_cache.update();
882  x_view.bind(lfs_cache);
883  x_view.read(xl);
884  x_view.unbind();
885 
886  // get Jacobian of geometry
887  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
888  JgeoIT = e.geometry().jacobianInverseTransposed(x);
889 
890  y = 0.0;
891 
892  // Loop over PowerLFS and calculate gradient for each child separately
893  for(unsigned int k = 0; k != TypeTree::degree(lfs); ++k)
894  {
895  // get local Jacobians/gradients of the shape functions
896  std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
897  lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
898 
899  Dune::FieldVector<RF,LBTraits::dimDomain> gradphi;
900  for (typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++)
901  {
902  gradphi = 0;
903  JgeoIT.umv(J[i][0], gradphi);
904 
905  y[k].axpy(xl[lfs.child(k).localIndex(i)], gradphi);
906  }
907  }
908  }
909 
910 
912  inline const typename Traits::GridViewType& getGridView () const
913  {
914  return pgfs->gridView();
915  }
916 
917  private:
920  typedef typename X::template ConstLocalView<LFSCache> XView;
921 
922  std::shared_ptr<GFS const> pgfs;
923  mutable LFS lfs;
924  mutable LFSCache lfs_cache;
925  mutable XView x_view;
926  mutable std::vector<RF> xl;
927  mutable std::vector<JT> J;
928  std::shared_ptr<const X> px; // FIXME: dummy pointer to make sure we take ownership of X
929  };
930 
944  template<typename Mat, typename RF, std::size_t size>
950  template<typename T>
951  static inline RF compute_derivative(const Mat& mat, const T& t, const unsigned int k)
952  {
953  // setting this to zero is just a test if the template specialization work
954  Dune::FieldVector<RF,size> grad_phi(0.0);
955  mat.umv(t,grad_phi);
956  return grad_phi[k];
957  // return 0.0;
958  }
959  };
960 
969  template<typename RF, std::size_t size>
970  struct SingleDerivativeComputationHelper<Dune::FieldMatrix<RF,size,size>,RF,size> {
976  template<typename T>
977  static inline RF compute_derivative(const Dune::FieldMatrix<RF,size,size>& mat, const T& t, const unsigned int k)
978  {
979  return mat[k]*t;
980  }
981  };
982 
993  template<typename RF, std::size_t size>
994  struct SingleDerivativeComputationHelper<Dune::DiagonalMatrix<RF,size>,RF,size> {
1000  template<typename T>
1001  static inline RF compute_derivative(const Dune::DiagonalMatrix<RF,size>& mat, const T& t, const unsigned int k)
1002  {
1003  return mat[k][k]*t[k];
1004  }
1005  };
1006 
1017  template<typename T, typename X>
1020  Dune::PDELab::GridFunctionTraits<
1021  typename T::Traits::GridViewType,
1022  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1023  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1024  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1025  VectorDiscreteGridFunctionDiv<T,X> >
1026  , public TypeTree::LeafNode
1027  {
1028  typedef T GFS;
1029 
1032  typename T::Traits::GridViewType,
1033  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1034  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1035  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1037  public :
1038  typedef typename BaseT::Traits Traits;
1039  typedef typename T::template Child<0>::Type ChildType;
1040  typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits;
1041 
1042  typedef typename LBTraits::RangeFieldType RF;
1043  typedef typename LBTraits::JacobianType JT;
1044 
1045  VectorDiscreteGridFunctionDiv(const GFS& gfs, const X& x_)
1046  : pgfs(stackobject_to_shared_ptr(gfs))
1047  , lfs(gfs)
1048  , lfs_cache(lfs)
1049  , x_view(x_)
1050  , xl(gfs.maxLocalSize())
1051  , J(gfs.maxLocalSize())
1052  {
1053  static_assert(LBTraits::dimDomain == TypeTree::StaticDegree<T>::value,
1054  "dimDomain and number of children has to be the same");
1055  }
1056 
1057  inline void evaluate(const typename Traits::ElementType& e,
1058  const typename Traits::DomainType& x,
1059  typename Traits::RangeType& y) const
1060  {
1061  // get and bind local functions space
1062  lfs.bind(e);
1063  lfs_cache.update();
1064  x_view.bind(lfs_cache);
1065  x_view.read(xl);
1066  x_view.unbind();
1067 
1068  // get Jacobian of geometry
1069  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1070  JgeoIT = e.geometry().jacobianInverseTransposed(x);
1071 
1072  const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1073  Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1074 
1075  y = 0.0;
1076 
1077  // loop over VectorGFS and calculate k-th derivative of k-th child
1078  for(unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1079 
1080  // get local Jacobians/gradients of the shape functions
1081  std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1082  lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1083 
1084  RF d_k_phi;
1085  for(typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1086  // compute k-th derivative of k-th child
1087  d_k_phi =
1089  typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1090  typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1091  N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1092  (JgeoIT,J[i][0],k);
1093 
1094  y += xl[lfs.child(k).localIndex(i)] * d_k_phi;
1095  }
1096  }
1097  }
1098 
1100  inline const typename Traits::GridViewType& getGridView() const
1101  {
1102  return pgfs->gridView();
1103  }
1104 
1105  private :
1108  typedef typename X::template ConstLocalView<LFSCache> XView;
1109 
1110  shared_ptr<GFS const> pgfs;
1111  mutable LFS lfs;
1112  mutable LFSCache lfs_cache;
1113  mutable XView x_view;
1114  mutable std::vector<RF> xl;
1115  mutable std::vector<JT> J;
1116  shared_ptr<const X> px;
1117  }; // end class VectorDiscreteGridFunctionDiv
1118 
1133  {
1134  typedef T GFS;
1135  public :
1136  VectorDiscreteGridFunctionCurl(const GFS& gfs, const X& x)
1137  {
1139  "Curl computation can only be done in two or three dimensions");
1140  }
1141  };
1142 
1153  template<typename T, typename X>
1156  Dune::PDELab::GridFunctionTraits<
1157  typename T::Traits::GridViewType,
1158  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1159  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1160  TypeTree::StaticDegree<T>::value,
1161  Dune::FieldVector<
1162  typename T::template Child<0>::Type::Traits::FiniteElementType
1163  ::Traits::LocalBasisType::Traits::RangeFieldType,
1164  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange
1165  TypeTree::StaticDegree<T>::value
1166  >
1167  >,
1168  VectorDiscreteGridFunctionCurl<T,X>
1169  >
1170  , public TypeTree::LeafNode
1171  {
1172  typedef T GFS;
1173 
1176  typename T::Traits::GridViewType,
1177  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1178  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1180  Dune::FieldVector<
1181  typename T::template Child<0>::Type::Traits::FiniteElementType
1182  ::Traits::LocalBasisType::Traits::RangeFieldType,
1183  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange
1184  TypeTree::StaticDegree<T>::value
1185  >
1186  >,
1188 
1189  public :
1190  typedef typename BaseT::Traits Traits;
1191  typedef typename T::template Child<0>::Type ChildType;
1192  typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits;
1193 
1194  typedef typename LBTraits::RangeFieldType RF;
1195  typedef typename LBTraits::JacobianType JT;
1196 
1197  VectorDiscreteGridFunctionCurl(const GFS& gfs, const X& x_)
1198  : pgfs(stackobject_to_shared_ptr(gfs))
1199  , lfs(gfs)
1200  , lfs_cache(lfs)
1201  , x_view(x_)
1202  , xl(gfs.maxLocalSize())
1203  , J(gfs.maxLocalSize())
1204  {
1205  static_assert(LBTraits::dimDomain == TypeTree::StaticDegree<T>::value,
1206  "dimDomain and number of children has to be the same");
1207  }
1208 
1209  inline void evaluate(const typename Traits::ElementType& e,
1210  const typename Traits::DomainType& x,
1211  typename Traits::RangeType& y) const
1212  {
1213  // get and bind local functions space
1214  lfs.bind(e);
1215  lfs_cache.update();
1216  x_view.bind(lfs_cache);
1217  x_view.read(xl);
1218  x_view.unbind();
1219 
1220  // get Jacobian of geometry
1221  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1222  JgeoIT = e.geometry().jacobianInverseTransposed(x);
1223 
1224  const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1225  Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1226 
1227  y = 0.0;
1228 
1229  // some handy variables for the curl in 3D
1230  int i1, i2;
1231 
1232  // loop over childs of VectorGFS
1233  for(unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1234 
1235  // get local Jacobians/gradients of the shape functions
1236  std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1237  lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1238 
1239  // pick out the right derivative and component for the curl computation
1240  i1 = (k+1)%3;
1241  i2 = (k+2)%3;
1242 
1243  RF d_k_phi;
1244  for(typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1245  // compute i2-th derivative of k-th child
1246  d_k_phi =
1248  typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1249  typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1250  N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1251  (JgeoIT,J[i][0],i2);
1252 
1253  y[i1] += xl[lfs.child(k).localIndex(i)] * d_k_phi;
1254 
1255  // compute i1-th derivative of k-th child
1256  d_k_phi =
1258  typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1259  typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1260  N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1261  (JgeoIT,J[i][0],i1);
1262 
1263  y[i2] -= xl[lfs.child(k).localIndex(i)] * d_k_phi;
1264  }
1265  }
1266  }
1267 
1269  inline const typename Traits::GridViewType& getGridView() const
1270  {
1271  return pgfs->gridView();
1272  }
1273 
1274  private :
1277  typedef typename X::template ConstLocalView<LFSCache> XView;
1278 
1279  shared_ptr<GFS const> pgfs;
1280  mutable LFS lfs;
1281  mutable LFSCache lfs_cache;
1282  mutable XView x_view;
1283  mutable std::vector<RF> xl;
1284  mutable std::vector<JT> J;
1285  shared_ptr<const X> px;
1286  }; // end class VectorDiscreteGridFunctionCurl (3D)
1287 
1298  template<typename T, typename X>
1301  Dune::PDELab::GridFunctionTraits<
1302  typename T::Traits::GridViewType,
1303  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1304  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1305  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1306  VectorDiscreteGridFunctionDiv<T,X> >
1307  , public TypeTree::LeafNode
1308  {
1309  typedef T GFS;
1310 
1313  typename T::Traits::GridViewType,
1314  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1315  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1316  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1318  public :
1319  typedef typename BaseT::Traits Traits;
1320  typedef typename T::template Child<0>::Type ChildType;
1321  typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits;
1322 
1323  typedef typename LBTraits::RangeFieldType RF;
1324  typedef typename LBTraits::JacobianType JT;
1325 
1326  VectorDiscreteGridFunctionCurl(const GFS& gfs, const X& x_)
1327  : pgfs(stackobject_to_shared_ptr(gfs))
1328  , lfs(gfs)
1329  , lfs_cache(lfs)
1330  , x_view(x_)
1331  , xl(gfs.maxLocalSize())
1332  , J(gfs.maxLocalSize())
1333  {
1334  static_assert(LBTraits::dimDomain == TypeTree::StaticDegree<T>::value,
1335  "dimDomain and number of children has to be the same");
1336  }
1337 
1338  inline void evaluate(const typename Traits::ElementType& e,
1339  const typename Traits::DomainType& x,
1340  typename Traits::RangeType& y) const
1341  {
1342  // get and bind local functions space
1343  lfs.bind(e);
1344  lfs_cache.update();
1345  x_view.bind(lfs_cache);
1346  x_view.read(xl);
1347  x_view.unbind();
1348 
1349  // get Jacobian of geometry
1350  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1351  JgeoIT = e.geometry().jacobianInverseTransposed(x);
1352 
1353  const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1354  Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1355 
1356  y = 0.0;
1357 
1358  // some handy variables for the curl computation in 2D
1359  RF sign = -1.0;
1360  int i2;
1361 
1362  // loop over childs of VectorGFS
1363  for(unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1364 
1365  // get local Jacobians/gradients of the shape functions
1366  std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1367  lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1368 
1369  RF d_k_phi;
1370  // pick out the right derivative
1371  i2 = 1-k;
1372  for(typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1373  // compute i2-th derivative of k-th child
1374  d_k_phi =
1376  typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1377  typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1378  N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1379  (JgeoIT,J[i][0],i2);
1380 
1381  y += sign * xl[lfs.child(k).localIndex(i)] * d_k_phi;
1382  }
1383  sign *= -1.0;
1384  }
1385  }
1386 
1388  inline const typename Traits::GridViewType& getGridView() const
1389  {
1390  return pgfs->gridView();
1391  }
1392 
1393  private :
1396  typedef typename X::template ConstLocalView<LFSCache> XView;
1397 
1398  shared_ptr<GFS const> pgfs;
1399  mutable LFS lfs;
1400  mutable LFSCache lfs_cache;
1401  mutable XView x_view;
1402  mutable std::vector<RF> xl;
1403  mutable std::vector<JT> J;
1404  shared_ptr<const X> px;
1405  }; // end class VectorDiscreteGridFunctionCurl (2D)
1406 
1408  } // namespace PDELab
1409 } // namespace Dune
1410 
1411 #endif // DUNE_PDELAB_GRIDFUNCTIONSPACE_GRIDFUNCTIONSPACEUTILITIES_HH
Helper class to compute a single derivative of scalar basis functions.
Definition: gridfunctionspaceutilities.hh:945
GridFunctionTraits< typename T::Traits::GridViewType, typename J2C::CurlField, J2C::dimCurl, typename J2C::Curl > Traits
Definition: gridfunctionspaceutilities.hh:222
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:912
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:858
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1197
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1323
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1040
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:875
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1039
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1195
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:114
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:148
T Traits
Export type traits.
Definition: function.hh:191
DiscreteGridFunctionPiola(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionPiola.
Definition: gridfunctionspaceutilities.hh:616
Helper class to calculate the Traits of DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:297
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:564
Definition: gridfunctionspaceutilities.hh:23
ChildType::Traits::FiniteElementType ::Traits::LocalBasisType::Traits::RangeType RT
Definition: gridfunctionspaceutilities.hh:730
DiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGradient.
Definition: gridfunctionspaceutilities.hh:519
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:117
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, const Remap &remap_)
construct
Definition: gridfunctionspaceutilities.hh:764
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:726
static RF compute_derivative(const Dune::FieldMatrix< RF, size, size > &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:977
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:626
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:528
VectorDiscreteGridFunctionDiv(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1045
convert a single component function space with a grid function representing the gradient ...
Definition: gridfunctionspaceutilities.hh:480
ChildType::Traits::FiniteElementType ::Traits::LocalBasisType::Traits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:728
DiscreteGridFunction(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:121
convert a grid function space and a coefficient vector into a grid function of the curl ...
Definition: gridfunctionspaceutilities.hh:199
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, std::size_t start=0)
construct
Definition: gridfunctionspaceutilities.hh:738
DiscreteGridFunctionGlobalCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGlobalCurl.
Definition: gridfunctionspaceutilities.hh:401
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:862
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1100
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:798
DiscreteGridFunction for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:690
a GridFunction maps x in DomainType to y in RangeType
Definition: function.hh:186
convert a single component function space with experimental global finite elements into a grid functi...
Definition: gridfunctionspaceutilities.hh:369
GV GridViewType
The type of the grid view the function lives on.
Definition: function.hh:114
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1269
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:169
set_entity_set_visitor(const ES &es_)
Definition: gridfunctionspaceutilities.hh:34
const ES & es
Definition: gridfunctionspaceutilities.hh:38
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1194
extract the curl of a function from the jacobian of that function
Definition: jacobiantocurl.hh:27
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1320
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1190
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:610
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x)
Definition: gridfunctionspaceutilities.hh:1136
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1057
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1038
Equivalent of DiscreteGridFunctionGradient for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:824
void leaf(const LeafGFS &leaf_gfs, TreePath tp)
Definition: gridfunctionspaceutilities.hh:29
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1321
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1324
DiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:229
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:863
Dune::FieldVector< GV::Grid::ctype, GV::dimension > DomainType
domain type in dim-size coordinates
Definition: function.hh:48
static RF compute_derivative(const Mat &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:951
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:450
static RF compute_derivative(const Dune::DiagonalMatrix< RF, size > &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:1001
Compute curl of vector-valued functions.
Definition: gridfunctionspaceutilities.hh:1132
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1388
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:241
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1209
const Entity & e
Definition: localfunctionspace.hh:111
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:413
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:653
VectorDiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:865
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1191
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:263
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:778
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1043
GridFunctionTraits< typename GFS::Traits::GridViewType, typename LBTraits::RangeFieldType, LBTraits::dimDomain, FieldVector< typename LBTraits::RangeFieldType, LBTraits::dimDomain > > Traits
Definition: gridfunctionspaceutilities.hh:506
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1319
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:860
DiscreteGridFunction(std::shared_ptr< const GFS > gfs, std::shared_ptr< const X > x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:136
Create a local function space from a global function space.
Definition: localfunctionspace.hh:670
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:859
convert a grid function space and a coefficient vector into a grid function
Definition: gridfunctionspaceutilities.hh:72
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1042
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1192
DiscreteGridFunctionCurlTraits< typename T::Traits::GridViewType, typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType, T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange > Traits
Definition: gridfunctionspaceutilities.hh:385
DiscreteGridFunction with Piola transformation.
Definition: gridfunctionspaceutilities.hh:586
traits class holding the function signature, same as in local function
Definition: function.hh:175
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1338
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:725
Compute divergence of vector-valued functions.
Definition: gridfunctionspaceutilities.hh:1018
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1326