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 
26 
27  //===============================================================
28  // output: convert grid function space to discrete grid function
29  //===============================================================
30 
31 
52  template<typename T, typename X>
54  : public TypeTree::LeafNode
56  GridFunctionTraits<
57  typename T::Traits::GridViewType,
58  typename BasisInterfaceSwitch<
59  typename FiniteElementInterfaceSwitch<
60  typename T::Traits::FiniteElementType
61  >::Basis
62  >::RangeField,
63  BasisInterfaceSwitch<
64  typename FiniteElementInterfaceSwitch<
65  typename T::Traits::FiniteElementType
66  >::Basis
67  >::dimRange,
68  typename BasisInterfaceSwitch<
69  typename FiniteElementInterfaceSwitch<
70  typename T::Traits::FiniteElementType
71  >::Basis
72  >::Range
73  >,
74  DiscreteGridFunction<T,X>
75  >
76  {
77  typedef T GFS;
78 
79  typedef typename Dune::BasisInterfaceSwitch<
80  typename FiniteElementInterfaceSwitch<
81  typename T::Traits::FiniteElementType
82  >::Basis
83  > BasisSwitch;
84  typedef GridFunctionInterface<
86  typename T::Traits::GridViewType,
87  typename BasisSwitch::RangeField,
88  BasisSwitch::dimRange,
89  typename BasisSwitch::Range
90  >,
92  > BaseT;
93 
94  public:
95  typedef typename BaseT::Traits Traits;
96 
102  DiscreteGridFunction (const GFS& gfs, const X& x_)
103  : pgfs(stackobject_to_shared_ptr(gfs))
104  , lfs(gfs)
105  , lfs_cache(lfs)
106  , x_view(x_)
107  , xl(gfs.maxLocalSize())
108  , yb(gfs.maxLocalSize())
109  {
110  }
111 
117  DiscreteGridFunction (std::shared_ptr<const GFS> gfs, std::shared_ptr<const X> x_)
118  : pgfs(gfs)
119  , lfs(*gfs)
120  , lfs_cache(lfs)
121  , x_view(*x_)
122  , xl(gfs->maxLocalSize())
123  , yb(gfs->maxLocalSize())
124  , px(x_) // FIXME: The LocalView should handle a shared_ptr correctly!
125  {
126  }
127 
128  // Evaluate
129  inline void evaluate (const typename Traits::ElementType& e,
130  const typename Traits::DomainType& x,
131  typename Traits::RangeType& y) const
132  {
133  typedef FiniteElementInterfaceSwitch<
135  > FESwitch;
136  lfs.bind(e);
137  lfs_cache.update();
138  x_view.bind(lfs_cache);
139  x_view.read(xl);
140  x_view.unbind();
141  FESwitch::basis(lfs.finiteElement()).evaluateFunction(x,yb);
142  y = 0;
143  for (unsigned int i=0; i<yb.size(); i++)
144  {
145  y.axpy(xl[i],yb[i]);
146  }
147  }
148 
150  inline const typename Traits::GridViewType& getGridView () const
151  {
152  return pgfs->gridView();
153  }
154 
155  private:
156 
159  typedef typename X::template ConstLocalView<LFSCache> XView;
160 
161  std::shared_ptr<GFS const> pgfs;
162  mutable LFS lfs;
163  mutable LFSCache lfs_cache;
164  mutable XView x_view;
165  mutable std::vector<typename Traits::RangeFieldType> xl;
166  mutable std::vector<typename Traits::RangeType> yb;
167  std::shared_ptr<const X> px; // FIXME: dummy pointer to make sure we take ownership of X
168  };
169 
179  template<typename T, typename X>
181  public GridFunctionInterface<
182  GridFunctionTraits<
183  typename T::Traits::GridViewType,
184  typename JacobianToCurl<typename T::Traits::FiniteElementType::
185  Traits::LocalBasisType::Traits::JacobianType>::CurlField,
186  JacobianToCurl<typename T::Traits::FiniteElementType::Traits::LocalBasisType::
187  Traits::JacobianType>::dimCurl,
188  typename JacobianToCurl<typename T::Traits::FiniteElementType::
189  Traits::LocalBasisType::Traits::JacobianType>::Curl
190  >,
191  DiscreteGridFunctionCurl<T,X>
192  >
193  {
194  typedef T GFS;
195  typedef typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
196  JacobianType Jacobian;
198 
199  public:
200  typedef GridFunctionTraits<
201  typename T::Traits::GridViewType,
202  typename J2C::CurlField, J2C::dimCurl, typename J2C::Curl
204 
210  DiscreteGridFunctionCurl(const GFS& gfs, const X& x_)
211  : pgfs(stackobject_to_shared_ptr(gfs))
212  , lfs(gfs)
213  , lfs_cache(lfs)
214  , x_view(x_)
215  , xl(gfs.maxLocalSize())
216  , jacobian(gfs.maxLocalSize())
217  , yb(gfs.maxLocalSize())
218  , px(stackobject_to_shared_ptr(x_))
219  {}
220 
221  // Evaluate
222  void evaluate (const typename Traits::ElementType& e,
223  const typename Traits::DomainType& x,
224  typename Traits::RangeType& y) const
225  {
226  static const J2C& j2C = J2C();
227 
228  lfs.bind();
229  lfs_cache.update();
230  x_view.bind(lfs_cache);
231  x_view.read(xl);
232  x_view.unbind();
233 
234  lfs.finiteElement().basis().evaluateJacobian(x,jacobian);
235 
236  y = 0;
237  for (std::size_t i=0; i < lfs.size(); i++) {
238  j2C(jacobian[i], yb);
239  y.axpy(xl[i], yb);
240  }
241  }
242 
244  const typename Traits::GridViewType& getGridView() const
245  { return pgfs->gridView(); }
246 
247  private:
249  BaseT;
252  typedef typename X::template ConstLocalView<LFSCache> XView;
253 
254  std::shared_ptr<GFS const> pgfs;
255  mutable LFS lfs;
256  mutable LFSCache lfs_cache;
257  mutable XView x_view;
258  mutable std::vector<typename Traits::RangeFieldType> xl;
259  mutable std::vector<Jacobian> jacobian;
260  mutable std::vector<typename Traits::RangeType> yb;
261  std::shared_ptr<const X> px; // FIXME: dummy pointer to make sure we take ownership of X
262  };
263 
265 
277  template<typename GV, typename RangeFieldType, int dimRangeOfBasis>
279  static_assert(AlwaysFalse<GV>::value,
280  "DiscreteGridFunctionCurl (and friends) work in 2D "
281  "and 3D only");
282  };
284 
290  template<typename GV, typename RangeFieldType>
292  : public GridFunctionTraits<GV,
293  RangeFieldType, 2,
294  FieldVector<RangeFieldType, 2> >
295  {
296  static_assert(GV::dimensionworld == 2,
297  "World dimension of grid must be 2 for the curl of a "
298  "scalar (1D) quantity");
299  };
301 
306  template<typename GV, typename RangeFieldType>
308  : public GridFunctionTraits<GV,
309  RangeFieldType, 1,
310  FieldVector<RangeFieldType, 1> >
311  {
312  static_assert(GV::dimensionworld == 2,
313  "World dimension of grid must be 2 for the curl of a"
314  "2D quantity");
315  };
317 
322  template<typename GV, typename RangeFieldType>
324  : public GridFunctionTraits<GV,
325  RangeFieldType, 3,
326  FieldVector<RangeFieldType, 3> >
327  {
328  static_assert(GV::dimensionworld == 3,
329  "World dimension of grid must be 3 for the curl of a"
330  "3D quantity");
331  };
332 
335 
349  template<typename T, typename X>
351  : public GridFunctionInterface<
352  DiscreteGridFunctionCurlTraits<
353  typename T::Traits::GridViewType,
354  typename T::Traits::FiniteElementType::Traits::
355  LocalBasisType::Traits::RangeFieldType,
356  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
357  dimRange>,
358  DiscreteGridFunctionGlobalCurl<T,X> >
359  {
360  public:
362  typename T::Traits::GridViewType,
363  typename T::Traits::FiniteElementType::Traits::
364  LocalBasisType::Traits::RangeFieldType,
365  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
366  dimRange> Traits;
367 
368  private:
369  typedef T GFS;
370  typedef GridFunctionInterface<
371  Traits,
373  typedef typename T::Traits::FiniteElementType::Traits::
374  LocalBasisType::Traits LBTraits;
375 
376  public:
382  DiscreteGridFunctionGlobalCurl (const GFS& gfs, const X& x_)
383  : pgfs(stackobject_to_shared_ptr(gfs))
384  , lfs(gfs)
385  , lfs_cache(lfs)
386  , x_view(x_)
387  , xl(gfs.maxLocalSize())
388  , J(gfs.maxLocalSize())
389  , px(stackobject_to_shared_ptr(x_))
390  {}
391 
392 
393  // Evaluate
394  inline void evaluate (const typename Traits::ElementType& e,
395  const typename Traits::DomainType& x,
396  typename Traits::RangeType& y) const
397  {
398  lfs.bind(e);
399  lfs_cache.update();
400  x_view.bind(lfs_cache);
401  x_view.read(xl);
402  x_view.unbind();
403 
404  lfs.finiteElement().localBasis().
405  evaluateJacobianGlobal(x,J,e.geometry());
406  y = 0;
407  for (unsigned int i=0; i<J.size(); i++)
408  // avoid a "case label value exceeds maximum value for type"
409  // warning: since dimRange is an anonymous enum, its type may
410  // contain only the values 0 and 1, resulting in a warning.
411  switch(unsigned(Traits::dimRange)) {
412  case 1:
413  y[0] += xl[i] * J[i][0][1];
414  y[1] += xl[i] * -J[i][0][0];
415  break;
416  case 2:
417  y[0] += xl[i]*(J[i][1][0] - J[i][0][1]);
418  break;
419  case 3:
420  y[0] += xl[i]*(J[i][2][1] - J[i][1][2]);
421  y[1] += xl[i]*(J[i][0][2] - J[i][2][0]);
422  y[2] += xl[i]*(J[i][1][0] - J[i][0][1]);
423  break;
424  default:
425  //how did that pass all the static asserts?
426  std::abort();
427  }
428  }
429 
431  inline const typename Traits::GridViewType& getGridView () const
432  {
433  return pgfs->gridView();
434  }
435 
436  private:
439  typedef typename X::template ConstLocalView<LFSCache> XView;
440 
441  std::shared_ptr<GFS const> pgfs;
442  mutable LFS lfs;
443  mutable LFSCache lfs_cache;
444  mutable XView x_view;
445  mutable std::vector<typename Traits::RangeFieldType> xl;
446  mutable std::vector<typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::JacobianType> J;
447  std::shared_ptr<const X> px; // FIXME: dummy pointer to make sure we take ownership of X
448  };
449 
452 
460  template<typename T, typename X>
462  : public GridFunctionInterface<
463  GridFunctionTraits<
464  typename T::Traits::GridViewType,
465  typename T::Traits::FiniteElementType::Traits::LocalBasisType
466  ::Traits::RangeFieldType,
467  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
468  ::dimDomain,
469  FieldVector<
470  typename T::Traits::FiniteElementType::Traits
471  ::LocalBasisType::Traits::RangeFieldType,
472  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
473  ::dimDomain> >,
474  DiscreteGridFunctionGradient<T,X> >
475  {
476  typedef T GFS;
477  typedef typename GFS::Traits::FiniteElementType::Traits::
478  LocalBasisType::Traits LBTraits;
479 
480  public:
481  typedef GridFunctionTraits<
482  typename GFS::Traits::GridViewType,
483  typename LBTraits::RangeFieldType,
484  LBTraits::dimDomain,
485  FieldVector<
486  typename LBTraits::RangeFieldType,
487  LBTraits::dimDomain> > Traits;
488 
489  private:
490  typedef GridFunctionInterface<
491  Traits,
493 
494  public:
500  DiscreteGridFunctionGradient (const GFS& gfs, const X& x_)
501  : pgfs(stackobject_to_shared_ptr(gfs))
502  , lfs(gfs)
503  , lfs_cache(lfs)
504  , x_view(x_)
505  , xl(lfs.size())
506  { }
507 
508  // Evaluate
509  inline void evaluate (const typename Traits::ElementType& e,
510  const typename Traits::DomainType& x,
511  typename Traits::RangeType& y) const
512  {
513  // get and bind local functions space
514  lfs.bind(e);
515  lfs_cache.update();
516  x_view.bind(lfs_cache);
517 
518  // get local coefficients
519  xl.resize(lfs.size());
520  x_view.read(xl);
521  x_view.unbind();
522 
523  // get Jacobian of geometry
524  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
525  JgeoIT = e.geometry().jacobianInverseTransposed(x);
526 
527  // get local Jacobians/gradients of the shape functions
528  std::vector<typename LBTraits::JacobianType> J(lfs.size());
529  lfs.finiteElement().localBasis().evaluateJacobian(x,J);
530 
531  typename Traits::RangeType gradphi;
532  y = 0;
533  for(unsigned int i = 0; i < lfs.size(); ++i) {
534  // compute global gradient of shape function i
535  gradphi = 0;
536  JgeoIT.umv(J[i][0], gradphi);
537 
538  // sum up global gradients, weighting them with the appropriate coeff
539  y.axpy(xl[i], gradphi);
540  }
541 
542  }
543 
545  inline const typename Traits::GridViewType& getGridView () const
546  {
547  return pgfs->gridView();
548  }
549 
550  private:
553  typedef typename X::template ConstLocalView<LFSCache> XView;
554 
555  std::shared_ptr<GFS const> pgfs;
556  mutable LFS lfs;
557  mutable LFSCache lfs_cache;
558  mutable XView x_view;
559  mutable std::vector<typename Traits::RangeFieldType> xl;
560  };
561 
566  template<typename T, typename X>
568  : public GridFunctionInterface<
569  GridFunctionTraits<
570  typename T::Traits::GridViewType,
571  typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
572  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
573  typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType
574  >,
575  DiscreteGridFunctionPiola<T,X>
576  >
577  {
578  typedef T GFS;
579 
580  typedef GridFunctionInterface<
582  typename T::Traits::GridViewType,
583  typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
584  T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
585  typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType
586  >,
588  > BaseT;
589 
590  public:
591  typedef typename BaseT::Traits Traits;
592 
597  DiscreteGridFunctionPiola (const GFS& gfs, const X& x_)
598  : pgfs(stackobject_to_shared_ptr(gfs))
599  , lfs(gfs)
600  , lfs_cache(lfs)
601  , x_view(x_)
602  , xl(pgfs->maxLocalSize())
603  , yb(pgfs->maxLocalSize())
604  {
605  }
606 
607  inline void evaluate (const typename Traits::ElementType& e,
608  const typename Traits::DomainType& x,
609  typename Traits::RangeType& y) const
610  {
611  // evaluate shape function on the reference element as before
612  lfs.bind(e);
613  lfs_cache.update();
614  x_view.bind(lfs_cache);
615  x_view.read(xl);
616  x_view.unbind();
617 
618  lfs.finiteElement().localBasis().evaluateFunction(x,yb);
619  typename Traits::RangeType yhat;
620  yhat = 0;
621  for (unsigned int i=0; i<yb.size(); i++)
622  yhat.axpy(xl[i],yb[i]);
623 
624  // apply Piola transformation
625  typename Traits::ElementType::Geometry::JacobianInverseTransposed
626  J = e.geometry().jacobianInverseTransposed(x);
627  J.invert();
628  y = 0;
629  J.umtv(yhat,y);
630  y /= J.determinant();
631  }
632 
634  inline const typename Traits::GridViewType& getGridView () const
635  {
636  return pgfs->gridView();
637  }
638 
639  private:
640 
643  typedef typename X::template ConstLocalView<LFSCache> XView;
644 
645  std::shared_ptr<GFS const> pgfs;
646  mutable LFS lfs;
647  mutable LFSCache lfs_cache;
648  mutable XView x_view;
649  mutable std::vector<typename Traits::RangeFieldType> xl;
650  mutable std::vector<typename Traits::RangeType> yb;
651 
652  };
653 
665 #ifdef __clang__
666  // clang is too stupid to correctly apply the constexpr qualifier of staticDegree in this context
668 #else
670 #endif
672  : public GridFunctionInterface<
673  GridFunctionTraits<
674  typename T::Traits::GridViewType,
675  typename T::template Child<0>::Type::Traits::FiniteElementType
676  ::Traits::LocalBasisType::Traits::RangeFieldType,
677  dimR,
678  Dune::FieldVector<
679  typename T::template Child<0>::Type::Traits::FiniteElementType
680  ::Traits::LocalBasisType::Traits::RangeFieldType,
681  dimR
682  >
683  >,
684  VectorDiscreteGridFunction<T,X>
685  >,
686  public TypeTree::LeafNode
687  {
688  typedef T GFS;
689 
690  typedef GridFunctionInterface<
692  typename T::Traits::GridViewType,
693  typename T::template Child<0>::Type::Traits::FiniteElementType
694  ::Traits::LocalBasisType::Traits::RangeFieldType,
695  dimR,
696  Dune::FieldVector<
697  typename T::template Child<0>::Type::Traits::FiniteElementType
698  ::Traits::LocalBasisType::Traits::RangeFieldType,
699  dimR
700  >
701  >,
703  > BaseT;
704 
705  public:
706  typedef typename BaseT::Traits Traits;
707  typedef typename T::template Child<0>::Type ChildType;
708  typedef typename ChildType::Traits::FiniteElementType
709  ::Traits::LocalBasisType::Traits::RangeFieldType RF;
710  typedef typename ChildType::Traits::FiniteElementType
711  ::Traits::LocalBasisType::Traits::RangeType RT;
712 
714 
719  VectorDiscreteGridFunction(const GFS& gfs, const X& x_,
720  std::size_t start = 0)
721  : pgfs(stackobject_to_shared_ptr(gfs))
722  , lfs(gfs)
723  , lfs_cache(lfs)
724  , x_view(x_)
725  , xl(gfs.maxLocalSize())
726  , yb(gfs.maxLocalSize())
727  {
728  for(std::size_t i = 0; i < dimR; ++i)
729  remap[i] = i + start;
730  }
731 
733 
744  template<class Remap>
745  VectorDiscreteGridFunction(const GFS& gfs, const X& x_,
746  const Remap &remap_)
747  : pgfs(stackobject_to_shared_ptr(gfs))
748  , lfs(gfs)
749  , lfs_cache(lfs)
750  , x_view(x_)
751  , xl(gfs.maxLocalSize())
752  , yb(gfs.maxLocalSize())
753  , px(stackobject_to_shared_ptr(x_))
754  {
755  for(std::size_t i = 0; i < dimR; ++i)
756  remap[i] = remap_[i];
757  }
758 
759  inline void evaluate (const typename Traits::ElementType& e,
760  const typename Traits::DomainType& x,
761  typename Traits::RangeType& y) const
762  {
763  lfs.bind(e);
764  lfs_cache.update();
765  x_view.bind(lfs_cache);
766  x_view.read(xl);
767  x_view.unbind();
768  for (unsigned int k=0; k < dimR; k++)
769  {
770  lfs.child(remap[k]).finiteElement().localBasis().
771  evaluateFunction(x,yb);
772  y[k] = 0.0;
773  for (unsigned int i=0; i<yb.size(); i++)
774  y[k] += xl[lfs.child(remap[k]).localIndex(i)]*yb[i];
775  }
776  }
777 
779  inline const typename Traits::GridViewType& getGridView () const
780  {
781  return pgfs->gridView();
782  }
783 
784  private:
787  typedef typename X::template ConstLocalView<LFSCache> XView;
788 
789  std::shared_ptr<GFS const> pgfs;
790  std::size_t remap[dimR];
791  mutable LFS lfs;
792  mutable LFSCache lfs_cache;
793  mutable XView x_view;
794  mutable std::vector<RF> xl;
795  mutable std::vector<RT> yb;
796  std::shared_ptr<const X> px; // FIXME: dummy pointer to make sure we take ownership of X
797  };
798 
804  template<typename T, typename X>
806  : public GridFunctionInterface<
807  GridFunctionTraits<
808  typename T::Traits::GridViewType,
809  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
810  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain,
811  TypeTree::StaticDegree<T>::value,
812  Dune::FieldMatrix<
813  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
814  TypeTree::StaticDegree<T>::value,
815  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain
816  >
817  >,
818  VectorDiscreteGridFunctionGradient<T,X>
819  >,
820  public TypeTree::LeafNode
821  {
822  typedef T GFS;
823 
824  typedef GridFunctionInterface<
826  typename T::Traits::GridViewType,
827  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
828  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain,
830  Dune::FieldMatrix<
831  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
832  TypeTree::StaticDegree<T>::value,
833  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain>
834  >,
836  > BaseT;
837 
838  public:
839  typedef typename BaseT::Traits Traits;
840  typedef typename T::template Child<0>::Type ChildType;
841  typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits;
842 
843  typedef typename LBTraits::RangeFieldType RF;
844  typedef typename LBTraits::JacobianType JT;
845 
846  VectorDiscreteGridFunctionGradient (const GFS& gfs, const X& x_)
847  : pgfs(stackobject_to_shared_ptr(gfs))
848  , lfs(gfs)
849  , lfs_cache(lfs)
850  , x_view(x_)
851  , xl(gfs.maxLocalSize())
852  , J(gfs.maxLocalSize())
853  {
854  }
855 
856  inline void evaluate(const typename Traits::ElementType& e,
857  const typename Traits::DomainType& x,
858  typename Traits::RangeType& y) const
859  {
860  // get and bind local functions space
861  lfs.bind(e);
862  lfs_cache.update();
863  x_view.bind(lfs_cache);
864  x_view.read(xl);
865  x_view.unbind();
866 
867  // get Jacobian of geometry
868  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
869  JgeoIT = e.geometry().jacobianInverseTransposed(x);
870 
871  y = 0.0;
872 
873  // Loop over PowerLFS and calculate gradient for each child separately
874  for(unsigned int k = 0; k != TypeTree::degree(lfs); ++k)
875  {
876  // get local Jacobians/gradients of the shape functions
877  std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
878  lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
879 
880  Dune::FieldVector<RF,LBTraits::dimDomain> gradphi;
881  for (typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++)
882  {
883  gradphi = 0;
884  JgeoIT.umv(J[i][0], gradphi);
885 
886  y[k].axpy(xl[lfs.child(k).localIndex(i)], gradphi);
887  }
888  }
889  }
890 
891 
893  inline const typename Traits::GridViewType& getGridView () const
894  {
895  return pgfs->gridView();
896  }
897 
898  private:
901  typedef typename X::template ConstLocalView<LFSCache> XView;
902 
903  std::shared_ptr<GFS const> pgfs;
904  mutable LFS lfs;
905  mutable LFSCache lfs_cache;
906  mutable XView x_view;
907  mutable std::vector<RF> xl;
908  mutable std::vector<JT> J;
909  std::shared_ptr<const X> px; // FIXME: dummy pointer to make sure we take ownership of X
910  };
911 
925  template<typename Mat, typename RF, std::size_t size>
931  template<typename T>
932  static inline RF compute_derivative(const Mat& mat, const T& t, const unsigned int k)
933  {
934  // setting this to zero is just a test if the template specialization work
935  Dune::FieldVector<RF,size> grad_phi(0.0);
936  mat.umv(t,grad_phi);
937  return grad_phi[k];
938  // return 0.0;
939  }
940  };
941 
950  template<typename RF, std::size_t size>
951  struct SingleDerivativeComputationHelper<Dune::FieldMatrix<RF,size,size>,RF,size> {
957  template<typename T>
958  static inline RF compute_derivative(const Dune::FieldMatrix<RF,size,size>& mat, const T& t, const unsigned int k)
959  {
960  return mat[k]*t;
961  }
962  };
963 
974  template<typename RF, std::size_t size>
975  struct SingleDerivativeComputationHelper<Dune::DiagonalMatrix<RF,size>,RF,size> {
981  template<typename T>
982  static inline RF compute_derivative(const Dune::DiagonalMatrix<RF,size>& mat, const T& t, const unsigned int k)
983  {
984  return mat[k][k]*t[k];
985  }
986  };
987 
998  template<typename T, typename X>
1001  Dune::PDELab::GridFunctionTraits<
1002  typename T::Traits::GridViewType,
1003  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1004  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1005  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1006  VectorDiscreteGridFunctionDiv<T,X> >
1007  , public TypeTree::LeafNode
1008  {
1009  typedef T GFS;
1010 
1013  typename T::Traits::GridViewType,
1014  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1015  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1016  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1018  public :
1019  typedef typename BaseT::Traits Traits;
1020  typedef typename T::template Child<0>::Type ChildType;
1021  typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits;
1022 
1023  typedef typename LBTraits::RangeFieldType RF;
1024  typedef typename LBTraits::JacobianType JT;
1025 
1026  VectorDiscreteGridFunctionDiv(const GFS& gfs, const X& x_)
1027  : pgfs(stackobject_to_shared_ptr(gfs))
1028  , lfs(gfs)
1029  , lfs_cache(lfs)
1030  , x_view(x_)
1031  , xl(gfs.maxLocalSize())
1032  , J(gfs.maxLocalSize())
1033  {
1034  static_assert(LBTraits::dimDomain == TypeTree::StaticDegree<T>::value,
1035  "dimDomain and number of children has to be the same");
1036  }
1037 
1038  inline void evaluate(const typename Traits::ElementType& e,
1039  const typename Traits::DomainType& x,
1040  typename Traits::RangeType& y) const
1041  {
1042  // get and bind local functions space
1043  lfs.bind(e);
1044  lfs_cache.update();
1045  x_view.bind(lfs_cache);
1046  x_view.read(xl);
1047  x_view.unbind();
1048 
1049  // get Jacobian of geometry
1050  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1051  JgeoIT = e.geometry().jacobianInverseTransposed(x);
1052 
1053  const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1054  Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1055 
1056  y = 0.0;
1057 
1058  // loop over VectorGFS and calculate k-th derivative of k-th child
1059  for(unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1060 
1061  // get local Jacobians/gradients of the shape functions
1062  std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1063  lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1064 
1065  RF d_k_phi;
1066  for(typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1067  // compute k-th derivative of k-th child
1068  d_k_phi =
1070  typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1071  typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1072  N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1073  (JgeoIT,J[i][0],k);
1074 
1075  y += xl[lfs.child(k).localIndex(i)] * d_k_phi;
1076  }
1077  }
1078  }
1079 
1081  inline const typename Traits::GridViewType& getGridView() const
1082  {
1083  return pgfs->gridView();
1084  }
1085 
1086  private :
1089  typedef typename X::template ConstLocalView<LFSCache> XView;
1090 
1091  shared_ptr<GFS const> pgfs;
1092  mutable LFS lfs;
1093  mutable LFSCache lfs_cache;
1094  mutable XView x_view;
1095  mutable std::vector<RF> xl;
1096  mutable std::vector<JT> J;
1097  shared_ptr<const X> px;
1098  }; // end class VectorDiscreteGridFunctionDiv
1099 
1114  {
1115  typedef T GFS;
1116  public :
1117  VectorDiscreteGridFunctionCurl(const GFS& gfs, const X& x)
1118  {
1120  "Curl computation can only be done in two or three dimensions");
1121  }
1122  };
1123 
1134  template<typename T, typename X>
1137  Dune::PDELab::GridFunctionTraits<
1138  typename T::Traits::GridViewType,
1139  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1140  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1141  TypeTree::StaticDegree<T>::value,
1142  Dune::FieldVector<
1143  typename T::template Child<0>::Type::Traits::FiniteElementType
1144  ::Traits::LocalBasisType::Traits::RangeFieldType,
1145  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange
1146  TypeTree::StaticDegree<T>::value
1147  >
1148  >,
1149  VectorDiscreteGridFunctionCurl<T,X>
1150  >
1151  , public TypeTree::LeafNode
1152  {
1153  typedef T GFS;
1154 
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,
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  >,
1169 
1170  public :
1171  typedef typename BaseT::Traits Traits;
1172  typedef typename T::template Child<0>::Type ChildType;
1173  typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits;
1174 
1175  typedef typename LBTraits::RangeFieldType RF;
1176  typedef typename LBTraits::JacobianType JT;
1177 
1178  VectorDiscreteGridFunctionCurl(const GFS& gfs, const X& x_)
1179  : pgfs(stackobject_to_shared_ptr(gfs))
1180  , lfs(gfs)
1181  , lfs_cache(lfs)
1182  , x_view(x_)
1183  , xl(gfs.maxLocalSize())
1184  , J(gfs.maxLocalSize())
1185  {
1186  static_assert(LBTraits::dimDomain == TypeTree::StaticDegree<T>::value,
1187  "dimDomain and number of children has to be the same");
1188  }
1189 
1190  inline void evaluate(const typename Traits::ElementType& e,
1191  const typename Traits::DomainType& x,
1192  typename Traits::RangeType& y) const
1193  {
1194  // get and bind local functions space
1195  lfs.bind(e);
1196  lfs_cache.update();
1197  x_view.bind(lfs_cache);
1198  x_view.read(xl);
1199  x_view.unbind();
1200 
1201  // get Jacobian of geometry
1202  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1203  JgeoIT = e.geometry().jacobianInverseTransposed(x);
1204 
1205  const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1206  Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1207 
1208  y = 0.0;
1209 
1210  // some handy variables for the curl in 3D
1211  int i1, i2;
1212 
1213  // loop over childs of VectorGFS
1214  for(unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1215 
1216  // get local Jacobians/gradients of the shape functions
1217  std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1218  lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1219 
1220  // pick out the right derivative and component for the curl computation
1221  i1 = (k+1)%3;
1222  i2 = (k+2)%3;
1223 
1224  RF d_k_phi;
1225  for(typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1226  // compute i2-th derivative of k-th child
1227  d_k_phi =
1229  typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1230  typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1231  N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1232  (JgeoIT,J[i][0],i2);
1233 
1234  y[i1] += xl[lfs.child(k).localIndex(i)] * d_k_phi;
1235 
1236  // compute i1-th derivative of k-th child
1237  d_k_phi =
1239  typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1240  typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1241  N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1242  (JgeoIT,J[i][0],i1);
1243 
1244  y[i2] -= xl[lfs.child(k).localIndex(i)] * d_k_phi;
1245  }
1246  }
1247  }
1248 
1250  inline const typename Traits::GridViewType& getGridView() const
1251  {
1252  return pgfs->gridView();
1253  }
1254 
1255  private :
1258  typedef typename X::template ConstLocalView<LFSCache> XView;
1259 
1260  shared_ptr<GFS const> pgfs;
1261  mutable LFS lfs;
1262  mutable LFSCache lfs_cache;
1263  mutable XView x_view;
1264  mutable std::vector<RF> xl;
1265  mutable std::vector<JT> J;
1266  shared_ptr<const X> px;
1267  }; // end class VectorDiscreteGridFunctionCurl (3D)
1268 
1279  template<typename T, typename X>
1282  Dune::PDELab::GridFunctionTraits<
1283  typename T::Traits::GridViewType,
1284  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1285  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1286  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1287  VectorDiscreteGridFunctionDiv<T,X> >
1288  , public TypeTree::LeafNode
1289  {
1290  typedef T GFS;
1291 
1294  typename T::Traits::GridViewType,
1295  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1296  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1297  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1299  public :
1300  typedef typename BaseT::Traits Traits;
1301  typedef typename T::template Child<0>::Type ChildType;
1302  typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits;
1303 
1304  typedef typename LBTraits::RangeFieldType RF;
1305  typedef typename LBTraits::JacobianType JT;
1306 
1307  VectorDiscreteGridFunctionCurl(const GFS& gfs, const X& x_)
1308  : pgfs(stackobject_to_shared_ptr(gfs))
1309  , lfs(gfs)
1310  , lfs_cache(lfs)
1311  , x_view(x_)
1312  , xl(gfs.maxLocalSize())
1313  , J(gfs.maxLocalSize())
1314  {
1315  static_assert(LBTraits::dimDomain == TypeTree::StaticDegree<T>::value,
1316  "dimDomain and number of children has to be the same");
1317  }
1318 
1319  inline void evaluate(const typename Traits::ElementType& e,
1320  const typename Traits::DomainType& x,
1321  typename Traits::RangeType& y) const
1322  {
1323  // get and bind local functions space
1324  lfs.bind(e);
1325  lfs_cache.update();
1326  x_view.bind(lfs_cache);
1327  x_view.read(xl);
1328  x_view.unbind();
1329 
1330  // get Jacobian of geometry
1331  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1332  JgeoIT = e.geometry().jacobianInverseTransposed(x);
1333 
1334  const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1335  Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1336 
1337  y = 0.0;
1338 
1339  // some handy variables for the curl computation in 2D
1340  RF sign = -1.0;
1341  int i2;
1342 
1343  // loop over childs of VectorGFS
1344  for(unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1345 
1346  // get local Jacobians/gradients of the shape functions
1347  std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1348  lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1349 
1350  RF d_k_phi;
1351  // pick out the right derivative
1352  i2 = 1-k;
1353  for(typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1354  // compute i2-th derivative of k-th child
1355  d_k_phi =
1357  typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1358  typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1359  N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1360  (JgeoIT,J[i][0],i2);
1361 
1362  y += sign * xl[lfs.child(k).localIndex(i)] * d_k_phi;
1363  }
1364  sign *= -1.0;
1365  }
1366  }
1367 
1369  inline const typename Traits::GridViewType& getGridView() const
1370  {
1371  return pgfs->gridView();
1372  }
1373 
1374  private :
1377  typedef typename X::template ConstLocalView<LFSCache> XView;
1378 
1379  shared_ptr<GFS const> pgfs;
1380  mutable LFS lfs;
1381  mutable LFSCache lfs_cache;
1382  mutable XView x_view;
1383  mutable std::vector<RF> xl;
1384  mutable std::vector<JT> J;
1385  shared_ptr<const X> px;
1386  }; // end class VectorDiscreteGridFunctionCurl (2D)
1387 
1389  } // namespace PDELab
1390 } // namespace Dune
1391 
1392 #endif // DUNE_PDELAB_GRIDFUNCTIONSPACE_GRIDFUNCTIONSPACEUTILITIES_HH
DiscreteGridFunction with Piola transformation.
Definition: gridfunctionspaceutilities.hh:567
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1300
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:759
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1020
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1172
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1021
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1369
static RF compute_derivative(const Dune::DiagonalMatrix< RF, size > &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:982
Equivalent of DiscreteGridFunctionGradient for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:805
DiscreteGridFunction(std::shared_ptr< const GFS > gfs, std::shared_ptr< const X > x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:117
Dune::FieldVector< GV::Grid::ctype, GV::dimension > DomainType
domain type in dim-size coordinates
Definition: function.hh:48
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:394
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1019
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:509
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:591
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1176
DiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:210
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1175
Compute curl of vector-valued functions.
Definition: gridfunctionspaceutilities.hh:1113
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:841
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:844
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:839
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:129
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1302
traits class holding the function signature, same as in local function
Definition: function.hh:175
ChildType::Traits::FiniteElementType ::Traits::LocalBasisType::Traits::RangeType RT
Definition: gridfunctionspaceutilities.hh:711
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1178
GridFunctionTraits< typename GFS::Traits::GridViewType, typename LBTraits::RangeFieldType, LBTraits::dimDomain, FieldVector< typename LBTraits::RangeFieldType, LBTraits::dimDomain > > Traits
Definition: gridfunctionspaceutilities.hh:487
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:607
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:856
DiscreteGridFunctionCurlTraits< typename T::Traits::GridViewType, typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType, T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange > Traits
Definition: gridfunctionspaceutilities.hh:366
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:779
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1081
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:244
DiscreteGridFunctionPiola(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionPiola.
Definition: gridfunctionspaceutilities.hh:597
Create a local function space from a global function space.
Definition: localfunctionspace.hh:670
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1173
convert a grid function space and a coefficient vector into a grid function
Definition: gridfunctionspaceutilities.hh:53
DiscreteGridFunction for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:671
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1038
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:431
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1307
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:843
static RF compute_derivative(const Dune::FieldMatrix< RF, size, size > &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:958
Compute divergence of vector-valued functions.
Definition: gridfunctionspaceutilities.hh:999
const Entity & e
Definition: localfunctionspace.hh:111
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:150
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
Helper class to compute a single derivative of scalar basis functions.
Definition: gridfunctionspaceutilities.hh:926
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:95
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1305
DiscreteGridFunction(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:102
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:545
DiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGradient.
Definition: gridfunctionspaceutilities.hh:500
ChildType::Traits::FiniteElementType ::Traits::LocalBasisType::Traits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:709
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:840
VectorDiscreteGridFunctionDiv(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1026
GridFunctionTraits< typename T::Traits::GridViewType, typename J2C::CurlField, J2C::dimCurl, typename J2C::Curl > Traits
Definition: gridfunctionspaceutilities.hh:203
T Traits
Export type traits.
Definition: function.hh:191
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:222
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1190
Helper class to calculate the Traits of DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:278
a GridFunction maps x in DomainType to y in RangeType
Definition: function.hh:186
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1304
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:706
extract the curl of a function from the jacobian of that function
Definition: jacobiantocurl.hh:27
convert a single component function space with a grid function representing the gradient ...
Definition: gridfunctionspaceutilities.hh:461
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1171
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:893
convert a grid function space and a coefficient vector into a grid function of the curl ...
Definition: gridfunctionspaceutilities.hh:180
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, const Remap &remap_)
construct
Definition: gridfunctionspaceutilities.hh:745
static RF compute_derivative(const Mat &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:932
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1319
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:707
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:117
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1023
DiscreteGridFunctionGlobalCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGlobalCurl.
Definition: gridfunctionspaceutilities.hh:382
GV GridViewType
The type of the grid view the function lives on.
Definition: function.hh:114
VectorDiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:846
void update()
Definition: lfsindexcache.hh:300
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1024
convert a single component function space with experimental global finite elements into a grid functi...
Definition: gridfunctionspaceutilities.hh:350
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, std::size_t start=0)
construct
Definition: gridfunctionspaceutilities.hh:719
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1250
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x)
Definition: gridfunctionspaceutilities.hh:1117
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1301
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:634