BasisFunction.h 13.8 KB
Newer Older
1
2
3
4
5
6
// ============================================================================
// ==                                                                        ==
// == AMDiS - Adaptive multidimensional simulations                          ==
// ==                                                                        ==
// ============================================================================
// ==                                                                        ==
7
// ==  TU Dresden                                                            ==
8
// ==                                                                        ==
9
10
11
// ==  Institut für Wissenschaftliches Rechnen                               ==
// ==  Zellescher Weg 12-14                                                  ==
// ==  01069 Dresden                                                         ==
12
13
14
15
// ==  germany                                                               ==
// ==                                                                        ==
// ============================================================================
// ==                                                                        ==
16
// ==  https://gforge.zih.tu-dresden.de/projects/amdis/                      ==
17
18
19
20
21
22
23
24
25
// ==                                                                        ==
// ============================================================================

/** \file BasisFunction.h */

#ifndef AMDIS_BASISFUNCTION_H
#define AMDIS_BASISFUNCTION_H

#include <string>
26
#include "AMDiS_fwd.h"
27
28
#include "Global.h"
#include "Boundary.h"
Thomas Witkowski's avatar
Thomas Witkowski committed
29
#include "MatrixVector.h"
30
31
32

namespace AMDiS {

33
  /// Function interface for evaluating basis functions.
34
35
36
  class BasFctType
  {
  public:
37
    BasFctType() {}
38

39
    virtual ~BasFctType() {}
40
41
42
43

    virtual double operator()(const DimVec<double>&) const = 0;
  };

44
  /// Function interface for evaluating gradients of basis functions. 
45
46
47
  class GrdBasFctType
  {
  public:
48
    GrdBasFctType() {}
49

50
    virtual ~GrdBasFctType() {}
51

52
    virtual void operator()(const DimVec<double>&, DimVec<double>&) const = 0;
53
54
  };
  
55
  /// Function interface for evaluating second derivative of basis functions.
56
57
58
  class D2BasFctType
  {
  public:
59
    D2BasFctType() {}
60

61
    virtual ~D2BasFctType() {}
62

63
    virtual void operator()(const DimVec<double>&, DimMat<double>&) const = 0;
64
65
  };
			    
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
  typedef BasFctType *BFptr;
  typedef GrdBasFctType *GBFptr;
  typedef D2BasFctType *DBFptr;

  /** \ingroup FEMSpace
   * \brief
   * Base class for finite element basis functions. In order to build up a
   * finite element space, we have to specify a set of local basis functions.
   * Together with the correspondig DOF administration and the underlying mesh,
   * the finite element space is given. 
   * This class holds the local basis functions and their derivatives of the
   * reference element. They are evaluated at barycentric coordinates, so they
   * can be used on every element of the mesh.  
   */
  class BasisFunction
  {  
  protected:
83
    /// Creates a BasisFunction object of given dim and degree 
Thomas Witkowski's avatar
Thomas Witkowski committed
84
    BasisFunction(std::string name, int dim, int degree);
85

86
    /// destructor
87
88
89
    virtual ~BasisFunction();

  public:
90
    /// compares two BasisFunction objects.
Thomas Witkowski's avatar
Thomas Witkowski committed
91
92
    virtual bool operator==(const BasisFunction& a) const 
    {
93
      return a.getName() == name;
94
    }
95

96
    /// Returns !(*this == b)
Thomas Witkowski's avatar
Thomas Witkowski committed
97
98
    inline bool operator!=(const BasisFunction& b) const 
    {
99
      return !(operator == (b));
100
    }
101

102
    /// Used by \ref getDOFIndices and \ref getVec
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
    virtual int* orderOfPositionIndices(const Element* el, GeoIndex position, 
					int positionIndex) const = 0;

