3 #ifndef DUNE_PDELAB_GRIDFUNCTIONSPACEUTILITIES_HH 4 #define DUNE_PDELAB_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" 52 template<
typename T,
typename X>
54 :
public TypeTree::LeafNode
57 typename T::Traits::GridViewType,
58 typename BasisInterfaceSwitch<
59 typename FiniteElementInterfaceSwitch<
60 typename T::Traits::FiniteElementType
64 typename FiniteElementInterfaceSwitch<
65 typename T::Traits::FiniteElementType
68 typename BasisInterfaceSwitch<
69 typename FiniteElementInterfaceSwitch<
70 typename T::Traits::FiniteElementType
74 DiscreteGridFunction<T,X>
79 typedef typename Dune::BasisInterfaceSwitch<
80 typename FiniteElementInterfaceSwitch<
81 typename T::Traits::FiniteElementType
86 typename T::Traits::GridViewType,
87 typename BasisSwitch::RangeField,
88 BasisSwitch::dimRange,
89 typename BasisSwitch::Range
103 : pgfs(stackobject_to_shared_ptr(gfs))
107 , xl(gfs.maxLocalSize())
108 , yb(gfs.maxLocalSize())
122 , xl(gfs->maxLocalSize())
123 , yb(gfs->maxLocalSize())
129 inline void evaluate (
const typename Traits::ElementType&
e,
130 const typename Traits::DomainType& x,
131 typename Traits::RangeType& y)
const 133 typedef FiniteElementInterfaceSwitch<
138 x_view.bind(lfs_cache);
141 FESwitch::basis(lfs.finiteElement()).evaluateFunction(x,yb);
143 for (
unsigned int i=0; i<yb.size(); i++)
152 return pgfs->gridView();
159 typedef typename X::template ConstLocalView<LFSCache> XView;
161 std::shared_ptr<GFS const> pgfs;
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;
179 template<
typename T,
typename X>
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
191 DiscreteGridFunctionCurl<T,X>
195 typedef typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
196 JacobianType Jacobian;
201 typename T::Traits::GridViewType,
202 typename J2C::CurlField, J2C::dimCurl,
typename J2C::Curl
211 : pgfs(stackobject_to_shared_ptr(gfs))
215 , xl(gfs.maxLocalSize())
216 , jacobian(gfs.maxLocalSize())
217 , yb(gfs.maxLocalSize())
218 , px(stackobject_to_shared_ptr(x_))
226 static const J2C& j2C = J2C();
230 x_view.bind(lfs_cache);
234 lfs.finiteElement().basis().evaluateJacobian(x,jacobian);
237 for (std::size_t i=0; i < lfs.size(); i++) {
238 j2C(jacobian[i], yb);
245 {
return pgfs->gridView(); }
252 typedef typename X::template ConstLocalView<LFSCache> XView;
254 std::shared_ptr<GFS const> pgfs;
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;
277 template<
typename GV,
typename RangeFieldType,
int dimRangeOfBasis>
280 "DiscreteGridFunctionCurl (and friends) work in 2D " 290 template<
typename GV,
typename RangeFieldType>
294 FieldVector<RangeFieldType, 2> >
296 static_assert(GV::dimensionworld == 2,
297 "World dimension of grid must be 2 for the curl of a " 298 "scalar (1D) quantity");
306 template<
typename GV,
typename RangeFieldType>
310 FieldVector<RangeFieldType, 1> >
312 static_assert(GV::dimensionworld == 2,
313 "World dimension of grid must be 2 for the curl of a" 322 template<
typename GV,
typename RangeFieldType>
326 FieldVector<RangeFieldType, 3> >
328 static_assert(GV::dimensionworld == 3,
329 "World dimension of grid must be 3 for the curl of a" 349 template<
typename T,
typename X>
352 DiscreteGridFunctionCurlTraits<
353 typename T::Traits::GridViewType,
354 typename T::Traits::FiniteElementType::Traits::
355 LocalBasisType::Traits::RangeFieldType,
356 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
358 DiscreteGridFunctionGlobalCurl<T,X> >
362 typename T::Traits::GridViewType,
363 typename T::Traits::FiniteElementType::Traits::
364 LocalBasisType::Traits::RangeFieldType,
365 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
373 typedef typename T::Traits::FiniteElementType::Traits::
374 LocalBasisType::Traits LBTraits;
383 : pgfs(stackobject_to_shared_ptr(gfs))
387 , xl(gfs.maxLocalSize())
388 , J(gfs.maxLocalSize())
389 , px(stackobject_to_shared_ptr(x_))
394 inline void evaluate (
const typename Traits::ElementType&
e,
395 const typename Traits::DomainType& x,
396 typename Traits::RangeType& y)
const 400 x_view.bind(lfs_cache);
404 lfs.finiteElement().localBasis().
405 evaluateJacobianGlobal(x,J,e.geometry());
407 for (
unsigned int i=0; i<J.size(); i++)
411 switch(
unsigned(Traits::dimRange)) {
413 y[0] += xl[i] * J[i][0][1];
414 y[1] += xl[i] * -J[i][0][0];
417 y[0] += xl[i]*(J[i][1][0] - J[i][0][1]);
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]);
433 return pgfs->gridView();
439 typedef typename X::template ConstLocalView<LFSCache> XView;
441 std::shared_ptr<GFS const> pgfs;
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;
460 template<
typename T,
typename X>
464 typename T::Traits::GridViewType,
465 typename T::Traits::FiniteElementType::Traits::LocalBasisType
466 ::Traits::RangeFieldType,
467 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
470 typename T::Traits::FiniteElementType::Traits
471 ::LocalBasisType::Traits::RangeFieldType,
472 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
474 DiscreteGridFunctionGradient<T,X> >
477 typedef typename GFS::Traits::FiniteElementType::Traits::
478 LocalBasisType::Traits LBTraits;
482 typename GFS::Traits::GridViewType,
483 typename LBTraits::RangeFieldType,
486 typename LBTraits::RangeFieldType,
501 : pgfs(stackobject_to_shared_ptr(gfs))
516 x_view.bind(lfs_cache);
519 xl.resize(lfs.size());
524 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
525 JgeoIT = e.geometry().jacobianInverseTransposed(x);
528 std::vector<typename LBTraits::JacobianType> J(lfs.size());
529 lfs.finiteElement().localBasis().evaluateJacobian(x,J);
533 for(
unsigned int i = 0; i < lfs.size(); ++i) {
536 JgeoIT.umv(J[i][0], gradphi);
539 y.axpy(xl[i], gradphi);
547 return pgfs->gridView();
553 typedef typename X::template ConstLocalView<LFSCache> XView;
555 std::shared_ptr<GFS const> pgfs;
557 mutable LFSCache lfs_cache;
558 mutable XView x_view;
559 mutable std::vector<typename Traits::RangeFieldType> xl;
566 template<
typename T,
typename X>
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
575 DiscreteGridFunctionPiola<T,X>
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
598 : pgfs(stackobject_to_shared_ptr(gfs))
602 , xl(pgfs->maxLocalSize())
603 , yb(pgfs->maxLocalSize())
607 inline void evaluate (
const typename Traits::ElementType&
e,
608 const typename Traits::DomainType& x,
609 typename Traits::RangeType& y)
const 614 x_view.bind(lfs_cache);
618 lfs.finiteElement().localBasis().evaluateFunction(x,yb);
619 typename Traits::RangeType yhat;
621 for (
unsigned int i=0; i<yb.size(); i++)
622 yhat.axpy(xl[i],yb[i]);
625 typename Traits::ElementType::Geometry::JacobianInverseTransposed
626 J = e.geometry().jacobianInverseTransposed(x);
630 y /= J.determinant();
636 return pgfs->gridView();
643 typedef typename X::template ConstLocalView<LFSCache> XView;
645 std::shared_ptr<GFS const> pgfs;
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;
665 template<
typename T,
typename X, std::
size_t dimR = T::CHILDREN>
669 typename T::Traits::GridViewType,
670 typename T::template Child<0>::Type::Traits::FiniteElementType
671 ::Traits::LocalBasisType::Traits::RangeFieldType,
674 typename T::template Child<0>::Type::Traits::FiniteElementType
675 ::Traits::LocalBasisType::Traits::RangeFieldType,
679 VectorDiscreteGridFunction<T,X>
681 public TypeTree::LeafNode
687 typename T::Traits::GridViewType,
688 typename T::template Child<0>::Type::Traits::FiniteElementType
689 ::Traits::LocalBasisType::Traits::RangeFieldType,
692 typename T::template Child<0>::Type::Traits::FiniteElementType
693 ::Traits::LocalBasisType::Traits::RangeFieldType,
703 typedef typename ChildType::Traits::FiniteElementType
704 ::Traits::LocalBasisType::Traits::RangeFieldType
RF;
705 typedef typename ChildType::Traits::FiniteElementType
706 ::Traits::LocalBasisType::Traits::RangeType
RT;
715 std::size_t start = 0)
716 : pgfs(stackobject_to_shared_ptr(gfs))
720 , xl(gfs.maxLocalSize())
721 , yb(gfs.maxLocalSize())
723 for(std::size_t i = 0; i < dimR; ++i)
724 remap[i] = i + start;
739 template<
class Remap>
742 : pgfs(stackobject_to_shared_ptr(gfs))
746 , xl(gfs.maxLocalSize())
747 , yb(gfs.maxLocalSize())
748 , px(stackobject_to_shared_ptr(x_))
750 for(std::size_t i = 0; i < dimR; ++i)
751 remap[i] = remap_[i];
754 inline void evaluate (
const typename Traits::ElementType&
e,
755 const typename Traits::DomainType& x,
756 typename Traits::RangeType& y)
const 760 x_view.bind(lfs_cache);
763 for (
unsigned int k=0; k < dimR; k++)
765 lfs.child(remap[k]).finiteElement().localBasis().
766 evaluateFunction(x,yb);
768 for (
unsigned int i=0; i<yb.size(); i++)
769 y[k] += xl[lfs.child(remap[k]).localIndex(i)]*yb[i];
776 return pgfs->gridView();
782 typedef typename X::template ConstLocalView<LFSCache> XView;
784 std::shared_ptr<GFS const> pgfs;
785 std::size_t remap[dimR];
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;
799 template<
typename T,
typename X>
803 typename T::Traits::GridViewType,
804 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
808 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
810 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain
813 VectorDiscreteGridFunctionGradient<T,X>
815 public TypeTree::LeafNode
821 typename T::Traits::GridViewType,
822 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
826 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
828 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain>
836 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
838 typedef typename LBTraits::RangeFieldType
RF;
839 typedef typename LBTraits::JacobianType
JT;
842 : pgfs(stackobject_to_shared_ptr(gfs))
846 , xl(gfs.maxLocalSize())
847 , J(gfs.maxLocalSize())
851 inline void evaluate(
const typename Traits::ElementType&
e,
852 const typename Traits::DomainType& x,
853 typename Traits::RangeType& y)
const 858 x_view.bind(lfs_cache);
863 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
864 JgeoIT = e.geometry().jacobianInverseTransposed(x);
869 for(
unsigned int k = 0; k != T::CHILDREN; ++k)
872 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
873 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
875 Dune::FieldVector<RF,LBTraits::dimDomain> gradphi;
876 for (
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++)
879 JgeoIT.umv(J[i][0], gradphi);
881 y[k].axpy(xl[lfs.child(k).localIndex(i)], gradphi);
890 return pgfs->gridView();
896 typedef typename X::template ConstLocalView<LFSCache> XView;
898 std::shared_ptr<GFS const> pgfs;
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;
920 template<
typename Mat,
typename RF, std::
size_t size>
930 Dune::FieldVector<RF,size> grad_phi(0.0);
945 template<
typename RF, std::
size_t size>
953 static inline RF
compute_derivative(
const Dune::FieldMatrix<RF,size,size>& mat,
const T& t,
const unsigned int k)
969 template<
typename RF, std::
size_t size>
977 static inline RF
compute_derivative(
const Dune::DiagonalMatrix<RF,size>& mat,
const T& t,
const unsigned int k)
979 return mat[k][k]*t[k];
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
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>,
1016 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
1018 typedef typename LBTraits::RangeFieldType
RF;
1019 typedef typename LBTraits::JacobianType
JT;
1022 : pgfs(stackobject_to_shared_ptr(gfs))
1026 , xl(gfs.maxLocalSize())
1027 , J(gfs.maxLocalSize())
1029 static_assert(LBTraits::dimDomain == T::CHILDREN,
1030 "dimDomain and number of children has to be the same");
1034 const typename Traits::DomainType& x,
1035 typename Traits::RangeType& y)
const 1040 x_view.bind(lfs_cache);
1045 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1046 JgeoIT = e.geometry().jacobianInverseTransposed(x);
1048 const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1049 Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1054 for(
unsigned int k=0; k != T::CHILDREN; ++k) {
1057 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1058 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1061 for(
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1065 typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1066 typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1067 N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1070 y += xl[lfs.child(k).localIndex(i)] * d_k_phi;
1078 return pgfs->gridView();
1084 typedef typename X::template ConstLocalView<LFSCache> XView;
1086 shared_ptr<GFS const> pgfs;
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;
1107 template<
typename T,
typename X, std::
size_t dimR = T::CHILDREN>
1115 "Curl computation can only be done in two or three dimensions");
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,
1138 typename T::template Child<0>::Type::Traits::FiniteElementType
1139 ::Traits::LocalBasisType::Traits::RangeFieldType,
1144 VectorDiscreteGridFunctionCurl<T,X>
1146 ,
public TypeTree::LeafNode
1152 typename T::Traits::GridViewType,
1153 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1157 typename T::template Child<0>::Type::Traits::FiniteElementType
1158 ::Traits::LocalBasisType::Traits::RangeFieldType,
1168 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
1170 typedef typename LBTraits::RangeFieldType
RF;
1171 typedef typename LBTraits::JacobianType
JT;
1174 : pgfs(stackobject_to_shared_ptr(gfs))
1178 , xl(gfs.maxLocalSize())
1179 , J(gfs.maxLocalSize())
1181 static_assert(LBTraits::dimDomain == T::CHILDREN,
1182 "dimDomain and number of children has to be the same");
1186 const typename Traits::DomainType& x,
1187 typename Traits::RangeType& y)
const 1192 x_view.bind(lfs_cache);
1197 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1198 JgeoIT = e.geometry().jacobianInverseTransposed(x);
1200 const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1201 Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1209 for(
unsigned int k=0; k != T::CHILDREN; ++k) {
1212 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1213 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1220 for(
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
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);
1229 y[i1] += xl[lfs.child(k).localIndex(i)] * 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);
1239 y[i2] -= xl[lfs.child(k).localIndex(i)] * d_k_phi;
1247 return pgfs->gridView();
1253 typedef typename X::template ConstLocalView<LFSCache> XView;
1255 shared_ptr<GFS const> pgfs;
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;
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
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>,
1297 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
1299 typedef typename LBTraits::RangeFieldType
RF;
1300 typedef typename LBTraits::JacobianType
JT;
1303 : pgfs(stackobject_to_shared_ptr(gfs))
1307 , xl(gfs.maxLocalSize())
1308 , J(gfs.maxLocalSize())
1310 static_assert(LBTraits::dimDomain == T::CHILDREN,
1311 "dimDomain and number of children has to be the same");
1315 const typename Traits::DomainType& x,
1316 typename Traits::RangeType& y)
const 1321 x_view.bind(lfs_cache);
1326 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1327 JgeoIT = e.geometry().jacobianInverseTransposed(x);
1329 const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1330 Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1339 for(
unsigned int k=0; k != T::CHILDREN; ++k) {
1342 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1343 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1348 for(
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
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);
1357 y += sign * xl[lfs.child(k).localIndex(i)] * d_k_phi;
1366 return pgfs->gridView();
1372 typedef typename X::template ConstLocalView<LFSCache> XView;
1374 shared_ptr<GFS const> pgfs;
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;
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
RangeFieldType RangeFieldType
Export type for range field.
Definition: function.hh:51
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
R RangeType
range type
Definition: function.hh:60
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