DOFVector.h 20.6 KB
Newer Older
1
2
3
4
5
6
// ============================================================================
// ==                                                                        ==
// == AMDiS - Adaptive multidimensional simulations                          ==
// ==                                                                        ==
// ============================================================================
// ==                                                                        ==
7
// ==  TU Dresden                                                            ==
8
// ==                                                                        ==
9
10
11
// ==  Institut fr 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
// ==                                                                        ==
// ============================================================================

/** \file DOFVector.h */

#ifndef AMDIS_DOFVECTOR_H 
#define AMDIS_DOFVECTOR_H 
 
25
26
27
28
#include <vector> 
#include <memory> 
#include <list> 

29
#include "AMDiS_fwd.h"
30
31
32
33
34
35
36
37
38
39
40
#include "FixVec.h"
#include "Global.h" 
#include "Flag.h" 
#include "RCNeighbourList.h" 
#include "DOFIterator.h"
#include "DOFIndexed.h"
#include "DOFContainer.h"
#include "Boundary.h"
#include "CreatorInterface.h"
#include "Serializable.h"
#include "DOFMatrix.h" 
Thomas Witkowski's avatar
Thomas Witkowski committed
41
#include "BasisFunction.h"
Thomas Witkowski's avatar
Thomas Witkowski committed
42
#include "FiniteElemSpace.h"
43

44
45
46
47
48
49
namespace AMDiS {

  template<typename T> 
  class DOFVectorBase : public DOFIndexed<T>
  {
  public:
Thomas Witkowski's avatar
Thomas Witkowski committed
50

51
52
    DOFVectorBase() 
      : feSpace(NULL),
53
	elementVector(3),
Thomas Witkowski's avatar
Thomas Witkowski committed
54
55
        boundaryManager(NULL),
        nBasFcts(0)
56
    {}
Thomas Witkowski's avatar
Thomas Witkowski committed
57
    
58
    DOFVectorBase(const FiniteElemSpace *f, std::string n);
59

Thomas Witkowski's avatar
Thomas Witkowski committed
60
    virtual ~DOFVectorBase();
61

Thomas Witkowski's avatar
Thomas Witkowski committed
62
63
64
65
66
    /** \brief 
     * For the given element, this function returns an array of all DOFs of this
     * DOFVector that are defined on this element.
     */
    virtual const T *getLocalVector(const Element *el, T* localVec) const;
67

Thomas Witkowski's avatar
Thomas Witkowski committed
68
69
70
71
    /** \brief
     * Evaluates the DOF vector at a set of quadrature points defined on the 
     * given element.
     */
72
73
    const T *getVecAtQPs(const ElInfo *elInfo, 
			 const Quadrature *quad,
74
			 const FastQuadrature *quadFast,
75
			 T *vecAtQPs) const;
76

Thomas Witkowski's avatar
Thomas Witkowski committed
77
78
79
80
81
82
    const T *getVecAtQPs(const ElInfo *smallElInfo, 
			 const ElInfo *largeElInfo,
			 const Quadrature *quad,
			 const FastQuadrature *quadFast,
			 T *vecAtQPs) const;

83
84
    const WorldVector<T> *getGrdAtQPs(const ElInfo *elInfo, 
				      const Quadrature *quad,
85
				      const FastQuadrature *quadFast,
86
				      WorldVector<T> *grdAtQPs) const;
87

Thomas Witkowski's avatar
Thomas Witkowski committed
88
89
90
91
92
93
    const WorldVector<T> *getGrdAtQPs(const ElInfo *smallElInfo, 
				      const ElInfo *largeElInfo,
				      const Quadrature *quad,
				      const FastQuadrature *quadFast,
				      WorldVector<T> *grdAtQPs) const;

94
95
    const WorldMatrix<T> *getD2AtQPs(const ElInfo *elInfo, 
				     const Quadrature *quad,
96
				     const FastQuadrature *quadFast,
97
				     WorldMatrix<T> *d2AtQPs) const;
98

Thomas Witkowski's avatar
Thomas Witkowski committed
99
    /// Returns the FE space the DOF vector is defined on.
Thomas Witkowski's avatar
Thomas Witkowski committed
100
101
    inline const FiniteElemSpace* getFESpace() const 
    {
102
      return feSpace;
103
    }
104

Thomas Witkowski's avatar
Thomas Witkowski committed
105
106
107
108
109
110
111
112
113
114
115
116
117
    /** \brief
     * Assembles the element vector for the given ellement and adds the
     * element matrix to the current DOF vector.
     */ 
    void assemble(T factor, ElInfo *elInfo,			    
		  const BoundaryType *bound, 
		  Operator *op = NULL);

