SubAssembler.h 6.83 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
#ifndef AMDIS_SUB_ASSEMBLER_H
#define AMDIS_SUB_ASSEMBLER_H

#include <map>
#include <vector>
#include "FixVec.h"
#include "Flag.h"
#include "MemoryManager.h"
#include "OpenMP.h"

namespace AMDiS {

  class Assembler;
  class Quadrature;
  class FastQuadrature;
  class Operator;
  class OperatorTerm;
  class BasisFunction;
  class FiniteElemSpace;
  class ElInfo;
  class ElementMatrix;
  class ElementVector;

  template<typename T> class DOFVectorBase;

  /**
   * \ingroup Assembler
   * 
   * \brief
   * Base class for SecondOrderAssembler, FirstOrderAssembler, 
   * ZeroOrderAssembler. The task of a SubAssembler is to assemble a list of 
   * terms of a special order and add their contributions to a DOFMatrix or a 
   * DOFVector. An Assembler can consist of up to four SubAssemblers: one 
   * SecondOrderAssembler for second order terms, one ZeroOrderAssembler for 
   * terms of order zero, and two FirstOrderAssemblers. One for terms with
   * derivatives of the basis functions connected to to row DOFs and one for 
   * those connected to the column DOFs.
   */
  class SubAssembler
  {
  public:
    MEMORY_MANAGED(SubAssembler);

    /** \brief
     * Creates a SubAssembler belonging to assembler for the terms of given 
     * order of Operator op. If order is equal to one, type spezifies what kind 
     * of FirstOrderType are to assemble. During construction of a SubAssembler
     * the needs and properties of the terms are considered.
     */
    SubAssembler(Operator *op,
		 Assembler *assembler,
		 Quadrature *quadrat,
		 int order, 
		 bool optimized,
		 FirstOrderType type = GRD_PHI);

    /** \brief
     * Destructor
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
60
    virtual ~SubAssembler() {}
61
62
63
64
65
66
67
68
69
70

    /** \brief
     * Calculates the element matrix for elInfo and adds it to mat. Memory
     * for mat must be provided by the caller.
     */
    virtual void calculateElementMatrix(const ElInfo *elInfo,
					ElementMatrix *mat) = 0;

    virtual void calculateElementMatrix(const ElInfo *rowElInfo,
					const ElInfo *colElInfo,
71
72
					const ElInfo *smallElInfo,
					const ElInfo *largeElInfo,
73
74
75
76
77
78
79
80
81
82
83
84
85
86
					ElementMatrix *mat) = 0;

    /** \brief
     * Calculates the element vector for elInfo and adds it to vec. Memory
     * for vec must be provided by the caller.
     */
    virtual void calculateElementVector(const ElInfo *elInfo,
					ElementVector *vec) = 0;

    /** \brief
     * Returns \ref terms
     */
    inline std::vector<OperatorTerm*> *getTerms() { 
      return &terms[omp_get_thread_num()]; 
Thomas Witkowski's avatar
Thomas Witkowski committed
87
    }
88
89
90
91
92
93

    /** \brief
     * Returns \ref quadrature.
     */ 
    inline Quadrature *getQuadrature() {
      return quadrature;
Thomas Witkowski's avatar
Thomas Witkowski committed
94
    }
95
96
97
98
99
100

    /** \brief
     * Sets \ref quadrature to q.
     */
    inline void setQuadrature(Quadrature* q) {
      quadrature = q;
Thomas Witkowski's avatar
Thomas Witkowski committed
101
    }
102
103
104
105
106
107
108
109
110
111
112
113
114
  
    /** \brief
     * Returns a vector with the world coordinates of the quadrature points
     * of \ref quadrature on the element of elInfo. 
     * Used by \ref OperatorTerm::initElement().
     */
    WorldVector<double>* getCoordsAtQPs(const ElInfo* elInfo, 
					Quadrature *quad = NULL);

    /** \brief
     * DOFVector dv evaluated at quadrature points.
     * Used by \ref OperatorTerm::initElement().
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
115
116
117
118
119
120
121
122
123
124
    double* getVectorAtQPs(DOFVectorBase<double>* dv, 
			   const ElInfo* elInfo,
			   Quadrature *quad = NULL);

    /** \brief
     *
     */
    double* getVectorAtQPs(DOFVectorBase<double>* dv, 
			   const ElInfo* smallElInfo,
			   const ElInfo* largeElInfo,
125
126
127
128
129
130
131
132
133
134
			   Quadrature *quad = NULL);

