dune-fem 2.8.0
Loading...
Searching...
No Matches
common/interpolate.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_SPACE_COMMON_INTERPOLATE_HH
2#define DUNE_FEM_SPACE_COMMON_INTERPOLATE_HH
3
4#include <algorithm>
5#include <type_traits>
6#include <utility>
7
8#include <dune/common/typetraits.hh>
9
10#include <dune/grid/common/partitionset.hh>
11#include <dune/grid/common/rangegenerators.hh>
12
14
21
22namespace Dune
23{
24
25 namespace Fem
26 {
27
28 // interpolate
29 // -----------
30
53 template< class GridFunction, class DiscreteFunction >
54 static inline void interpolate ( const GridFunction &u, DiscreteFunction &v )
55 {
56 // just call interpolate for the all partition
57 interpolate( u, v, Partitions::all );
58 }
59
60 template< class Function, class DiscreteFunction, unsigned int partitions >
61 static inline std::enable_if_t< !std::is_convertible< Function, HasLocalFunction >::value >
62 interpolate ( const Function &u, DiscreteFunction &v, PartitionSet< partitions > ps )
63 {
64 typedef typename DiscreteFunction :: DiscreteFunctionSpaceType :: GridPartType GridPartType;
65 typedef GridFunctionAdapter< Function, GridPartType > GridFunctionType;
66
67 GridFunctionType uGrid( "uGrid", u, v.space().gridPart() );
68
69 interpolate( uGrid, v, ps );
70 }
71
72 template< class GridFunction, 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 )
75 {
79 interpolation( v.space() );
80
81 // iterate over selected partition
82 for( const auto entity : elements( v.gridPart(), ps ) )
83 {
84 // initialize u to entity
85 auto uGuard = bindGuard( uLocal, entity );
86
87 // bind v to entity
88 auto vGuard = bindGuard( vLocal, entity );
89
90 // bind interpolation to entity
91 auto iGuard = bindGuard( interpolation, entity );
92
93 // perform local interpolation
94 interpolation( uLocal, vLocal );
95 }
96 }
97
98
99
100 namespace Impl
101 {
102
103 template< class Entity, class FunctionSpace, class Weight >
104 struct WeightLocalFunction
105 {
106 typedef Entity EntityType;
107 typedef FunctionSpace FunctionSpaceType;
108
109 typedef typename FunctionSpaceType::DomainFieldType DomainFieldType;
110 typedef typename FunctionSpaceType::RangeFieldType RangeFieldType;
111
112 typedef typename FunctionSpaceType::DomainType DomainType;
113 typedef typename FunctionSpaceType::RangeType RangeType;
114 typedef typename FunctionSpaceType::JacobianRangeType JacobianRangeType;
115 typedef typename FunctionSpaceType::HessianRangeType HessianRangeType;
116
117 typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
118
119 static constexpr int dimDomain = FunctionSpaceType::dimDomain;
120 static constexpr int dimRange = FunctionSpaceType::dimRange;
121
122 explicit WeightLocalFunction ( Weight &weight, int order ) : weight_( weight ), order_(order) {}
123
124 void bind ( const EntityType &entity ) { entity_ = entity; weight_.setEntity( entity ); }
125 void unbind () {}
126
127 const EntityType entity() const { return entity_; }
128
129 template< class Point >
130 void evaluate ( const Point &x, RangeType &value ) const
131 {
132 const RangeFieldType weight = weight_( coordinate( x ) );
133 for( int i = 0; i < dimRange; ++i )
134 value[ i ] = weight;
135 }
136
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 >
140 {
141 for( const auto &qp : quadrature )
142 evaluate( qp, values[ qp.index() ] );
143 }
144
145 int order() const
146 { return order_; }
147 private:
148 Weight &weight_;
149 int order_;
150 Entity entity_;
151 };
152
153 } // namespace Impl
154
155
156
157 // interpolate with weights
158 // ------------------------
159
160 template< class GridFunction, 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 >
163 {
164 DiscreteFunction w( v );
165 interpolate( u, v, std::forward< Weight >( weight ), w );
166 }
167
168 template< class GridFunction, 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 >
171 {
172 interpolate( gridFunctionAdapter( u, v.gridPart() ), v, std::forward< Weight >( weight ) );
173 }
174
175 template< class GridFunction, class DiscreteFunction, class Weight >
176 inline static auto interpolate ( const GridFunction &u, DiscreteFunction &v, Weight &&weight, DiscreteFunction &w )
177 -> std::enable_if_t< std::is_base_of< HasLocalFunction, GridFunction >::value && Capabilities::hasInterpolation< typename DiscreteFunction::DiscreteFunctionSpaceType >::v,
178 void_t< decltype( std::declval< Weight >().setEntity( std::declval< const typename DiscreteFunction::DiscreteFunctionSpaceType::EntityType & >() ) ) > >
179 {
180 typedef typename DiscreteFunction::DiscreteFunctionSpaceType::EntityType EntityType;
181
182 const auto &space = w.space();
183 Impl::WeightLocalFunction< EntityType, std::remove_reference_t< typename DiscreteFunction::FunctionSpaceType >, Weight > localWeight( weight, w.order() );
185 interpolate( u, v, [ &interpolation, &localWeight ] ( const EntityType &entity, AddLocalContribution< DiscreteFunction > &w ) {
186 auto weightGuard = bindGuard( localWeight, entity );
187 auto iGuard = bindGuard( interpolation, entity );
188 interpolation( localWeight, w );
189 }, w );
190 }
191
192 template< class GridFunction, class DiscreteFunction, class Weight >
193 inline static auto interpolate ( const GridFunction &u, DiscreteFunction &v, Weight &&weight, DiscreteFunction &w )
194 -> std::enable_if_t< std::is_base_of< HasLocalFunction, GridFunction >::value && Capabilities::hasInterpolation< typename DiscreteFunction::DiscreteFunctionSpaceType >::v,
195 void_t< decltype( std::declval< Weight >()( std::declval< typename DiscreteFunction::DiscreteFunctionSpaceType::EntityType & >(), std::declval< AddLocalContribution< DiscreteFunction > & >() ) ) > >
196 {
197 typedef typename DiscreteFunction::DofType DofType;
198
199 v.clear();
200 w.clear();
201
202 {
203 ConstLocalFunction< GridFunction > uLocal( u );
204 AddLocalContribution< DiscreteFunction > vLocal( v ), wLocal( w );
205 LocalInterpolation< typename DiscreteFunction::DiscreteFunctionSpaceType >
206 interpolation( v.space() );
207
208 for( const auto &entity : v.space() )
209 {
210 auto uGuard = bindGuard( uLocal, entity );
211 auto vGuard = bindGuard( vLocal, entity );
212 auto wGuard = bindGuard( wLocal, entity );
213 auto iGuard = bindGuard( interpolation, entity );
214
215 // interpolate u and store in v
216 interpolation( uLocal, vLocal );
217
218 // evaluate DoF-wise weight
219 weight( entity, wLocal );
220
221 // multiply interpolated values by weight
222 std::transform( vLocal.begin(), vLocal.end(), wLocal.begin(), vLocal.begin(), std::multiplies< DofType >() );
223 }
224 } // ensure the local contributions go out of scope, here (communication)
225
226 std::transform( v.dbegin(), v.dend(), w.dbegin(), v.dbegin(), [] ( DofType v, DofType w ) {
227 // weights are non-negative, so cancellation cannot occur
228 return (w > DofType( 0 ) ? v / w : v);
229 } );
230 }
231
232 } // namespace Fem
233
234} // namespace Dune
235
236#endif // #ifndef DUNE_FEM_SPACE_COMMON_INTERPOLATE_HH
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