    void assemble2(T factor, 
		   ElInfo *mainElInfo, ElInfo *auxElInfo,
		   ElInfo *smallElInfo, ElInfo *largeElInfo,
		   const BoundaryType *bound, 
		   Operator *op = NULL);
118
119

    void addElementVector(T sign,
120
			  const ElementVector& elVec, 
121
			  const BoundaryType *bound,
122
			  ElInfo *elInfo,
123
			  bool add = true); 
Thomas Witkowski's avatar
Thomas Witkowski committed
124
125
126
127
128
129
130

    /* \brief
     * That function must be called after the matrix assembling has been finished. 
     * This makes it possible to start some cleanup or matrix data compressing 
     * procedures.
     */
    void finishAssembling();
131
132
133
134
135
136
137
138
 
    inline void addOperator(Operator* op, 
			    double *factor = NULL,
			    double *estFactor = NULL) 
    {
      operators.push_back(op);
      operatorFactor.push_back(factor);
      operatorEstFactor.push_back(estFactor);
139
    }
140

Thomas Witkowski's avatar
Thomas Witkowski committed
141
142
    inline std::vector<double*>::iterator getOperatorFactorBegin() 
    {
143
      return operatorFactor.begin();
144
    }
145

Thomas Witkowski's avatar
Thomas Witkowski committed
146
147
    inline std::vector<double*>::iterator getOperatorFactorEnd() 
    {
148
      return operatorFactor.end();
149
    }
150

Thomas Witkowski's avatar
Thomas Witkowski committed
151
152
    inline std::vector<double*>::iterator getOperatorEstFactorBegin() 
    {
153
      return operatorEstFactor.begin();
154
    }
155

Thomas Witkowski's avatar
Thomas Witkowski committed
156
157
    inline std::vector<double*>::iterator getOperatorEstFactorEnd() 
    {
158
      return operatorEstFactor.end();
159
    }
160
161


Thomas Witkowski's avatar
Thomas Witkowski committed
162
163
    inline std::vector<Operator*>::iterator getOperatorsBegin() 
    {
164
      return operators.begin();
165
    }
166

Thomas Witkowski's avatar
Thomas Witkowski committed
167
168
    inline std::vector<Operator*>::iterator getOperatorsEnd() 
    {
169
      return operators.end();
170
    }
171
172
173
174
175
176
177
178
179
180
181

    Flag getAssembleFlag();

    /** \brief
     * Evaluates \f[ u_h(x(\lambda)) = \sum_{i=0}^{m-1} vec[ind[i]] * 
     * \varphi^i(\lambda) \f] where \f$ \varphi^i \f$ is the i-th basis function,
     * \f$ x(\lambda) \f$ are the world coordinates of lambda and
     * \f$ m \f$ is the number of basis functions
     */
    T evalUh(const DimVec<double>& lambda, DegreeOfFreedom* ind);

Thomas Witkowski's avatar
Thomas Witkowski committed
182
183
    inline std::vector<Operator*>& getOperators() 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
184
      return operators; 
185
    }
Thomas Witkowski's avatar
Thomas Witkowski committed
186

Thomas Witkowski's avatar
Thomas Witkowski committed
187
188
    inline std::vector<double*>& getOperatorFactor() 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
189
      return operatorFactor; 
190
    }
191

Thomas Witkowski's avatar
Thomas Witkowski committed
192
193
    inline std::vector<double*>& getOperatorEstFactor() 
    { 
194
      return operatorEstFactor; 
195
    }
196

Thomas Witkowski's avatar
Thomas Witkowski committed
197
    /// Returns \ref name
Thomas Witkowski's avatar
Thomas Witkowski committed
198
    inline std::string getName() const 
Thomas Witkowski's avatar
Thomas Witkowski committed
199
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
200
      return name; 
201
    } 
202

Thomas Witkowski's avatar
Thomas Witkowski committed
203
204
    inline BoundaryManager* getBoundaryManager() const 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
205
      return boundaryManager; 
206
    }
207

