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" 52 template<
typename T,
typename X>
56 typename T::Traits::GridViewType,
57 typename BasisInterfaceSwitch<
58 typename FiniteElementInterfaceSwitch<
59 typename T::Traits::FiniteElementType
63 typename FiniteElementInterfaceSwitch<
64 typename T::Traits::FiniteElementType
67 typename BasisInterfaceSwitch<
68 typename FiniteElementInterfaceSwitch<
69 typename T::Traits::FiniteElementType
73 DiscreteGridFunction<T,X>
78 typedef typename Dune::BasisInterfaceSwitch<
79 typename FiniteElementInterfaceSwitch<
80 typename T::Traits::FiniteElementType
85 typename T::Traits::GridViewType,
86 typename BasisSwitch::RangeField,
87 BasisSwitch::dimRange,
88 typename BasisSwitch::Range
102 : pgfs(stackobject_to_shared_ptr(gfs))
106 , xl(gfs.maxLocalSize())
107 , yb(gfs.maxLocalSize())
121 , xl(gfs->maxLocalSize())
122 , yb(gfs->maxLocalSize())
128 inline void evaluate (
const typename Traits::ElementType&
e,
129 const typename Traits::DomainType& x,
130 typename Traits::RangeType& y)
const 132 typedef FiniteElementInterfaceSwitch<
137 x_view.bind(lfs_cache);
140 FESwitch::basis(lfs.finiteElement()).evaluateFunction(x,yb);
142 for (
unsigned int i=0; i<yb.size(); i++)
151 return pgfs->gridView();
158 typedef typename X::template ConstLocalView<LFSCache> XView;
160 std::shared_ptr<GFS const> pgfs;
162 mutable LFSCache lfs_cache;
163 mutable XView x_view;
164 mutable std::vector<typename Traits::RangeFieldType> xl;
165 mutable std::vector<typename Traits::RangeType> yb;
166 std::shared_ptr<const X> px;
178 template<
typename T,
typename X>
182 typename T::Traits::GridViewType,
183 typename JacobianToCurl<typename T::Traits::FiniteElementType::
184 Traits::LocalBasisType::Traits::JacobianType>::CurlField,
185 JacobianToCurl<typename T::Traits::FiniteElementType::Traits::LocalBasisType::
186 Traits::JacobianType>::dimCurl,
187 typename JacobianToCurl<typename T::Traits::FiniteElementType::
188 Traits::LocalBasisType::Traits::JacobianType>::Curl
190 DiscreteGridFunctionCurl<T,X>
194 typedef typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
195 JacobianType Jacobian;
200 typename T::Traits::GridViewType,
201 typename J2C::CurlField, J2C::dimCurl,
typename J2C::Curl
210 : pgfs(stackobject_to_shared_ptr(gfs))
214 , xl(gfs.maxLocalSize())
215 , jacobian(gfs.maxLocalSize())
216 , yb(gfs.maxLocalSize())
217 , px(stackobject_to_shared_ptr(x_))
225 static const J2C& j2C = J2C();
229 x_view.bind(lfs_cache);
233 lfs.finiteElement().basis().evaluateJacobian(x,jacobian);
236 for (std::size_t i=0; i < lfs.size(); i++) {
237 j2C(jacobian[i], yb);
244 {
return pgfs->gridView(); }
251 typedef typename X::template ConstLocalView<LFSCache> XView;
253 std::shared_ptr<GFS const> pgfs;
255 mutable LFSCache lfs_cache;
256 mutable XView x_view;
257 mutable std::vector<typename Traits::RangeFieldType> xl;
258 mutable std::vector<Jacobian> jacobian;
259 mutable std::vector<typename Traits::RangeType> yb;
260 std::shared_ptr<const X> px;
276 template<
typename GV,
typename RangeFieldType,
int dimRangeOfBasis>
279 "DiscreteGridFunctionCurl (and friends) work in 2D " 289 template<
typename GV,
typename RangeFieldType>
293 FieldVector<RangeFieldType, 2> >
295 static_assert(GV::dimensionworld == 2,
296 "World dimension of grid must be 2 for the curl of a " 297 "scalar (1D) quantity");
305 template<
typename GV,
typename RangeFieldType>
309 FieldVector<RangeFieldType, 1> >
311 static_assert(GV::dimensionworld == 2,
312 "World dimension of grid must be 2 for the curl of a" 321 template<
typename GV,
typename RangeFieldType>
325 FieldVector<RangeFieldType, 3> >
327 static_assert(GV::dimensionworld == 3,
328 "World dimension of grid must be 3 for the curl of a" 348 template<
typename T,
typename X>
351 DiscreteGridFunctionCurlTraits<
352 typename T::Traits::GridViewType,
353 typename T::Traits::FiniteElementType::Traits::
354 LocalBasisType::Traits::RangeFieldType,
355 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
357 DiscreteGridFunctionGlobalCurl<T,X> >
361 typename T::Traits::GridViewType,
362 typename T::Traits::FiniteElementType::Traits::
363 LocalBasisType::Traits::RangeFieldType,
364 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::
372 typedef typename T::Traits::FiniteElementType::Traits::
373 LocalBasisType::Traits LBTraits;
382 : pgfs(stackobject_to_shared_ptr(gfs))
386 , xl(gfs.maxLocalSize())
387 , J(gfs.maxLocalSize())
388 , px(stackobject_to_shared_ptr(x_))
393 inline void evaluate (
const typename Traits::ElementType&
e,
394 const typename Traits::DomainType& x,
395 typename Traits::RangeType& y)
const 399 x_view.bind(lfs_cache);
403 lfs.finiteElement().localBasis().
404 evaluateJacobianGlobal(x,J,e.geometry());
406 for (
unsigned int i=0; i<J.size(); i++)
410 switch(
unsigned(Traits::dimRange)) {
412 y[0] += xl[i] * J[i][0][1];
413 y[1] += xl[i] * -J[i][0][0];
416 y[0] += xl[i]*(J[i][1][0] - J[i][0][1]);
419 y[0] += xl[i]*(J[i][2][1] - J[i][1][2]);
420 y[1] += xl[i]*(J[i][0][2] - J[i][2][0]);
421 y[2] += xl[i]*(J[i][1][0] - J[i][0][1]);
432 return pgfs->gridView();
438 typedef typename X::template ConstLocalView<LFSCache> XView;
440 std::shared_ptr<GFS const> pgfs;
442 mutable LFSCache lfs_cache;
443 mutable XView x_view;
444 mutable std::vector<typename Traits::RangeFieldType> xl;
445 mutable std::vector<typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::JacobianType> J;
446 std::shared_ptr<const X> px;
459 template<
typename T,
typename X>
463 typename T::Traits::GridViewType,
464 typename T::Traits::FiniteElementType::Traits::LocalBasisType
465 ::Traits::RangeFieldType,
466 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
469 typename T::Traits::FiniteElementType::Traits
470 ::LocalBasisType::Traits::RangeFieldType,
471 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
473 DiscreteGridFunctionGradient<T,X> >
476 typedef typename GFS::Traits::FiniteElementType::Traits::
477 LocalBasisType::Traits LBTraits;
481 typename GFS::Traits::GridViewType,
482 typename LBTraits::RangeFieldType,
485 typename LBTraits::RangeFieldType,
500 : pgfs(stackobject_to_shared_ptr(gfs))
515 x_view.bind(lfs_cache);
518 xl.resize(lfs.size());
523 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
524 JgeoIT = e.geometry().jacobianInverseTransposed(x);
527 std::vector<typename LBTraits::JacobianType> J(lfs.size());
528 lfs.finiteElement().localBasis().evaluateJacobian(x,J);
532 for(
unsigned int i = 0; i < lfs.size(); ++i) {
535 JgeoIT.umv(J[i][0], gradphi);
538 y.axpy(xl[i], gradphi);
546 return pgfs->gridView();
552 typedef typename X::template ConstLocalView<LFSCache> XView;
554 std::shared_ptr<GFS const> pgfs;
556 mutable LFSCache lfs_cache;
557 mutable XView x_view;
558 mutable std::vector<typename Traits::RangeFieldType> xl;
565 template<
typename T,
typename X>
569 typename T::Traits::GridViewType,
570 typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
571 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
572 typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType
574 DiscreteGridFunctionPiola<T,X>
581 typename T::Traits::GridViewType,
582 typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
583 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
584 typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType
597 : pgfs(stackobject_to_shared_ptr(gfs))
601 , xl(pgfs->maxLocalSize())
602 , yb(pgfs->maxLocalSize())
606 inline void evaluate (
const typename Traits::ElementType&
e,
607 const typename Traits::DomainType& x,
608 typename Traits::RangeType& y)
const 613 x_view.bind(lfs_cache);
617 lfs.finiteElement().localBasis().evaluateFunction(x,yb);
618 typename Traits::RangeType yhat;
620 for (
unsigned int i=0; i<yb.size(); i++)
621 yhat.axpy(xl[i],yb[i]);
624 typename Traits::ElementType::Geometry::JacobianInverseTransposed
625 J = e.geometry().jacobianInverseTransposed(x);
629 y /= J.determinant();
635 return pgfs->gridView();
642 typedef typename X::template ConstLocalView<LFSCache> XView;
644 std::shared_ptr<GFS const> pgfs;
646 mutable LFSCache lfs_cache;
647 mutable XView x_view;
648 mutable std::vector<typename Traits::RangeFieldType> xl;
649 mutable std::vector<typename Traits::RangeType> yb;
673 typename T::Traits::GridViewType,
674 typename T::template Child<0>::Type::Traits::FiniteElementType
675 ::Traits::LocalBasisType::Traits::RangeFieldType,
678 typename T::template Child<0>::Type::Traits::FiniteElementType
679 ::Traits::LocalBasisType::Traits::RangeFieldType,
683 VectorDiscreteGridFunction<T,X>
690 typename T::Traits::GridViewType,
691 typename T::template Child<0>::Type::Traits::FiniteElementType
692 ::Traits::LocalBasisType::Traits::RangeFieldType,
695 typename T::template Child<0>::Type::Traits::FiniteElementType
696 ::Traits::LocalBasisType::Traits::RangeFieldType,
706 typedef typename ChildType::Traits::FiniteElementType
707 ::Traits::LocalBasisType::Traits::RangeFieldType
RF;
708 typedef typename ChildType::Traits::FiniteElementType
709 ::Traits::LocalBasisType::Traits::RangeType
RT;
718 std::size_t start = 0)
719 : pgfs(stackobject_to_shared_ptr(gfs))
723 , xl(gfs.maxLocalSize())
724 , yb(gfs.maxLocalSize())
726 for(std::size_t i = 0; i < dimR; ++i)
727 remap[i] = i + start;
742 template<
class Remap>
745 : pgfs(stackobject_to_shared_ptr(gfs))
749 , xl(gfs.maxLocalSize())
750 , yb(gfs.maxLocalSize())
751 , px(stackobject_to_shared_ptr(x_))
753 for(std::size_t i = 0; i < dimR; ++i)
754 remap[i] = remap_[i];
757 inline void evaluate (
const typename Traits::ElementType&
e,
758 const typename Traits::DomainType& x,
759 typename Traits::RangeType& y)
const 763 x_view.bind(lfs_cache);
766 for (
unsigned int k=0; k < dimR; k++)
768 lfs.child(remap[k]).finiteElement().localBasis().
769 evaluateFunction(x,yb);
771 for (
unsigned int i=0; i<yb.size(); i++)
772 y[k] += xl[lfs.child(remap[k]).localIndex(i)]*yb[i];
779 return pgfs->gridView();
785 typedef typename X::template ConstLocalView<LFSCache> XView;
787 std::shared_ptr<GFS const> pgfs;
788 std::size_t remap[dimR];
790 mutable LFSCache lfs_cache;
791 mutable XView x_view;
792 mutable std::vector<RF> xl;
793 mutable std::vector<RT> yb;
794 std::shared_ptr<const X> px;
802 template<
typename T,
typename X>
806 typename T::Traits::GridViewType,
807 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
809 TypeTree::StaticDegree<T>::value,
811 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
812 TypeTree::StaticDegree<T>::value,
813 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain
816 VectorDiscreteGridFunctionGradient<T,X>
823 typename T::Traits::GridViewType,
824 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
828 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
829 TypeTree::StaticDegree<T>::value,
830 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain>
838 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
840 typedef typename LBTraits::RangeFieldType
RF;
841 typedef typename LBTraits::JacobianType
JT;
844 : pgfs(stackobject_to_shared_ptr(gfs))
848 , xl(gfs.maxLocalSize())
849 , J(gfs.maxLocalSize())
853 inline void evaluate(
const typename Traits::ElementType&
e,
854 const typename Traits::DomainType& x,
855 typename Traits::RangeType& y)
const 860 x_view.bind(lfs_cache);
865 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
866 JgeoIT = e.geometry().jacobianInverseTransposed(x);
871 for(
unsigned int k = 0; k != TypeTree::degree(lfs); ++k)
874 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
875 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
877 Dune::FieldVector<RF,LBTraits::dimDomain> gradphi;
878 for (
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++)
881 JgeoIT.umv(J[i][0], gradphi);
883 y[k].axpy(xl[lfs.child(k).localIndex(i)], gradphi);
892 return pgfs->gridView();
898 typedef typename X::template ConstLocalView<LFSCache> XView;
900 std::shared_ptr<GFS const> pgfs;
902 mutable LFSCache lfs_cache;
903 mutable XView x_view;
904 mutable std::vector<RF> xl;
905 mutable std::vector<JT> J;
906 std::shared_ptr<const X> px;
922 template<
typename Mat,
typename RF, std::
size_t size>
932 Dune::FieldVector<RF,size> grad_phi(0.0);
947 template<
typename RF, std::
size_t size>
955 static inline RF
compute_derivative(
const Dune::FieldMatrix<RF,size,size>& mat,
const T& t,
const unsigned int k)
971 template<
typename RF, std::
size_t size>
979 static inline RF
compute_derivative(
const Dune::DiagonalMatrix<RF,size>& mat,
const T& t,
const unsigned int k)
981 return mat[k][k]*t[k];
995 template<
typename T,
typename X>
998 Dune::PDELab::GridFunctionTraits<
999 typename T::Traits::GridViewType,
1000 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1001 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1002 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1003 VectorDiscreteGridFunctionDiv<T,X> >
1009 typename T::Traits::GridViewType,
1010 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1011 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1012 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1017 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
1019 typedef typename LBTraits::RangeFieldType
RF;
1020 typedef typename LBTraits::JacobianType
JT;
1023 : pgfs(stackobject_to_shared_ptr(gfs))
1027 , xl(gfs.maxLocalSize())
1028 , J(gfs.maxLocalSize())
1031 "dimDomain and number of children has to be the same");
1035 const typename Traits::DomainType& x,
1036 typename Traits::RangeType& y)
const 1041 x_view.bind(lfs_cache);
1046 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1047 JgeoIT = e.geometry().jacobianInverseTransposed(x);
1049 const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1050 Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1055 for(
unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1058 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1059 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1062 for(
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1066 typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1067 typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1068 N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1071 y += xl[lfs.child(k).localIndex(i)] * d_k_phi;
1079 return pgfs->gridView();
1085 typedef typename X::template ConstLocalView<LFSCache> XView;
1087 shared_ptr<GFS const> pgfs;
1089 mutable LFSCache lfs_cache;
1090 mutable XView x_view;
1091 mutable std::vector<RF> xl;
1092 mutable std::vector<JT> J;
1093 shared_ptr<const X> px;
1116 "Curl computation can only be done in two or three dimensions");
1130 template<
typename T,
typename X>
1133 Dune::PDELab::GridFunctionTraits<
1134 typename T::Traits::GridViewType,
1135 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1137 TypeTree::StaticDegree<T>::value,
1139 typename T::template Child<0>::Type::Traits::FiniteElementType
1140 ::Traits::LocalBasisType::Traits::RangeFieldType,
1142 TypeTree::StaticDegree<T>::value
1145 VectorDiscreteGridFunctionCurl<T,X>
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,
1160 TypeTree::StaticDegree<T>::value
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 == TypeTree::StaticDegree<T>::value,
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 != TypeTree::degree(lfs); ++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> >
1288 typename T::Traits::GridViewType,
1289 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1290 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1291 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1296 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
1298 typedef typename LBTraits::RangeFieldType
RF;
1299 typedef typename LBTraits::JacobianType
JT;
1302 : pgfs(stackobject_to_shared_ptr(gfs))
1306 , xl(gfs.maxLocalSize())
1307 , J(gfs.maxLocalSize())
1310 "dimDomain and number of children has to be the same");
1314 const typename Traits::DomainType& x,
1315 typename Traits::RangeType& y)
const 1320 x_view.bind(lfs_cache);
1325 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1326 JgeoIT = e.geometry().jacobianInverseTransposed(x);
1328 const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1329 Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1338 for(
unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1341 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1342 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1347 for(
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1351 typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1352 typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1353 N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1354 (JgeoIT,J[i][0],i2);
1356 y += sign * xl[lfs.child(k).localIndex(i)] * d_k_phi;
1365 return pgfs->gridView();
1371 typedef typename X::template ConstLocalView<LFSCache> XView;
1373 shared_ptr<GFS const> pgfs;
1375 mutable LFSCache lfs_cache;
1376 mutable XView x_view;
1377 mutable std::vector<RF> xl;
1378 mutable std::vector<JT> J;
1379 shared_ptr<const X> px;
1386 #endif // DUNE_PDELAB_GRIDFUNCTIONSPACE_GRIDFUNCTIONSPACEUTILITIES_HH GridFunctionTraits< typename T::Traits::GridViewType, typename J2C::CurlField, J2C::dimCurl, typename J2C::Curl > Traits
Definition: gridfunctionspaceutilities.hh:202
DiscreteGridFunction with Piola transformation.
Definition: gridfunctionspaceutilities.hh:566
const Entity & e
Definition: localfunctionspace.hh:111
Helper class to calculate the Traits of DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:277
static RF compute_derivative(const Dune::DiagonalMatrix< RF, size > &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:979
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1301
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1173
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1166
leaf of a function tree
Definition: function.hh:298
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:221
DiscreteGridFunctionGlobalCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGlobalCurl.
Definition: gridfunctionspaceutilities.hh:381
convert a single component function space with a grid function representing the gradient ...
Definition: gridfunctionspaceutilities.hh:460
Dune::FieldVector< GV::Grid::ctype, GV::dimension > DomainType
domain type in dim-size coordinates
Definition: function.hh:49
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1313
DiscreteGridFunction(std::shared_ptr< const GFS > gfs, std::shared_ptr< const X > x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:116
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:149
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1168
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1020
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1295
DiscreteGridFunction for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:670
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:393
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:118
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1034
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:705
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1077
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1298
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:840
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:704
void update()
Definition: lfsindexcache.hh:300
GV GridViewType
The type of the grid view the function lives on.
Definition: function.hh:115
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:94
Compute divergence of vector-valued functions.
Definition: gridfunctionspaceutilities.hh:996
VectorDiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:843
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:838
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:590
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:128
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:777
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:606
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1015
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1185
static RF compute_derivative(const Dune::FieldMatrix< RF, size, size > &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:955
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1167
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
DiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGradient.
Definition: gridfunctionspaceutilities.hh:499
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1170
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1245
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:544
Helper class to compute a single derivative of scalar basis functions.
Definition: gridfunctionspaceutilities.hh:923
DiscreteGridFunctionPiola(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionPiola.
Definition: gridfunctionspaceutilities.hh:596
DiscreteGridFunction(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:101
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:836
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:853
T Traits
Export type traits.
Definition: function.hh:192
R RangeType
range type
Definition: function.hh:61
VectorDiscreteGridFunctionDiv(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1022
Equivalent of DiscreteGridFunctionGradient for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:803
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:837
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1296
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1294
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:890
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:430
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:841
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1017
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:757
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1171
RangeFieldType RangeFieldType
Export type for range field.
Definition: function.hh:52
Compute curl of vector-valued functions.
Definition: gridfunctionspaceutilities.hh:1109
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, std::size_t start=0)
construct
Definition: gridfunctionspaceutilities.hh:717
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x)
Definition: gridfunctionspaceutilities.hh:1113
static RF compute_derivative(const Mat &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:929
convert a grid function space and a coefficient vector into a grid function of the curl ...
Definition: gridfunctionspaceutilities.hh:179
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1363
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1019
extract the curl of a function from the jacobian of that function
Definition: jacobiantocurl.hh:27
DiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:209
DiscreteGridFunctionCurlTraits< typename T::Traits::GridViewType, typename T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType, T::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange > Traits
Definition: gridfunctionspaceutilities.hh:365
ChildType::Traits::FiniteElementType ::Traits::LocalBasisType::Traits::RangeType RT
Definition: gridfunctionspaceutilities.hh:709
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:508
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:243
GridFunctionTraits< typename GFS::Traits::GridViewType, typename LBTraits::RangeFieldType, LBTraits::dimDomain, FieldVector< typename LBTraits::RangeFieldType, LBTraits::dimDomain > > Traits
Definition: gridfunctionspaceutilities.hh:486
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1299
convert a single component function space with experimental global finite elements into a grid functi...
Definition: gridfunctionspaceutilities.hh:349
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
Create a local function space from a global function space.
Definition: localfunctionspace.hh:670
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:633
convert a grid function space and a coefficient vector into a grid function
Definition: gridfunctionspaceutilities.hh:53
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, const Remap &remap_)
construct
Definition: gridfunctionspaceutilities.hh:743
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1016
ChildType::Traits::FiniteElementType ::Traits::LocalBasisType::Traits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:707
traits class holding the function signature, same as in local function
Definition: function.hh:176