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>
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;
674 typename T::Traits::GridViewType,
675 typename T::template Child<0>::Type::Traits::FiniteElementType
676 ::Traits::LocalBasisType::Traits::RangeFieldType,
679 typename T::template Child<0>::Type::Traits::FiniteElementType
680 ::Traits::LocalBasisType::Traits::RangeFieldType,
684 VectorDiscreteGridFunction<T,X>
686 public TypeTree::LeafNode
692 typename T::Traits::GridViewType,
693 typename T::template Child<0>::Type::Traits::FiniteElementType
694 ::Traits::LocalBasisType::Traits::RangeFieldType,
697 typename T::template Child<0>::Type::Traits::FiniteElementType
698 ::Traits::LocalBasisType::Traits::RangeFieldType,
708 typedef typename ChildType::Traits::FiniteElementType
709 ::Traits::LocalBasisType::Traits::RangeFieldType
RF;
710 typedef typename ChildType::Traits::FiniteElementType
711 ::Traits::LocalBasisType::Traits::RangeType
RT;
720 std::size_t start = 0)
721 : pgfs(stackobject_to_shared_ptr(gfs))
725 , xl(gfs.maxLocalSize())
726 , yb(gfs.maxLocalSize())
728 for(std::size_t i = 0; i < dimR; ++i)
729 remap[i] = i + start;
744 template<
class Remap>
747 : pgfs(stackobject_to_shared_ptr(gfs))
751 , xl(gfs.maxLocalSize())
752 , yb(gfs.maxLocalSize())
753 , px(stackobject_to_shared_ptr(x_))
755 for(std::size_t i = 0; i < dimR; ++i)
756 remap[i] = remap_[i];
759 inline void evaluate (
const typename Traits::ElementType&
e,
760 const typename Traits::DomainType& x,
761 typename Traits::RangeType& y)
const 765 x_view.bind(lfs_cache);
768 for (
unsigned int k=0; k < dimR; k++)
770 lfs.child(remap[k]).finiteElement().localBasis().
771 evaluateFunction(x,yb);
773 for (
unsigned int i=0; i<yb.size(); i++)
774 y[k] += xl[lfs.child(remap[k]).localIndex(i)]*yb[i];
781 return pgfs->gridView();
787 typedef typename X::template ConstLocalView<LFSCache> XView;
789 std::shared_ptr<GFS const> pgfs;
790 std::size_t remap[dimR];
792 mutable LFSCache lfs_cache;
793 mutable XView x_view;
794 mutable std::vector<RF> xl;
795 mutable std::vector<RT> yb;
796 std::shared_ptr<const X> px;
804 template<
typename T,
typename X>
808 typename T::Traits::GridViewType,
809 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
811 TypeTree::StaticDegree<T>::value,
813 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
814 TypeTree::StaticDegree<T>::value,
815 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain
818 VectorDiscreteGridFunctionGradient<T,X>
820 public TypeTree::LeafNode
826 typename T::Traits::GridViewType,
827 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
831 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
832 TypeTree::StaticDegree<T>::value,
833 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimDomain>
841 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
843 typedef typename LBTraits::RangeFieldType
RF;
844 typedef typename LBTraits::JacobianType
JT;
847 : pgfs(stackobject_to_shared_ptr(gfs))
851 , xl(gfs.maxLocalSize())
852 , J(gfs.maxLocalSize())
856 inline void evaluate(
const typename Traits::ElementType&
e,
857 const typename Traits::DomainType& x,
858 typename Traits::RangeType& y)
const 863 x_view.bind(lfs_cache);
868 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
869 JgeoIT = e.geometry().jacobianInverseTransposed(x);
874 for(
unsigned int k = 0; k != TypeTree::degree(lfs); ++k)
877 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
878 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
880 Dune::FieldVector<RF,LBTraits::dimDomain> gradphi;
881 for (
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++)
884 JgeoIT.umv(J[i][0], gradphi);
886 y[k].axpy(xl[lfs.child(k).localIndex(i)], gradphi);
895 return pgfs->gridView();
901 typedef typename X::template ConstLocalView<LFSCache> XView;
903 std::shared_ptr<GFS const> pgfs;
905 mutable LFSCache lfs_cache;
906 mutable XView x_view;
907 mutable std::vector<RF> xl;
908 mutable std::vector<JT> J;
909 std::shared_ptr<const X> px;
925 template<
typename Mat,
typename RF, std::
size_t size>
935 Dune::FieldVector<RF,size> grad_phi(0.0);
950 template<
typename RF, std::
size_t size>
958 static inline RF
compute_derivative(
const Dune::FieldMatrix<RF,size,size>& mat,
const T& t,
const unsigned int k)
974 template<
typename RF, std::
size_t size>
982 static inline RF
compute_derivative(
const Dune::DiagonalMatrix<RF,size>& mat,
const T& t,
const unsigned int k)
984 return mat[k][k]*t[k];
998 template<
typename T,
typename X>
1001 Dune::PDELab::GridFunctionTraits<
1002 typename T::Traits::GridViewType,
1003 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1004 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1005 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1006 VectorDiscreteGridFunctionDiv<T,X> >
1007 ,
public TypeTree::LeafNode
1013 typename T::Traits::GridViewType,
1014 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1015 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1016 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1021 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
1023 typedef typename LBTraits::RangeFieldType
RF;
1024 typedef typename LBTraits::JacobianType
JT;
1027 : pgfs(stackobject_to_shared_ptr(gfs))
1031 , xl(gfs.maxLocalSize())
1032 , J(gfs.maxLocalSize())
1035 "dimDomain and number of children has to be the same");
1039 const typename Traits::DomainType& x,
1040 typename Traits::RangeType& y)
const 1045 x_view.bind(lfs_cache);
1050 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1051 JgeoIT = e.geometry().jacobianInverseTransposed(x);
1053 const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1054 Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1059 for(
unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1062 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1063 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1066 for(
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1070 typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1071 typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1072 N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1075 y += xl[lfs.child(k).localIndex(i)] * d_k_phi;
1083 return pgfs->gridView();
1089 typedef typename X::template ConstLocalView<LFSCache> XView;
1091 shared_ptr<GFS const> pgfs;
1093 mutable LFSCache lfs_cache;
1094 mutable XView x_view;
1095 mutable std::vector<RF> xl;
1096 mutable std::vector<JT> J;
1097 shared_ptr<const X> px;
1120 "Curl computation can only be done in two or three dimensions");
1134 template<
typename T,
typename X>
1137 Dune::PDELab::GridFunctionTraits<
1138 typename T::Traits::GridViewType,
1139 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1141 TypeTree::StaticDegree<T>::value,
1143 typename T::template Child<0>::Type::Traits::FiniteElementType
1144 ::Traits::LocalBasisType::Traits::RangeFieldType,
1146 TypeTree::StaticDegree<T>::value
1149 VectorDiscreteGridFunctionCurl<T,X>
1151 ,
public TypeTree::LeafNode
1157 typename T::Traits::GridViewType,
1158 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1162 typename T::template Child<0>::Type::Traits::FiniteElementType
1163 ::Traits::LocalBasisType::Traits::RangeFieldType,
1165 TypeTree::StaticDegree<T>::value
1173 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
1175 typedef typename LBTraits::RangeFieldType
RF;
1176 typedef typename LBTraits::JacobianType
JT;
1179 : pgfs(stackobject_to_shared_ptr(gfs))
1183 , xl(gfs.maxLocalSize())
1184 , J(gfs.maxLocalSize())
1186 static_assert(LBTraits::dimDomain == TypeTree::StaticDegree<T>::value,
1187 "dimDomain and number of children has to be the same");
1191 const typename Traits::DomainType& x,
1192 typename Traits::RangeType& y)
const 1197 x_view.bind(lfs_cache);
1202 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1203 JgeoIT = e.geometry().jacobianInverseTransposed(x);
1205 const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1206 Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1214 for(
unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1217 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1218 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1225 for(
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1229 typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1230 typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1231 N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1232 (JgeoIT,J[i][0],i2);
1234 y[i1] += xl[lfs.child(k).localIndex(i)] * d_k_phi;
1239 typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1240 typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1241 N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1242 (JgeoIT,J[i][0],i1);
1244 y[i2] -= xl[lfs.child(k).localIndex(i)] * d_k_phi;
1252 return pgfs->gridView();
1258 typedef typename X::template ConstLocalView<LFSCache> XView;
1260 shared_ptr<GFS const> pgfs;
1262 mutable LFSCache lfs_cache;
1263 mutable XView x_view;
1264 mutable std::vector<RF> xl;
1265 mutable std::vector<JT> J;
1266 shared_ptr<const X> px;
1279 template<
typename T,
typename X>
1282 Dune::PDELab::GridFunctionTraits<
1283 typename T::Traits::GridViewType,
1284 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1285 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1286 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1287 VectorDiscreteGridFunctionDiv<T,X> >
1288 ,
public TypeTree::LeafNode
1294 typename T::Traits::GridViewType,
1295 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeFieldType,
1296 T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::dimRange,
1297 typename T::template Child<0>::Type::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType>,
1302 typedef typename ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits
LBTraits;
1304 typedef typename LBTraits::RangeFieldType
RF;
1305 typedef typename LBTraits::JacobianType
JT;
1308 : pgfs(stackobject_to_shared_ptr(gfs))
1312 , xl(gfs.maxLocalSize())
1313 , J(gfs.maxLocalSize())
1316 "dimDomain and number of children has to be the same");
1320 const typename Traits::DomainType& x,
1321 typename Traits::RangeType& y)
const 1326 x_view.bind(lfs_cache);
1331 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
1332 JgeoIT = e.geometry().jacobianInverseTransposed(x);
1334 const typename Traits::ElementType::Geometry::JacobianInverseTransposed::size_type N =
1335 Traits::ElementType::Geometry::JacobianInverseTransposed::rows;
1344 for(
unsigned int k=0; k != TypeTree::degree(lfs); ++k) {
1347 std::vector<typename LBTraits::JacobianType> J(lfs.child(k).size());
1348 lfs.child(k).finiteElement().localBasis().evaluateJacobian(x,J);
1353 for(
typename LFS::Traits::SizeType i=0; i<lfs.child(k).size(); i++) {
1357 typename Traits::ElementType::Geometry::JacobianInverseTransposed,
1358 typename Traits::ElementType::Geometry::JacobianInverseTransposed::field_type,
1359 N>::template compute_derivative<typename LBTraits::JacobianType::row_type>
1360 (JgeoIT,J[i][0],i2);
1362 y += sign * xl[lfs.child(k).localIndex(i)] * d_k_phi;
1371 return pgfs->gridView();
1377 typedef typename X::template ConstLocalView<LFSCache> XView;
1379 shared_ptr<GFS const> pgfs;
1381 mutable LFSCache lfs_cache;
1382 mutable XView x_view;
1383 mutable std::vector<RF> xl;
1384 mutable std::vector<JT> J;
1385 shared_ptr<const X> px;
1392 #endif // DUNE_PDELAB_GRIDFUNCTIONSPACE_GRIDFUNCTIONSPACEUTILITIES_HH DiscreteGridFunction with Piola transformation.
Definition: gridfunctionspaceutilities.hh:567
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1300
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:759
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1020
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1172
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1021
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1369
static RF compute_derivative(const Dune::DiagonalMatrix< RF, size > &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:982
Equivalent of DiscreteGridFunctionGradient for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:805
DiscreteGridFunction(std::shared_ptr< const GFS > gfs, std::shared_ptr< const X > x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:117
Dune::FieldVector< GV::Grid::ctype, GV::dimension > DomainType
domain type in dim-size coordinates
Definition: function.hh:48
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:394
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1019
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:509
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:591
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1176
DiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:210
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1175
Compute curl of vector-valued functions.
Definition: gridfunctionspaceutilities.hh:1113
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:841
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:844
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:839
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:129
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1302
traits class holding the function signature, same as in local function
Definition: function.hh:175
ChildType::Traits::FiniteElementType ::Traits::LocalBasisType::Traits::RangeType RT
Definition: gridfunctionspaceutilities.hh:711
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1178
GridFunctionTraits< typename GFS::Traits::GridViewType, typename LBTraits::RangeFieldType, LBTraits::dimDomain, FieldVector< typename LBTraits::RangeFieldType, LBTraits::dimDomain > > Traits
Definition: gridfunctionspaceutilities.hh:487
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:607
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:856
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
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:779
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1081
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:244
DiscreteGridFunctionPiola(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionPiola.
Definition: gridfunctionspaceutilities.hh:597
Create a local function space from a global function space.
Definition: localfunctionspace.hh:670
ChildType::Traits::FiniteElementType::Traits::LocalBasisType::Traits LBTraits
Definition: gridfunctionspaceutilities.hh:1173
convert a grid function space and a coefficient vector into a grid function
Definition: gridfunctionspaceutilities.hh:53
DiscreteGridFunction for vector-valued functions.
Definition: gridfunctionspaceutilities.hh:671
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1038
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:431
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1307
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:843
static RF compute_derivative(const Dune::FieldMatrix< RF, size, size > &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:958
Compute divergence of vector-valued functions.
Definition: gridfunctionspaceutilities.hh:999
const Entity & e
Definition: localfunctionspace.hh:111
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:150
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
Helper class to compute a single derivative of scalar basis functions.
Definition: gridfunctionspaceutilities.hh:926
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:95
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1305
DiscreteGridFunction(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunction.
Definition: gridfunctionspaceutilities.hh:102
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:545
DiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGradient.
Definition: gridfunctionspaceutilities.hh:500
ChildType::Traits::FiniteElementType ::Traits::LocalBasisType::Traits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:709
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:840
VectorDiscreteGridFunctionDiv(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:1026
GridFunctionTraits< typename T::Traits::GridViewType, typename J2C::CurlField, J2C::dimCurl, typename J2C::Curl > Traits
Definition: gridfunctionspaceutilities.hh:203
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:222
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1190
Helper class to calculate the Traits of DiscreteGridFunctionCurl.
Definition: gridfunctionspaceutilities.hh:278
a GridFunction maps x in DomainType to y in RangeType
Definition: function.hh:186
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1304
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:706
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
convert a single component function space with a grid function representing the gradient ...
Definition: gridfunctionspaceutilities.hh:461
BaseT::Traits Traits
Definition: gridfunctionspaceutilities.hh:1171
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:893
convert a grid function space and a coefficient vector into a grid function of the curl ...
Definition: gridfunctionspaceutilities.hh:180
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, const Remap &remap_)
construct
Definition: gridfunctionspaceutilities.hh:745
static RF compute_derivative(const Mat &mat, const T &t, const unsigned int k)
Definition: gridfunctionspaceutilities.hh:932
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: gridfunctionspaceutilities.hh:1319
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:707
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:117
LBTraits::RangeFieldType RF
Definition: gridfunctionspaceutilities.hh:1023
DiscreteGridFunctionGlobalCurl(const GFS &gfs, const X &x_)
Construct a DiscreteGridFunctionGlobalCurl.
Definition: gridfunctionspaceutilities.hh:382
GV GridViewType
The type of the grid view the function lives on.
Definition: function.hh:114
VectorDiscreteGridFunctionGradient(const GFS &gfs, const X &x_)
Definition: gridfunctionspaceutilities.hh:846
void update()
Definition: lfsindexcache.hh:300
LBTraits::JacobianType JT
Definition: gridfunctionspaceutilities.hh:1024
convert a single component function space with experimental global finite elements into a grid functi...
Definition: gridfunctionspaceutilities.hh:350
VectorDiscreteGridFunction(const GFS &gfs, const X &x_, std::size_t start=0)
construct
Definition: gridfunctionspaceutilities.hh:719
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:1250
VectorDiscreteGridFunctionCurl(const GFS &gfs, const X &x)
Definition: gridfunctionspaceutilities.hh:1117
T::template Child< 0 >::Type ChildType
Definition: gridfunctionspaceutilities.hh:1301
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: gridfunctionspaceutilities.hh:634
R RangeType
range type
Definition: function.hh:60