Thomas Witkowski's avatar
Thomas Witkowski committed
208
209
    inline void setBoundaryManager(BoundaryManager *bm) 
    {
210
      boundaryManager = bm;
211
    }
212

Thomas Witkowski's avatar
Thomas Witkowski committed
213
#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
214
    inline void setRankDofs(std::map<DegreeOfFreedom, bool> &dofmap)
Thomas Witkowski's avatar
Thomas Witkowski committed
215
    {
216
      //      rankDofs = &dofmap;
Thomas Witkowski's avatar
Thomas Witkowski committed
217
218
219
220
    }

    inline bool isRankDof(DegreeOfFreedom dof)
    {
221
222
      TEST_EXIT_DBG(rankDofs)("No rank dofs set!\n");
      return (*rankDofs)[dof];
Thomas Witkowski's avatar
Thomas Witkowski committed
223
224
225
    }
#endif

226
  protected:
Thomas Witkowski's avatar
Thomas Witkowski committed
227
    ///
228
229
    const FiniteElemSpace *feSpace;

Thomas Witkowski's avatar
Thomas Witkowski committed
230
    ///
231
    std::string name;
232

Thomas Witkowski's avatar
Thomas Witkowski committed
233
    ///
234
    ElementVector elementVector;
235

Thomas Witkowski's avatar
Thomas Witkowski committed
236
    ///
237
    std::vector<Operator*> operators;
238

Thomas Witkowski's avatar
Thomas Witkowski committed
239
    ///
240
    std::vector<double*> operatorFactor;
241

Thomas Witkowski's avatar
Thomas Witkowski committed
242
    ///
243
    std::vector<double*> operatorEstFactor;
244

Thomas Witkowski's avatar
Thomas Witkowski committed
245
    ///
246
247
    BoundaryManager *boundaryManager;

Thomas Witkowski's avatar
Thomas Witkowski committed
248
    /// Number of basis functions of the used finite element space.
Thomas Witkowski's avatar
Thomas Witkowski committed
249
    int nBasFcts;
250

Thomas Witkowski's avatar
Thomas Witkowski committed
251
    /// Are used to store temporary local dofs of an element. Thread safe.
252
253
    std::vector<DegreeOfFreedom*> localIndices;

Thomas Witkowski's avatar
Thomas Witkowski committed
254
    /// Are used to store temporary local values of an element. Thread safe.
Thomas Witkowski's avatar
Thomas Witkowski committed
255
256
    std::vector<T*> localVectors;

Thomas Witkowski's avatar
Thomas Witkowski committed
257
    /// Temporarly used in \ref getGrdAtQPs. Thread safe.
Thomas Witkowski's avatar
Thomas Witkowski committed
258
259
    std::vector<DimVec<double>*> grdTmp;

Thomas Witkowski's avatar
Thomas Witkowski committed
260
    /// Temporarly used in \ref getGrdAtQPs. Thread safe.
261
262
    std::vector<DimVec<double>*> grdPhis;

Thomas Witkowski's avatar
Thomas Witkowski committed
263
    /// Temporarly used in \ref getD2AtQPs. Thread safe.
264
    std::vector<DimMat<double>*> D2Phis;
Thomas Witkowski's avatar
Thomas Witkowski committed
265

Thomas Witkowski's avatar
Thomas Witkowski committed
266
    /// Dimension of the mesh this DOFVectorBase belongs to
Thomas Witkowski's avatar
Thomas Witkowski committed
267
    int dim;
Thomas Witkowski's avatar
Thomas Witkowski committed
268
269

#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
270
  public:
271
    std::map<DegreeOfFreedom, bool> *rankDofs;
Thomas Witkowski's avatar
Thomas Witkowski committed
272
#endif
Thomas Witkowski's avatar
Thomas Witkowski committed
273
  };
274

