DOFVector.h 21.7 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
57
58
59
60
    {
#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
      rankDofs = NULL;
#endif
    }
Thomas Witkowski's avatar
Thomas Witkowski committed
61
    
62
    DOFVectorBase(const FiniteElemSpace *f, std::string n);
63

Thomas Witkowski's avatar
Thomas Witkowski committed
64
    virtual ~DOFVectorBase();
65

66
67
68
    /// Creates temporary used data structures.
    void createTempData();

Thomas Witkowski's avatar
Thomas Witkowski committed
69
70
71
72
    /** \brief 
     * For the given element, this function returns an array of all DOFs of this
     * DOFVector that are defined on this element.
     */
73
    virtual const T *getLocalVector(const Element *el, T* localVec) const;
74

Thomas Witkowski's avatar
Thomas Witkowski committed
75
76
77
78
    /** \brief
     * Evaluates the DOF vector at a set of quadrature points defined on the 
     * given element.
     */
79
80
    const T *getVecAtQPs(const ElInfo *elInfo, 
			 const Quadrature *quad,
81
			 const FastQuadrature *quadFast,
82
			 T *vecAtQPs) const;
83

Thomas Witkowski's avatar
Thomas Witkowski committed
84
85
86
87
88
89
    const T *getVecAtQPs(const ElInfo *smallElInfo, 
			 const ElInfo *largeElInfo,
			 const Quadrature *quad,
			 const FastQuadrature *quadFast,
			 T *vecAtQPs) const;

90
91
    const WorldVector<T> *getGrdAtQPs(const ElInfo *elInfo, 
				      const Quadrature *quad,
92
				      const FastQuadrature *quadFast,
93
				      WorldVector<T> *grdAtQPs) const;
94

Thomas Witkowski's avatar
Thomas Witkowski committed
95
96
97
98
99
100
    const WorldVector<T> *getGrdAtQPs(const ElInfo *smallElInfo, 
				      const ElInfo *largeElInfo,
				      const Quadrature *quad,
				      const FastQuadrature *quadFast,
				      WorldVector<T> *grdAtQPs) const;

101
102
    const WorldMatrix<T> *getD2AtQPs(const ElInfo *elInfo, 
				     const Quadrature *quad,
103
				     const FastQuadrature *quadFast,
104
				     WorldMatrix<T> *d2AtQPs) const;
105

Thomas Witkowski's avatar
Thomas Witkowski committed
106
    /// Returns the FE space the DOF vector is defined on.
107
    inline const FiniteElemSpace* getFeSpace() const 
Thomas Witkowski's avatar
Thomas Witkowski committed
108
    {
109
      return feSpace;
110
    }
111

