3 #ifndef DUNE_PDELAB_GRIDFUNCTIONSPACEUTILITIES_HH
4 #define DUNE_PDELAB_GRIDFUNCTIONSPACEUTILITIES_HH
9 #include<dune/common/exceptions.hh>
10 #include <dune/common/fvector.hh>
11 #include <dune/common/static_assert.hh>
13 #include <dune/localfunctions/common/interfaceswitch.hh>
15 #include"../common/function.hh"
53 template<
typename T,
typename X>
55 :
public TypeTree::LeafNode
58 typename T::Traits::GridViewType,
59 typename BasisInterfaceSwitch<
60 typename FiniteElementInterfaceSwitch<
61 typename T::Traits::FiniteElementType
65 typename FiniteElementInterfaceSwitch<
66 typename T::Traits::FiniteElementType
69 typename BasisInterfaceSwitch<
70 typename FiniteElementInterfaceSwitch<
71 typename T::Traits::FiniteElementType
75 DiscreteGridFunction<T,X>
80 typedef typename Dune::BasisInterfaceSwitch<
81 typename FiniteElementInterfaceSwitch<
82 typename T::Traits::FiniteElementType
87 typename T::Traits::GridViewType,
88 typename BasisSwitch::RangeField,
89 BasisSwitch::dimRange,
90 typename BasisSwitch::Range
104 : pgfs(stackobject_to_shared_ptr(gfs))
108 , xl(gfs.maxLocalSize())
109 , yb(gfs.maxLocalSize())
123 , xl(gfs->maxLocalSize())
124 , yb(gfs->maxLocalSize())
130 inline void evaluate (
const typename Traits::ElementType&
e,
131 const typename Traits::DomainType& x,
132 typename Traits::RangeType& y)
const
134 typedef FiniteElementInterfaceSwitch<
139 x_view.bind(lfs_cache);
142 FESwitch::basis(lfs.finiteElement()).evaluateFunction(x,yb);
144 for (
unsigned int i=0; i<yb.size(); i++)
153 return pgfs->gridView();
160 typedef typename X::template ConstLocalView<LFSCache> XView;
162 shared_ptr<GFS const> pgfs;
164 mutable LFSCache lfs_cache;
165 mutable XView x_view;
166 mutable std::vector<typename Traits::RangeFieldType> xl;
167 mutable std::vector<typename Traits::RangeType> yb;
168 shared_ptr<const X> px;
180 template<
typename T,
typename X>
184 typename T::Traits::GridViewType,
185 typename JacobianToCurl<typename T::Traits::FiniteElementType::
186 Traits::LocalBasisType::Traits::JacobianType>::CurlField,
187 JacobianToCurl<typename T::Traits::FiniteElementType::Traits::LocalBasisType::
188 Traits::JacobianType>::dimCurl,
189 typename JacobianToCurl<typename T::Traits::FiniteElementType::
190 Traits::LocalBasisType::Traits::JacobianType>::Curl
192 DiscreteGridFunctionCurl<T,X>
196 typedef typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
197 JacobianType Jacobian;
202 typename T::Traits::GridViewType,
203 typename J2C::CurlField, J2C::dimCurl,
typename J2C::Curl
212 : pgfs(stackobject_to_shared_ptr(gfs))
216 , xl(gfs.maxLocalSize())
217 , jacobian(gfs.maxLocalSize())
218 , yb(gfs.maxLocalSize())
219 , px(stackobject_to_shared_ptr(x_))
227 static const J2C& j2C =
J2C();
231 x_view.bind(lfs_cache);
235 lfs.finiteElement().basis().evaluateJacobian(x,jacobian);
238 for (std::size_t i=0; i < lfs.
size(); i++) {
239 j2C(jacobian[i], yb);
246 {
return pgfs->gridView(); }
253 typedef typename X::template ConstLocalView<LFSCache> XView;
255 shared_ptr<GFS const> pgfs;
257 mutable LFSCache lfs_cache;
258 mutable XView x_view;
259 mutable std::vector<typename Traits::RangeFieldType> xl;
260 mutable std::vector<Jacobian> jacobian;
261 mutable std::vector<typename Traits::RangeType> yb;
262 shared_ptr<const X> px;
278 template<
typename GV,
typename RangeFieldType,
int dimRangeOfBasis>
281 "DiscreteGridFunctionCurl (and friends) work in 2D "
291 template<
typename GV,
typename RangeFieldType>
295 FieldVector<RangeFieldType, 2> >
298 "World dimension of grid must be 2 for the curl of a "
299 "scalar (1D) quantity");
307 template<
typename GV,
typename RangeFieldType>
311 FieldVector<RangeFieldType, 1> >
314 "World dimension of grid must be 2 for the curl of a"
323 template<
typename GV,
typename RangeFieldType>
327 FieldVector<RangeFieldType, 3> >
330 "World dimension of grid must be 3 for the curl of a"
350 template<
typename T,
typename X>
353 DiscreteGridFunctionCurlTraits<
354 typename T::Traits::GridViewType,
355 typename T::Traits::FiniteElementType::Traits::
356 LocalBasisType::Traits::RangeFieldType,
357 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
359 DiscreteGridFunctionGlobalCurl<T,X> >
363 typename T::Traits::GridViewType,
364 typename T::Traits::FiniteElementType::Traits::
365 LocalBasisType::Traits::RangeFieldType,
366 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
374 typedef typename T::Traits::FiniteElementType::Traits::
375 LocalBasisType::Traits LBTraits;
384 : pgfs(stackobject_to_shared_ptr(gfs))
388 , xl(gfs.maxLocalSize())
389 , J(gfs.maxLocalSize())
390 , px(stackobject_to_shared_ptr(x_))
395 inline void evaluate (
const typename Traits::ElementType&
e,
396 const typename Traits::DomainType& x,
397 typename Traits::RangeType& y)
const
401 x_view.bind(lfs_cache);
405 lfs.finiteElement().localBasis().
406 evaluateJacobianGlobal(x,J,e.geometry());
408 for (
unsigned int i=0; i<J.size(); i++)
412 switch(
unsigned(Traits::dimRange)) {
414 y[0] += xl[i] * J[i][0][1];
415 y[1] += xl[i] * -J[i][0][0];
418 y[0] += xl[i]*(J[i][1][0] - J[i][0][1]);
421 y[0] += xl[i]*(J[i][2][1] - J[i][1][2]);
422 y[1] += xl[i]*(J[i][0][2] - J[i][2][0]);
423 y[2] += xl[i]*(J[i][1][0] - J[i][0][1]);
434 return pgfs->gridView();
440 typedef typename X::template ConstLocalView<LFSCache> XView;
442 shared_ptr<GFS const> pgfs;
444 mutable LFSCache lfs_cache;
445 mutable XView x_view;
446 mutable std::vector<typename Traits::RangeFieldType> xl;
447 mutable std::vector<typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::JacobianType> J;
448 shared_ptr<const X> px;
461 template<
typename T,
typename X>
465 typename T::Traits::GridViewType,
466 typename T::Traits::FiniteElementType::Traits::LocalBasisType
467 ::Traits::RangeFieldType,
468 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
471 typename T::Traits::FiniteElementType::Traits
472 ::LocalBasisType::Traits::RangeFieldType,
473 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
475 DiscreteGridFunctionGradient<T,X> >
478 typedef typename GFS::Traits::FiniteElementType::Traits::
479 LocalBasisType::Traits LBTraits;
483 typename GFS::Traits::GridViewType,
484 typename LBTraits::RangeFieldType,
487 typename LBTraits::RangeFieldType,
502 : pgfs(stackobject_to_shared_ptr(gfs))
517 x_view.bind(lfs_cache);
520 xl.resize(lfs.
size());
525 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
526 JgeoIT = e.geometry().jacobianInverseTransposed(x);
529 std::vector<typename LBTraits::JacobianType> J(lfs.
size());
530 lfs.finiteElement().localBasis().evaluateJacobian(x,J);
534 for(
unsigned int i = 0; i < lfs.
size(); ++i) {
537 JgeoIT.umv(J[i][0], gradphi);
540 y.axpy(xl[i], gradphi);
548 return pgfs->gridView();
554 typedef typename X::template ConstLocalView<LFSCache> XView;
556 shared_ptr<GFS const> pgfs;
558 mutable LFSCache lfs_cache;
559 mutable XView x_view;
560 mutable std::vector<typename Traits::RangeFieldType> xl;
567 template<
typename T,
typename X>
571 typename T::Traits::GridViewType,
572 typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
573 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
574 typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType
576 DiscreteGridFunctionPiola<T,X>
583 typename T::Traits::GridViewType,
584 typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
585 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
586 typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType
599 : pgfs(stackobject_to_shared_ptr(gfs))
603 , xl(pgfs->maxLocalSize())
604 , yb(pgfs->maxLocalSize())
608 inline void evaluate (
const typename Traits::ElementType&
e,
609 const typename Traits::DomainType& x,
610 typename Traits::RangeType& y)
const
615 x_view.bind(lfs_cache);
619 lfs.finiteElement().localBasis().evaluateFunction(x,yb);
620 typename Traits::RangeType yhat;
622 for (
unsigned int i=0; i<yb.size(); i++)
623 yhat.axpy(xl[i],yb[i]);
626 typename Traits::ElementType::Geometry::JacobianInverseTransposed
627 J = e.geometry().jacobianInverseTransposed(x);
631 y /= J.determinant();
637 return pgfs->gridView();
644 typedef typename X::template ConstLocalView<LFSCache> XView;
646 shared_ptr<GFS const> pgfs;
648 mutable LFSCache lfs_cache;
649 mutable XView x_view;
650 mutable std::vector<typename Traits::RangeFieldType> xl;
651 mutable std::vector<typename Traits::RangeType> yb;
666 template<
typename T,
typename X, std::
size_t dimR = T::CHILDREN>
670 typename T::Traits::GridViewType,
671 typename T::template Child<0>::Type::Traits::FiniteElementType
672 ::Traits::LocalBasisType::Traits::RangeFieldType,
675 typename T::template Child<0>::Type::Traits::FiniteElementType
676 ::Traits::LocalBasisType::Traits::RangeFieldType,
680 VectorDiscreteGridFunction<T,X>
682 public TypeTree::LeafNode
688 typename T::Traits::GridViewType,
689 typename T::template Child<0>::Type::Traits::FiniteElementType
690 ::Traits::LocalBasisType::Traits::RangeFieldType,
693 typename T::template Child<0>::Type::Traits::FiniteElementType
694 ::Traits::LocalBasisType::Traits::RangeFieldType,
704 typedef typename ChildType::Traits::FiniteElementType
705 ::Traits::LocalBasisType::Traits::RangeFieldType
RF;
706 typedef typename ChildType::Traits::FiniteElementType
707 ::Traits::LocalBasisType::Traits::RangeType
RT;
716 std::size_t start = 0)
717 : pgfs(stackobject_to_shared_ptr(gfs))
721 , xl(gfs.maxLocalSize())
722 , yb(gfs.maxLocalSize())
724 for(std::size_t i = 0; i < dimR; ++i)
725 remap[i] = i + start;
740 template<
class Remap>
743 : pgfs(stackobject_to_shared_ptr(gfs))
747 , xl(gfs.maxLocalSize())
748 , yb(gfs.maxLocalSize())
749 , px(stackobject_to_shared_ptr(x_))
751 for(std::size_t i = 0; i < dimR; ++i)
752 remap[i] = remap_[i];
755 inline void evaluate (
const typename Traits::ElementType&
e,
756 const typename Traits::DomainType& x,
757 typename Traits::RangeType& y)
const
761 x_view.bind(lfs_cache);
764 for (
unsigned int k=0; k < dimR; k++)
766 lfs.child(remap[k]).finiteElement().localBasis().
767 evaluateFunction(x,yb);
769 for (
unsigned int i=0; i<yb.size(); i++)
770 y[k] += xl[lfs.child(remap[k]).localIndex(i)]*yb[i];
777 return pgfs->gridView();
783 typedef typename X::template ConstLocalView<LFSCache> XView;
785 shared_ptr<GFS const> pgfs;
786 std::size_t remap[dimR];
788 mutable LFSCache lfs_cache;
789 mutable XView x_view;
790 mutable std::vector<RF> xl;
791 mutable std::vector<RT> yb;
792 shared_ptr<const X> px;
800 template<
typename T,
typename X>
804 typename T::Traits::GridViewType,
805 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
809 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
811 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain
814 VectorDiscreteGridFunctionGradient<T,X>
816 public TypeTree::LeafNode
822 typename T::Traits::GridViewType,
823 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
827 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
829 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain>
837 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
839 typedef typename LBTraits::RangeFieldType
RF;
840 typedef typename LBTraits::JacobianType
JT;
843 : pgfs(stackobject_to_shared_ptr(gfs))
847 , xl(gfs.maxLocalSize())
848 , J(gfs.maxLocalSize())
852 inline void evaluate(
const typename Traits::ElementType&
e,
853 const typename Traits::DomainType& x,
854 typename Traits::RangeType& y)
const
859 x_view.bind(lfs_cache);
864 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
865 JgeoIT = e.geometry().jacobianInverseTransposed(x);
870 for(
unsigned int k = 0; k != T::CHILDREN; ++k)
873 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
874 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
876 Dune::FieldVector<RF,LBTraits::dimDomain> gradphi;
877 for (
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++)
880 JgeoIT.umv(J[i][0], gradphi);
882 y[k].axpy(xl[lfs.child(k).localIndex(i)], gradphi);
891 return pgfs->gridView();
897 typedef typename X::template ConstLocalView<LFSCache> XView;
899 shared_ptr<GFS const> pgfs;
901 mutable LFSCache lfs_cache;
902 mutable XView x_view;
903 mutable std::vector<RF> xl;
904 mutable std::vector<JT> J;
905 shared_ptr<const X> px;
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:432
DiscreteGridFunction(shared_ptr< const GFS > gfs, shared_ptr< const X > x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:118
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:835
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:775
DiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGradient.
Definition: gridfunctionspaceutilities.hh:501
VectorDiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:842
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:755
convert a single component function space with experimental global finite elements into a grid functi...
Definition: gridfunctionspaceutilities.hh:351
dune_static_assert(AlwaysFalse< GV >::value,"DiscreteGridFunctionCurl (and friends) work in 2D ""and 3D only")
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:223
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:836
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:118
void update()
Definition: lfsindexcache.hh:300
Equivalent of DiscreteGridFunctionGradient for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:801
DiscreteGridFunctionGlobalCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGlobalCurl.
Definition: gridfunctionspaceutilities.hh:383
Helper class to calculate the Traits of DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:279
RangeFieldType RangeFieldType
Export type for range field.
Definition: function.hh:52
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:705
Dune::FieldVector< GV::Grid::ctype, GV::dimension > DomainType
domain type in dim-size coordinates
Definition: function.hh:49
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:635
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, std::size_t start=0)
construct
Definition: gridfunctionspaceutilities.hh:715
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:130
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:510
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:889
extract the curl of a function from the jacobian of that function
Definition: jacobiantocurl.hh:28
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:546
a GridFunction maps x in DomainType to y in RangeType
Definition: function.hh:187
R RangeType
range type
Definition: function.hh:61
DiscreteGridFunction for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:667
DiscreteGridFunction with Piola transformation.
Definition: gridfunctionspaceutilities.hh:568
DiscreteGridFunctionCurlTraits< typename T::Traits::GridViewType, typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType, T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange > Traits
Definition: gridfunctionspaceutilities.hh:367
Traits::IndexContainer::size_type size() const
get current size
Definition: localfunctionspace.hh:206
GridFunctionTraits< typename T::Traits::GridViewType, typename J2C::CurlField, J2C::dimCurl, typename J2C::Curl > Traits
Definition: gridfunctionspaceutilities.hh:204
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:96
static const unsigned int value
Definition: gridfunctionspace/tags.hh:175
Create a local function space from a global function space.
Definition: localfunctionspace.hh:702
DiscreteGridFunctionPiola(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionPiola.
Definition: gridfunctionspaceutilities.hh:598
traits class holding the function signature, same as in local function
Definition: function.hh:176
convert a grid function space and a coefficient vector into a grid function of the curl ...
Definition: gridfunctionspaceutilities.hh:181
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:151
T Traits
Export type traits.
Definition: function.hh:192
DiscreteGridFunction(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:103
DiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:211
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:592
GV GridViewType
The type of the grid view the function lives on.
Definition: function.hh:115
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:702
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType RT
Definition: gridfunctionspaceutilities.hh:707
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:703
const E & e
Definition: interpolate.hh:172
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, const Remap &remap_)
construct
Definition: gridfunctionspaceutilities.hh:741
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:839
convert a single component function space with a grid function representing the gradient ...
Definition: gridfunctionspaceutilities.hh:462
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:395
convert a grid function space and a coefficient vector into a grid function
Definition: gridfunctionspaceutilities.hh:54
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:840
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:608
GridFunctionTraits< typename GFS::Traits::GridViewType, typename LBTraits::RangeFieldType, LBTraits::dimDomain, FieldVector< typename LBTraits::RangeFieldType, LBTraits::dimDomain > > Traits
Definition: gridfunctionspaceutilities.hh:488
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:245
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:852
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:837