Thomas Witkowski's avatar
Thomas Witkowski committed
275
  /// Specifies which operation should be done after coarsening
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
  typedef enum{
    NO_OPERATION = 0,   
    COARSE_RESTRICT = 1,
    COARSE_INTERPOL = 2 
  } CoarsenOperation;
 

  /** \ingroup DOFAdministration 
   * \brief
   * The DOFs described above are just integers that can be used as indices into 
   * vectors and matrices. During refinement and coarsening of the mesh, the 
   * number of used DOFs, the meaning of one integer index, and even the total 
   * range of DOFs change. To be able to handle these changes automatically for 
   * all vectors, which are indexed by the DOFs, special data structures are 
   * used which contain such vector data. Lists of these structures are kept in 
   * DOFAdmin, so that all vectors in the lists can be resized together with the
   * range of DOFs. During refinement and coarsening of elements, values can be
   * interpolated automatically to new DOFs, and restricted from old DOFs.
   */
  template<typename T> 
  class DOFVector : public DOFVectorBase<T>, public Serializable
  {  
  public:
    /** \ingroup DOFAdministration
     * \brief
     * Enables the access of DOFVector<T>::Iterator. Alias for DOFIterator<T>
     */
    class Iterator : public DOFIterator<T> {
    public:
      Iterator(DOFIndexed<T> *c, DOFIteratorType type)
	: DOFIterator<T>(c, type)
307
      {}
308
309
310

      Iterator(DOFAdmin *admin, DOFIndexed<T> *c, DOFIteratorType type)
	: DOFIterator<T>(admin, c, type)
311
      {}
312
313
314
315
    };

    class Creator : public CreatorInterface<DOFVector<T> > {
    public:
316
317
318
      Creator(FiniteElemSpace *feSpace_) 
        : feSpace(feSpace_) 
      {}
319

Thomas Witkowski's avatar
Thomas Witkowski committed
320
321
      DOFVector<T> *create() 
      { 
322
	return new DOFVector<T>(feSpace, ""); 
323
      }
324

Thomas Witkowski's avatar
Thomas Witkowski committed
325
326
      void free(DOFVector<T> *vec) 
      { 
327
	delete vec; 
328
      }
329
330
331
332
333
334

    private:
      FiniteElemSpace *feSpace;
    };

  public:
Thomas Witkowski's avatar
Thomas Witkowski committed
335
    /// Empty constructor. No initialization!
336
337
338
    DOFVector() 
      : DOFVectorBase<T>(), 
	coarsenOperation(NO_OPERATION)
339
    {}
340

Thomas Witkowski's avatar
Thomas Witkowski committed
341
    /// Constructs a DOFVector with name n belonging to FiniteElemSpace f
342
    DOFVector(const FiniteElemSpace* f, std::string n); 
343

Thomas Witkowski's avatar
Thomas Witkowski committed
344
    /// Initialization.
345
    void init(const FiniteElemSpace* f, std::string n);
346

Thomas Witkowski's avatar
Thomas Witkowski committed
347
    /// Copy Constructor
Thomas Witkowski's avatar
Thomas Witkowski committed
348
349
    DOFVector(const DOFVector& rhs) 
    {
Thomas Witkowski's avatar
Thomas Witkowski committed
350
      *this = rhs;   
351
352
353
      this->name = rhs.name + "copy";
      if (this->feSpace && this->feSpace->getAdmin())
	(dynamic_cast<DOFAdmin*>(this->feSpace->getAdmin()))->addDOFIndexed(this);
354
    }
355

Thomas Witkowski's avatar
Thomas Witkowski committed
356
    /// Destructor
357
358
    virtual ~DOFVector();

Thomas Witkowski's avatar
Thomas Witkowski committed
359
    /// Returns iterator to the begin of \ref vec
Thomas Witkowski's avatar
Thomas Witkowski committed
360
361
    typename std::vector<T>::iterator begin() 
    { 
362
      return vec.begin(); 
363
    }
364

Thomas Witkowski's avatar
Thomas Witkowski committed
365
    /// Returns iterator to the end of \ref vec
Thomas Witkowski's avatar
Thomas Witkowski committed
366
367
    typename std::vector<T>::iterator end() 
    { 
368
      return vec.end(); 
369
    }
370
371
372
373
374
375
  
    /** \brief
     * Used by DOFAdmin to compress this DOFVector. Implementation of
     * DOFIndexedBase::compress()
     */
    virtual void compressDOFIndexed(int first, int last,
376
				    std::vector<DegreeOfFreedom> &newDof);
Thomas Witkowski's avatar
Thomas Witkowski committed
377
378

    /// Sets \ref coarsenOperation to op
Thomas Witkowski's avatar
Thomas Witkowski committed
379
380
    inline void setCoarsenOperation(CoarsenOperation op) 
    { 
381
      coarsenOperation = op; 
382
    }
383

Thomas Witkowski's avatar
Thomas Witkowski committed
384
    /// Returns \ref coarsenOperation
Thomas Witkowski's avatar
Thomas Witkowski committed
385
386
    inline CoarsenOperation getCoarsenOperation() 
    { 
387
      return coarsenOperation; 
388
    }
389

Thomas Witkowski's avatar
Thomas Witkowski committed
390
    /// Restriction after coarsening. Implemented for DOFVector<double>
391
    inline void coarseRestrict(RCNeighbourList&, int) {}
392

Thomas Witkowski's avatar
Thomas Witkowski committed
393
    /// Interpolation after refinement.
394
    inline void refineInterpol(RCNeighbourList&, int) {}
395

Thomas Witkowski's avatar
Thomas Witkowski committed
396
    /// Returns \ref vec
Thomas Witkowski's avatar
Thomas Witkowski committed
397
398
    std::vector<T>& getVector() 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
399
      return vec;
400
    }
