ProblemStat.h 21.2 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
/******************************************************************************
 *
 * AMDiS - Adaptive multidimensional simulations
 *
 * Copyright (C) 2013 Dresden University of Technology. All Rights Reserved.
 * Web: https://fusionforge.zih.tu-dresden.de/projects/amdis
 *
 * Authors: 
 * Simon Vey, Thomas Witkowski, Andreas Naumann, Simon Praetorius, et al.
 *
 * This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
 * WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
 *
 *
 * This file is part of AMDiS
 *
 * See also license.opensource.txt in the distribution.
 * 
 ******************************************************************************/
20
21


22

23
/** \file ProblemStat.h */
24

25
26
#ifndef AMDIS_PROBLEM_STAT_H
#define AMDIS_PROBLEM_STAT_H
27

Thomas Witkowski's avatar
Thomas Witkowski committed
28
29
#include <vector>
#include <list>
30
#include "AMDiS_fwd.h"
31
#include "ProblemStatBase.h"
32
#include "Initfile.h"
33
34
35
#include "Boundary.h"
#include "MatrixVector.h"
#include "StandardProblemIteration.h"
36
#include "io/ElementFileWriter.h"
Thomas Witkowski's avatar
Thomas Witkowski committed
37
#include "ComponentTraverseInfo.h"
38
#include "AbstractFunction.h"
39
#include "solver/SolverMatrix.h"
40
#include "SystemVector.h"
41
42
43

