dune-fem  2.4.1-rc
shapefunctionset/simple.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_SHAPEFUNCTIONSET_SIMPLE_HH
2 #define DUNE_FEM_SHAPEFUNCTIONSET_SIMPLE_HH
3 
4 // C++ includes
5 #include <cstddef>
6 #include <vector>
7 
9 
10 namespace Dune
11 {
12 
13  namespace Fem
14  {
15 
16  // AbstractShapeFunction
17  // ---------------------
18 
19  template< class FunctionSpace >
21  {
23 
24  public:
26 
31 
32  virtual ~AbstractShapeFunction () {}
33 
34  virtual void evaluate ( const DomainType &x, RangeType &value ) const = 0;
35 
36  virtual void jacobian ( const DomainType &x, JacobianRangeType &jacobian ) const = 0;
37 
38  virtual void hessian ( const DomainType &x, HessianRangeType &hessian ) const = 0;
39 
40  const ThisType *clone () const = 0;
41  };
42 
43 
44 
45  // SimpleShapeFunctionSet
46  // ----------------------
47 
48  template< class ShapeFunction >
50  {
52 
53  public:
54  typedef ShapeFunction ShapeFunctionType;
55 
56  typedef typename ShapeFunction::FunctionSpaceType FunctionSpaceType;
61 
62  template< class Factory >
63  explicit SimpleShapeFunctionSet ( const Factory &factory );
64 
65  SimpleShapeFunctionSet ( const ThisType &other );
66 
67  const ThisType &operator= ( const ThisType &other );
68 
70 
71  int order () const { return order_; }
72 
73  // Shape Function Set Interface Methods
74  std::size_t size () const { return shapeFunctions_.size(); }
75 
76  template< class Point, class Functor >
77  void evaluateEach ( const Point &x, Functor functor ) const;
78 
79  template< class Point, class Functor >
80  void jacobianEach ( const Point &x, Functor functor ) const;
81 
82  template< class Point, class Functor >
83  void hessianEach ( const Point &x, Functor functor ) const;
84 
85  protected:
86  std::vector< const ShapeFunctionType * > shapeFunctions_;
87  int order_;
88  };
89 
90 
91 
92  // Implementation of SimpleShapeFunctionSet
93  // ----------------------------------------
94 
95  template< class ShapeFunction >
96  template< class Factory >
98  ::SimpleShapeFunctionSet ( const Factory &factory )
99  {
100  const std::size_t numShapeFunctions = factory.numShapeFunctions();
101  shapeFunctions_.resize( numShapeFunctions );
102  for( std::size_t i = 0; i < numShapeFunctions; ++i )
103  shapeFunctions_[ i ] = factory.createShapeFunction( i );
104  order_ = factory.order();
105  }
106 
107  template< class ShapeFunction >
110  {
111  *this = other;
112  }
113 
114  template< class ShapeFunction >
117  {
118  if( this == &other )
119  return *this;
120 
121  for( std::size_t i = 0; i < size(); ++i )
122  delete shapeFunctions_[ i ];
123 
124  const std::size_t numShapeFunctions = other.size();
125  shapeFunctions_.resize( numShapeFunctions );
126  for( std::size_t i = 0; i < numShapeFunctions; ++i )
127  shapeFunctions_[ i ] = other.shapeFunctions_[ i ]->clone();
128 
129  order_ = other.order_;
130  return *this;
131  }
132 
133 
134  template< class ShapeFunction >
136  {
137  for( std::size_t i = 0; i < size(); ++i )
138  delete shapeFunctions_[ i ];
139  }
140 
141 
142  template< class ShapeFunction >
143  template< class Point, class Functor >
145  ::evaluateEach ( const Point &x, Functor functor ) const
146  {
147  for( std::size_t i = 0; i < size(); ++i )
148  {
149  RangeType value;
150  shapeFunctions_[ i ]->evaluate( coordinate( x ), value );
151  functor( i, value );
152  }
153  }
154 
155 
156  template< class ShapeFunction >
157  template< class Point, class Functor >
159  ::jacobianEach ( const Point &x, Functor functor ) const
160  {
161  for( std::size_t i = 0; i < size(); ++i )
162  {
164  shapeFunctions_[ i ]->jacobian( coordinate( x ), jacobian );
165  functor( i, jacobian );
166  }
167  }
168 
169 
170  template< class ShapeFunction >
171  template< class Point, class Functor >
173  ::hessianEach ( const Point &x, Functor functor ) const
174  {
175  for( std::size_t i = 0; i < size(); ++i )
176  {
178  shapeFunctions_[ i ]->hessian( coordinate( x ), hessian );
179  functor( i, hessian );
180  }
181  }
182 
183  } // namespace Fem
184 
185 } // namespace Dune
186 
187 #endif // #ifndef DUNE_FEM_SHAPEFUNCTIONSET_SIMPLE_HH
std::size_t size() const
Definition: shapefunctionset/simple.hh:74
FunctionSpaceType::RangeType RangeType
Definition: shapefunctionset/simple.hh:58
Definition: shapefunctionset/simple.hh:20
void evaluateEach(const Point &x, Functor functor) const
Definition: shapefunctionset/simple.hh:145
A vector valued function space.
Definition: functionspace.hh:16
void hessianEach(const Point &x, Functor functor) const
Definition: shapefunctionset/simple.hh:173
virtual void jacobian(const DomainType &x, JacobianRangeType &jacobian) const =0
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::LinearMappingType JacobianRangeType
Intrinsic type used for the jacobian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:74
FieldVector< FieldMatrix< RangeFieldType, dimDomain, dimDomain >, dimRange > HessianRangeType
Intrinsic type used for the hessian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:78
FunctionSpaceType::JacobianRangeType JacobianRangeType
Definition: shapefunctionset/simple.hh:59
std::vector< const ShapeFunctionType * > shapeFunctions_
Definition: shapefunctionset/simple.hh:86
FunctionSpaceType::JacobianRangeType JacobianRangeType
Definition: shapefunctionset/simple.hh:29
virtual void hessian(const DomainType &x, HessianRangeType &hessian) const =0
Definition: shapefunctionset/simple.hh:49
FunctionSpaceType::HessianRangeType HessianRangeType
Definition: shapefunctionset/simple.hh:60
FunctionSpaceType::DomainType DomainType
Definition: shapefunctionset/simple.hh:57
const ThisType * clone() const =0
Definition: coordinate.hh:4
void jacobianEach(const Point &x, Functor functor) const
Definition: shapefunctionset/simple.hh:159
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::DomainType DomainType
Type of domain vector (using type of domain field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:66
FunctionSpaceType::HessianRangeType HessianRangeType
Definition: shapefunctionset/simple.hh:30
~SimpleShapeFunctionSet()
Definition: shapefunctionset/simple.hh:135
ShapeFunction ShapeFunctionType
Definition: shapefunctionset/simple.hh:54
FunctionSpace FunctionSpaceType
Definition: shapefunctionset/simple.hh:25
FunctionSpaceType::RangeType RangeType
Definition: shapefunctionset/simple.hh:28
virtual ~AbstractShapeFunction()
Definition: shapefunctionset/simple.hh:32
int order_
Definition: shapefunctionset/simple.hh:87
int order() const
Definition: shapefunctionset/simple.hh:71
SimpleShapeFunctionSet(const Factory &factory)
Definition: shapefunctionset/simple.hh:98
virtual void evaluate(const DomainType &x, RangeType &value) const =0
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
FunctionSpaceType::DomainType DomainType
Definition: shapefunctionset/simple.hh:27
ShapeFunction::FunctionSpaceType FunctionSpaceType
Definition: shapefunctionset/simple.hh:56
static const Point & coordinate(const Point &x)
Definition: coordinate.hh:11
const ThisType & operator=(const ThisType &other)
Definition: shapefunctionset/simple.hh:116