1 #ifndef DUNE_FEM_VTXPROJECTION_HH 2 #define DUNE_FEM_VTXPROJECTION_HH 16 template <
class Gr
idPartType>
18 typedef FieldVector<double,GridPartType::dimension>
DomainType;
19 template <
class EntityType>
21 volume_ = en.geometry().volume();
31 template<
class Function,
class DiscreteFunction,
class Weight >
35 typedef typename Function::LocalFunctionType LocalFunctionType;
37 typedef typename DiscreteFunction::DiscreteFunctionSpaceType DiscreteFunctionSpaceType;
38 typedef typename DiscreteFunction::LocalFunctionType LocalDiscreteFunctionType;
40 typedef typename DiscreteFunctionSpaceType::GridPartType GridPartType;
41 typedef typename DiscreteFunctionSpaceType::LagrangePointSetType LagrangePointSetType;
43 typedef typename GridPartType::IntersectionIteratorType IntersectionIteratorType;
44 typedef typename GridPartType::IntersectionType IntersectionType;
45 typedef typename GridPartType::template Codim< 0 >::EntityType EntityType;
46 typedef typename GridPartType::template Codim< 0 >::GeometryType GeometryType;;
48 typedef typename LagrangePointSetType::template Codim< 1 >::SubEntityIteratorType FaceDofIteratorType;
50 typedef typename FunctionSpaceType::RangeFieldType RangeFieldType;
51 typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
52 typedef typename FunctionSpaceType::RangeType RangeType;
53 typedef typename GeometryType::LocalCoordinate LocalCoordinateType;
55 const unsigned int dimRange = FunctionSpaceType::dimRange;
56 const DiscreteFunctionSpaceType &space = u.space();
59 DiscreteFunction w (
"weight", space );
62 for (
const auto& entity : space )
64 weight.setEntity( entity );
66 LocalDiscreteFunctionType lw = w.localFunction( entity );
67 LocalDiscreteFunctionType lu = u.localFunction( entity );
69 const LocalFunctionType lf = f.localFunction( entity );
70 const LagrangePointSetType &lagrangePointSet = space.lagrangePointSet( entity );
72 const unsigned int numPoints = lagrangePointSet.nop();
73 for(
unsigned int pt = 0; pt < numPoints; ++pt )
76 lf.
evaluate( lagrangePointSet[ pt ], val );
78 RealType wght = weight( lagrangePointSet.point( pt ) );
91 typedef typename DiscreteFunction::DofIteratorType DofIteratorType;
93 const DofIteratorType udend = u.dend();
94 DofIteratorType udit = u.dbegin();
95 DofIteratorType wdit = w.dbegin();
96 for( ; udit != udend; ++udit, ++wdit )
100 if ( weight > 1e-12 )
108 const GridPartType &gridPart = space.gridPart();
109 for(
const auto& entity : space )
111 const LagrangePointSetType &lagrangePointSet = space.lagrangePointSet( entity );
113 const IntersectionIteratorType iend = gridPart.iend( entity );
114 for( IntersectionIteratorType iit = gridPart.ibegin( entity ); iit != iend; ++iit )
116 const IntersectionType &intersection = *iit;
118 if( intersection.neighbor() )
121 EntityType neighbor =
make_entity( intersection.outside() );
124 if( entity.level() > neighbor.level() )
126 const int indexInInside = intersection.indexInInside();
128 typedef typename IntersectionType::LocalGeometry LocalGeometryType;
129 const LocalGeometryType &geoIn = intersection.geometryInInside();
130 const LocalGeometryType &geoOut = intersection.geometryInOutside();
132 LocalDiscreteFunctionType uIn = u.localFunction( entity );
133 LocalDiscreteFunctionType uOut = u.localFunction( neighbor );
135 const FaceDofIteratorType fdend = lagrangePointSet.template endSubEntity< 1 >( indexInInside );
136 FaceDofIteratorType fdit = lagrangePointSet.template beginSubEntity< 1 >( indexInInside );
137 for( ; fdit != fdend; ++fdit )
139 const LocalCoordinateType &xIn = lagrangePointSet.point( *fdit );
140 const LocalCoordinateType xOut = geoOut.global( geoIn.local( xIn ) );
143 uOut.evaluate( xOut, val );
163 template <
typename DType,
typename RType >
171 typedef typename Dune::FieldTraits< RangeFieldType >::real_type
RealType;
172 typedef typename RType :: DiscreteFunctionSpaceType :: GridPartType
GridPartType;
177 template <
class WeightType>
179 WeightType& weight)
const 184 void operator() (
const DomainType& f, RangeType& discFunc)
const 197 #endif // #ifndef DUNE_FEM_VTXPROJECTION_HH static void project(const Function &f, DiscreteFunction &u, Weight &weight)
Definition: vtxprojection.hh:32
DType DomainType
Definition: vtxprojection.hh:167
FunctionSpaceImp FunctionSpaceType
type of function space this function belongs to
Definition: function.hh:56
double operator()(const DomainType &point)
Definition: vtxprojection.hh:23
specialize with 'false' if grid part has no underlying dune grid (default=true)
Definition: gridpart/common/capabilities.hh:17
VtxProjection()
Constructor.
Definition: vtxprojection.hh:174
void setEntity(const EntityType &en)
Definition: vtxprojection.hh:20
Definition: vtxprojection.hh:17
Dune::FieldTraits< RangeFieldType >::real_type RealType
Definition: vtxprojection.hh:171
Double abs(const Double &a)
Definition: double.hh:860
abstract operator
Definition: operator.hh:25
Dune::EntityPointer< Grid, Implementation >::Entity make_entity(const Dune::EntityPointer< Grid, Implementation > &entityPointer)
Definition: compatibility.hh:23
Definition: coordinate.hh:4
FieldVector< double, GridPartType::dimension > DomainType
Definition: vtxprojection.hh:18
RType::RangeFieldType RangeFieldType
Definition: vtxprojection.hh:170
RType RangeType
Definition: vtxprojection.hh:168
RType::DiscreteFunctionSpaceType::GridPartType GridPartType
Definition: vtxprojection.hh:172
Definition: vtxprojection.hh:29
Abstract class representing a function.
Definition: function.hh:43
double volume_
Definition: vtxprojection.hh:26
DomainType::RangeFieldType DomainFieldType
Definition: vtxprojection.hh:169
static const Point & coordinate(const Point &x)
Definition: coordinate.hh:11
The Projection class which average discontinuous function in the Lagrangepoints.
Definition: vtxprojection.hh:164
void evaluate(const DomainType &x, RangeType &value) const
evaluate the function
Definition: function.hh:111