3 #ifndef DUNE_PDELAB_GRIDFUNCTIONSPACE_GRIDFUNCTIONSPACEUTILITIES_HH 4 #define DUNE_PDELAB_GRIDFUNCTIONSPACE_GRIDFUNCTIONSPACEUTILITIES_HH 9 #include<dune/common/exceptions.hh> 10 #include <dune/common/fvector.hh> 12 #include <dune/localfunctions/common/interfaceswitch.hh> 14 #include"../common/function.hh" 22 template <
typename ES>
24 :
public TypeTree::TreeVisitor
25 ,
public TypeTree::DynamicTraversal
28 template<
typename LeafGFS,
typename TreePath>
29 void leaf(
const LeafGFS& leaf_gfs, TreePath tp)
71 template<
typename T,
typename X>
73 :
public TypeTree::LeafNode
76 typename T::Traits::GridViewType,
77 typename BasisInterfaceSwitch<
78 typename FiniteElementInterfaceSwitch<
79 typename T::Traits::FiniteElementType
83 typename FiniteElementInterfaceSwitch<
84 typename T::Traits::FiniteElementType
87 typename BasisInterfaceSwitch<
88 typename FiniteElementInterfaceSwitch<
89 typename T::Traits::FiniteElementType
93 DiscreteGridFunction<T,X>
98 typedef typename Dune::BasisInterfaceSwitch<
99 typename FiniteElementInterfaceSwitch<
100 typename T::Traits::FiniteElementType
105 typename T::Traits::GridViewType,
106 typename BasisSwitch::RangeField,
107 BasisSwitch::dimRange,
108 typename BasisSwitch::Range
122 : pgfs(stackobject_to_shared_ptr(gfs))
126 , xl(gfs.maxLocalSize())
127 , yb(gfs.maxLocalSize())
141 , xl(gfs->maxLocalSize())
142 , yb(gfs->maxLocalSize())
148 inline void evaluate (
const typename Traits::ElementType&
e,
149 const typename Traits::DomainType& x,
150 typename Traits::RangeType& y)
const 152 typedef FiniteElementInterfaceSwitch<
157 x_view.bind(lfs_cache);
160 FESwitch::basis(lfs.finiteElement()).evaluateFunction(x,yb);
162 for (
unsigned int i=0; i<yb.size(); i++)
171 return pgfs->gridView();
178 typedef typename X::template ConstLocalView<LFSCache> XView;
180 std::shared_ptr<GFS const> pgfs;
182 mutable LFSCache lfs_cache;
183 mutable XView x_view;
184 mutable std::vector<typename Traits::RangeFieldType> xl;
185 mutable std::vector<typename Traits::RangeType> yb;
186 std::shared_ptr<const X> px;
198 template<
typename T,
typename X>
202 typename T::Traits::GridViewType,
203 typename JacobianToCurl<typename T::Traits::FiniteElementType::
204 Traits::LocalBasisType::Traits::JacobianType>::CurlField,
205 JacobianToCurl<typename T::Traits::FiniteElementType::Traits::LocalBasisType::
206 Traits::JacobianType>::dimCurl,
207 typename JacobianToCurl<typename T::Traits::FiniteElementType::
208 Traits::LocalBasisType::Traits::JacobianType>::Curl
210 DiscreteGridFunctionCurl<T,X>
214 typedef typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
215 JacobianType Jacobian;
220 typename T::Traits::GridViewType,
221 typename J2C::CurlField, J2C::dimCurl,
typename J2C::Curl
230 : pgfs(stackobject_to_shared_ptr(gfs))
234 , xl(gfs.maxLocalSize())
235 , jacobian(gfs.maxLocalSize())
236 , yb(gfs.maxLocalSize())
237 , px(stackobject_to_shared_ptr(x_))
245 static const J2C& j2C = J2C();
249 x_view.bind(lfs_cache);
253 lfs.finiteElement().basis().evaluateJacobian(x,jacobian);
256 for (std::size_t i=0; i < lfs.size(); i++) {
257 j2C(jacobian[i], yb);
264 {
return pgfs->gridView(); }
271 typedef typename X::template ConstLocalView<LFSCache> XView;
273 std::shared_ptr<GFS const> pgfs;
275 mutable LFSCache lfs_cache;
276 mutable XView x_view;
277 mutable std::vector<typename Traits::RangeFieldType> xl;
278 mutable std::vector<Jacobian> jacobian;
279 mutable std::vector<typename Traits::RangeType> yb;
280 std::shared_ptr<const X> px;
296 template<
typename GV,
typename RangeFieldType,
int dimRangeOfBasis>
299 "DiscreteGridFunctionCurl (and friends) work in 2D " 309 template<
typename GV,
typename RangeFieldType>
313 FieldVector<RangeFieldType, 2> >
315 static_assert(GV::dimensionworld == 2,
316 "World dimension of grid must be 2 for the curl of a " 317 "scalar (1D) quantity");
325 template<
typename GV,
typename RangeFieldType>
329 FieldVector<RangeFieldType, 1> >
331 static_assert(GV::dimensionworld == 2,
332 "World dimension of grid must be 2 for the curl of a" 341 template<
typename GV,
typename RangeFieldType>
345 FieldVector<RangeFieldType, 3> >
347 static_assert(GV::dimensionworld == 3,
348 "World dimension of grid must be 3 for the curl of a" 368 template<
typename T,
typename X>
371 DiscreteGridFunctionCurlTraits<
372 typename T::Traits::GridViewType,
373 typename T::Traits::FiniteElementType::Traits::
374 LocalBasisType::Traits::RangeFieldType,
375 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
377 DiscreteGridFunctionGlobalCurl<T,X> >
381 typename T::Traits::GridViewType,
382 typename T::Traits::FiniteElementType::Traits::
383 LocalBasisType::Traits::RangeFieldType,
384 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
392 typedef typename T::Traits::FiniteElementType::Traits::
393 LocalBasisType::Traits LBTraits;
402 : pgfs(stackobject_to_shared_ptr(gfs))
406 , xl(gfs.maxLocalSize())
407 , J(gfs.maxLocalSize())
408 , px(stackobject_to_shared_ptr(x_))
413 inline void evaluate (
const typename Traits::ElementType&
e,
414 const typename Traits::DomainType& x,
415 typename Traits::RangeType& y)
const 419 x_view.bind(lfs_cache);
423 lfs.finiteElement().localBasis().
424 evaluateJacobianGlobal(x,J,e.geometry());
426 for (
unsigned int i=0; i<J.size(); i++)
430 switch(
unsigned(Traits::dimRange)) {
432 y[0] += xl[i] * J[i][0][1];
433 y[1] += xl[i] * -J[i][0][0];
436 y[0] += xl[i]*(J[i][1][0] - J[i][0][1]);
439 y[0] += xl[i]*(J[i][2][1] - J[i][1][2]);
440 y[1] += xl[i]*(J[i][0][2] - J[i][2][0]);
441 y[2] += xl[i]*(J[i][1][0] - J[i][0][1]);
452 return pgfs->gridView();
458 typedef typename X::template ConstLocalView<LFSCache> XView;
460 std::shared_ptr<GFS const> pgfs;
462 mutable LFSCache lfs_cache;
463 mutable XView x_view;
464 mutable std::vector<typename Traits::RangeFieldType> xl;
465 mutable std::vector<typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::JacobianType> J;
466 std::shared_ptr<const X> px;
479 template<
typename T,
typename X>
483 typename T::Traits::GridViewType,
484 typename T::Traits::FiniteElementType::Traits::LocalBasisType
485 ::Traits::RangeFieldType,
486 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
489 typename T::Traits::FiniteElementType::Traits
490 ::LocalBasisType::Traits::RangeFieldType,
491 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
493 DiscreteGridFunctionGradient<T,X> >
496 typedef typename GFS::Traits::FiniteElementType::Traits::
497 LocalBasisType::Traits LBTraits;
501 typename GFS::Traits::GridViewType,
502 typename LBTraits::RangeFieldType,
505 typename LBTraits::RangeFieldType,
520 : pgfs(stackobject_to_shared_ptr(gfs))
535 x_view.bind(lfs_cache);
538 xl.resize(lfs.size());
543 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
544 JgeoIT = e.geometry().jacobianInverseTransposed(x);
547 std::vector<typename LBTraits::JacobianType> J(lfs.size());
548 lfs.finiteElement().localBasis().evaluateJacobian(x,J);
552 for(
unsigned int i = 0; i < lfs.size(); ++i) {
555 JgeoIT.umv(J[i][0], gradphi);
558 y.axpy(xl[i], gradphi);
566 return pgfs->gridView();
572 typedef typename X::template ConstLocalView<LFSCache> XView;
574 std::shared_ptr<GFS const> pgfs;
576 mutable LFSCache lfs_cache;
577 mutable XView x_view;
578 mutable std::vector<typename Traits::RangeFieldType> xl;
585 template<
typename T,
typename X>
589 typename T::Traits::GridViewType,
590 typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
591 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
592 typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType
594 DiscreteGridFunctionPiola<T,X>
601 typename T::Traits::GridViewType,
602 typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
603 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
604 typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType
617 : pgfs(stackobject_to_shared_ptr(gfs))
621 , xl(pgfs->maxLocalSize())
622 , yb(pgfs->maxLocalSize())
626 inline void evaluate (
const typename Traits::ElementType&
e,
627 const typename Traits::DomainType& x,
628 typename Traits::RangeType& y)
const 633 x_view.bind(lfs_cache);
637 lfs.finiteElement().localBasis().evaluateFunction(x,yb);
638 typename Traits::RangeType yhat;
640 for (
unsigned int i=0; i<yb.size(); i++)
641 yhat.axpy(xl[i],yb[i]);
644 typename Traits::ElementType::Geometry::JacobianInverseTransposed
645 J = e.geometry().jacobianInverseTransposed(x);
649 y /= J.determinant();
655 return pgfs->gridView();
662 typedef typename X::template ConstLocalView<LFSCache> XView;
664 std::shared_ptr<GFS const> pgfs;
666 mutable LFSCache lfs_cache;
667 mutable XView x_view;
668 mutable std::vector<typename Traits::RangeFieldType> xl;
669 mutable std::vector<typename Traits::RangeType> yb;
693 typename T::Traits::GridViewType,
694 typename T::template Child<0>::Type::Traits::FiniteElementType
695 ::Traits::LocalBasisType::Traits::RangeFieldType,
698 typename T::template Child<0>::Type::Traits::FiniteElementType
699 ::Traits::LocalBasisType::Traits::RangeFieldType,
703 VectorDiscreteGridFunction<T,X>
705 public TypeTree::LeafNode
711 typename T::Traits::GridViewType,
712 typename T::template Child<0>::Type::Traits::FiniteElementType
713 ::Traits::LocalBasisType::Traits::RangeFieldType,
716 typename T::template Child<0>::Type::Traits::FiniteElementType
717 ::Traits::LocalBasisType::Traits::RangeFieldType,
727 typedef typename ChildType::Traits::FiniteElementType
728 ::Traits::LocalBasisType::Traits::RangeFieldType
RF;
729 typedef typename ChildType::Traits::FiniteElementType
730 ::Traits::LocalBasisType::Traits::RangeType
RT;
739 std::size_t start = 0)
740 : pgfs(stackobject_to_shared_ptr(gfs))
744 , xl(gfs.maxLocalSize())
745 , yb(gfs.maxLocalSize())
747 for(std::size_t i = 0; i < dimR; ++i)
748 remap[i] = i + start;
763 template<
class Remap>
766 : pgfs(stackobject_to_shared_ptr(gfs))
770 , xl(gfs.maxLocalSize())
771 , yb(gfs.maxLocalSize())
772 , px(stackobject_to_shared_ptr(x_))
774 for(std::size_t i = 0; i < dimR; ++i)
775 remap[i] = remap_[i];
778 inline void evaluate (
const typename Traits::ElementType&
e,
779 const typename Traits::DomainType& x,
780 typename Traits::RangeType& y)
const 784 x_view.bind(lfs_cache);
787 for (
unsigned int k=0; k < dimR; k++)
789 lfs.child(remap[k]).finiteElement().localBasis().
790 evaluateFunction(x,yb);
792 for (
unsigned int i=0; i<yb.size(); i++)
793 y[k] += xl[lfs.child(remap[k]).localIndex(i)]*yb[i];
800 return pgfs->gridView();
806 typedef typename X::template ConstLocalView<LFSCache> XView;
808 std::shared_ptr<GFS const> pgfs;
809 std::size_t remap[dimR];
811 mutable LFSCache lfs_cache;
812 mutable XView x_view;
813 mutable std::vector<RF> xl;
814 mutable std::vector<RT> yb;
815 std::shared_ptr<const X> px;
823 template<
typename T,
typename X>
827 typename T::Traits::GridViewType,
828 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
830 TypeTree::StaticDegree<T>::value,
832 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
833 TypeTree::StaticDegree<T>::value,
834 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain
837 VectorDiscreteGridFunctionGradient<T,X>
839 public TypeTree::LeafNode
845 typename T::Traits::GridViewType,
846 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
850 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
851 TypeTree::StaticDegree<T>::value,
852 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain>
860 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
862 typedef typename LBTraits::RangeFieldType
RF;
863 typedef typename LBTraits::JacobianType
JT;
866 : pgfs(stackobject_to_shared_ptr(gfs))
870 , xl(gfs.maxLocalSize())
871 , J(gfs.maxLocalSize())
875 inline void evaluate(
const typename Traits::ElementType&
e,
876 const typename Traits::DomainType& x,
877 typename Traits::RangeType& y)
const 882 x_view.bind(lfs_cache);
887 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
888 JgeoIT = e.geometry().jacobianInverseTransposed(x);
893 for(
unsigned int k = 0; k != TypeTree::degree(lfs); ++k)
896 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
897 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
899 Dune::FieldVector<RF,LBTraits::dimDomain> gradphi;
900 for (
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++)
903 JgeoIT.umv(J[i][0], gradphi);
905 y[k].axpy(xl[lfs.child(k).localIndex(i)], gradphi);
914 return pgfs->gridView();
920 typedef typename X::template ConstLocalView<LFSCache> XView;
922 std::shared_ptr<GFS const> pgfs;
924 mutable LFSCache lfs_cache;
925 mutable XView x_view;
926 mutable std::vector<RF> xl;
927 mutable std::vector<JT> J;
928 std::shared_ptr<const X> px;
944 template<
typename Mat,
typename RF, std::
size_t size>
954 Dune::FieldVector<RF,size> grad_phi(0.0);
969 template<
typename RF, std::
size_t size>
977 static inline RF
compute_derivative(
const Dune::FieldMatrix<RF,size,size>& mat,
const T& t,
const unsigned int k)
993 template<
typename RF, std::
size_t size>
1000 template<
typename T>
1001 static inline RF
compute_derivative(
const Dune::DiagonalMatrix<RF,size>& mat,
const T& t,
const unsigned int k)
1003 return mat[k][k]*t[k];
1017 template<
typename T,
typename X>
1020 Dune::PDELab::GridFunctionTraits<
1021 typename T::Traits::GridViewType,
1022 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1023 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1024 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1025 VectorDiscreteGridFunctionDiv<T,X> >
1026 ,
public TypeTree::LeafNode
1032 typename T::Traits::GridViewType,
1033 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1034 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1035 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1040 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
1042 typedef typename LBTraits::RangeFieldType
RF;
1043 typedef typename LBTraits::JacobianType
JT;
1046 : pgfs(stackobject_to_shared_ptr(gfs))
1050 , xl(gfs.maxLocalSize())
1051 , J(gfs.maxLocalSize())
1054 "dimDomain and number of children has to be the same");
1058 const typename Traits::DomainType& x,
1059 typename Traits::RangeType& y)
const 1064 x_view.bind(lfs_cache);
1069 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1070 JgeoIT = e.geometry().jacobianInverseTransposed(x);
1072 const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1073 Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1078 for(
unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1081 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1082 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1085 for(
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1089 typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1090 typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1091 N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1094 y += xl[lfs.child(k).localIndex(i)] * d_k_phi;
1102 return pgfs->gridView();
1108 typedef typename X::template ConstLocalView<LFSCache> XView;
1110 shared_ptr<GFS const> pgfs;
1112 mutable LFSCache lfs_cache;
1113 mutable XView x_view;
1114 mutable std::vector<RF> xl;
1115 mutable std::vector<JT> J;
1116 shared_ptr<const X> px;
1139 "Curl computation can only be done in two or three dimensions");
1153 template<
typename T,
typename X>
1156 Dune::PDELab::GridFunctionTraits<
1157 typename T::Traits::GridViewType,
1158 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1160 TypeTree::StaticDegree<T>::value,
1162 typename T::template Child<0>::Type::Traits::FiniteElementType
1163 ::Traits::LocalBasisType::Traits::RangeFieldType,
1165 TypeTree::StaticDegree<T>::value
1168 VectorDiscreteGridFunctionCurl<T,X>
1170 ,
public TypeTree::LeafNode
1176 typename T::Traits::GridViewType,
1177 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1181 typename T::template Child<0>::Type::Traits::FiniteElementType
1182 ::Traits::LocalBasisType::Traits::RangeFieldType,
1184 TypeTree::StaticDegree<T>::value
1192 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
1194 typedef typename LBTraits::RangeFieldType
RF;
1195 typedef typename LBTraits::JacobianType
JT;
1198 : pgfs(stackobject_to_shared_ptr(gfs))
1202 , xl(gfs.maxLocalSize())
1203 , J(gfs.maxLocalSize())
1205 static_assert(LBTraits::dimDomain == TypeTree::StaticDegree<T>::value,
1206 "dimDomain and number of children has to be the same");
1210 const typename Traits::DomainType& x,
1211 typename Traits::RangeType& y)
const 1216 x_view.bind(lfs_cache);
1221 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1222 JgeoIT = e.geometry().jacobianInverseTransposed(x);
1224 const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1225 Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1233 for(
unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1236 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1237 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1244 for(
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1248 typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1249 typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1250 N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1251 (JgeoIT,J[i][0],i2);
1253 y[i1] += xl[lfs.child(k).localIndex(i)] * d_k_phi;
1258 typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1259 typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1260 N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1261 (JgeoIT,J[i][0],i1);
1263 y[i2] -= xl[lfs.child(k).localIndex(i)] * d_k_phi;
1271 return pgfs->gridView();
1277 typedef typename X::template ConstLocalView<LFSCache> XView;
1279 shared_ptr<GFS const> pgfs;
1281 mutable LFSCache lfs_cache;
1282 mutable XView x_view;
1283 mutable std::vector<RF> xl;
1284 mutable std::vector<JT> J;
1285 shared_ptr<const X> px;
1298 template<
typename T,
typename X>
1301 Dune::PDELab::GridFunctionTraits<
1302 typename T::Traits::GridViewType,
1303 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1304 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1305 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1306 VectorDiscreteGridFunctionDiv<T,X> >
1307 ,
public TypeTree::LeafNode
1313 typename T::Traits::GridViewType,
1314 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1315 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1316 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1321 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
1323 typedef typename LBTraits::RangeFieldType
RF;
1324 typedef typename LBTraits::JacobianType
JT;
1327 : pgfs(stackobject_to_shared_ptr(gfs))
1331 , xl(gfs.maxLocalSize())
1332 , J(gfs.maxLocalSize())
1335 "dimDomain and number of children has to be the same");
1339 const typename Traits::DomainType& x,
1340 typename Traits::RangeType& y)
const 1345 x_view.bind(lfs_cache);
1350 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1351 JgeoIT = e.geometry().jacobianInverseTransposed(x);
1353 const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1354 Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1363 for(
unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1366 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1367 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1372 for(
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1376 typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1377 typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1378 N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1379 (JgeoIT,J[i][0],i2);
1381 y += sign * xl[lfs.child(k).localIndex(i)] * d_k_phi;
1390 return pgfs->gridView();
1396 typedef typename X::template ConstLocalView<LFSCache> XView;
1398 shared_ptr<GFS const> pgfs;
1400 mutable LFSCache lfs_cache;
1401 mutable XView x_view;
1402 mutable std::vector<RF> xl;
1403 mutable std::vector<JT> J;
1404 shared_ptr<const X> px;
1411 #endif // DUNE_PDELAB_GRIDFUNCTIONSPACE_GRIDFUNCTIONSPACEUTILITIES_HH Helper class to compute a single derivative of scalar basis functions.
Definition: gridfunctionspaceutilities.hh:945
GridFunctionTraits< typename T::Traits::GridViewType, typename J2C::CurlField, J2C::dimCurl, typename J2C::Curl > Traits
Definition: gridfunctionspaceutilities.hh:222
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:912
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:858
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1197
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1323
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1040
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:875
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1039
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1195
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:114
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:148
T Traits
Export type traits.
Definition: function.hh:191
DiscreteGridFunctionPiola(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionPiola.
Definition: gridfunctionspaceutilities.hh:616
Helper class to calculate the Traits of DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:297
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:564
Definition: gridfunctionspaceutilities.hh:23
ChildType::Traits::FiniteElementType ::Traits::LocalBasisType::Traits::RangeType RT
Definition: gridfunctionspaceutilities.hh:730
DiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGradient.
Definition: gridfunctionspaceutilities.hh:519
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:117
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, const Remap &remap_)
construct
Definition: gridfunctionspaceutilities.hh:764
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:726
static RF compute_derivative(const Dune::FieldMatrix< RF, size, size > &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:977
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:626
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:528
VectorDiscreteGridFunctionDiv(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1045
convert a single component function space with a grid function representing the gradient ...
Definition: gridfunctionspaceutilities.hh:480
ChildType::Traits::FiniteElementType ::Traits::LocalBasisType::Traits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:728
DiscreteGridFunction(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:121
convert a grid function space and a coefficient vector into a grid function of the curl ...
Definition: gridfunctionspaceutilities.hh:199
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, std::size_t start=0)
construct
Definition: gridfunctionspaceutilities.hh:738
DiscreteGridFunctionGlobalCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGlobalCurl.
Definition: gridfunctionspaceutilities.hh:401
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:862
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1100
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:798
DiscreteGridFunction for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:690
a GridFunction maps x in DomainType to y in RangeType
Definition: function.hh:186
convert a single component function space with experimental global finite elements into a grid functi...
Definition: gridfunctionspaceutilities.hh:369
GV GridViewType
The type of the grid view the function lives on.
Definition: function.hh:114
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1269
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:169
set_entity_set_visitor(const ES &es_)
Definition: gridfunctionspaceutilities.hh:34
const ES & es
Definition: gridfunctionspaceutilities.hh:38
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1194
RangeFieldType RangeFieldType
Export type for range field.
Definition: function.hh:51
extract the curl of a function from the jacobian of that function
Definition: jacobiantocurl.hh:27
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1320
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1190
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:610
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x)
Definition: gridfunctionspaceutilities.hh:1136
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1057
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1038
Equivalent of DiscreteGridFunctionGradient for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:824
void leaf(const LeafGFS &leaf_gfs, TreePath tp)
Definition: gridfunctionspaceutilities.hh:29
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1321
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1324
DiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:229
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:863
Dune::FieldVector< GV::Grid::ctype, GV::dimension > DomainType
domain type in dim-size coordinates
Definition: function.hh:48
static RF compute_derivative(const Mat &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:951
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:450
static RF compute_derivative(const Dune::DiagonalMatrix< RF, size > &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:1001
Compute curl of vector-valued functions.
Definition: gridfunctionspaceutilities.hh:1132
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1388
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:241
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1209
const Entity & e
Definition: localfunctionspace.hh:111
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:413
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:653
VectorDiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:865
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1191
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:263
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:778
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1043
GridFunctionTraits< typename GFS::Traits::GridViewType, typename LBTraits::RangeFieldType, LBTraits::dimDomain, FieldVector< typename LBTraits::RangeFieldType, LBTraits::dimDomain > > Traits
Definition: gridfunctionspaceutilities.hh:506
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1319
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:860
DiscreteGridFunction(std::shared_ptr< const GFS > gfs, std::shared_ptr< const X > x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:136
Create a local function space from a global function space.
Definition: localfunctionspace.hh:670
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:859
convert a grid function space and a coefficient vector into a grid function
Definition: gridfunctionspaceutilities.hh:72
R RangeType
range type
Definition: function.hh:60
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1042
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1192
DiscreteGridFunctionCurlTraits< typename T::Traits::GridViewType, typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType, T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange > Traits
Definition: gridfunctionspaceutilities.hh:385
DiscreteGridFunction with Piola transformation.
Definition: gridfunctionspaceutilities.hh:586
traits class holding the function signature, same as in local function
Definition: function.hh:175
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1338
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:725
Compute divergence of vector-valued functions.
Definition: gridfunctionspaceutilities.hh:1018
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1326