    /** \brief
     * Pointer to a function which connects the set of local basis functions
     * with its global DOFs.
     * getDOFIndices(el, admin, dof) returns a pointer to a const vector of 
     * length \ref nBasFcts where the i-th entry is the index of the DOF 
     * associated to the i-th basis function; arguments are the actual element 
     * el and the DOF admin admin of the corresponding finite element space 
     * (these indices depend on all defined DOF admins and thus on the 
     * corresponding admin); if the last argument dof is NULL, getDOFndices 
     * has to provide memory for storing this vector, which is overwritten on the
     * next call of getDOFIndices; if dof is not NULL, dof is a pointer to a 
     * vector which has to be filled;   
     */
    virtual const DegreeOfFreedom* getDOFIndices(const Element*,
						 const DOFAdmin&, 
						 DegreeOfFreedom*) const = 0;

    /** \brief
Thomas Witkowski's avatar
Thomas Witkowski committed
124
125
126
127
128
129
130
131
132
     * The second argument 'bound' has to be a pointer to a vector which has 
     * to be filled. Its length is \ref nBasFcts (the number of basis functions
     * in the used finite element space). After calling this function, the i-th 
     * entry of the array is the boundary type of the i-th basis function of this
     * element.
     * 
     * This function needs boundary information within the ElInfo object; thus, 
     * all routines using this function on the elements need the FILL_BOUND 
     * flag during mesh traversal;
133
     */
134
    virtual void getBound(const ElInfo*, BoundaryType *) const {}
135

Thomas Witkowski's avatar
Thomas Witkowski committed
136
    /// Returns \ref degree of BasisFunction
Thomas Witkowski's avatar
Thomas Witkowski committed
137
138
    inline const int getDegree() const 
    { 
139
      return degree; 
140
    }
141

Thomas Witkowski's avatar
Thomas Witkowski committed
142
    /// Returns \ref dim of BasisFunction
Thomas Witkowski's avatar
Thomas Witkowski committed
143
144
    inline const int getDim() const 
    { 
145
      return dim; 
146
    }
147

Thomas Witkowski's avatar
Thomas Witkowski committed
148
    /// Returns \ref nBasFcts which is the number of local basis functions
Thomas Witkowski's avatar
Thomas Witkowski committed
149
150
    inline const int getNumber() const 
    { 
151
      return nBasFcts; 
152
    }
153

Thomas Witkowski's avatar
Thomas Witkowski committed
154
    /// Returns \ref name of BasisFunction
Thomas Witkowski's avatar
Thomas Witkowski committed
155
    inline std::string getName() const 
Thomas Witkowski's avatar
Thomas Witkowski committed
156
    { 
157
      return name; 
158
    }
159

Thomas Witkowski's avatar
Thomas Witkowski committed
160
    /// Returns \ref nDOF[i]
161
162
    int getNumberOfDOFs(int i) const;

Thomas Witkowski's avatar
Thomas Witkowski committed
163
    /// Returns \ref nDOF
Thomas Witkowski's avatar
Thomas Witkowski committed
164
165
    inline DimVec<int>* getNumberOfDOFs() const 
    { 
166
      return nDOF; 
167
    }
168

Thomas Witkowski's avatar
Thomas Witkowski committed
169
    /// Initialisation of the \ref nDOF vector. Must be implemented by sub classes
170
171
    virtual void setNDOF() = 0;

Thomas Witkowski's avatar
Thomas Witkowski committed
172
    /// Returns the barycentric coordinates of the i-th basis function.
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
    virtual DimVec<double> *getCoords(int i) const = 0;