namespace AMDiS {

44
45
46
47
48
49
  struct OperatorPos 
  {
    int row, col;
    double *factor, *estFactor;
  };

50

51
52
53
54
55
56
57
  /**
   * \ingroup Problem 
   *
   * \brief
   * This class defines the stationary problem definition in sequential
   * computations. For parallel computations, see \ref ParallelProblemStatBase.
   **/
58
  class ProblemStatSeq : public ProblemStatBase
59
  {
60
  protected:
61
    // Defines a mapping type from dof indices to world coordinates.
62
    typedef std::map<DegreeOfFreedom, WorldVector<double> > DofToCoord;
63
64

    // Defines a mapping type from dof indices to world coordinates.
65
    typedef std::map<WorldVector<double>, DegreeOfFreedom> CoordToDof;
66

67
  public:
Thomas Witkowski's avatar
Thomas Witkowski committed
68
    /// Constructor
69
    ProblemStatSeq(std::string nameStr,
70
		   ProblemIterationInterface *problemIteration = NULL);
71

Thomas Witkowski's avatar
Thomas Witkowski committed
72
    /// Destructor
Thomas Witkowski's avatar
Thomas Witkowski committed
73
    virtual ~ProblemStatSeq();
74

Thomas Witkowski's avatar
Thomas Witkowski committed
75
    /// Initialisation of the problem.
76
    virtual void initialize(Flag initFlag,
77
			    ProblemStatSeq *adoptProblem = NULL,
78
79
			    Flag adoptFlag = INIT_NOTHING);

Thomas Witkowski's avatar
Thomas Witkowski committed
80
    /// Used in \ref initialize().
81
82
    virtual void createMesh();

83
84
85
    /// Used in \ref initialize().
    virtual void createRefCoarseManager();

Thomas Witkowski's avatar
Thomas Witkowski committed
86
    /// Used in \ref initialize().
87
    virtual void createFeSpace(DOFAdmin *admin);
88

Thomas Witkowski's avatar
Thomas Witkowski committed
89
    /// Used in \ref initialize().
90
91
    virtual void createMatricesAndVectors();

Thomas Witkowski's avatar
Thomas Witkowski committed
92
    /// Used in \ref initialize().
93
94
    virtual void createSolver();

Thomas Witkowski's avatar
Thomas Witkowski committed
95
    /// Used in \ref initialize().
96
97
    virtual void createEstimator();

Thomas Witkowski's avatar
Thomas Witkowski committed
98
    /// Used in \ref initialize().
99
100
    virtual void createMarker();

Thomas Witkowski's avatar
Thomas Witkowski committed
101
    /// Used in \ref initialize().
102
103
    virtual void createFileWriter();

104
105
    /// Used in \ref initialize(). This function is deprecated and should not 
    /// be used anymore. There is no guarantee that it will work in the future.
106
107
    virtual void doOtherStuff();

108
109
    /// Implementation of ProblemStatBase::solve(). Deligates the solving
    /// of problems system to \ref solver.
110
111
    void solve(AdaptInfo *adaptInfo,
	       bool createMatrixData = true,
112
	       bool storeMatrixData = false) override;
113

114
115
    /// Implementation of ProblemStatBase::estimate(). Deligates the estimation
    /// to \ref estimator.
116
    void estimate(AdaptInfo *adaptInfo) override;
117

118
119
    /// Implementation of ProblemStatBase::markElements().
    /// Deligated to \ref adapt.
120
    Flag markElements(AdaptInfo *adaptInfo) override;
121

122
123
    /// Implementation of ProblemStatBase::refineMesh(). Deligated to the
    /// RefinementManager of \ref adapt.
124
    Flag refineMesh(AdaptInfo *adaptInfo) override;
125

126
127
    /// Implementation of ProblemStatBase::coarsenMesh(). Deligated to the
    /// CoarseningManager of \ref adapt.
128
    Flag coarsenMesh(AdaptInfo *adaptInfo) override;
129

130
131
    /// Implementation of ProblemStatBase::buildBeforeRefine().
    /// Does nothing here.
132
    void buildBeforeRefine(AdaptInfo *adaptInfo, Flag) override {}
133

134
135
    /// Implementation of ProblemStatBase::buildBeforeCoarsen().
    /// Does nothing here.
136
    void buildBeforeCoarsen(AdaptInfo *adaptInfo, Flag) override {}
137

138
139
140
    /// Implementation of ProblemStatBase::buildAfterCoarsen().
    /// Assembles \ref A and \ref rhs. With the last two parameters, assembling
    /// can be restricted to matrices or vectors only.
141
    void buildAfterCoarsen(AdaptInfo *adaptInfo, Flag flag,
142
				   bool assembleMatrix = true,
143
				   bool assembleVector = true) override;
144

145
146
    bool dualMeshTraverseRequired();

147
148
    void dualAssemble(AdaptInfo *adaptInfo, Flag flag, 
		      bool asmMatrix = true, bool asmVector = true);
149

150
    /// Implementation of ProblemStatBase::getNumComponents(), TODO: Wrong!!
151
152
    virtual int getNumComponents() 
    { 
153
      return nComponents; 
Thomas Witkowski's avatar
Thomas Witkowski committed
154
    }
155
156
157
158
159
    
    virtual int getNumAddComponents()
    {
      return nAddComponents;
    }
160

161
    /// Writes output files. TODO: Make obsolete.
162
163
    void writeFiles(AdaptInfo *adaptInfo, bool force);

164
165
166
    /// Writes output files.
    void writeFiles(AdaptInfo &adaptInfo, bool force);

Thomas Witkowski's avatar
Thomas Witkowski committed
167
    /// Interpolates fct to \ref solution.
168
    void interpolInitialSolution(std::vector<AbstractFunction<double, WorldVector<double> >*> *fct);
169

170
    /// Adds an operator to \ref A.
171
    void addMatrixOperator(Operator *op, int i, int j,
172
			   double *factor = NULL, double *estFactor = NULL);
173

174
    /// Adds an operator to \ref A.
175
    void addMatrixOperator(Operator &op, int i, int j,
176
			   double *factor = NULL, double *estFactor = NULL);
177

178
    /// Adds an operator to \ref rhs.
179
    void addVectorOperator(Operator *op, int i,
180
			   double *factor = NULL, double *estFactor = NULL);
181

182
    /// Adds an operator to \ref rhs.
183
    void addVectorOperator(Operator &op, int i,
184
			   double *factor = NULL, double *estFactor = NULL);
185

186
187
    /// Adds a Dirichlet boundary condition, where the rhs is given by an 
    /// abstract function.
188
    virtual void addDirichletBC(BoundaryType type, int row, int col,
189
190
				AbstractFunction<double, WorldVector<double> > *b);

191
192
    /// Adds a Dirichlet boundary condition, where the rhs is given by a DOF
    /// vector.
193
194
195
    virtual void addDirichletBC(BoundaryType type, int row, int col,
				DOFVector<double> *vec);

196
197
    /// Adds a Neumann boundary condition, where the rhs is given by an
    /// abstract function.
198
199
200
    virtual void addNeumannBC(BoundaryType type, int row, int col, 
			      AbstractFunction<double, WorldVector<double> > *n);

201
202
203
204
205
206
    /// Adds a Neumann boundary condition, where the rhs is given by an DOF
    /// vector.
    virtual void addNeumannBC(BoundaryType type, int row, int col, 
			      DOFVector<double> *n);

    /// Adds Robin boundary condition.
207
208
209
210
    virtual void addRobinBC(BoundaryType type, int row, int col, 
			    AbstractFunction<double, WorldVector<double> > *n,
			    AbstractFunction<double, WorldVector<double> > *r);

211
212
213
214
215
216
217
218
219
220
    /// Adds Robin boundary condition.
    virtual void addRobinBC(BoundaryType type, int row, int col, 
			    DOFVector<double> *n,
			    DOFVector<double> *r);

    /// Adds Robin boundary condition.
    virtual void addRobinBC(BoundaryType type, int row, int col, 
			    Operator *n,
			    Operator *r);

221
    /// Adds a periodic boundary condition.
222
223
    virtual void addPeriodicBC(BoundaryType type, int row, int col);

224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
    /// add boundary operator to matrix side
    virtual void addBoundaryMatrixOperator(BoundaryType type, 
          Operator *op, int row, int col);

    virtual void addBoundaryMatrixOperator(BoundaryType type, 
          Operator &op, int row, int col)
    {
      addBoundaryMatrixOperator(type, &op, row, col);
    }

    /// add boundary operator to rhs (vector) side
    virtual void addBoundaryVectorOperator(BoundaryType type, 
          Operator *op, int row);

    virtual void addBoundaryVectorOperator(BoundaryType type, 
          Operator &op, int row)
    {
      addBoundaryVectorOperator(type, &op, row);
    }

244
245
    /// This function assembles a DOFMatrix and a DOFVector for the case,
    /// the meshes from row and col FE-space are equal.
246
    void assembleOnOneMesh(const FiniteElemSpace *feSpace, 
Thomas Witkowski's avatar
Thomas Witkowski committed
247
			   Flag assembleFlag,
Thomas Witkowski's avatar
Thomas Witkowski committed
248
249
			   DOFMatrix *matrix, DOFVector<double> *vector);

Thomas Witkowski's avatar
Thomas Witkowski committed
250
251

    ///
252
253
254
255
    void assembleBoundaryConditions(DOFVector<double> *rhs,
				    DOFVector<double> *solution,
				    Mesh *mesh,
				    Flag assembleFlag);
256
257

 
258
259
260
261
    /** \name getting methods
     * \{ 
     */

Thomas Witkowski's avatar
Thomas Witkowski committed
262
    /// Returns \ref solution.
263
264
    inline SystemVector* getSolution() 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
265
266
      return solution; 
    }
267

268
    inline DOFVector<double>* getSolution(int i)
269
270
271
272
    {
      return solution->getDOFVector(i);
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
273
    /// Returns \ref rhs.
274
    inline SystemVector* getRhs() 
275
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
276
277
      return rhs; 
    }
278

279
    inline DOFVector<double>* getRhsVector(int i = 0)
280
281
282
283
    {
      return rhs->getDOFVector(i);
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
284
    /// Returns \ref systemMatrix.
285
286
    inline Matrix<DOFMatrix*> *getSystemMatrix() 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
287
288
      return systemMatrix; 
    }
289

Thomas Witkowski's avatar
Thomas Witkowski committed
290
    /// Returns a pointer to the corresponding DOFMatrix.
291
292
    inline DOFMatrix* getSystemMatrix(int row, int col) 
    {
Thomas Witkowski's avatar
Thomas Witkowski committed
293
294
295
296
      return (*systemMatrix)[row][col];
    }