401

Thomas Witkowski's avatar
Thomas Witkowski committed
402
    /// Returns size of \ref vec
Thomas Witkowski's avatar
Thomas Witkowski committed
403
404
    inline int getSize() const 
    { 
405
      return vec.size();
406
    } 
407

Thomas Witkowski's avatar
Thomas Witkowski committed
408
    /// Returns used size of the vector.
Thomas Witkowski's avatar
Thomas Witkowski committed
409
410
    inline int getUsedSize() const 
    { 
411
      return this->feSpace->getAdmin()->getUsedSize(); 
412
    }
413

Thomas Witkowski's avatar
Thomas Witkowski committed
414
    /// Resizes \ref vec to n
Thomas Witkowski's avatar
Thomas Witkowski committed
415
416
    inline void resize(int n) 
    { 
417
418
      FUNCNAME("DOFVector<T>::resize()");
      TEST_EXIT_DBG((n >= 0))("Can't resize DOFVector to negative size\n"); 
419
      vec.resize(n);
420
    } 
421

Thomas Witkowski's avatar
Thomas Witkowski committed
422
    /// Resizes \ref vec to n and inits new values with init
Thomas Witkowski's avatar
Thomas Witkowski committed
423
424
    inline void resize(int n, T init) 
    { 
425
426
      FUNCNAME("DOFVector<T>::resize()");
      TEST_EXIT_DBG((n >= 0))("Can't resize DOFVector to negative size\n"); 
427
      vec.resize(n, init);
428
    } 
429

Thomas Witkowski's avatar
Thomas Witkowski committed
430
    /// Returns \ref vec[i]
Thomas Witkowski's avatar
Thomas Witkowski committed
431
432
    inline const T& operator[](DegreeOfFreedom i) const 
    {
433
      FUNCNAME("DOFVector<T>::operator[]");
434
435
      TEST_EXIT_DBG((i>= 0) && (i < static_cast<int>(vec.size())))
	("Illegal vector index %d.\n",i);
436
      return vec[i];
437
    } 
438

Thomas Witkowski's avatar
Thomas Witkowski committed
439
    /// Returns \ref vec[i]
Thomas Witkowski's avatar
Thomas Witkowski committed
440
441
    inline T& operator[](DegreeOfFreedom i) 
    {
442
      FUNCNAME("DOFVector<T>::operator[]");
443
444
      TEST_EXIT_DBG((i >= 0) && (i < static_cast<int>(vec.size())))
	("Illegal vector index %d.\n",i);
445
      return vec[i];
446
    }
Thomas Witkowski's avatar
Thomas Witkowski committed
447
 
Thomas Witkowski's avatar
Thomas Witkowski committed
448
    /// Calculates Integral of this DOFVector
449
450
    double Int(Quadrature* q = NULL) const;

Thomas Witkowski's avatar
Thomas Witkowski committed
451
    /// Calculates L1 norm of this DOFVector
452
    double L1Norm(Quadrature* q = NULL) const;
Thomas Witkowski's avatar
Thomas Witkowski committed
453
 
Thomas Witkowski's avatar
Thomas Witkowski committed
454
    /// Calculates L2 norm of this DOFVector
Thomas Witkowski's avatar
Thomas Witkowski committed
455
456
    inline double L2Norm(Quadrature* q = NULL) const 
    {
457
      return sqrt(L2NormSquare());
458
    }
459

Thomas Witkowski's avatar
Thomas Witkowski committed
460
    /// Calculates square of L2 norm of this DOFVector