    /** \brief
     * Returns a pointer to a const vector with interpolation coefficients of the
     * function f; if indices is a pointer to NULL, the coefficient for all 
     * basis functions are calculated and the i-th entry in the vector is the 
     * coefficient of the i-th basis function; if indices is non NULL, only the 
     * coefficients for a subset of the local basis functions has to be 
     * calculated; n is the number of those basis functions, indices[0], . . . 
     * , indices[n-1] are the local indices of the basis functions where the
     * coefficients have to be calculated, and the i-th entry in the return 
     * vector is then the coefficient of the indices[i]-th basis function; coeff 
     * may be a pointer to a vector which has to be filled 
     * (compare the dof argument of \ref getDOFIndices());
     * such a function usually needs vertex coordinate information; thus, all 
     * routines using this function on the elements need the FILL COORDS flag 
     * during mesh traversal.
     * Must be implemented by sub classes.
     */
    virtual const double* interpol(const ElInfo *el_info, 
				   int n, const int *indices, 
				   AbstractFunction<double, WorldVector<double> > *f,
				   double *coeff) = 0;

197
    /// WorldVector<double> valued interpol function.
198
199
200
201
202
    virtual const WorldVector<double>* interpol(const ElInfo *el_info, int no, 
						const int *b_no,
						AbstractFunction<WorldVector<double>,
						             WorldVector<double> > *f, 
						WorldVector<double> *vec) = 0;
203

204
    /// Returns the i-th local basis function
Thomas Witkowski's avatar
Thomas Witkowski committed
205
206
    inline BasFctType *getPhi(int i) const 
    { 
207
      return (*phi)[i]; 
208
    }
209

210
    /// Returns the gradient of the i-th local basis function
Thomas Witkowski's avatar
Thomas Witkowski committed
211
212
    inline GrdBasFctType *getGrdPhi(int i) const 
    { 
213
      return (*grdPhi)[i]; 
214
    }
215

216
    /// Returns the second derivative of the i-th local basis function
Thomas Witkowski's avatar
Thomas Witkowski committed
217
218
    inline D2BasFctType *getD2Phi(int i) const 
    { 
219
      return (*d2Phi)[i]; 
220
    }
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243

    /** \brief
     * Approximates the L2 scalar products of a given function with all basis 
     * functions by numerical quadrature and adds the corresponding values to a 
     * DOF vector;
     * f is a pointer for the evaluation of the given function in world 
     * coordinates x and returns the value of that function at x; if f is a NULL
     *  pointer, nothing is done;
     * fh is the DOF vector where at the i-th entry the approximation of the L2 
     * scalar product of the given function with the i-th global basis function 
     * of fh->feSpace is stored;
     * quad is the quadrature for the approximation of the integral on each leaf 
     * element of fh->feSpace->mesh; if quad is a NULL pointer, a default 
     * quadrature which is exact of degree 2*fh->feSpace->basFcts->degree-2 is 
     * used.
     * The integrals are approximated by looping over all leaf elements, 
     * computing the approximations to the element contributions and adding 
     * these values to the vector fh by add element vec().
     * The vector fh is not initialized with 0.0; only the new contributions are 
     * added
     */
    virtual void l2ScpFctBas(Quadrature*,
			     AbstractFunction<double, WorldVector<double> >* /*f*/,
244
245
			     DOFVector<double>* /*fh*/)
    {}
246

247
    /// WorldVector<double> valued l2ScpFctBas function
248
249
    virtual void l2ScpFctBas(Quadrature* ,
			     AbstractFunction<WorldVector<double>, WorldVector<double> >* /*f*/,
250
251
			     DOFVector<WorldVector<double> >* /*fh*/) 
    {}
252
253


254
255
    /// Interpolates a DOFIndexed<double> after refinement
    virtual void refineInter(DOFIndexed<double> *, RCNeighbourList*, int)
256
    {}
257

258
259
    /// Interpolates a DOFIndexed<double> after coarsening
    virtual void coarseInter(DOFIndexed<double> *, RCNeighbourList*, int)
260
    {}
261

262
263
    /// Restricts a DOFIndexed<double> after coarsening
    virtual void coarseRestr(DOFIndexed<double> *, RCNeighbourList*, int)
264
    {}
265

266
267
    /// Interpolates a DOFVector<WorldVector<double> > after refinement
    virtual void refineInter(DOFVector<WorldVector<double> >*, RCNeighbourList*, int)
268
    {}
269

270
271
    /// Interpolates a DOFVector<WorldVector<double> > after coarsening
    virtual void coarseInter(DOFVector<WorldVector<double> >*, RCNeighbourList*, int)
272
    {}
273

274
275
    /// Restricts a DOFVector<WorldVector<double> > after coarsening
    virtual void coarseRestr(DOFVector<WorldVector<double> >*, RCNeighbourList*, int)
276
    {}
277

278
    /// Returns local dof indices of the element for the given fe space.
Thomas Witkowski's avatar
Thomas Witkowski committed
279
280
281
    virtual const DegreeOfFreedom *getLocalIndices(const Element *el,
						   const DOFAdmin *admin,
						   DegreeOfFreedom *dofPtr) const
282
283
    {
      return NULL;
284
    }
285

286
287
288
289
290
291
292
293
294
295
296
    inline void getLocalIndices(const Element *el,
				const DOFAdmin *admin,
				std::vector<DegreeOfFreedom> &indices) const
    {
      FUNCNAME("BasisFunction::getLocalIndices()");
      
      TEST_EXIT_DBG(static_cast<int>(indices.size()) >= nBasFcts)
	("Index vector is too small!\n");

      getLocalIndices(el, admin, &(indices[0]));
    }
Thomas Witkowski's avatar
Thomas Witkowski committed
297

Thomas Witkowski's avatar
Thomas Witkowski committed
298
299
300
    virtual void getLocalDofPtrVec(const Element *el, 
				   const DOFAdmin *admin,
				   std::vector<const DegreeOfFreedom*>& vec) const
301
    {}
Thomas Witkowski's avatar
Thomas Witkowski committed
302

303
304
305
306
307
308

