dune-pdelab  2.4.1
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_GRIDFUNCTIONSPACEUTILITIES_HH
4 #define DUNE_PDELAB_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  template<typename T, typename X, std::size_t dimR = T::CHILDREN>
667  : public GridFunctionInterface<
668  GridFunctionTraits<
669  typename T::Traits::GridViewType,
670  typename T::template Child<0>::Type::Traits::FiniteElementType
671  ::Traits::LocalBasisType::Traits::RangeFieldType,
672  dimR,
673  Dune::FieldVector<
674  typename T::template Child<0>::Type::Traits::FiniteElementType
675  ::Traits::LocalBasisType::Traits::RangeFieldType,
676  dimR
677  >
678  >,
679  VectorDiscreteGridFunction<T,X>
680  >,
681  public TypeTree::LeafNode
682  {
683  typedef T GFS;
684 
685  typedef GridFunctionInterface<
687  typename T::Traits::GridViewType,
688  typename T::template Child<0>::Type::Traits::FiniteElementType
689  ::Traits::LocalBasisType::Traits::RangeFieldType,
690  dimR,
691  Dune::FieldVector<
692  typename T::template Child<0>::Type::Traits::FiniteElementType
693  ::Traits::LocalBasisType::Traits::RangeFieldType,
694  dimR
695  >
696  >,
698  > BaseT;
699 
700  public:
701  typedef typename BaseT::Traits Traits;
702  typedef typename T::template Child<0>::Type ChildType;
703  typedef typename ChildType::Traits::FiniteElementType
704  ::Traits::LocalBasisType::Traits::RangeFieldType RF;
705  typedef typename ChildType::Traits::FiniteElementType
706  ::Traits::LocalBasisType::Traits::RangeType RT;
707 
709 
714  VectorDiscreteGridFunction(const GFS& gfs, const X& x_,
715  std::size_t start = 0)
716  : pgfs(stackobject_to_shared_ptr(gfs))
717  , lfs(gfs)
718  , lfs_cache(lfs)
719  , x_view(x_)
720  , xl(gfs.maxLocalSize())
721  , yb(gfs.maxLocalSize())
722  {
723  for(std::size_t i = 0; i < dimR; ++i)
724  remap[i] = i + start;
725  }
726 
728 
739  template<class Remap>
740  VectorDiscreteGridFunction(const GFS& gfs, const X& x_,
741  const Remap &remap_)
742  : pgfs(stackobject_to_shared_ptr(gfs))
743  , lfs(gfs)
744  , lfs_cache(lfs)
745  , x_view(x_)
746  , xl(gfs.maxLocalSize())
747  , yb(gfs.maxLocalSize())
748  , px(stackobject_to_shared_ptr(x_))
749  {
750  for(std::size_t i = 0; i < dimR; ++i)
751  remap[i] = remap_[i];
752  }
753 
754  inline void evaluate (const typename Traits::ElementType& e,
755  const typename Traits::DomainType& x,
756  typename Traits::RangeType& y) const
757  {
758  lfs.bind(e);
759  lfs_cache.update();
760  x_view.bind(lfs_cache);
761  x_view.read(xl);
762  x_view.unbind();
763  for (unsigned int k=0; k < dimR; k++)
764  {
765  lfs.child(remap[k]).finiteElement().localBasis().
766  evaluateFunction(x,yb);
767  y[k] = 0.0;
768  for (unsigned int i=0; i<yb.size(); i++)
769  y[k] += xl[lfs.child(remap[k]).localIndex(i)]*yb[i];
770  }
771  }
772 
774  inline const typename Traits::GridViewType& getGridView () const
775  {
776  return pgfs->gridView();
777  }
778 
779  private:
782  typedef typename X::template ConstLocalView<LFSCache> XView;
783 
784  std::shared_ptr<GFS const> pgfs;
785  std::size_t remap[dimR];
786  mutable LFS lfs;
787  mutable LFSCache lfs_cache;
788  mutable XView x_view;
789  mutable std::vector<RF> xl;
790  mutable std::vector<RT> yb;
791  std::shared_ptr<const X> px; // FIXME: dummy pointer to make sure we take ownership of X
792  };
793 
799  template<typename T, typename X>
801  : public GridFunctionInterface<
802  GridFunctionTraits<
803  typename T::Traits::GridViewType,
804  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
805  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain,
806  T::CHILDREN,
807  Dune::FieldMatrix<
808  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
809  T::CHILDREN,
810  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain
811  >
812  >,
813  VectorDiscreteGridFunctionGradient<T,X>
814  >,
815  public TypeTree::LeafNode
816  {
817  typedef T GFS;
818 
819  typedef GridFunctionInterface<
821  typename T::Traits::GridViewType,
822  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
823  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain,
824  T::CHILDREN,
825  Dune::FieldMatrix<
826  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
827  T::CHILDREN,
828  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain>
829  >,
831  > BaseT;
832 
833  public:
834  typedef typename BaseT::Traits Traits;
835  typedef typename T::template Child<0>::Type ChildType;
836  typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits;
837 
838  typedef typename LBTraits::RangeFieldType RF;
839  typedef typename LBTraits::JacobianType JT;
840 
841  VectorDiscreteGridFunctionGradient (const GFS& gfs, const X& x_)
842  : pgfs(stackobject_to_shared_ptr(gfs))
843  , lfs(gfs)
844  , lfs_cache(lfs)
845  , x_view(x_)
846  , xl(gfs.maxLocalSize())
847  , J(gfs.maxLocalSize())
848  {
849  }
850 
851  inline void evaluate(const typename Traits::ElementType& e,
852  const typename Traits::DomainType& x,
853  typename Traits::RangeType& y) const
854  {
855  // get and bind local functions space
856  lfs.bind(e);
857  lfs_cache.update();
858  x_view.bind(lfs_cache);
859  x_view.read(xl);
860  x_view.unbind();
861 
862  // get Jacobian of geometry
863  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
864  JgeoIT = e.geometry().jacobianInverseTransposed(x);
865 
866  y = 0.0;
867 
868  // Loop over PowerLFS and calculate gradient for each child separately
869  for(unsigned int k = 0; k != T::CHILDREN; ++k)
870  {
871  // get local Jacobians/gradients of the shape functions
872  std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
873  lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
874 
875  Dune::FieldVector<RF,LBTraits::dimDomain> gradphi;
876  for (typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++)
877  {
878  gradphi = 0;
879  JgeoIT.umv(J[i][0], gradphi);
880 
881  y[k].axpy(xl[lfs.child(k).localIndex(i)], gradphi);
882  }
883  }
884  }
885 
886 
888  inline const typename Traits::GridViewType& getGridView () const
889  {
890  return pgfs->gridView();
891  }
892 
893  private:
896  typedef typename X::template ConstLocalView<LFSCache> XView;
897 
898  std::shared_ptr<GFS const> pgfs;
899  mutable LFS lfs;
900  mutable LFSCache lfs_cache;
901  mutable XView x_view;
902  mutable std::vector<RF> xl;
903  mutable std::vector<JT> J;
904  std::shared_ptr<const X> px; // FIXME: dummy pointer to make sure we take ownership of X
905  };
906 
920  template<typename Mat, typename RF, std::size_t size>
926  template<typename T>
927  static inline RF compute_derivative(const Mat& mat, const T& t, const unsigned int k)
928  {
929  // setting this to zero is just a test if the template specialization work
930  Dune::FieldVector<RF,size> grad_phi(0.0);
931  mat.umv(t,grad_phi);
932  return grad_phi[k];
933  // return 0.0;
934  }
935  };
936 
945  template<typename RF, std::size_t size>
946  struct SingleDerivativeComputationHelper<Dune::FieldMatrix<RF,size,size>,RF,size> {
952  template<typename T>
953  static inline RF compute_derivative(const Dune::FieldMatrix<RF,size,size>& mat, const T& t, const unsigned int k)
954  {
955  return mat[k]*t;
956  }
957  };
958 
969  template<typename RF, std::size_t size>
970  struct SingleDerivativeComputationHelper<Dune::DiagonalMatrix<RF,size>,RF,size> {
976  template<typename T>
977  static inline RF compute_derivative(const Dune::DiagonalMatrix<RF,size>& mat, const T& t, const unsigned int k)
978  {
979  return mat[k][k]*t[k];
980  }
981  };
982 
993  template<typename T, typename X>
996  Dune::PDELab::GridFunctionTraits<
997  typename T::Traits::GridViewType,
998  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
999  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1000  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1001  VectorDiscreteGridFunctionDiv<T,X> >
1002  , public TypeTree::LeafNode
1003  {
1004  typedef T GFS;
1005 
1008  typename T::Traits::GridViewType,
1009  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1010  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1011  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1013  public :
1014  typedef typename BaseT::Traits Traits;
1015  typedef typename T::template Child<0>::Type ChildType;
1016  typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits;
1017 
1018  typedef typename LBTraits::RangeFieldType RF;
1019  typedef typename LBTraits::JacobianType JT;
1020 
1021  VectorDiscreteGridFunctionDiv(const GFS& gfs, const X& x_)
1022  : pgfs(stackobject_to_shared_ptr(gfs))
1023  , lfs(gfs)
1024  , lfs_cache(lfs)
1025  , x_view(x_)
1026  , xl(gfs.maxLocalSize())
1027  , J(gfs.maxLocalSize())
1028  {
1029  static_assert(LBTraits::dimDomain == T::CHILDREN,
1030  "dimDomain and number of children has to be the same");
1031  }
1032 
1033  inline void evaluate(const typename Traits::ElementType& e,
1034  const typename Traits::DomainType& x,
1035  typename Traits::RangeType& y) const
1036  {
1037  // get and bind local functions space
1038  lfs.bind(e);
1039  lfs_cache.update();
1040  x_view.bind(lfs_cache);
1041  x_view.read(xl);
1042  x_view.unbind();
1043 
1044  // get Jacobian of geometry
1045  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1046  JgeoIT = e.geometry().jacobianInverseTransposed(x);
1047 
1048  const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1049  Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1050 
1051  y = 0.0;
1052 
1053  // loop over VectorGFS and calculate k-th derivative of k-th child
1054  for(unsigned int k=0; k != T::CHILDREN; ++k) {
1055 
1056  // get local Jacobians/gradients of the shape functions
1057  std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1058  lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1059 
1060  RF d_k_phi;
1061  for(typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1062  // compute k-th derivative of k-th child
1063  d_k_phi =
1065  typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1066  typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1067  N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1068  (JgeoIT,J[i][0],k);
1069 
1070  y += xl[lfs.child(k).localIndex(i)] * d_k_phi;
1071  }
1072  }
1073  }
1074 
1076  inline const typename Traits::GridViewType& getGridView() const
1077  {
1078  return pgfs->gridView();
1079  }
1080 
1081  private :
1084  typedef typename X::template ConstLocalView<LFSCache> XView;
1085 
1086  shared_ptr<GFS const> pgfs;
1087  mutable LFS lfs;
1088  mutable LFSCache lfs_cache;
1089  mutable XView x_view;
1090  mutable std::vector<RF> xl;
1091  mutable std::vector<JT> J;
1092  shared_ptr<const X> px;
1093  }; // end class VectorDiscreteGridFunctionDiv
1094 
1107  template<typename T, typename X, std::size_t dimR = T::CHILDREN>
1109  {
1110  typedef T GFS;
1111  public :
1112  VectorDiscreteGridFunctionCurl(const GFS& gfs, const X& x)
1113  {
1115  "Curl computation can only be done in two or three dimensions");
1116  }
1117  };
1118 
1129  template<typename T, typename X>
1132  Dune::PDELab::GridFunctionTraits<
1133  typename T::Traits::GridViewType,
1134  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1135  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1136  T::CHILDREN,
1137  Dune::FieldVector<
1138  typename T::template Child<0>::Type::Traits::FiniteElementType
1139  ::Traits::LocalBasisType::Traits::RangeFieldType,
1140  //T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange
1141  T::CHILDREN
1142  >
1143  >,
1144  VectorDiscreteGridFunctionCurl<T,X>
1145  >
1146  , public TypeTree::LeafNode
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,
1155  T::CHILDREN,
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  T::CHILDREN
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 == T::CHILDREN,
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 != T::CHILDREN; ++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  , public TypeTree::LeafNode
1284  {
1285  typedef T GFS;
1286 
1289  typename T::Traits::GridViewType,
1290  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1291  T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1292  typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1294  public :
1295  typedef typename BaseT::Traits Traits;
1296  typedef typename T::template Child<0>::Type ChildType;
1297  typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits;
1298 
1299  typedef typename LBTraits::RangeFieldType RF;
1300  typedef typename LBTraits::JacobianType JT;
1301 
1302  VectorDiscreteGridFunctionCurl(const GFS& gfs, const X& x_)
1303  : pgfs(stackobject_to_shared_ptr(gfs))
1304  , lfs(gfs)
1305  , lfs_cache(lfs)
1306  , x_view(x_)
1307  , xl(gfs.maxLocalSize())
1308  , J(gfs.maxLocalSize())
1309  {
1310  static_assert(LBTraits::dimDomain == T::CHILDREN,
1311  "dimDomain and number of children has to be the same");
1312  }
1313 
1314  inline void evaluate(const typename Traits::ElementType& e,
1315  const typename Traits::DomainType& x,
1316  typename Traits::RangeType& y) const
1317  {
1318  // get and bind local functions space
1319  lfs.bind(e);
1320  lfs_cache.update();
1321  x_view.bind(lfs_cache);
1322  x_view.read(xl);
1323  x_view.unbind();
1324 
1325  // get Jacobian of geometry
1326  const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1327  JgeoIT = e.geometry().jacobianInverseTransposed(x);
1328 
1329  const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1330  Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1331 
1332  y = 0.0;
1333 
1334  // some handy variables for the curl computation in 2D
1335  RF sign = -1.0;
1336  int i2;
1337 
1338  // loop over childs of VectorGFS
1339  for(unsigned int k=0; k != T::CHILDREN; ++k) {
1340 
1341  // get local Jacobians/gradients of the shape functions
1342  std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1343  lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1344 
1345  RF d_k_phi;
1346  // pick out the right derivative
1347  i2 = 1-k;
1348  for(typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1349  // compute i2-th derivative of k-th child
1350  d_k_phi =
1352  typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1353  typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1354  N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1355  (JgeoIT,J[i][0],i2);
1356 
1357  y += sign * xl[lfs.child(k).localIndex(i)] * d_k_phi;
1358  }
1359  sign *= -1.0;
1360  }
1361  }
1362 
1364  inline const typename Traits::GridViewType& getGridView() const
1365  {
1366  return pgfs->gridView();
1367  }
1368 
1369  private :
1372  typedef typename X::template ConstLocalView<LFSCache> XView;
1373 
1374  shared_ptr<GFS const> pgfs;
1375  mutable LFS lfs;
1376  mutable LFSCache lfs_cache;
1377  mutable XView x_view;
1378  mutable std::vector<RF> xl;
1379  mutable std::vector<JT> J;
1380  shared_ptr<const X> px;
1381  }; // end class VectorDiscreteGridFunctionCurl (2D)
1382 
1384  } // namespace PDELab
1385 } // namespace Dune
1386 
1387 #endif
DiscreteGridFunction(std::shared_ptr< const GFS > gfs, std::shared_ptr< const X > x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:117
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:222
convert a single component function space with a grid function representing the gradient ...
Definition: gridfunctionspaceutilities.hh:461
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, const Remap &remap_)
construct
Definition: gridfunctionspaceutilities.hh:740
DiscreteGridFunctionGlobalCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGlobalCurl.
Definition: gridfunctionspaceutilities.hh:382
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1295
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:851
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:834
convert a grid function space and a coefficient vector into a grid function
Definition: gridfunctionspaceutilities.hh:53
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1364
static RF compute_derivative(const Dune::DiagonalMatrix< RF, size > &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:977
static const unsigned int value
Definition: gridfunctionspace/tags.hh:177
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:545
GridFunctionTraits< typename GFS::Traits::GridViewType, typename LBTraits::RangeFieldType, LBTraits::dimDomain, FieldVector< typename LBTraits::RangeFieldType, LBTraits::dimDomain > > Traits
Definition: gridfunctionspaceutilities.hh:487
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:835
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1168
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1296
DiscreteGridFunctionPiola(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionPiola.
Definition: gridfunctionspaceutilities.hh:597
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:702
DiscreteGridFunction for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:666
DiscreteGridFunction with Piola transformation.
Definition: gridfunctionspaceutilities.hh:567
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1302
extract the curl of a function from the jacobian of that function
Definition: jacobiantocurl.hh:27
convert a single component function space with experimental global finite elements into a grid functi...
Definition: gridfunctionspaceutilities.hh:350
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:836
GV GridViewType
The type of the grid view the function lives on.
Definition: function.hh:114
Equivalent of DiscreteGridFunctionGradient for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:800
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1173
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:509
Definition: adaptivity.hh:27
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType RT
Definition: gridfunctionspaceutilities.hh:706
DiscreteGridFunction(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:102
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:117
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1015
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:95
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1033
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:888
DiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:210
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:754
Compute curl of vector-valued functions.
Definition: gridfunctionspaceutilities.hh:1108
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:394
Dune::FieldVector< GV::Grid::ctype, GV::dimension > DomainType
domain type in dim-size coordinates
Definition: function.hh:48
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1019
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:634
Create a local function space from a global function space.
Definition: localfunctionspace.hh:670
Helper class to calculate the Traits of DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:278
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:129
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x)
Definition: gridfunctionspaceutilities.hh:1112
GridFunctionTraits< typename T::Traits::GridViewType, typename J2C::CurlField, J2C::dimCurl, typename J2C::Curl > Traits
Definition: gridfunctionspaceutilities.hh:203
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1300
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:774
DiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGradient.
Definition: gridfunctionspaceutilities.hh:500
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:838
Compute divergence of vector-valued functions.
Definition: gridfunctionspaceutilities.hh:994
static RF compute_derivative(const Dune::FieldMatrix< RF, size, size > &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:953
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:704
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1314
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1167
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1170
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1297
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1018
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1299
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, std::size_t start=0)
construct
Definition: gridfunctionspaceutilities.hh:714
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
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1185
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1016
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:150
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:431
VectorDiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:841
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1014
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:607
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:591
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1076
void update()
Definition: lfsindexcache.hh:300
convert a grid function space and a coefficient vector into a grid function of the curl ...
Definition: gridfunctionspaceutilities.hh:180
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:701
a GridFunction maps x in DomainType to y in RangeType
Definition: function.hh:186
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:244
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1166
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:839
Helper class to compute a single derivative of scalar basis functions.
Definition: gridfunctionspaceutilities.hh:921
VectorDiscreteGridFunctionDiv(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1021
traits class holding the function signature, same as in local function
Definition: function.hh:175
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1171
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1245
static RF compute_derivative(const Mat &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:927