Thomas Witkowski's avatar
Thomas Witkowski committed
112
113
114
115
116
117
118
119
120
121
122
123
124
    /** \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);
125
126

    void addElementVector(T sign,
127
			  const ElementVector& elVec, 
128
			  const BoundaryType *bound,
129
			  ElInfo *elInfo,
130
			  bool add = true); 
Thomas Witkowski's avatar
Thomas Witkowski committed
131
132
133
134
135
136
137

    /* \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();
138
139
140
141
142
143
144
145
 
    inline void addOperator(Operator* op, 
			    double *factor = NULL,
			    double *estFactor = NULL) 
    {
      operators.push_back(op);
      operatorFactor.push_back(factor);
      operatorEstFactor.push_back(estFactor);
146
    }
147

Thomas Witkowski's avatar
Thomas Witkowski committed
148
149
    inline std::vector<double*>::iterator getOperatorFactorBegin() 
    {
150
      return operatorFactor.begin();
151
    }
152

Thomas Witkowski's avatar
Thomas Witkowski committed
153
154
    inline std::vector<double*>::iterator getOperatorFactorEnd() 
    {
155
      return operatorFactor.end();
156
    }
157

Thomas Witkowski's avatar
Thomas Witkowski committed
158
159
    inline std::vector<double*>::iterator getOperatorEstFactorBegin() 
    {
160
      return operatorEstFactor.begin();
161
    }
162

Thomas Witkowski's avatar
Thomas Witkowski committed
163
164
    inline std::vector<double*>::iterator getOperatorEstFactorEnd() 
    {
165
      return operatorEstFactor.end();
166
    }
167
168


Thomas Witkowski's avatar
Thomas Witkowski committed
169
170
    inline std::vector<Operator*>::iterator getOperatorsBegin() 
    {
171
      return operators.begin();
172
    }
173

Thomas Witkowski's avatar
Thomas Witkowski committed
174
175
    inline std::vector<Operator*>::iterator getOperatorsEnd() 
    {
176
      return operators.end();
177
    }
178
179
180
181
182
183
184
185
186
187
188

    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
189
190
    inline std::vector<Operator*>& getOperators() 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
191
      return operators; 
192
    }
Thomas Witkowski's avatar
Thomas Witkowski committed
193

Thomas Witkowski's avatar
Thomas Witkowski committed
194
195
    inline std::vector<double*>& getOperatorFactor() 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
196
      return operatorFactor; 
197
    }
198

Thomas Witkowski's avatar
Thomas Witkowski committed
199
200
    inline std::vector<double*>& getOperatorEstFactor() 
    { 
201
      return operatorEstFactor; 
202
    }
203

Thomas Witkowski's avatar
Thomas Witkowski committed
204
    /// Returns \ref name
Thomas Witkowski's avatar
Thomas Witkowski committed
205
    inline std::string getName() const 
Thomas Witkowski's avatar
Thomas Witkowski committed
206
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
207
      return name; 
208
    } 
209

210
211
212
213
214
    inline void setName(std::string n)
    {
      name=n;
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
215
216
    inline BoundaryManager* getBoundaryManager() const 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
217
      return boundaryManager; 
218
    }
219

Thomas Witkowski's avatar
Thomas Witkowski committed
220
221
    inline void setBoundaryManager(BoundaryManager *bm) 
    {
222
      boundaryManager = bm;
223
    }
224

Thomas Witkowski's avatar
Thomas Witkowski committed
225
#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
226
    inline void setRankDofs(std::map<DegreeOfFreedom, bool> &dofmap)
Thomas Witkowski's avatar
Thomas Witkowski committed
227
    {
228
      rankDofs = &dofmap;
Thomas Witkowski's avatar
Thomas Witkowski committed
229
230
231
232
    }

    inline bool isRankDof(DegreeOfFreedom dof)
    {
233
      TEST_EXIT_DBG(rankDofs)("No rank dofs set!\n");
234

235
      return (*rankDofs)[dof];
Thomas Witkowski's avatar
Thomas Witkowski committed
236
237
238
    }
#endif

239
  protected:
Thomas Witkowski's avatar
Thomas Witkowski committed
240
    ///
241
242
    const FiniteElemSpace *feSpace;

Thomas Witkowski's avatar
Thomas Witkowski committed
243
    ///
244
    std::string name;
245

Thomas Witkowski's avatar
Thomas Witkowski committed
246
    ///
247
    ElementVector elementVector;
248

Thomas Witkowski's avatar
Thomas Witkowski committed
249
    ///
250
    std::vector<Operator*> operators;
251

Thomas Witkowski's avatar
Thomas Witkowski committed
252
    ///
253
    std::vector<double*> operatorFactor;
254

Thomas Witkowski's avatar
Thomas Witkowski committed
255
    ///
256
    std::vector<double*> operatorEstFactor;
257

Thomas Witkowski's avatar
Thomas Witkowski committed
258
    ///
259
260
    BoundaryManager *boundaryManager;

Thomas Witkowski's avatar
Thomas Witkowski committed
261
    /// Number of basis functions of the used finite element space.
Thomas Witkowski's avatar
Thomas Witkowski committed
262
    int nBasFcts;
263

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

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

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

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

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

Thomas Witkowski's avatar
Thomas Witkowski committed
279
    /// Dimension of the mesh this DOFVectorBase belongs to
Thomas Witkowski's avatar
Thomas Witkowski committed
280
    int dim;
Thomas Witkowski's avatar
Thomas Witkowski committed
281
282

#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
283
    std::map<DegreeOfFreedom, bool> *rankDofs;
Thomas Witkowski's avatar
Thomas Witkowski committed
284
#endif
Thomas Witkowski's avatar
Thomas Witkowski committed
285
  };
286

Thomas Witkowski's avatar
Thomas Witkowski committed
287
  /// Specifies which operation should be done after coarsening
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
  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)
319
      {}
320
321
322

      Iterator(DOFAdmin *admin, DOFIndexed<T> *c, DOFIteratorType type)
	: DOFIterator<T>(admin, c, type)
323
      {}
324
325
326
327
    };

    class Creator : public CreatorInterface<DOFVector<T> > {
    public:
328
329
330
      Creator(FiniteElemSpace *feSpace_) 
        : feSpace(feSpace_) 
      {}
331

Thomas Witkowski's avatar
Thomas Witkowski committed
332
333
      DOFVector<T> *create() 
      { 
334
	return new DOFVector<T>(feSpace, ""); 
335
      }
336

Thomas Witkowski's avatar
Thomas Witkowski committed
337
338
      void free(DOFVector<T> *vec) 
      { 
339
	delete vec; 
340
      }
341
342
343
344
345
346

    private:
      FiniteElemSpace *feSpace;
    };

  public:
Thomas Witkowski's avatar
Thomas Witkowski committed
347
    /// Empty constructor. No initialization!
348
349
350
    DOFVector() 
      : DOFVectorBase<T>(), 
	coarsenOperation(NO_OPERATION)
351
    {}
352

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

Thomas Witkowski's avatar
Thomas Witkowski committed
356
    /// Initialization.
357
    void init(const FiniteElemSpace* f, std::string n);
358

Thomas Witkowski's avatar
Thomas Witkowski committed
359
    /// Copy Constructor
Thomas Witkowski's avatar
Thomas Witkowski committed
360
361
    DOFVector(const DOFVector& rhs) 
    {
Thomas Witkowski's avatar
Thomas Witkowski committed
362
      *this = rhs;   
363
364
365
      this->name = rhs.name + "copy";
      if (this->feSpace && this->feSpace->getAdmin())
	(dynamic_cast<DOFAdmin*>(this->feSpace->getAdmin()))->addDOFIndexed(this);
366
367

      this->createTempData();
368
    }
369

Thomas Witkowski's avatar
Thomas Witkowski committed
370
    /// Destructor
371
372
    virtual ~DOFVector();

Thomas Witkowski's avatar
Thomas Witkowski committed
373
    /// Returns iterator to the begin of \ref vec
Thomas Witkowski's avatar
Thomas Witkowski committed
374
375
    typename std::vector<T>::iterator begin() 
    { 
376
      return vec.begin(); 
377
    }
378

Thomas Witkowski's avatar
Thomas Witkowski committed
379
    /// Returns iterator to the end of \ref vec
Thomas Witkowski's avatar
Thomas Witkowski committed
380
381
    typename std::vector<T>::iterator end() 
    { 
382
      return vec.end(); 
383
    }
384
385
386
387
388
389
  
    /** \brief
     * Used by DOFAdmin to compress this DOFVector. Implementation of
     * DOFIndexedBase::compress()
     */
    virtual void compressDOFIndexed(int first, int last,
390
				    std::vector<DegreeOfFreedom> &newDof);
