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