dune-fem  2.4.1-rc
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 <vector>
5 #include <type_traits>
6 
7 #include <dune/grid/common/partitionset.hh>
8 #include <dune/grid/common/rangegenerators.hh>
9 
12 
13 namespace Dune
14 {
15 
16  namespace Fem
17  {
40  template< class GridFunction, class DiscreteFunction >
41  static inline void interpolate ( const GridFunction &u, DiscreteFunction &v )
42  {
43  // just call interpolate for the all partition
44  interpolate( u, v, Partitions::all );
45  }
46 
47  // interpolate implementations
48  // ---------------------------
49  template< class GridFunction, class DiscreteFunction, unsigned int partitions >
50  static inline void interpolate ( const GridFunction &u, DiscreteFunction &v, PartitionSet< partitions > ps )
51  {
52  const bool hasLocalFunction = std::is_convertible< GridFunction, HasLocalFunction >::value;
53  interpolate( u, v, ps, std::integral_constant< bool, hasLocalFunction >() );
54  }
55 
56  template< class Function, class DiscreteFunction, unsigned int partitions >
57  static inline void interpolate ( const Function &u, DiscreteFunction &v, PartitionSet< partitions > ps, std::integral_constant< bool, false > )
58  {
59  typedef typename DiscreteFunction :: DiscreteFunctionSpaceType :: GridPartType GridPartType;
60  typedef GridFunctionAdapter< Function, GridPartType > GridFunctionType;
61 
62  GridFunctionType uGrid( "uGrid", u, v.space().gridPart() );
63 
64  interpolate( uGrid, v, ps );
65  }
66 
67  template< class GridFunction, class DiscreteFunction, unsigned int partitions >
68  static inline void interpolate ( const GridFunction &u, DiscreteFunction &v, PartitionSet< partitions > ps, std::integral_constant< bool, true > )
69  {
70 
71  // obtain grid view
72  const auto& gv = static_cast< typename DiscreteFunction::GridPartType::GridViewType >( v.gridPart() );
73 
74  // reserve memory for local dof vector
75  std::vector< typename DiscreteFunction::RangeFieldType > ldv;
76  ldv.reserve( v.space().blockMapper().maxNumDofs() * DiscreteFunction::DiscreteFunctionSpaceType::localBlockSize );
77 
78  typename GridFunction::LocalFunctionType uLocal( u );
79 
80  // iterate over selected partition
81  for( const auto entity : elements( gv, ps ) )
82  {
83  // obtain local interpolation
84  const auto interpolation = v.space().interpolation( entity );
85 
86  // resize local dof vector
87  ldv.resize( v.space().basisFunctionSet( entity ).size() );
88 
89  // interpolate u locally
90  uLocal.init( entity );
91  interpolation( uLocal, ldv );
92 
93  // write local dofs into v
94  v.setLocalDofs( entity, ldv );
95  }
96  }
97 
98  } // namespace Fem
99 
100 } // namespace Dune
101 
102 #endif // #ifndef DUNE_FEM_SPACE_COMMON_INTERPOLATE_HH
Definition: coordinate.hh:4
GridFunctionAdapter provides local functions for a Function.
Definition: gridfunctionadapter.hh:36
Abstract class representing a function.
Definition: function.hh:43
static void interpolate(const GridFunction &u, DiscreteFunction &v)
perform native interpolation of a discrete function space
Definition: interpolate.hh:41