461
462
    double L2NormSquare(Quadrature* q = NULL) const;

Thomas Witkowski's avatar
Thomas Witkowski committed
463
    /// Calculates H1 norm of this DOFVector
Thomas Witkowski's avatar
Thomas Witkowski committed
464
465
    inline double H1Norm(Quadrature* q = NULL) const 
    {
466
      return sqrt(H1NormSquare());
467
    }
468

Thomas Witkowski's avatar
Thomas Witkowski committed
469
    /// Calculates square of H1 norm of this DOFVector
470
471
    double H1NormSquare(Quadrature* q = NULL) const;  

Thomas Witkowski's avatar
Thomas Witkowski committed
472
    /// Calculates euclidian norm of this DOFVector
473
474
    double nrm2() const; 

Thomas Witkowski's avatar
Thomas Witkowski committed
475
    /// Returns square of the euclidian norm.
476
477
    double squareNrm2() const;

Thomas Witkowski's avatar
Thomas Witkowski committed
478
    /// Calculates l2 norm of this DOFVector
Thomas Witkowski's avatar
Thomas Witkowski committed
479
480
    inline double l2norm() const 
    { 
481
      return nrm2();
482
    }
483

Thomas Witkowski's avatar
Thomas Witkowski committed
484
    /// Calculates the absolute sum of this DOFVector
485
486
    T asum() const; 

Thomas Witkowski's avatar
Thomas Witkowski committed
487
    /// Calculates the l1 norm of this DOFVector
Thomas Witkowski's avatar
Thomas Witkowski committed
488
489
    inline double l1norm() const 
    { 
490
      return asum();
491
    } 
492

Thomas Witkowski's avatar
Thomas Witkowski committed
493
    /// Calculates doublewell of this DOFVector
494
    double DoubleWell(Quadrature* q = NULL) const;
Thomas Witkowski's avatar
Thomas Witkowski committed
495
 
Thomas Witkowski's avatar
Thomas Witkowski committed
496
    /// Calculates the sum of this DOFVector
497
    T sum() const; 
Thomas Witkowski's avatar
Thomas Witkowski committed
498
 
Thomas Witkowski's avatar
Thomas Witkowski committed
499
    /// Sets \ref vec[i] = val, i=0 , ... , size
500
501
    void set(T val); 

Thomas Witkowski's avatar
Thomas Witkowski committed
502
    /// Assignment operator for setting complete vector to a certain value d
Thomas Witkowski's avatar
Thomas Witkowski committed
503
504
    inline DOFVector<T>& operator=(T d) 
    {
505
506
      set(d); 
      return *this;
507
    } 
508

Thomas Witkowski's avatar
Thomas Witkowski committed
509
    /// Assignment operator between two vectors
510
    DOFVector<T>& operator=(const DOFVector<T>&); 
511

Thomas Witkowski's avatar
Thomas Witkowski committed
512
    /// vec[i] = v.vec[i]
513
514
    void copy(const DOFVector<T>& v); 

Thomas Witkowski's avatar
Thomas Witkowski committed
515
    /// Returns minimum of DOFVector
516
517
    T min() const; 

Thomas Witkowski's avatar
Thomas Witkowski committed
518
    /// Returns maximum of DOFVector
Thomas Witkowski's avatar
Thomas Witkowski committed
519
520
    T max() const;

521
    /// Returns absolute maximum of DOFVector
Thomas Witkowski's avatar
Thomas Witkowski committed
522
    T absMax() const;
523

524
    /// Used by interpol while mesh traversal
525
    void interpolFct(ElInfo* elinfo);
526

527
    /// Prints \ref vec to stdout
528
529
    void print() const; 

530
    ///
531
532
    int calcMemoryUsage() const;

533
534
535
536
537
538
539
540
    /** \brief
     * Computes the coefficients of the interpolant of the function fct and
     * stores these in the DOFVector
     */
    void interpol(AbstractFunction<T, WorldVector<double> > *fct);

    void interpol(DOFVector<T> *v, double factor);

541
    /// Writes the data of the DOFVector to an output stream.
Thomas Witkowski's avatar
Thomas Witkowski committed
542
543
    void serialize(std::ostream &out) 
    {
544
545
546
      unsigned int size = vec.size();
      out.write(reinterpret_cast<const char*>(&size), sizeof(unsigned int));
      out.write(reinterpret_cast<const char*>(&(vec[0])), size * sizeof(T));
547
    }
