1#ifndef DUNE_FEM_SPACE_COMMON_INTERPOLATE_HH
2#define DUNE_FEM_SPACE_COMMON_INTERPOLATE_HH
8#include <dune/common/typetraits.hh>
10#include <dune/grid/common/partitionset.hh>
11#include <dune/grid/common/rangegenerators.hh>
53 template<
class Gr
idFunction,
class DiscreteFunction >
54 static inline void interpolate (
const GridFunction &u, DiscreteFunction &v )
60 template<
class Function,
class DiscreteFunction,
unsigned int partitions >
61 static inline std::enable_if_t< !std::is_convertible< Function, HasLocalFunction >::value >
64 typedef typename DiscreteFunction :: DiscreteFunctionSpaceType :: GridPartType GridPartType;
67 GridFunctionType uGrid(
"uGrid", u, v.space().gridPart() );
72 template<
class Gr
idFunction,
class DiscreteFunction,
unsigned int partitions >
74 interpolate (
const GridFunction &u, DiscreteFunction &v, PartitionSet< partitions > ps )
79 interpolation( v.space() );
82 for(
const auto entity : elements( v.gridPart(), ps ) )
85 auto uGuard =
bindGuard( uLocal, entity );
88 auto vGuard =
bindGuard( vLocal, entity );
91 auto iGuard =
bindGuard( interpolation, entity );
94 interpolation( uLocal, vLocal );
103 template<
class Entity,
class FunctionSpace,
class Weight >
104 struct WeightLocalFunction
106 typedef Entity EntityType;
107 typedef FunctionSpace FunctionSpaceType;
117 typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
122 explicit WeightLocalFunction ( Weight &weight,
int order ) : weight_( weight ), order_(order) {}
124 void bind (
const EntityType &entity ) { entity_ = entity; weight_.setEntity( entity ); }
127 const EntityType entity()
const {
return entity_; }
129 template<
class Po
int >
130 void evaluate (
const Point &x, RangeType &value )
const
132 const RangeFieldType weight = weight_(
coordinate( x ) );
133 for(
int i = 0; i < dimRange; ++i )
137 template<
class Quadrature,
class Values >
138 auto evaluateQuadrature (
const Quadrature &quadrature, Values &values )
const
139 -> std::enable_if_t< std::is_same<
decltype( values[ 0 ] ), RangeType & >::value >
141 for(
const auto &qp : quadrature )
142 evaluate( qp, values[ qp.index() ] );
160 template<
class Gr
idFunction,
class DiscreteFunction,
class Weight >
161 inline static auto interpolate (
const GridFunction &u, DiscreteFunction &v, Weight &&weight )
162 -> std::enable_if_t< !std::is_const< std::remove_reference_t< Weight > >::value >
164 DiscreteFunction w( v );
165 interpolate( u, v, std::forward< Weight >( weight ), w );
168 template<
class Gr
idFunction,
class DiscreteFunction,
class Weight >
169 inline static auto interpolate (
const GridFunction &u, DiscreteFunction &v, Weight &&weight )
170 -> std::enable_if_t< !std::is_const< std::remove_reference_t< Weight > >::value && !std::is_base_of< HasLocalFunction, GridFunction >::value >
175 template<
class Gr
idFunction,
class DiscreteFunction,
class Weight >
176 inline static auto interpolate (
const GridFunction &u, DiscreteFunction &v, Weight &&weight, DiscreteFunction &w )
178 void_t< decltype( std::declval< Weight >().setEntity( std::declval< const typename DiscreteFunction::DiscreteFunctionSpaceType::EntityType & >() ) ) > >
180 typedef typename DiscreteFunction::DiscreteFunctionSpaceType::EntityType EntityType;
182 const auto &space = w.space();
183 Impl::WeightLocalFunction< EntityType, std::remove_reference_t< typename DiscreteFunction::FunctionSpaceType >, Weight > localWeight( weight, w.order() );
186 auto weightGuard =
bindGuard( localWeight, entity );
187 auto iGuard =
bindGuard( interpolation, entity );
188 interpolation( localWeight, w );
192 template<
class Gr
idFunction,
class DiscreteFunction,
class Weight >
193 inline static auto interpolate (
const GridFunction &u, DiscreteFunction &v, Weight &&weight, DiscreteFunction &w )
195 void_t< decltype( std::declval< Weight >()( std::declval< typename DiscreteFunction::DiscreteFunctionSpaceType::EntityType & >(), std::declval< AddLocalContribution< DiscreteFunction > & >() ) ) > >
197 typedef typename DiscreteFunction::DofType DofType;
203 ConstLocalFunction< GridFunction > uLocal( u );
204 AddLocalContribution< DiscreteFunction > vLocal( v ), wLocal( w );
205 LocalInterpolation< typename DiscreteFunction::DiscreteFunctionSpaceType >
206 interpolation( v.space() );
208 for(
const auto &entity : v.space() )
210 auto uGuard =
bindGuard( uLocal, entity );
211 auto vGuard =
bindGuard( vLocal, entity );
212 auto wGuard =
bindGuard( wLocal, entity );
213 auto iGuard =
bindGuard( interpolation, entity );
216 interpolation( uLocal, vLocal );
219 weight( entity, wLocal );
222 std::transform( vLocal.begin(), vLocal.end(), wLocal.begin(), vLocal.begin(), std::multiplies< DofType >() );
226 std::transform( v.dbegin(), v.dend(), w.dbegin(), v.dbegin(), [] ( DofType v, DofType w ) {
228 return (w > DofType( 0 ) ? v / w : v);
static GridFunctionAdapter< Function, GridPart > gridFunctionAdapter(std::string name, const Function &function, const GridPart &gridPart, unsigned int order)
convert a function to a grid function
Definition: gridfunctionadapter.hh:384
static void interpolate(const GridFunction &u, DiscreteFunction &v)
perform native interpolation of a discrete function space
Definition: common/interpolate.hh:54
Definition: bindguard.hh:11
typename Impl::ConstLocalFunction< GridFunction >::Type ConstLocalFunction
Definition: const.hh:604
static const Point & coordinate(const Point &x)
Definition: coordinate.hh:14
static auto bindGuard(Object &object, Args &&... args) -> std::enable_if_t< isBindable< Object, Args... >::value, BindGuard< Object > >
Definition: bindguard.hh:67
Definition: common/localcontribution.hh:14
Abstract class representing a function.
Definition: common/function.hh:50
BasicGridFunctionAdapter provides local functions for a Function.
Definition: gridfunctionadapter.hh:79
determine whether a discrete function space provides a (local) interpolation
Definition: space/common/capabilities.hh:165
static const bool v
Definition: space/common/capabilities.hh:166
ExplicitFieldVector< FieldMatrix< RangeFieldType, dimDomain, dimDomain >, dimRange > HessianRangeType
Intrinsic type used for the hessian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:79
FunctionSpaceTraits::DomainFieldType DomainFieldType
Intrinsic type used for values in the domain field (usually a double)
Definition: functionspaceinterface.hh:60
FunctionSpaceTraits::RangeType RangeType
Type of range vector (using type of range field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:71
@ dimDomain
dimension of domain vector space
Definition: functionspaceinterface.hh:46
@ dimRange
dimension of range vector space
Definition: functionspaceinterface.hh:48
FunctionSpaceTraits::LinearMappingType JacobianRangeType
Intrinsic type used for the jacobian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:75
FunctionSpaceTraits::DomainType DomainType
Type of domain vector (using type of domain field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:67
FunctionSpaceTraits::RangeFieldType RangeFieldType
Intrinsic type used for values in the range field (usually a double)
Definition: functionspaceinterface.hh:63
Definition: common/localinterpolation.hh:22