Thomas Witkowski's avatar
Thomas Witkowski committed
391
392

    /// Sets \ref coarsenOperation to op
Thomas Witkowski's avatar
Thomas Witkowski committed
393
394
    inline void setCoarsenOperation(CoarsenOperation op) 
    { 
395
      coarsenOperation = op; 
396
    }
397

Thomas Witkowski's avatar
Thomas Witkowski committed
398
    /// Returns \ref coarsenOperation
Thomas Witkowski's avatar
Thomas Witkowski committed
399
400
    inline CoarsenOperation getCoarsenOperation() 
    { 
401
      return coarsenOperation; 
402
    }
403

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

Thomas Witkowski's avatar
Thomas Witkowski committed
407
    /// Interpolation after refinement.
408
    inline void refineInterpol(RCNeighbourList&, int) {}
409

Thomas Witkowski's avatar
Thomas Witkowski committed
410
    /// Returns \ref vec
Thomas Witkowski's avatar
Thomas Witkowski committed
411
412
    std::vector<T>& getVector() 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
413
      return vec;
414
    }
415

Thomas Witkowski's avatar
Thomas Witkowski committed
416
    /// Returns size of \ref vec
Thomas Witkowski's avatar
Thomas Witkowski committed
417
418
    inline int getSize() const 
    { 
419
      return vec.size();
420
    } 