    /** \brief
     * Gradients of DOFVector dv evaluated at quadrature points.
     * Used by \ref OperatorTerm::initElement().
     */
    WorldVector<double>* getGradientsAtQPs(DOFVectorBase<double>* dv, 
					   const ElInfo* elInfo,
					   Quadrature *quad = NULL);

Thomas Witkowski's avatar
Thomas Witkowski committed
135
136
137
138
139
140
141
142
    /** \brief
     *
     */
    WorldVector<double>* getGradientsAtQPs(DOFVectorBase<double>* dv, 
					   const ElInfo* smallElInfo,
					   const ElInfo* largeElInfo,
					   Quadrature *quad = NULL);

143
144
145
146
147
    /** \brief
     * Called once for each ElInfo when \ref calculateElementMatrix() or
     * \ref calculateElementVector() is called for the first time for this
     * Element.
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
148
149
    virtual void initElement(const ElInfo *smallElInfo,
			     const ElInfo *largeElInfo = NULL,
150
151
152
153
154
155
156
			     Quadrature *quad = NULL);

    /** \brief
     * Returns \ref psiFast.
     */
    const FastQuadrature *getPsiFast() const { 
      return psiFast; 
Thomas Witkowski's avatar
Thomas Witkowski committed
157
    }
158
159
160
161
162
163

    /** \brief
     * Returns \ref phiFast.
     */
    const FastQuadrature *getPhiFast() const { 
      return phiFast; 
Thomas Witkowski's avatar
Thomas Witkowski committed
164
    }
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281

  protected:
    /** \brief
     * Updates \ref psiFast and \ref phiFast.
     */
    FastQuadrature *updateFastQuadrature(FastQuadrature *quadFast,
					 const BasisFunction *psi,
					 Flag updateFlag);
  
  protected:
    /** \brief
     * Problem dimension
     */
    int dim;

    /** \brief
     * Number of rows of the element matrix and length of the element
     * vector. Is equal to the number of row basis functions
     */
    int nRow;

    /** \brief
     * Number of columns of the element matrix. Is equal to the number
     * of column basis functions
     */
    int nCol;

    /** \brief
     * Used for \ref getVectorAtQPs().
     */
    class ValuesAtQPs {
    public:
      Vector<double> values;
      bool valid;
    };

    /** \brief
     * Used for \ref getGradientsAtQPs().
     */
    class GradientsAtQPs {
    public:
      Vector<WorldVector<double> > values;
      bool valid;
    };

    /** \brief
     * Used for \ref getVectorAtQPs().
     */
    std::map<const DOFVectorBase<double>*, ValuesAtQPs*> valuesAtQPs;

    /** \brief
     * Used for \ref getGradientsAtQPs().
     */
    std::map<const DOFVectorBase<double>*, GradientsAtQPs*> gradientsAtQPs;

    /** \brief
     * Set and updated by \ref initElement() for each ElInfo. 
     * coordsAtQPs[i] points to the coordinates of the i-th quadrature point.
     */
    WorldVector<double> *coordsAtQPs;

    /** \brief
     * Used for \ref getCoordsAtQPs().
     */
    bool coordsValid;

    /** \brief
     * Used for \ref getCoordsAtQP(). Stores the number of allocated WorldVectors.
     */
    int coordsNumAllocated;

    /** \brief
     * Needed Quadrature. Constructed in the constructor of SubAssembler
     */
    Quadrature *quadrature;

    /** \brief
     * FastQuadrature for row basis functions
     */
    FastQuadrature *psiFast;

    /** \brief
     * FastQuadrature for column basis functions
     */
    FastQuadrature *phiFast;

    /** \brief
     * Corresponding Assembler
     */
    Assembler* owner;

    /** \brief
     * Flag that specifies whether the element matrix is symmetric.
     */
    bool symmetric;

    /** \brief
     * List of all terms with a contribution to this SubAssembler
     */
    std::vector< std::vector<OperatorTerm*> > terms;

    /** \brief
     *
     */
    bool opt;

    /** \brief
     *
     */
    bool firstCall;

    friend class Assembler;
  };

}

#endif // AMDSIS_SUB_ASSEMBLER_H