    /// Returns mesh of given component
297
    inline Mesh* getMesh(int comp = 0) 
298
    {
299
      FUNCNAME("ProblemStatSeq::getMesh()");
Thomas Witkowski's avatar
Thomas Witkowski committed
300
      TEST_EXIT(comp < static_cast<int>(componentMeshes.size()) && comp >= 0)
301
	("invalid component number\n");
Thomas Witkowski's avatar
Thomas Witkowski committed
302
      return componentMeshes[comp]; 
Thomas Witkowski's avatar
Thomas Witkowski committed
303
    }
304

Thomas Witkowski's avatar
Thomas Witkowski committed
305
    /// Returns \ref meshes
306
    inline std::vector<Mesh*> getMeshes() 
307
    {
Thomas Witkowski's avatar
Thomas Witkowski committed
308
309
      return meshes; 
    }
310

Thomas Witkowski's avatar
Thomas Witkowski committed
311
    /// Returns \ref feSpace_.
312
    inline const FiniteElemSpace* getFeSpace(int comp = 0) 
313
    { 
314
      FUNCNAME("ProblemStatSeq::getFeSpace()");
315
      TEST_EXIT(comp < static_cast<int>(componentSpaces.size()) && comp >= 0)
316
	("invalid component number: %d\n", comp);
317
      return componentSpaces[comp]; 
Thomas Witkowski's avatar
Thomas Witkowski committed
318
    }
319

Thomas Witkowski's avatar
Thomas Witkowski committed
320
    /// Returns \ref feSpaces.
321
    inline std::vector<const FiniteElemSpace*>& getFeSpaces() 
322
    { 
323
      return feSpaces; 
Thomas Witkowski's avatar
Thomas Witkowski committed
324
    }
325

Thomas Witkowski's avatar
Thomas Witkowski committed
326
    /// Returns \ref componentSpaces;
327
    inline std::vector<const FiniteElemSpace*>& getComponentSpaces() 
328
    {
329
      return componentSpaces;
330
331
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
332
    /// Returns \ref estimator.
333
    inline std::vector<Estimator*> getEstimators() 
334
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
335
336
      return estimator; 
    }
337

Thomas Witkowski's avatar
Thomas Witkowski committed
338
    /// Returns \ref estimator.
339
    inline Estimator* getEstimator(int comp = 0) 
340
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
341
342
      return estimator[comp]; 
    }