421

Thomas Witkowski's avatar
Thomas Witkowski committed
422
    /// Returns used size of the vector.
Thomas Witkowski's avatar
Thomas Witkowski committed
423
424
    inline int getUsedSize() const 
    { 
425
      return this->feSpace->getAdmin()->getUsedSize(); 
426
    }
427

Thomas Witkowski's avatar
Thomas Witkowski committed
428
    /// Resizes \ref vec to n
Thomas Witkowski's avatar
Thomas Witkowski committed
429
430
    inline void resize(int n) 
    { 
431
      FUNCNAME("DOFVector<T>::resize()");
432
      TEST_EXIT_DBG(n >= 0)("Can't resize DOFVector to negative size\n"); 
433
      vec.resize(n);
434
    } 
435

Thomas Witkowski's avatar
Thomas Witkowski committed
436
    /// Resizes \ref vec to n and inits new values with init
Thomas Witkowski's avatar
Thomas Witkowski committed
437
438
    inline void resize(int n, T init) 
    { 
439
      FUNCNAME("DOFVector<T>::resize()");
440
      TEST_EXIT_DBG(n >= 0)("Can't resize DOFVector to negative size\n"); 
441
      vec.resize(n, init);
442
    } 
443

Thomas Witkowski's avatar
Thomas Witkowski committed
444
    /// Returns \ref vec[i]
Thomas Witkowski's avatar
Thomas Witkowski committed
445
446
    inline const T& operator[](DegreeOfFreedom i) const 
    {
447
      FUNCNAME("DOFVector<T>::operator[]");
448
      TEST_EXIT_DBG(i >= 0 && i < static_cast<int>(vec.size()))
449
	("Illegal vector index %d.\n",i);
450
      return vec[i];
451
    } 
452

Thomas Witkowski's avatar
Thomas Witkowski committed
453
    /// Returns \ref vec[i]
Thomas Witkowski's avatar
Thomas Witkowski committed
454
455
    inline T& operator[](DegreeOfFreedom i) 
    {
456
      FUNCNAME("DOFVector<T>::operator[]");
457
458
459
460

      TEST_EXIT_DBG(i >= 0 && i < static_cast<int>(vec.size())) 
 	("Illegal vector index %d.\n",i); 

461
      return vec[i];
462
    }
Thomas Witkowski's avatar
Thomas Witkowski committed
463
 
Thomas Witkowski's avatar
Thomas Witkowski committed
464
    /// Calculates Integral of this DOFVector
465
466
    double Int(Quadrature* q = NULL) const;

Thomas Witkowski's avatar
Thomas Witkowski committed
467
    /// Calculates L1 norm of this DOFVector
468
    double L1Norm(Quadrature* q = NULL) const;
Thomas Witkowski's avatar
Thomas Witkowski committed
469
 
Thomas Witkowski's avatar
Thomas Witkowski committed
470
    /// Calculates L2 norm of this DOFVector
Thomas Witkowski's avatar
Thomas Witkowski committed
471
472
    inline double L2Norm(Quadrature* q = NULL) const 
    {
473
      return sqrt(L2NormSquare());
474
    }
475

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

Thomas Witkowski's avatar
Thomas Witkowski committed
479
    /// Calculates H1 norm of this DOFVector
Thomas Witkowski's avatar
Thomas Witkowski committed
480
481
    inline double H1Norm(Quadrature* q = NULL) const 
    {
482
      return sqrt(H1NormSquare());
483
    }
484

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

Thomas Witkowski's avatar
Thomas Witkowski committed
488
    /// Calculates euclidian norm of this DOFVector
489
490
    double nrm2() const; 

Thomas Witkowski's avatar
Thomas Witkowski committed
491
    /// Returns square of the euclidian norm.
492
493
    double squareNrm2() const;

Thomas Witkowski's avatar
Thomas Witkowski committed
494
    /// Calculates l2 norm of this DOFVector
Thomas Witkowski's avatar
Thomas Witkowski committed
495
496
    inline double l2norm() const 
    { 
497
      return nrm2();
498
    }
499

Thomas Witkowski's avatar
Thomas Witkowski committed
500
    /// Calculates the absolute sum of this DOFVector
501
502
    T asum() const; 

