dune-fem  2.4.1-rc
lagrange/shapefunctionset.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_SPACE_LAGRANGE_SHAPEFUNCTIONSET_HH
2 #define DUNE_FEM_SPACE_LAGRANGE_SHAPEFUNCTIONSET_HH
3 
4 #include <cassert>
5 #include <cstdlib>
6 
7 #include <dune/common/fvector.hh>
8 #include <dune/common/nullptr.hh>
9 #include <dune/common/forloop.hh>
10 
11 #include <dune/geometry/genericgeometry/topologytypes.hh>
12 #include <dune/geometry/type.hh>
13 
16 
17 #include "genericbasefunctions.hh"
18 #include "genericlagrangepoints.hh"
19 
20 /*
21  @file
22  @brief Shape function set for Lagrange space
23  @author Christoph Gersbacher
24 */
25 
26 namespace Dune
27 {
28 
29  namespace Fem
30  {
31 
32  // LagrangeShapeFunctionInterface
33  // ------------------------------
34 
40  template< class FunctionSpace >
42  {
44 
45  static_assert( (FunctionSpace::dimRange == 1), "FunctionSpace must be scalar." );
46 
47  public:
49 
54 
56 
57  virtual void evaluate ( const DomainType &x, RangeType &value ) const = 0;
58  virtual void jacobian ( const DomainType &x, JacobianRangeType &jacobian ) const = 0;
59  virtual void hessian ( const DomainType &x, HessianRangeType &hessian ) const = 0;
60 
61  virtual int order () const = 0;
62 
63  virtual const ThisType *clone () const = 0;
64  };
65 
66 
67 
68  // LagrangeShapeFunction
69  // ---------------------
70 
79  template< class FunctionSpace, class GeometryType, unsigned int polOrder >
81  : public LagrangeShapeFunctionInterface< FunctionSpace >
82  {
85 
86  static const int dimension = FunctionSpace::dimDomain;
87 
88  public:
91 
92  typedef typename BaseType::DomainType DomainType;
93  typedef typename BaseType::RangeType RangeType;
96 
97  explicit LagrangeShapeFunction ( const GenericBaseFunctionType &genericShapeFunction )
98  : genericShapeFunction_( genericShapeFunction )
99  {}
100 
101  virtual void evaluate ( const DomainType &x, RangeType &value ) const;
102  virtual void jacobian ( const DomainType &x, JacobianRangeType &jacobian ) const;
103  virtual void hessian ( const DomainType &x, HessianRangeType &hessian ) const;
104 
105  virtual int order () const { return polOrder; }
106 
107  virtual const BaseType *clone () const { return new ThisType( *this ); }
108 
109  protected:
111  };
112 
113 
114 
115  // LagrangeShapeFunctionFactory
116  // ----------------------------
117 
125  template< class FunctionSpace, int maxPolOrder >
127  {
128  static const int dimension = FunctionSpace::dimDomain;
129 
130  public:
132 
133  private:
134  template< class Topology >
135  struct Switch;
136 
137  public:
138  explicit LagrangeShapeFunctionFactory ( const Dune::GeometryType &type, const int order = maxPolOrder )
139  : topologyId_( type.id() ),
140  order_( order )
141  {}
142 
143  int order () const;
144 
145  std::size_t numShapeFunctions () const;
146 
147  ShapeFunctionType *createShapeFunction( std::size_t i ) const;
148 
149  private:
150  const unsigned int topologyId_;
151  const int order_;
152  };
153 
154 
155 
156  // LagrangeShapeFunctionSet
157  // ------------------------
158 
165  template< class FunctionSpace, int polOrder >
167  : public SimpleShapeFunctionSet<
168  typename LagrangeShapeFunctionFactory< FunctionSpace, polOrder >::ShapeFunctionType
169  >
170  {
171  static_assert( (FunctionSpace::dimRange == 1), "FunctionSpace must be scalar." );
172 
175 
176  public:
177  LagrangeShapeFunctionSet ( const Dune::GeometryType &type )
178  : BaseType( ShapeFunctionFactoryType( type ) )
179  {}
180  };
181 
182 
183 
184  // Implementation of LagrangeShapeFunction
185  // ---------------------------------------
186 
187  template< class FunctionSpace, class GeometryType, unsigned int polOrder >
189  ::evaluate ( const DomainType &x, RangeType &value ) const
190  {
191  FieldVector< int, 0 > diffVariable;
192  return genericShapeFunction_.evaluate( diffVariable, x, value );
193  }
194 
195 
196  template< class FunctionSpace, class GeometryType, unsigned int polOrder >
199  {
200  FieldVector< int, 1 > diffVariable;
201  RangeType tmp;
202 
203  int &i = diffVariable[ 0 ];
204  for( i = 0; i < dimension; ++i )
205  {
206  genericShapeFunction_.evaluate( diffVariable, x, tmp );
207  jacobian[ 0 ][ i ] = tmp[ 0 ];
208  }
209  }
210 
211 
212  template< class FunctionSpace, class GeometryType, unsigned int polOrder >
215  {
216  FieldVector< int, 2 > diffVariable;
217  RangeType tmp;
218 
219  int &i = diffVariable[ 0 ];
220  for( i = 0; i < dimension; ++i )
221  {
222  // we use symmetrized evaluation of the hessian, since calling
223  // evaluate is in general quite expensive
224  int &j = diffVariable[ 1 ];
225  for( j = 0; j < i; ++j )
226  {
227  genericShapeFunction_.evaluate( diffVariable, x, tmp );
228  hessian[ 0 ][ i ][ j ] = hessian[ 0 ][ j ][ i ] = tmp[ 0 ];
229  }
230 
231  assert( j == i );
232  genericShapeFunction_.evaluate( diffVariable, x, tmp );
233  hessian[ 0 ][ i ][ i ] = tmp[ 0 ];
234  }
235  }
236 
237 
238 
239  // LagrangeShapeFunctionFactory::Switch
240  // ------------------------------------
241 
242  template< class FunctionSpace, int maxPolOrder >
243  template< class Topology >
244  struct LagrangeShapeFunctionFactory< FunctionSpace, maxPolOrder >::Switch
245  {
246  // get generic geometry type
247  static const unsigned int topologyId = Topology::id;
249  ::GenericGeometryType GenericGeometryType;
250 
251  template <int polOrd>
252  struct CheckOrder
253  {
254  // type of scalar shape function for current polynomial order
257 
258  static void apply ( const int order, std::size_t &size )
259  {
260  if( order == polOrd )
261  {
263  }
264  }
265 
266  static void apply ( const std::size_t &i, const int order, ShapeFunctionType *&shapeFunction )
267  {
268  if( order == polOrd )
269  {
270  shapeFunction = new ShapeFunctionImpl( GenericBaseFunctionType( i ) );
271  }
272  }
273  };
274 
275  static void apply ( const int order, std::size_t &size )
276  {
277  // loop over all possible polynomial order to find correct size
278  Dune::ForLoop< CheckOrder, 0, maxPolOrder >::apply( order, size );
279  }
280 
281  static void apply ( const std::size_t &i, const int order, ShapeFunctionType *&shapeFunction )
282  {
283  // loop over all possible polynomial order to create correct shape function
284  Dune::ForLoop< CheckOrder, 0, maxPolOrder >::apply( i, order, shapeFunction );
285  }
286  };
287 
288 
289 
290  // Implementation of LagrangeShapeFunctionFactory
291  // ----------------------------------------------
292 
293  template< class FunctionSpace, int polOrder >
295 
296 
297  template< class FunctionSpace, int polOrder >
299  ::order () const
300  {
301  return order_;
302  }
303 
304 
305  template< class FunctionSpace, int polOrder >
308  {
309  std::size_t numShapeFunctions;
310  GenericGeometry::IfTopology< Switch, dimension >::apply( topologyId_, order_, numShapeFunctions );
311  return numShapeFunctions;
312  }
313 
314 
315  template< class FunctionSpace, int polOrder >
318  ::createShapeFunction( const std::size_t i ) const
319  {
320  ShapeFunctionType *shapeFunction( nullptr );
321  GenericGeometry::IfTopology< Switch, dimension >::apply( topologyId_, i, order_, shapeFunction );
322  assert( shapeFunction );
323  return shapeFunction;
324  }
325 
326 
327 
328  // Extern Template Instantiations
329  // ------------------------------
330 
337 
344 
351 
352  } // namespace Fem
353 
354 } // namespace Dune
355 
356 #endif // #ifndef DUNE_FEM_SPACE_LAGRANGE_SHAPEFUNCTIONSET_HH
virtual void jacobian(const DomainType &x, JacobianRangeType &jacobian) const =0
abstract base class for Lagrange shape functions
Definition: lagrange/shapefunctionset.hh:41
virtual void evaluate(const DomainType &x, RangeType &value) const =0
Definition: genericlagrangepoints.hh:21
LagrangeShapeFunctionFactory(const Dune::GeometryType &type, const int order=maxPolOrder)
Definition: lagrange/shapefunctionset.hh:138
A vector valued function space.
Definition: functionspace.hh:16
virtual int order() const
Definition: lagrange/shapefunctionset.hh:105
virtual void hessian(const DomainType &x, HessianRangeType &hessian) const
Definition: lagrange/shapefunctionset.hh:214
virtual void evaluate(const DomainType &x, RangeType &value) const
Definition: lagrange/shapefunctionset.hh:189
VectorSpaceTraits< DomainField, RangeField, dimD, dimR >::LinearMappingType JacobianRangeType
Intrinsic type used for the jacobian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:74
virtual const ThisType * clone() const =0
FieldVector< FieldMatrix< RangeFieldType, dimDomain, dimDomain >, dimRange > HessianRangeType
Intrinsic type used for the hessian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:78
virtual void hessian(const DomainType &x, HessianRangeType &hessian) const =0
dimension of domain vector space
Definition: functionspaceinterface.hh:45
ShapeFunctionType * createShapeFunction(std::size_t i) const
Definition: lagrange/shapefunctionset.hh:318
static void apply(const int order, std::size_t &size)
Definition: lagrange/shapefunctionset.hh:258
Definition: lagrange/shapefunctionset.hh:252
dimension of range vector space
Definition: functionspaceinterface.hh:47
virtual void jacobian(const DomainType &x, JacobianRangeType &jacobian) const
Definition: lagrange/shapefunctionset.hh:198
ShapeFunctionImpl::GenericBaseFunctionType GenericBaseFunctionType
Definition: lagrange/shapefunctionset.hh:256
Definition: shapefunctionset/simple.hh:49
int order() const
Definition: lagrange/shapefunctionset.hh:299
LagrangeShapeFunction< FunctionSpace, GenericGeometryType, polOrd > ShapeFunctionImpl
Definition: lagrange/shapefunctionset.hh:255
factory class
Definition: lagrange/shapefunctionset.hh:126
FunctionSpaceType::HessianRangeType HessianRangeType
Definition: lagrange/shapefunctionset.hh:53
BaseType::DomainType DomainType
Definition: lagrange/shapefunctionset.hh:92
Definition: coordinate.hh:4
Definition: genericgeometry.hh:172
BaseType::RangeType RangeType
Definition: lagrange/shapefunctionset.hh:93
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
BaseType::JacobianRangeType JacobianRangeType
Definition: lagrange/shapefunctionset.hh:94
FunctionSpaceType::JacobianRangeType JacobianRangeType
Definition: lagrange/shapefunctionset.hh:52
LagrangeShapeFunctionInterface< FunctionSpace > ShapeFunctionType
Definition: lagrange/shapefunctionset.hh:131
static void apply(const std::size_t &i, const int order, ShapeFunctionType *&shapeFunction)
Definition: lagrange/shapefunctionset.hh:266
GenericBaseFunctionType genericShapeFunction_
Definition: lagrange/shapefunctionset.hh:110
LagrangeShapeFunction(const GenericBaseFunctionType &genericShapeFunction)
Definition: lagrange/shapefunctionset.hh:97
Lagrange shape function set.
Definition: lagrange/shapefunctionset.hh:166
virtual const BaseType * clone() const
Definition: lagrange/shapefunctionset.hh:107
GenericLagrangeBaseFunction< FunctionSpace, GeometryType, polOrder > GenericBaseFunctionType
Definition: lagrange/shapefunctionset.hh:90
implementation of Lagrange shape function using generic Lagrange shape functions
Definition: lagrange/shapefunctionset.hh:80
LagrangeShapeFunctionSet(const Dune::GeometryType &type)
Definition: lagrange/shapefunctionset.hh:177
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
BaseType::HessianRangeType HessianRangeType
Definition: lagrange/shapefunctionset.hh:95
std::size_t numShapeFunctions() const
Definition: lagrange/shapefunctionset.hh:307
FunctionSpaceType::RangeType RangeType
Definition: lagrange/shapefunctionset.hh:51
virtual ~LagrangeShapeFunctionInterface()
Definition: lagrange/shapefunctionset.hh:55
FunctionSpace FunctionSpaceType
Definition: lagrange/shapefunctionset.hh:45
FunctionSpaceType::DomainType DomainType
Definition: lagrange/shapefunctionset.hh:50