343

Thomas Witkowski's avatar
Thomas Witkowski committed
344
    /// Returns \ref refinementManager.
345
    inline RefinementManager* getRefinementManager(int comp = 0) 
346
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
347
348
      return refinementManager; 
    }
349

Thomas Witkowski's avatar
Thomas Witkowski committed
350
    /// Returns \ref refinementManager.
351
    inline CoarseningManager* getCoarseningManager(int comp = 0) 
352
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
353
354
      return coarseningManager; 
    }
355

Thomas Witkowski's avatar
Thomas Witkowski committed
356
    /// Returns \ref solver.
357
    inline LinearSolverInterface* getSolver() 
358
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
359
360
      return solver; 
    }
361

Thomas Witkowski's avatar
Thomas Witkowski committed
362
    /// Returns \ref marker.
363
    inline Marker *getMarker(int comp = 0) 
364
    { 
365
      return marker[comp]; 
Thomas Witkowski's avatar
Thomas Witkowski committed
366
    }
367

Thomas Witkowski's avatar
Thomas Witkowski committed
368
    /// Returns \ref marker.
369
    inline std::vector<Marker*> getMarkers() 
370
    { 
371
      return marker; 
Thomas Witkowski's avatar
Thomas Witkowski committed
372
    }
373

Thomas Witkowski's avatar
Thomas Witkowski committed
374
    /// Returns the name of the problem
375
    inline std::string getName() override
376
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
377
378
      return name; 
    }
379

Praetorius, Simon's avatar
Praetorius, Simon committed
380
    /// Returns the name of the problem
381
    inline std::string getComponentName(int comp = 0)
Praetorius, Simon's avatar
Praetorius, Simon committed
382
383
384
385
386
387
    {
      TEST_EXIT(comp < static_cast<int>(componentNames.size()) && comp >= 0)
	("invalid component number\n");
      return componentNames[comp];
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
388
    /// Returns \ref useGetBound.
389
390
    inline bool getBoundUsed() 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
391
392
      return useGetBound; 
    }
393

394
395
396
397
398
399
    /// Returns \ref info.
    int getInfo()
    {
      return info;
    }

400
401
402
403
404
405
    /// Returns \ref deserialized;
    bool isDeserialized()
    {
      return deserialized;
    }

406
407
408
409
410
411
    /** \} */

    /** \name setting methods
     * \{ 
     */

Thomas Witkowski's avatar
Thomas Witkowski committed
412
    /// Sets \ref estimator.
413
    inline void setEstimator(std::vector<Estimator*> est) 
414
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
415
416
      estimator = est; 
    }
417

Thomas Witkowski's avatar
Thomas Witkowski committed
418
    /// Sets the FE space for the given component.
