1#ifndef DUNE_FEM_SPACE_COMMON_INTERPOLATE_HH
2#define DUNE_FEM_SPACE_COMMON_INTERPOLATE_HH
10#include <dune/grid/common/partitionset.hh>
11#include <dune/grid/common/rangegenerators.hh>
13#include <dune/fem/common/bindguard.hh>
15#include <dune/fem/function/common/discretefunction.hh>
16#include <dune/fem/function/common/gridfunctionadapter.hh>
17#include <dune/fem/function/common/localcontribution.hh>
18#include <dune/fem/function/localfunction/const.hh>
19#include <dune/fem/space/common/capabilities.hh>
20#include <dune/fem/space/common/localinterpolation.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;
65 typedef GridFunctionAdapter< Function, GridPartType > GridFunctionType;
67 GridFunctionType uGrid(
"uGrid", u, v.space().gridPart() );
72 template<
class Gr
idFunction,
class DiscreteFunction,
unsigned int partitions >
73 static inline std::enable_if_t< std::is_convertible< GridFunction, HasLocalFunction >::value && Capabilities::hasInterpolation< typename DiscreteFunction::DiscreteFunctionSpaceType >::v >
74 interpolate (
const GridFunction &u, DiscreteFunction &v, PartitionSet< partitions > ps )
76 ConstLocalFunction< GridFunction > uLocal( u );
77 LocalContribution< DiscreteFunction, Assembly::Set > vLocal( v,
false );
78 LocalInterpolation< typename DiscreteFunction::DiscreteFunctionSpaceType >
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;
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 )
136 template<
class Po
int >
137 void jacobian (
const Point &x, JacobianRangeType &value )
const
139 const RangeFieldType weight = weight_( coordinate( x ) );
140 for(
int i = 0; i < dimRange; ++i )
144 template<
class Quadrature,
class Values >
145 auto evaluateQuadrature (
const Quadrature &quadrature, Values &values )
const
146 -> std::enable_if_t< std::is_same<
decltype( values[ 0 ] ), RangeType & >::value >
148 for(
const auto &qp : quadrature )
149 evaluate( qp, values[ qp.index() ] );
167 template<
class Gr
idFunction,
class DiscreteFunction,
class Weight >
168 inline static auto interpolate (
const GridFunction &u, DiscreteFunction &v, Weight &&weight )
169 -> std::enable_if_t< !std::is_const< std::remove_reference_t< Weight > >::value >
171 DiscreteFunction w( v );
172 interpolate( u, v, std::forward< Weight >( weight ), w );
175 template<
class Gr
idFunction,
class DiscreteFunction,
class Weight >
176 inline static auto interpolate (
const GridFunction &u, DiscreteFunction &v, Weight &&weight )
177 -> std::enable_if_t< !std::is_const< std::remove_reference_t< Weight > >::value && !std::is_base_of< HasLocalFunction, GridFunction >::value >
182 template<
class Gr
idFunction,
class DiscreteFunction,
class Weight >
183 inline static auto interpolate (
const GridFunction &u, DiscreteFunction &v, Weight &&weight, DiscreteFunction &w )
184 -> std::enable_if_t< std::is_base_of< HasLocalFunction, GridFunction >::value && Capabilities::hasInterpolation< typename DiscreteFunction::DiscreteFunctionSpaceType >::v,
185 void_t< decltype( std::declval< Weight >().setEntity( std::declval< const typename DiscreteFunction::DiscreteFunctionSpaceType::EntityType & >() ) ) > >
187 typedef typename DiscreteFunction::DiscreteFunctionSpaceType::EntityType EntityType;
189 const auto &space = w.space();
190 Impl::WeightLocalFunction< EntityType, std::remove_reference_t< typename DiscreteFunction::FunctionSpaceType >, Weight > localWeight( weight, w.order() );
191 LocalInterpolation< typename DiscreteFunction::DiscreteFunctionSpaceType > interpolation( space );
192 interpolate( u, v, [ &interpolation, &localWeight ] (
const EntityType &entity, AddLocalContribution< DiscreteFunction > &w ) {
193 auto weightGuard = bindGuard( localWeight, entity );
194 auto iGuard = bindGuard( interpolation, entity );
195 interpolation( localWeight, w );
199 template<
class Gr
idFunction,
class DiscreteFunction,
class Weight >
200 inline static auto interpolate (
const GridFunction &u, DiscreteFunction &v, Weight &&weight, DiscreteFunction &w )
201 -> std::enable_if_t< std::is_base_of< HasLocalFunction, GridFunction >::value && Capabilities::hasInterpolation< typename DiscreteFunction::DiscreteFunctionSpaceType >::v,
202 void_t< decltype( std::declval< Weight >()( std::declval< typename DiscreteFunction::DiscreteFunctionSpaceType::EntityType & >(), std::declval< AddLocalContribution< DiscreteFunction > & >() ) ) > >
204 typedef typename DiscreteFunction::DofType DofType;
210 ConstLocalFunction< GridFunction > uLocal( u );
211 AddLocalContribution< DiscreteFunction > vLocal( v ), wLocal( w );
212 LocalInterpolation< typename DiscreteFunction::DiscreteFunctionSpaceType >
213 interpolation( v.space() );
215 for(
const auto &entity : v.space() )
217 auto uGuard = bindGuard( uLocal, entity );
218 auto vGuard = bindGuard( vLocal, entity );
219 auto wGuard = bindGuard( wLocal, entity );
220 auto iGuard = bindGuard( interpolation, entity );
223 interpolation( uLocal, vLocal );
226 weight( entity, wLocal );
229 std::transform( vLocal.begin(), vLocal.end(), wLocal.begin(), vLocal.begin(), std::multiplies< DofType >() );
233 std::transform( v.dbegin(), v.dend(), w.dbegin(), v.dbegin(), [] ( DofType v, DofType w ) {
235 return (w > DofType( 0 ) ? v / w : v);
Definition: explicitfieldvector.hh:75
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
A vector valued function space.
Definition: functionspace.hh:60
static void interpolate(const GridFunction &u, DiscreteFunction &v)
perform native interpolation of a discrete function space
Definition: interpolate.hh:54
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
constexpr All all
PartitionSet for all partitions.
Definition: partitionset.hh:295
Dune namespace.
Definition: alignedallocator.hh:13
A set of PartitionType values.
Definition: partitionset.hh:137
Traits for type conversions and type information.