Thomas Witkowski's avatar
Thomas Witkowski committed
503
    /// Calculates the l1 norm of this DOFVector
Thomas Witkowski's avatar
Thomas Witkowski committed
504
505
    inline double l1norm() const 
    { 
506
      return asum();
507
    } 
508

Thomas Witkowski's avatar
Thomas Witkowski committed
509
    /// Calculates doublewell of this DOFVector
510
    double DoubleWell(Quadrature* q = NULL) const;
Thomas Witkowski's avatar
Thomas Witkowski committed
511
 
Thomas Witkowski's avatar
Thomas Witkowski committed
512
    /// Calculates the sum of this DOFVector
513
    T sum() const; 
Thomas Witkowski's avatar
Thomas Witkowski committed
514
 
Thomas Witkowski's avatar
Thomas Witkowski committed
515
    /// Sets \ref vec[i] = val, i=0 , ... , size
516
517
    void set(T val); 

Thomas Witkowski's avatar
Thomas Witkowski committed
518
    /// Assignment operator for setting complete vector to a certain value d
Thomas Witkowski's avatar
Thomas Witkowski committed
519
520
    inline DOFVector<T>& operator=(T d) 
    {
521
522
      set(d); 
      return *this;
523
    } 
524

Thomas Witkowski's avatar
Thomas Witkowski committed
525
    /// Assignment operator between two vectors
526
    DOFVector<T>& operator=(const DOFVector<T>&); 
527

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

531
    /// Returns minimum of DOFVector.
532
533
    T min() const; 

534
    /// Returns maximum of DOFVector.
Thomas Witkowski's avatar
Thomas Witkowski committed
535
536
    T max() const;

537
    /// Returns absolute maximum of DOFVector.
Thomas Witkowski's avatar
Thomas Witkowski committed
538
    T absMax() const;
539

540
541
542
543
    /// Returns the average value of the DOFVector.
    T average() const;

    /// Used by interpol while mesh traversal.
544
    void interpolFct(ElInfo* elinfo);
545

546
    /// Prints \ref vec to stdout.
547
548
    void print() const; 

549
    ///
550
551
    int calcMemoryUsage() const;

552
553
554
555
556
557
    /** \brief
     * Computes the coefficients of the interpolant of the function fct and
     * stores these in the DOFVector
     */
    void interpol(AbstractFunction<T, WorldVector<double> > *fct);

558
    void interpol(DOFVector<T> *v, double factor=1.0);
559

560
    /// Writes the data of the DOFVector to an output stream.
Thomas Witkowski's avatar
Thomas Witkowski committed
561
562
    void serialize(std::ostream &out) 
    {
563
564
565
      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));
566
    }
567

568
    /// Reads data of an DOFVector from an input stream.
Thomas Witkowski's avatar
Thomas Witkowski committed
569
570
    void deserialize(std::istream &in) 
    {
571
572
573
574
      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));
575
    }
576
577
578
579
580

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

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

581
    DOFVector<WorldVector<T> >* getRecoveryGradient(DOFVector<WorldVector<T> >*) const;
582
583

  protected: 
584
    /// Data container
585
    std::vector<T> vec; 
586
 
587
    /// Specifies what operation should be performed after coarsening
588
589
    CoarsenOperation coarsenOperation;

590
    /// Used by \ref interpol
591
592
    AbstractFunction<T, WorldVector<double> > *interFct;

593
    /// Used for mesh traversal
594
595
596
    static DOFVector<T> *traverseVector;
  }; 

Thomas Witkowski's avatar
Thomas Witkowski committed
597

598
599
600
601
602
603
  template<>
  void DOFVector<double>::refineInterpol(RCNeighbourList&, int); 

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

Thomas Witkowski's avatar
Thomas Witkowski committed
604
605
  inline double min(const DOFVector<double>& v) 
  {
Thomas Witkowski's avatar
Thomas Witkowski committed
606
    return v.min();
607
  } 
Thomas Witkowski's avatar
Thomas Witkowski committed
608

Thomas Witkowski's avatar
Thomas Witkowski committed
609
610
  inline double max(const DOFVector<double>& v) 
  {
Thomas Witkowski's avatar
Thomas Witkowski committed
611
    return v.max();
612
  }
613
614
615
616
617
618
619
620
621
622
623
624
625
626


  /** \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
     */