419
    inline void setFeSpace(const FiniteElemSpace *feSpace, int comp = 0) 
420
    {
421
      feSpaces[comp] = feSpace;
Thomas Witkowski's avatar
Thomas Witkowski committed
422
    }
423

424
425
426
427
428
    void setFeSpaces(std::vector<FiniteElemSpace const*> feSpaces_)
    {
      feSpaces = feSpaces_;
    }

429
    void setComponentSpace(int comp, const FiniteElemSpace *feSpace)
430
431
432
    {
      if (static_cast<int>(componentSpaces.size()) < nComponents)
        componentSpaces.resize(nComponents);
433
      TEST_EXIT(comp >= 0 && comp < nComponents + nAddComponents)
434
435
436
437
438
        ("Component number not in feasable range!");

      componentSpaces[comp] = feSpace;
    }
    
Thomas Witkowski's avatar
Thomas Witkowski committed
439
    /// Sets \ref estimator.
440
    inline void setEstimator(Estimator* est, int comp = 0) 
441
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
442
443
      estimator[comp] = est; 
    }
444

Thomas Witkowski's avatar
Thomas Witkowski committed
445
    /// Sets \ref marker.
446
    inline void setMarker(Marker* mark, int comp = 0) 
447
    { 
448
      marker[comp] = mark; 
Thomas Witkowski's avatar
Thomas Witkowski committed
449
    }
450

Thomas Witkowski's avatar
Thomas Witkowski committed
451
    /// Sets \ref solver.
452
    inline void setSolver(LinearSolverInterface* sol) 
453
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
454
455
      solver = sol; 
    }
456