548

549
    /// Reads data of an DOFVector from an input stream.
Thomas Witkowski's avatar
Thomas Witkowski committed
550
551
    void deserialize(std::istream &in) 
    {
552
553
554
555
      unsigned int size;
      in.read(reinterpret_cast<char*>(&size), sizeof(unsigned int));
      vec.resize(size);
      in.read(reinterpret_cast<char*>(&(vec[0])), size * sizeof(T));
556
    }
557
558
559
560
561

    DOFVector<WorldVector<T> > *getGradient(DOFVector<WorldVector<T> >*) const;

    WorldVector<DOFVector<T>*> *getGradient(WorldVector<DOFVector<T>*> *grad) const;

562
    DOFVector<WorldVector<T> >* getRecoveryGradient(DOFVector<WorldVector<T> >*) const;
563
564

  protected: 
565
    /// Data container
566
    std::vector<T> vec; 
567
 
568
    /// Specifies what operation should be performed after coarsening
569
570
    CoarsenOperation coarsenOperation;

571
    /// Used by \ref interpol
572
573
    AbstractFunction<T, WorldVector<double> > *interFct;

574
    /// Used for mesh traversal
575
576
577
    static DOFVector<T> *traverseVector;
  }; 

Thomas Witkowski's avatar
Thomas Witkowski committed
578

579
580
581
582
583
584
  template<>
  void DOFVector<double>::refineInterpol(RCNeighbourList&, int); 

  template<>
  void DOFVector<double>::coarseRestrict(RCNeighbourList&, int); 

Thomas Witkowski's avatar
Thomas Witkowski committed
585
586
  inline double min(const DOFVector<double>& v) 
  {
Thomas Witkowski's avatar
Thomas Witkowski committed
587
    return v.min();
588
  } 
Thomas Witkowski's avatar
Thomas Witkowski committed
589

Thomas Witkowski's avatar
Thomas Witkowski committed
590
591
  inline double max(const DOFVector<double>& v) 
  {
Thomas Witkowski's avatar
Thomas Witkowski committed
592
    return v.max();
593
  }