627
    DOFVectorDOF(const FiniteElemSpace* feSpace_, std::string name_)
628
629
630
      : DOFVector<DegreeOfFreedom>(feSpace_, name_)
    {
      feSpace->getAdmin()->addDOFContainer(this);
631
    }
632
  
633
    /// Deregisters itself at DOFAdmin.
Thomas Witkowski's avatar
Thomas Witkowski committed
634
635
    ~DOFVectorDOF() 
    {
636
      feSpace->getAdmin()->removeDOFContainer(this);
637
    }
638
639
640
641
642

    /** \brief
     * Implements DOFContainer::operator[]() by calling 
     * DOFVEctor<DegreeOfFreedom>::operator[]()
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
643
644
    DegreeOfFreedom& operator[](DegreeOfFreedom i) 
    {
645
      return DOFVector<DegreeOfFreedom>::operator[](i);
646
    }
647

Thomas Witkowski's avatar
Thomas Witkowski committed
648
649
    const DegreeOfFreedom& operator[](DegreeOfFreedom i) const 
    {
650
651
652
      return DOFVector<DegreeOfFreedom>::operator[](i);
    }

653
    /// Implements DOFIndexedBase::getSize()
Thomas Witkowski's avatar
Thomas Witkowski committed
654
655
    int getSize() const 
    {
656
      return DOFVector<DegreeOfFreedom>::getSize();
657
    }
658

659
    /// Implements DOFIndexedBase::resize()
Thomas Witkowski's avatar
Thomas Witkowski committed
660
661
    void resize(int size) 
    {
662
663
664
665
666
667
668
669
670
      DOFVector<DegreeOfFreedom>::resize(size);
    }

    void freeDOFContent(DegreeOfFreedom dof);

  protected:
    DOFVectorDOF();
  };

Thomas Witkowski's avatar
Thomas Witkowski committed
671

672
  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
673
674
  double norm(DOFVector<T> *vec) 
  {
675
    return vec->nrm2();
676
  }
677
678

  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
679
680
  double L2Norm(DOFVector<T> *vec) 
  {
681
    return vec->L2Norm();
682
  }
683
684

  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
685
686
  double H1Norm(DOFVector<T> *vec) 
  {
687
    return vec->H1Norm();
688
  }
689
690

  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
691
692
  void print(DOFVector<T> *vec) 
{
693
    vec->print();
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
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740

  // 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
741
742
  inline void scal(T a, DOFVector<T>& y) 
  {
743
744
    y *= a;
  }
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760

  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);

761
762
763
764
765
766
767
768
769
  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);


770
  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
771
772
  inline void set(DOFVector<T>& vec, T d) 
  {
773
    vec.set(d);
774
  }
775
776

  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
777
778
  inline void setValue(DOFVector<T>& vec, T d) 
  {
779
    vec.set(d);
780
  }
781
782

  template<typename T>
Thomas Witkowski's avatar
Thomas Witkowski committed
783
784
  inline int size(DOFVector<T> *vec) 
  {
785
    return vec->getUsedSize();
786
  }
787

Thomas Witkowski's avatar
Thomas Witkowski committed
788
789
790
791
792
793
794
795
796
797
  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());
  }

798
799
  WorldVector<DOFVector<double>*> *transform(DOFVector<WorldVector<double> > *vec,
					     WorldVector<DOFVector<double>*> *result);
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819

  /** \brief
   * Computes the integral of a function that includes two different DOFVectors. This
   * function works also for the case that the DOFVectors are defined on two different
   * meshes.
   */
  double integrate(const DOFVector<double> &vec1,
		   const DOFVector<double> &vec2,
		   BinaryAbstractFunction<double, double, double> *fct);

  /// Computes the integral of a function with two DOFVectors defined on the same mesh.
  double intSingleMesh(const DOFVector<double> &vec1,
		       const DOFVector<double> &vec2,
		       BinaryAbstractFunction<double, double, double> *fct);

  /// Computes the integral of a function with two DOFVectors defined on different meshes.
  double intDualMesh(const DOFVector<double> &vec1,
		     const DOFVector<double> &vec2,
		     BinaryAbstractFunction<double, double, double> *fct);

820
821
822
823
824
}

#include "DOFVector.hh"

#endif  // !_DOFVECTOR_H_