dune-fem  2.4.1-rc
rannacherturek/localinterpolation.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_SPACE_RANNACHERTUREK_LOCALINTERPOLATION_HH
2 #define DUNE_FEM_SPACE_RANNACHERTUREK_LOCALINTERPOLATION_HH
3 
4 #include <cstddef>
5 
6 #include <utility>
7 #include <vector>
8 
9 #include <dune/common/fvector.hh>
10 
11 
12 namespace Dune
13 {
14 
15  namespace Fem
16  {
17 
18  // RannacherTurekLocalInterpolation
19  // --------------------------------
20 
21  template< class BasisFunctionSet, class LocalInterpolation >
23  {
25 
26  public:
28  typedef LocalInterpolation LocalInterpolationType;
29 
30  private:
32 
33  typedef typename FunctionSpaceType::RangeType RangeType;
34  typedef typename FunctionSpaceType::RangeFieldType RangeFieldType;
35  static const int dimRange = FunctionSpaceType::dimRange;
36 
37  typedef std::size_t size_type;
38 
39  template< class LocalFunction >
40  struct LocalFunctionWrapper;
41 
42  public:
43  explicit RannacherTurekLocalInterpolation ( const BasisFunctionSetType &basisFunctionSet,
44  const LocalInterpolationType &localInterpolation = LocalInterpolationType() )
45  : basisFunctionSet_( basisFunctionSet ),
46  localInterpolation_( localInterpolation )
47  {}
48 
49  RannacherTurekLocalInterpolation ( const ThisType & ) = default;
50 
51  RannacherTurekLocalInterpolation ( ThisType &&other )
52  : basisFunctionSet_( std::move( other.basisFunctionSet_ ) ),
53  localInterpolation_( std::move( other.localInterpolation_ ) )
54  {}
55 
56  RannacherTurekLocalInterpolation &operator= ( const ThisType & ) = default;
57 
59  {
60  basisFunctionSet_ = std::move( other.basisFunctionSet_ );
61  localInterpolation_ = std::move( other.localInterpolation_ );
62  return *this;
63  }
64 
65  BasisFunctionSetType basisFunctionSet () const
66  {
67  return basisFunctionSet_;
68  }
69 
70  template< class LocalFunction, class LocalDofVector >
71  void operator() ( const LocalFunction &localFunction, LocalDofVector &localDofVector ) const
72  {
73  apply( localFunction, localDofVector );
74  }
75 
76  template< class LocalFunction, class LocalDofVector >
77  void apply ( const LocalFunction &localFunction, LocalDofVector &localDofVector ) const;
78 
79  protected:
80  const LocalInterpolationType &localInterpolation () const { return localInterpolation_; }
81 
82  private:
83  BasisFunctionSet basisFunctionSet_;
84  LocalInterpolationType localInterpolation_;
85  };
86 
87 
88 
89  // Implementation of RannacherTurekLocalInterpolation::LocalFunctionWrapper
90  // ------------------------------------------------------------------------
91 
92  template< class LocalInterpolation, class RangeVector >
93  template< class LocalFunction >
94  struct RannacherTurekLocalInterpolation< LocalInterpolation, RangeVector >::LocalFunctionWrapper
95  {
96  LocalFunctionWrapper ( const LocalFunction &localFunction, size_type component )
97  : localFunction_( localFunction ),
98  component_( component )
99  {}
100 
101  template< class RangeFieldType >
102  void evaluate ( const typename LocalFunction::DomainType &x, Dune::FieldVector< RangeFieldType, 1 > &y ) const
103  {
104  typename LocalFunction::RangeType z;
105  localFunction().evaluate( x, z );
106  y = z[ component() ];
107  }
108 
109  private:
110  const LocalFunction &localFunction () const { return localFunction_; }
111 
112  size_type component () const { return component_; }
113 
114  const LocalFunction &localFunction_;
115  size_type component_;
116  };
117 
118 
119 
120  // Implementation of RannacherTurekLocalInterpolation
121  // --------------------------------------------------
122 
123 
124  template< class LocalInterpolation, class RangeVector >
125  template< class LocalFunction, class LocalDofVector >
127  ::apply ( const LocalFunction &localFunction, LocalDofVector &localDofVector ) const
128  {
129  std::vector< RangeFieldType > phi;
130 
131  for( int k = 0; k < dimRange; ++k )
132  {
133  LocalFunctionWrapper< LocalFunction > localFunctionWrapper( localFunction, k );
134  localInterpolation().interpolate( localFunctionWrapper, phi );
135 
136  const size_type size = phi.size();
137  for( size_type i = 0; i < size; ++i )
138  localDofVector[ i*dimRange + k ] = phi[ i ];
139  }
140  }
141 
142  } // namespace Fem
143 
144 } // namespace Dune
145 
146 #endif // #ifndef DUNE_FEM_SPACE_RANNACHERTUREK_LOCALINTERPOLATION_HH
LocalInterpolation LocalInterpolationType
Definition: rannacherturek/localinterpolation.hh:28
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::RangeFieldType RangeFieldType
Intrinsic type used for values in the range field (usually a double)
Definition: functionspaceinterface.hh:62
A vector valued function space.
Definition: functionspace.hh:16
BasisFunctionSet BasisFunctionSetType
Definition: rannacherturek/localinterpolation.hh:27
dimension of range vector space
Definition: functionspaceinterface.hh:47
interface for local functions
Definition: localfunction.hh:41
RannacherTurekLocalInterpolation & operator=(const ThisType &)=default
RannacherTurekLocalInterpolation(const BasisFunctionSetType &basisFunctionSet, const LocalInterpolationType &localInterpolation=LocalInterpolationType())
Definition: rannacherturek/localinterpolation.hh:43
Definition: coordinate.hh:4
const LocalInterpolationType & localInterpolation() const
Definition: rannacherturek/localinterpolation.hh:80
STL namespace.
RannacherTurekLocalInterpolation(ThisType &&other)
Definition: rannacherturek/localinterpolation.hh:51
void move(ArrayInterface< T > &array, const unsigned int oldOffset, const unsigned int newOffset, const unsigned int length)
Definition: array_inline.hh:38
void evaluate(const PointType &x, RangeType &ret) const
evaluate the local function
Definition: localfunction.hh:300
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::RangeType RangeType
Type of range vector (using type of range field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:70
void apply(const LocalFunction &localFunction, LocalDofVector &localDofVector) const
Definition: rannacherturek/localinterpolation.hh:127
FunctionSpaceType::DomainType DomainType
type of domain vectors, i.e., type of coordinates
Definition: localfunction.hh:69
Definition: rannacherturek/localinterpolation.hh:22
BasisFunctionSetType basisFunctionSet() const
Definition: rannacherturek/localinterpolation.hh:65
Definition: basisfunctionset/basisfunctionset.hh:31
void operator()(const LocalFunction &localFunction, LocalDofVector &localDofVector) const
Definition: rannacherturek/localinterpolation.hh:71
FunctionSpaceType::RangeType RangeType
type of range vectors, i.e., type of function values
Definition: localfunction.hh:71