594
595
596
597
598
599
600
601
602
603
604
605
606
607


  /** \ingroup DOFAdministration
   * \brief
   * A DOFVector that stores DOF indices.
   */
  class DOFVectorDOF : public DOFVector<DegreeOfFreedom>,
		       public DOFContainer
  {
  public:  
    /** \brief
     * Calls constructor of DOFVector<DegreeOfFreedom> and registers itself
     * as DOFContainer at DOFAdmin
     */
608
    DOFVectorDOF(const FiniteElemSpace* feSpace_, std::string name_)
609
610
611
      : DOFVector<DegreeOfFreedom>(feSpace_, name_)
    {
      feSpace->getAdmin()->addDOFContainer(this);
612
    }
613
  
614
    /// Deregisters itself at DOFAdmin.
Thomas Witkowski's avatar
Thomas Witkowski committed
615
616
    ~DOFVectorDOF() 
    {
617
      feSpace->getAdmin()->removeDOFContainer(this);
618
    }
619
620
621
622
623

    /** \brief
     * Implements DOFContainer::operator[]() by calling 
     * DOFVEctor<DegreeOfFreedom>::operator[]()
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
624
625
    DegreeOfFreedom& operator[](DegreeOfFreedom i) 
    {
626
      return DOFVector<DegreeOfFreedom>::operator[](i);
627
    }
628

Thomas Witkowski's avatar
Thomas Witkowski committed
629
630
    const DegreeOfFreedom& operator[](DegreeOfFreedom i) const 
    {
631
632
633
      return DOFVector<DegreeOfFreedom>::operator[](i);
    }

634
    /// Implements DOFIndexedBase::getSize()
Thomas Witkowski's avatar
Thomas Witkowski committed
635
636
    int getSize() const 
    {
637
      return DOFVector<DegreeOfFreedom>::getSize();
638
    }
639

640
    /// Implements DOFIndexedBase::resize()
Thomas Witkowski's avatar
Thomas Witkowski committed
641
642
    void resize(int size) 
    {
643
644
645
646
647
648
649
650
651
      DOFVector<DegreeOfFreedom>::resize(size);
    }

    void freeDOFContent(DegreeOfFreedom dof);

  protected:
    DOFVectorDOF();
  };

Thomas Witkowski's avatar
Thomas Witkowski committed
652

653
  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
654
655
  double norm(DOFVector<T> *vec) 
  {
656
    return vec->nrm2();
657
  }
658
659

  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
660
661
  double L2Norm(DOFVector<T> *vec) 
  {
662
    return vec->L2Norm();
663
  }
664
665

  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
666
667
  double H1Norm(DOFVector<T> *vec) 
  {
668
    return vec->H1Norm();
669
  }
670
671

  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
672
673
  void print(DOFVector<T> *vec) 
{
674
    vec->print();
675
  }
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721

  // point wise multiplication
  template<typename T>
  const DOFVector<T>& operator*=(DOFVector<T>& x, const DOFVector<T>& y);

  // multiplication with scalar
  template<typename T>
  const DOFVector<T>& operator*=(DOFVector<T>& x, T scal);

  // scalar product
  template<typename T>
  T operator*(DOFVector<T>& x, DOFVector<T>& y);

  // addition
  template<typename T>
  const DOFVector<T>& operator+=(DOFVector<T>& x, const DOFVector<T>& y);

  // subtraction
  template<typename T>
  const DOFVector<T>& operator-=(DOFVector<T>& x, const DOFVector<T>& y);

  template<typename T>
  const DOFVector<T>& operator*(const DOFVector<T>& v, double d);

  template<typename T>
  const DOFVector<T>& operator*(double d, const DOFVector<T>& v);

  template<typename T>
  const DOFVector<T>& operator+(const DOFVector<T>&v1 , const DOFVector<T>& v2);

  // y = a*x + y
  template<typename T>
  void axpy(double a,const DOFVector<T>& x, DOFVector<T>& y);

  // matrix vector product
  template<typename T>
  void mv(MatrixTranspose transpose, 
	  const DOFMatrix &a, 
	  const DOFVector<T> &x,
	  DOFVector<T> &result,
	  bool add = false); 

  template<typename T>
  void xpay(double a,const DOFVector<T>& x,DOFVector<T>& y);

  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
722
723
  inline void scal(T a, DOFVector<T>& y) 
  {
724
725
    y *= a;
  }
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741

  template<typename T>
  inline const DOFVector<T>& mult(double scal,
				  const DOFVector<T>& v,
				  DOFVector<T>& result);

  template<typename T>
  inline const DOFVector<T>& add(const DOFVector<T>& v,
				 double scal,
				 DOFVector<T>& result);

  template<typename T>
  inline const DOFVector<T>& add(const DOFVector<T>& v1,
				 const DOFVector<T>& v2,
				 DOFVector<T>& result);

742
743
744
745
746
747
748
749
750
  template<typename T>
  inline const DOFVector<T>& mod(const DOFVector<T>& v,
				 DOFVector<T>& result);

  template<typename T>
  inline const DOFVector<T>& Tanh(const DOFVector<T>& v,
				  DOFVector<T>& result);


751
  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
752
753
  inline void set(DOFVector<T>& vec, T d) 
  {
754
    vec.set(d);
755
  }
756
757

  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
758
759
  inline void setValue(DOFVector<T>& vec, T d) 
  {
760
    vec.set(d);
761
  }
762
763

  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
764
765
  inline int size(DOFVector<T> *vec) 
  {
766
    return vec->getUsedSize();
767
  }
768

Thomas Witkowski's avatar
Thomas Witkowski committed
769
770
771
772
773
774
775
776
777
778
  template<typename T>
  inline void checkFeSpace(const FiniteElemSpace* feSpace, const std::vector<T>& vec)
  {
    TEST_EXIT_DBG(feSpace)("feSpace is NULL\n");
    TEST_EXIT_DBG(feSpace->getAdmin())("admin is NULL\n");
    TEST_EXIT_DBG(static_cast<int>(vec.size()) >= feSpace->getAdmin()->getUsedSize())
      ("size = %d too small: admin->sizeUsed = %d\n", vec.size(),
       feSpace->getAdmin()->getUsedSize());
  }

779
780
781
782
783
784
785
  WorldVector<DOFVector<double>*> *transform(DOFVector<WorldVector<double> > *vec,
					     WorldVector<DOFVector<double>*> *result);
}

#include "DOFVector.hh"

#endif  // !_DOFVECTOR_H_