    /** \brief
     * Evaluates elements value at barycentric coordinates lambda with local 
     * coefficient vector uh.
     */
    double evalUh(const DimVec<double>& lambda, const double* uh) const;
309

310
311
312
313
314
315
316
    /** \brief
     * Evaluates elements value at barycentric coordinates lambda with local 
     * coefficient vector uh. If val is not NULL the result will be stored 
     * there, otherwise a pointer to a static local variable is returned which 
     * will be overwritten after the next call.
     */
    const WorldVector<double>& evalUh(const DimVec<double>& lambda, 
317
318
				      const WorldVector<double>* uh, 
				      WorldVector<double>* val) const;
319

320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
    /** \brief
     * Evaluates the gradient at barycentric coordinates lambda. Lambda is the
     * Jacobian of the barycentric coordinates. uh is the local coefficient
     * vector. If val is not NULL the result will be stored 
     * there, otherwise a pointer to a static local variable is returned which 
     * will be overwritten after the next call.
     */
    const WorldVector<double>& evalGrdUh(const DimVec<double>& lambda,
					 const DimVec<WorldVector<double> >& Lambda,
					 const double* uh,  WorldVector<double>* val) const;

    /** \brief
     * Evaluates the second derivative at barycentric coordinates lambda. 
     * Lambda is the Jacobian of the barycentric coordinates. uh is the local 
     * coefficient vector. If val is not NULL the result will be stored 
     * there, otherwise a pointer to a static local variable is returned which 
     * will be overwritten after the next call.
     */
    const WorldMatrix<double>& evalD2Uh(const DimVec<double>& lambda,
					const DimVec<WorldVector<double> >& Lambda,
					const double* uh, WorldMatrix<double>* val) const;

  protected:
343
    /// Textual description
344
    std::string name;     
345

346
    /// Number of basisfunctions on one Element                 
347
348
    int nBasFcts;

349
    /// Maximal degree of the basis functions                 
350
351
    int degree;

352
    /// Dimension of the basis functions                  
353
354
    int dim;

355
    /// Dimension of the world.
Thomas Witkowski's avatar
Thomas Witkowski committed
356
357
    int dow;

358
    /// Number of DOFs at the different positions                  
359
360
    DimVec<int> *nDOF;

361
    /// Vector of the local functions
362
    std::vector<BasFctType*> *phi;
363

364
    /// Vector of gradients
365
    std::vector<GrdBasFctType*> *grdPhi;
366

367
    /// Vector of second derivatives
368
    std::vector<D2BasFctType*> *d2Phi;
Thomas Witkowski's avatar
Thomas Witkowski committed
369
370
371
372
373

    /** \brief
     * Is used by function evalGrdUh. To make it thread safe, we need a
     * temporary DimVec for every thread.
     */
374
    std::vector<DimVec<double>* > grdTmpVec1;
Thomas Witkowski's avatar
Thomas Witkowski committed
375
376
377
378
379

    /** \brief
     * Is used by function evalGrdUh. To make it thread safe, we need a
     * temporary DimVec for every thread.
     */
380
    std::vector<DimVec<double>* > grdTmpVec2;
381
382
383
384
385
  };

}

#endif  // AMDIS_BASISFUNCTION_H