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
151
    /// Returns nr of components \ref nComponents
    virtual int getNumComponents()
152
    { 
153
      return nComponents; 
Thomas Witkowski's avatar
Thomas Witkowski committed
154
    }
155
    
156
    /// Returns nr of additional components \ref nAddComponents
157
158
159
160
    virtual int getNumAddComponents()
    {
      return nAddComponents;
    }
161

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

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

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

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

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

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

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

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

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

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

202
203
204
205
206
207
    /// 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.
208
209
210
211
    virtual void addRobinBC(BoundaryType type, int row, int col, 
			    AbstractFunction<double, WorldVector<double> > *n,
			    AbstractFunction<double, WorldVector<double> > *r);

212
213
214
215
216
217
218
219
220
221
    /// 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);

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

225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
    /// 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);
    }

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

Thomas Witkowski's avatar
Thomas Witkowski committed
251
252

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

      componentMeshes[comp] = mesh;
    }

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

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

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

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

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

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

Thomas Witkowski's avatar
Thomas Witkowski committed
555
  protected:
556
557
    /// 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
558
559
    void computeError(AdaptInfo *adaptInfo);

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

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

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

575
576
577
    /// 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
578
579
580
    int nMeshes;

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

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

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

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

592
593
594
    /// 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
595
596
597
    ComponentTraverseInfo traverseInfo;
    
    /// Responsible for element marking.
598
    std::vector<Marker*> marker;
599

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

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

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

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

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

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

618
619
620
621
622
    /// 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.
623
    std::vector<std::vector<bool> > assembleMatrixOnlyOnce;
624

625
626
627
628
    /// 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.
629
    std::vector<std::vector<bool> > assembledMatrix;
630

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

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

637
638
639
    /// 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
640
    RefinementManager *refinementManager;
641

642
643
644
    /// 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
645
    CoarseningManager *coarseningManager;
646
  
Thomas Witkowski's avatar
Thomas Witkowski committed
647
648
    /// Info level.
    int info;
649

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

654
655
656
    /// 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.
657
    std::vector<AbstractFunction<double, WorldVector<double> >*> exactSolutionFcts;
658

659
660
    /// If true, the error is not estimated but computed from the exact solution
    /// defined by \ref exactSolutionFcts.
661
    bool computeExactError;
662
    
663
664
665
666
    /// 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.
667
    bool boundaryConditionSet;
668

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

673
    std::map<Operator*, std::vector<OperatorPos> > operators;
674
    
675
    /// time needed to solve the linear system
676
    double solutionTime;
677
678
    
    /// time needed to assemble the linear system
679
    double buildTime;
680
681
    
    template <class> friend class detail::CouplingProblemStat;
682
  };
683
684
  
  namespace detail
685
  {
686
687
688
    template <class ProblemStatType>
    struct ProblemStat : public ProblemStatType, 
			 public StandardProblemIteration
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
732
      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;
733
#endif
734
735
736
}

#endif