2 #ifndef DUNE_PDELAB_LOCALOPERATOR_DARCYFEM_HH 3 #define DUNE_PDELAB_LOCALOPERATOR_DARCYFEM_HH 5 #include <dune/common/fvector.hh> 6 #include <dune/geometry/referenceelements.hh> 22 template<
typename P,
typename T,
typename X>
25 Dune::PDELab::GridFunctionTraits<
26 typename T::Traits::GridViewType,
27 typename T::Traits::FiniteElementType::Traits::LocalBasisType
28 ::Traits::RangeFieldType,
29 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
32 typename T::Traits::FiniteElementType::Traits
33 ::LocalBasisType::Traits::RangeFieldType,
34 T::Traits::FiniteElementType::Traits::LocalBasisType::Traits
36 DarcyVelocityFromHeadFEM<P,T,X> >
39 using LBTraits =
typename GFS::Traits::FiniteElementType::
40 Traits::LocalBasisType::Traits;
44 using LView =
typename X::template LocalView<LFSCache>;
48 typename GFS::Traits::GridViewType,
49 typename LBTraits::RangeFieldType,
52 typename LBTraits::RangeFieldType,
53 LBTraits::dimDomain> >;
67 : pgfs(stackobject_to_shared_ptr(gfs))
68 , pxg(stackobject_to_shared_ptr(x_))
69 , pp(stackobject_to_shared_ptr(p))
85 std::vector<typename Traits::RangeFieldType> xl(lfs.
size());
86 lview.bind(lfs_cache);
91 auto geo = e.geometry();
94 const typename Traits::ElementType::Geometry::JacobianInverseTransposed
95 JgeoIT(geo.jacobianInverseTransposed(x));
98 std::vector<typename LBTraits::JacobianType> J(lfs.
size());
99 lfs.finiteElement().localBasis().evaluateJacobian(x,J);
103 for(
unsigned int i = 0; i < lfs.
size(); ++i) {
106 JgeoIT.umv(J[i][0], gradphi);
109 minusgrad.axpy(-xl[i], gradphi);
114 typename P::Traits::PermTensorType A(pp->A(e,inside_cell_center_local));
121 return pgfs->gridView();
125 std::shared_ptr<const GFS> pgfs;
126 std::shared_ptr<X> pxg;
127 std::shared_ptr<const P> pp;
133 #endif // DUNE_PDELAB_LOCALOPERATOR_DARCYFEM_HH
const Entity & e
Definition: localfunctionspace.hh:111
leaf of a function tree
Definition: function.hh:298
Dune::FieldVector< GV::Grid::ctype, GV::dimension > DomainType
domain type in dim-size coordinates
Definition: function.hh:49
Traits::IndexContainer::size_type size() const
get current size
Definition: localfunctionspace.hh:206
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:118
void update()
Definition: lfsindexcache.hh:300
GV GridViewType
The type of the grid view the function lives on.
Definition: function.hh:115
const Traits::GridViewType & getGridView() const
get a reference to the GridView
Definition: darcyfem.hh:119
const P & p
Definition: constraints.hh:147
DarcyVelocityFromHeadFEM(const P &p, const GFS &gfs, X &x_)
Construct a DarcyVelocityFromHeadFEM.
Definition: darcyfem.hh:66
R RangeType
range type
Definition: function.hh:61
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Definition: darcyfem.hh:76
ReferenceElementWrapper< ReferenceElement< typename Geometry::ctype, Geometry::mydimension > > referenceElement(const Geometry &geo)
Returns the reference element for the given geometry.
Definition: referenceelements.hh:144
Dune::PDELab::GridFunctionTraits< typename GFS::Traits::GridViewType, typename LBTraits::RangeFieldType, LBTraits::dimDomain, Dune::FieldVector< typename LBTraits::RangeFieldType, LBTraits::dimDomain > > Traits
Definition: darcyfem.hh:53
Provide Darcy velocity as a vector-valued grid function.
Definition: darcyfem.hh:23
traits class holding the function signature, same as in local function
Definition: function.hh:176