457
458
459
460
461
    void setSolverMatrix(SolverMatrix<Matrix<DOFMatrix*> >& solverMatrix_)
    {
      solverMatrix.setMatrix(*solverMatrix_.getOriginalMat());
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
462
    ///
463
    inline void setAssembleMatrixOnlyOnce(int i = 0, int j = 0, bool value = true) 
464
    {
Thomas Witkowski's avatar
Thomas Witkowski committed
465
      assembleMatrixOnlyOnce[i][j] = value;
466
467
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
468
    ///
469
    void setExactSolutionFct(AbstractFunction<double, WorldVector<double> > *fct,
470
471
			     int component) 
    {
472
473
474
      exactSolutionFcts[component] = fct;
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
475
    ///
476
    AbstractFunction<double, WorldVector<double> > *getExactSolutionFct(int i = 0) 
477
    {
478
      return exactSolutionFcts[i];
479
480
    }

481
    ///
482
    std::vector< AbstractFunction<double, WorldVector<double> >* > getExactSolutionFcts() 
483
484
485
486
    {
      return exactSolutionFcts;
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
487
    ///
488
489
    void setComputeExactError(bool v) 
    {
490
491
492
      computeExactError = v;
    }

493
494
495
496
497
498
    /// Sets \ref writeAsmInfo;
    void setWriteAsmInfo(bool b)
    {
      writeAsmInfo = b;
    }

499
    void setMeshes(std::vector<Mesh*> meshes_)
500
501
502
503
504
505
506
507
508
    {
      meshes = meshes_;
      nMeshes = static_cast<int>(meshes.size());
    }

    void setComponentMesh(int comp, Mesh* mesh)
    {
      if (static_cast<int>(componentMeshes.size()) < nComponents)
        componentMeshes.resize(nComponents);
509
      TEST_EXIT(comp >= 0 && comp < nComponents + nAddComponents)
510
511
512
513
514
515
516
517
518
519
520
521
522
523
        ("Component number not in feasable range!");

      componentMeshes[comp] = mesh;
    }

    void setRefinementManager(RefinementManager *ref)
    {
      refinementManager = ref;
    }

    void setCoarseningManager(CoarseningManager *coarse)
    {
      coarseningManager = coarse;
    }
524
525
    /** \} */

526
527
    /// Outputs the mesh of the given component, but the values are taken from
    /// the residual error estimator. 
528
    void writeResidualMesh(int comp, AdaptInfo *adaptInfo, std::string name);
529

Thomas Witkowski's avatar
Thomas Witkowski committed
530
    /// Function that implements the serialization procedure.
531
    void serialize(std::ostream &out) override;
532

Thomas Witkowski's avatar
Thomas Witkowski committed
533
    /// Function that implements the deserialization procedure.
534
    void deserialize(std::istream &in) override;
Thomas Witkowski's avatar
Thomas Witkowski committed
535
536

    /// Returns \ref fileWriters.
537
    std::vector<FileWriterInterface*>& getFileWriterList() 
538
    {
Thomas Witkowski's avatar
Thomas Witkowski committed
539
540
      return fileWriters;
    }
541
    
542
    /// Returns \ref solutionTime.
543
544
545
546
547
    double getSolutionTime()
    {
      return solutionTime;
    }
    
548
    /// Returns \ref buildTime.
549
550
551
552
    double getBuildTime()
    {
      return buildTime;
    }
553

Thomas Witkowski's avatar
Thomas Witkowski committed
554
  protected:
555
556
    /// If the exact solution is known, the problem can compute the exact
    /// error instead of the error estimation. This is done in this function.
Thomas Witkowski's avatar
Thomas Witkowski committed
557
558
    void computeError(AdaptInfo *adaptInfo);

559
  protected:
Thomas Witkowski's avatar
Thomas Witkowski committed
560
561
    
    /// Name of this problem.
562
    std::string name;
563

Thomas Witkowski's avatar
Thomas Witkowski committed
564
    /// Number of problem components
565
    int nComponents;
566
567
568
    
    /// Number of additional components
    int nAddComponents;
569

570
571
    /// Stores the names for all components. Is used for naming the solution
    /// vectors, \ref solution.
572
    std::vector<std::string> componentNames;
573

574
575
576
    /// Number of problem meshes. If all components are defined on the same mesh,
    /// this number is 1. Otherwise, this variable is the number of different meshes
    /// within the problem.
Thomas Witkowski's avatar
Thomas Witkowski committed
577
578
579
    int nMeshes;

    /// FE spaces of this problem.
580
    std::vector<const FiniteElemSpace*> feSpaces;
581

Thomas Witkowski's avatar
Thomas Witkowski committed
582
    /// Meshes of this problem.
583
    std::vector<Mesh*> meshes;
584

Thomas Witkowski's avatar
Thomas Witkowski committed
585
    /// Pointer to the fe spaces for the different problem components
586
    std::vector<const FiniteElemSpace*> componentSpaces;
587

Thomas Witkowski's avatar
Thomas Witkowski committed
588
    /// Pointer to the meshes for the different problem components
589
    std::vector<Mesh*> componentMeshes;
590

591
592
593
    /// Stores information about which meshes must be traversed to assemble the
    /// specific components. I.e., it was implemented to make use of different
    /// meshes for different components.
Thomas Witkowski's avatar
Thomas Witkowski committed
594
595
596
    ComponentTraverseInfo traverseInfo;
    
    /// Responsible for element marking.
597
    std::vector<Marker*> marker;
598

Thomas Witkowski's avatar
Thomas Witkowski committed
599
    /// Estimator of this problem. Used in \ref estimate().
600
    std::vector<Estimator*> estimator;
601

Thomas Witkowski's avatar
Thomas Witkowski committed
602
    /// Linear solver of this problem. Used in \ref solve().
603
    LinearSolverInterface *solver;
604

Thomas Witkowski's avatar
Thomas Witkowski committed
605
606
    /// System vector  storing the calculated solution of the problem.
    SystemVector *solution;
607

Thomas Witkowski's avatar
Thomas Witkowski committed
608
609
    /// System vector for the right hand side 
    SystemVector *rhs;
610

Thomas Witkowski's avatar
Thomas Witkowski committed
611
612
    /// System matrix
    Matrix<DOFMatrix*> *systemMatrix;
613

614
615
616
    /// Composed system matrix
    SolverMatrix<Matrix<DOFMatrix*> > solverMatrix;

617
618
619
620
621
    /// Some DOFMatrices of the systemMatrix may be assembled only once (for
    /// example if they are independent of the time or older solutions). If
    /// [i][j] of this field is set to true, the corresponding DOFMatrix will
    /// be assembled only once. All other matrices will be assembled at every
    /// time step.
622
    std::vector<std::vector<bool> > assembleMatrixOnlyOnce;
623

624
625
626
627
    /// If [i][j] of this field is set to true, the corresponding DOFMatrix of
    /// the systemMatrix has been assembled at least once. This field is used
    /// to determine, if assembling of a matrix can be ommitted, if it is set
    /// to be assembled only once.
628
    std::vector<std::vector<bool> > assembledMatrix;
629

Thomas Witkowski's avatar
Thomas Witkowski committed
630
631
    /// Determines whether domain boundaries should be considered at assembling.
    bool useGetBound;
632

Thomas Witkowski's avatar
Thomas Witkowski committed
633
    /// Writes the meshes and solution after the adaption loop.
634
    std::vector<FileWriterInterface*> fileWriters;
635

636
637
638
    /// All actions of mesh refinement are performed by refinementManager.
    /// If new refinement algorithms should be realized, one has to override
    /// RefinementManager and give one instance of it to AdaptStationary.
Thomas Witkowski's avatar
Thomas Witkowski committed
639
    RefinementManager *refinementManager;
640

641
642
643
    /// All actions of mesh coarsening are performed by coarseningManager.
    /// If new coarsening algorithms should be realized, one has to override
    /// CoarseningManager and give one instance of it to AdaptStationary.
Thomas Witkowski's avatar
Thomas Witkowski committed
644
    CoarseningManager *coarseningManager;
645
  
Thomas Witkowski's avatar
Thomas Witkowski committed
646
647
    /// Info level.
    int info;
648

649
650
    /// If true, the stationary problem was deserialized from some serialization
    /// file.
651
652
    bool deserialized;

653
654
655
    /// This vectors stores pointers to functions defining the exact solution of
    /// the problem. This may be used to compute the real error of the computed
    /// solution.
656
    std::vector<AbstractFunction<double, WorldVector<double> >*> exactSolutionFcts;
657

658
659
    /// If true, the error is not estimated but computed from the exact solution
    /// defined by \ref exactSolutionFcts.
660
    bool computeExactError;
661
    
662
663
664
665
    /// If at least on boundary condition is set, this variable is true. It is 
    /// used to ensure that no operators are added after boundary condition were
    /// set. If this would happen, boundary conditions could set wrong on off 
    /// diagonal matrices.
666
    bool boundaryConditionSet;
667

668
669
    /// If true, AMDiS prints information about the assembling procedure to 
    /// the screen.
670
    bool writeAsmInfo;
671

672
    std::map<Operator*, std::vector<OperatorPos> > operators;
673
    
674
    /// time needed to solve the linear system
675
    double solutionTime;
676
677
    
    /// time needed to assemble the linear system
678
    double buildTime;
679
680
    
    template <class> friend class detail::CouplingProblemStat;
681
  };
682
683
  
  namespace detail
684
  {
685
686
687
    template <class ProblemStatType>
    struct ProblemStat : public ProblemStatType, 
			 public StandardProblemIteration
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
722
723
724
725
726
727
728
729
730
731
      using ProblemStatType::getName;
      
      /// Constructor
      ProblemStat(std::string nameStr,
		  ProblemIterationInterface *problemIteration = NULL)
	: ProblemStatType(nameStr, problemIteration),
	  StandardProblemIteration(this)
      { }

      /// Returns number of managed problems
      // implements StandardProblemIteration::getNumProblems()
      virtual int getNumProblems() override
      { 
	return 1; 
      }

      /// Returns the problem with the given number. If only one problem
      /// is managed by this master problem, the number hasn't to be given.
      // implements StandardProblemIteration::getProblem(int)
      virtual ProblemStatBase *getProblem(int number = 0) override
      { 
	return this; 
      }

      /// Determines the execution order of the single adaption steps. If adapt is
      /// true, mesh adaption will be performed. This allows to avoid mesh adaption,
      /// e.g. in timestep adaption loops of timestep adaptive strategies.
      // implements StandardProblemIteration::oneIteration(AdaptInfo*, Flag)
      virtual Flag oneIteration(AdaptInfo *adaptInfo, Flag toDo = FULL_ITERATION) override
      {	
	for (int i = 0; i < ProblemStatType::getNumComponents(); i++)
	  if (adaptInfo->spaceToleranceReached(i))
	    adaptInfo->allowRefinement(false, i);
	  else
	    adaptInfo->allowRefinement(true, i);	

	return StandardProblemIteration::oneIteration(adaptInfo, toDo);
      }
    };
  }
  
#ifndef HAVE_PARALLEL_DOMAIN_AMDIS
  typedef detail::ProblemStat<ProblemStatSeq> ProblemStat;
732
#endif
733
734
735
}

#endif