AdaptInfo.h 18.1 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
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
// ==                                                                        ==
// ============================================================================

/** \file AdaptInfo.h */

#ifndef AMDIS_ADAPTINFO_H
#define AMDIS_ADAPTINFO_H

#include "MatrixVector.h"
#include "Parameters.h"
#include "Serializable.h"

namespace AMDiS {

  /**
   * \ingroup Adaption
   * 
   * \brief
   * Holds adapt parameters and infos about the problem. Base class
   * for AdaptInfoScal and AdaptInfoVec.
   */
  class AdaptInfo : public Serializable
  {
  protected:
    /** \brief
     * Stores adapt infos for a scalar problem or for one component of a 
     * vector valued problem.
     */
    class ScalContent {
    public:
47
      /// Constructor.
Thomas Witkowski's avatar
Thomas Witkowski committed
48
      ScalContent(std::string prefix) 
49
50
51
52
53
54
55
56
57
58
	: est_sum(0.0),
	  est_t_sum(0.0),
	  est_max(0.0),
	  est_t_max(0.0),
	  spaceTolerance(1.0),
	  timeTolerance(1.0),
	  timeErrLow(1.0),
	  coarsenAllowed(0),
	  refinementAllowed(1),
	  refineBisections(1),
59
  	  coarseBisections(1)
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
      {
	double totalTol = 1.0;
	double relSpaceErr = 1.0;
	double relTimeErr = 0.5;
	double timeTheta1 = 1.0;
	double timeTheta2 = 0.3;

	GET_PARAMETER(0, prefix + "->tolerance", "%f", &totalTol);
	GET_PARAMETER(0, prefix + "->rel space error", "%f", &relSpaceErr);
	GET_PARAMETER(0, prefix + "->rel time error", "%f", &relTimeErr);
	GET_PARAMETER(0, prefix + "->time theta 1", "%f", &timeTheta1);
	GET_PARAMETER(0, prefix + "->time theta 2", "%f", &timeTheta2);
	GET_PARAMETER(0, prefix + "->coarsen allowed", "%d", &coarsenAllowed);
	GET_PARAMETER(0, prefix + "->refinement allowed", "%d", &refinementAllowed);
	GET_PARAMETER(0, prefix + "->refine bisections", "%d", &refineBisections);
	GET_PARAMETER(0, prefix + "->coarsen bisections", "%d", &coarseBisections);

	spaceTolerance = totalTol * relSpaceErr;
	timeTolerance = totalTol * relTimeErr * timeTheta1;
	timeErrLow = totalTol * relTimeErr * timeTheta2;
80
      }
81

82
      /// Sum of all error estimates
83
84
      double est_sum;

85
      /// Sum of all time error estimates
86
87
      double est_t_sum;

88
      /// maximal local error estimate.
89
90
      double est_max;

91
      /// Maximum of all time error estimates
92
93
      double est_t_max;

94
      /// Tolerance for the (absolute or relative) error
95
96
      double spaceTolerance;

97
      /// Time tolerance.
98
99
      double timeTolerance;

100
      /// Lower bound for the time error.
101
102
      double timeErrLow;

103
      /// true if coarsening is allowed, false otherwise.
104
105
      int coarsenAllowed;

106
      /// true if refinement is allowed, false otherwise.
107
108
109
110
111
112
113
      int refinementAllowed;

      /** \brief
       * parameter to tell the marking strategy how many bisections should be 
       * performed when an element is marked for refinement; usually the value is
       * 1 or DIM
       */
114
      int refineBisections;
115
116
117
118
119
120

      /** \brief
       * parameter to tell the marking strategy how many bisections should
       * be undone when an element is marked for coarsening; usually the value is 
       * 1 or DIM
       */                          
121
      int coarseBisections;    
122
123
124
    };

  public:
125
    /// Constructor.
Thomas Witkowski's avatar
Thomas Witkowski committed
126
    AdaptInfo(std::string name_, int size = 1) 
127
128
129
130
131
132
133
134
135
136
137
138
139
140
      : name(name_), 
	spaceIteration(-1),
	maxSpaceIteration(-1),
	timestepIteration(0),
	maxTimestepIteration(30),
	timeIteration(0),
	maxTimeIteration(30),
	time(0.0),
	startTime(0.0),
	endTime(1.0),
	timestep(0.0),
	minTimestep(0.0),
	maxTimestep(1.0),
	timestepNumber(0),
Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
141
	nTimesteps(0),
142
143
144
145
	solverIterations(0),
	maxSolverIterations(0),
	solverTolerance(1e-8),
	solverResidual(0.0),
146
147
        scalContents(size),
	isDeserialized_(false)
148
149
150
151
152
153
154
155
156
157
158
159
    {
      GET_PARAMETER(0, name_ + "->start time", "%f", &startTime);
      time = startTime;
      GET_PARAMETER(0, name_ + "->timestep", "%f", &timestep);
      GET_PARAMETER(0, name_ + "->end time", "%f", &endTime);
      GET_PARAMETER(0, name_ + "->max iteration", "%d", &maxSpaceIteration);
      GET_PARAMETER(0, name_ + "->max timestep iteration", "%d", &maxTimestepIteration);
      GET_PARAMETER(0, name_ + "->max time iteration", "%d", &maxTimeIteration);

      GET_PARAMETER(0, name_ + "->min timestep", "%f", &minTimestep);
      GET_PARAMETER(0, name_ + "->max timestep", "%f", &maxTimestep);

Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
160
161
      GET_PARAMETER(0, name_ + "->number of timesteps", "%d", &nTimesteps);

162
163
164
165
166
167
      if (size == 1) {
	scalContents[0] = new ScalContent(name); 
      } else {
	char number[5];
	for (int i = 0; i < size; i++) {
	  sprintf(number, "[%d]", i);
Thomas Witkowski's avatar
Thomas Witkowski committed
168
	  scalContents[i] = new ScalContent(name + std::string(number)); 
169
170
	}
      }
171
    }
172

173
    /// Destructor.
Thomas Witkowski's avatar
Thomas Witkowski committed
174
175
176
    virtual ~AdaptInfo() 
    {
      for (unsigned int i = 0;  i < scalContents.size(); i++)
177
	delete scalContents[i];
178
    }
179
180
181
182
183
184
185
186
187
188
189

    inline void reset() 
    {
      spaceIteration = -1;
      timestepIteration = 0;
      timeIteration = 0;
      time = 0.0;
      timestep = 0.0;
      timestepNumber = 0;
      solverIterations = 0;
      solverResidual = 0.0;
190
191

      GET_PARAMETER(0, name + "->timestep", "%f", &timestep);
192
    }
193

194
    /// Returns whether space tolerance is reached.
195
196
    virtual bool spaceToleranceReached() 
    {
Thomas Witkowski's avatar
Thomas Witkowski committed
197
198
199
200
      for (unsigned int i = 0; i < scalContents.size(); i++) {
	std::cout << "est_sum:" <<scalContents[i]->est_sum 
		  << " spaceTol: " << scalContents[i]->spaceTolerance 
		  << std::endl;
201
	if (!(scalContents[i]->est_sum < scalContents[i]->spaceTolerance))
202
	  return false;
203
      }
204

205
      return true;
206
    }
207

208
    /// Returns whether space tolerance of component i is reached.
209
210
    virtual bool spaceToleranceReached(int i) 
    {
211
      if (!(scalContents[i]->est_sum < scalContents[i]->spaceTolerance))
212
	return false;
213
      else
214
	return true;
215
    }
216

217
    /// Returns whether time tolerance is reached.
218
219
    virtual bool timeToleranceReached() 
    {
Thomas Witkowski's avatar
Thomas Witkowski committed
220
      for (unsigned int i = 0; i < scalContents.size(); i++)
221
	if (!(scalContents[i]->est_t_sum < scalContents[i]->timeTolerance))
222
	  return false;
223

224
      return true;
225
    }
226

227
    /// Returns whether time tolerance of component i is reached.
228
229
    virtual bool timeToleranceReached(int i) 
    {
230
      if (!(scalContents[i]->est_t_sum < scalContents[i]->timeTolerance))
231
	return false;
232
      else
233
	return true;
234
    }
235

236
    /// Returns whether time error is under its lower bound.
237
238
    virtual bool timeErrorLow() 
    {
Thomas Witkowski's avatar
Thomas Witkowski committed
239
      for (unsigned int i = 0; i < scalContents.size(); i++)
240
	if (!(scalContents[i]->est_t_sum < scalContents[i]->timeErrLow))
241
	  return false;
242

243
      return true;
244
245
    }

246
    /// Print debug information about time error and its bound.
247
248
    void printTimeErrorLowInfo() 
    {
Thomas Witkowski's avatar
Thomas Witkowski committed
249
      for (unsigned int i = 0; i < scalContents.size(); i++)
250
251
252
	std::cout << "    Time error estimate = " << scalContents[i]->est_t_sum 
		  << "    Time error bound = " << scalContents[i]->timeErrLow << "\n";
    }
253

254
    /// Returns \ref spaceIteration.
255
256
    inline int getSpaceIteration() 
    { 
257
      return spaceIteration; 
258
    }
259

260
    /// Sets \ref spaceIteration.
261
262
    inline void setSpaceIteration(int it) 
    { 
263
      spaceIteration = it; 
264
    }
265
  
266
    /// Returns \ref maxSpaceIteration.
267
268
    inline int getMaxSpaceIteration() 
    { 
269
      return maxSpaceIteration;
270
    }
271

272
    /// Sets \ref maxSpaceIteration.
273
274
    inline void setMaxSpaceIteration(int it) 
    { 
275
      maxSpaceIteration = it; 
276
    }
277
  
278
    /// Increments \ref spaceIteration by 1;
279
280
    inline void incSpaceIteration() 
    { 
281
      spaceIteration++; 
282
    }
283

284
    /// Sets \ref timestepIteration.
285
286
    inline void setTimestepIteration(int it) 
    { 
287
      timestepIteration = it; 
288
    }
289
  
290
    /// Returns \ref timestepIteration.
291
292
    inline int getTimestepIteration() 
    { 
293
      return timestepIteration; 
294
    }
295

296
    /// Increments \ref timestepIteration by 1;
297
298
    inline void incTimestepIteration() 
    { 
299
      timestepIteration++; 
300
    }
301

302
    /// Returns \ref maxTimestepIteration.
303
304
    inline int getMaxTimestepIteration() 
    { 
305
      return maxTimestepIteration; 
306
    }
307

308
    /// Sets \ref maxTimestepIteration.
309
310
    inline void setMaxTimestepIteration(int it) 
    { 
311
      maxTimestepIteration = it; 
312
    }
313
  
314
    /// Sets \ref timeIteration.
315
316
    inline void setTimeIteration(int it) 
    { 
317
      timeIteration = it; 
318
    }
319
  
320
    /// Returns \ref timeIteration.
321
322
    inline int getTimeIteration() 
    { 
323
      return timeIteration; 
324
    }
325

326
    /// Increments \ref timesIteration by 1;
327
328
    inline void incTimeIteration() 
    { 
329
      timeIteration++; 
330
    }
331

332
    /// Returns \ref maxTimeIteration.
333
334
    inline int getMaxTimeIteration() 
    { 
335
      return maxTimeIteration; 
336
    }
337

338
    /// Sets \ref maxTimeIteration.
339
340
    inline void setMaxTimeIteration(int it) 
    { 
341
      maxTimeIteration = it; 
342
    }
343
  
344
    /// Returns \ref timestepNumber.
345
346
    inline int getTimestepNumber() 
    { 
347
      return timestepNumber; 
348
    }
349

Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
350
    /// Returns \ref nTimesteps.
351
352
    inline int getNumberOfTimesteps() 
    {
Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
353
354
355
      return nTimesteps;
    }

356
    /// Increments \ref timestepNumber by 1;
357
358
    inline void incTimestepNumber() 
    { 
359
      timestepNumber++; 
360
    }
361

362
    /// Sets \ref est_sum.
363
364
    inline void setEstSum(double e, int index) 
    {
365
      scalContents[index]->est_sum = e;
366
    }
367

368
    /// Sets \ref est_max.
369
370
    inline void setEstMax(double e, int index) 
    {
371
      scalContents[index]->est_max = e;
372
    }
373

374
    /// Sets \ref est_max.
375
376
    inline void setTimeEstMax(double e, int index) 
    {
377
      scalContents[index]->est_t_max = e;
378
    }
379

380
    /// Sets \ref est_t_sum.
381
382
    inline void setTimeEstSum(double e, int index) 
    {
383
      scalContents[index]->est_t_sum = e;
384
    }
385

386
    /// Returns \ref est_sum.
387
388
    inline double getEstSum(int index) 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
389
390
391
392
393
      FUNCNAME("AdaptInfo::getEstSum()");

      TEST_EXIT_DBG(static_cast<unsigned int>(index) < scalContents.size())
	("Wrong index for adaptInfo!\n");

394
      return scalContents[index]->est_sum; 
395
    }
396

397
    /// Returns \ref est_t_sum.
398
399
    inline double getEstTSum(int index) 
    { 
400
      return scalContents[index]->est_t_sum; 
401
    }
402

403
    /// Returns \ref est_max.
404
405
    inline double getEstMax(int index) 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
406
407
408
409
410
      FUNCNAME("AdaptInfo::getEstSum()");

      TEST_EXIT_DBG(static_cast<unsigned int>(index) < scalContents.size())
	("Wrong index for adaptInfo!\n");

411
      return scalContents[index]->est_max; 
412
    }
413

414
    /// Returns \ref est_max.
415
416
    inline double getTimeEstMax(int index) 
    { 
417
      return scalContents[index]->est_t_max; 
418
    }
419

420
    /// Returns \ref est_t_sum.
421
422
    inline double getTimeEstSum(int index) 
    { 
423
      return scalContents[index]->est_t_sum; 
424
    }
425

426
    /// Returns \ref spaceTolerance.
427
428
    inline double getSpaceTolerance(int index) 
    { 
429
      return scalContents[index]->spaceTolerance; 
430
    }  
431

432
    /// Sets \ref spaceTolerance.
433
434
    inline void setSpaceTolerance(int index, double tol) 
    { 
435
      scalContents[index]->spaceTolerance = tol; 
436
    }  
437

438
    /// Returns \ref timeTolerance.
439
440
    inline double getTimeTolerance(int index) 
    { 
441
      return scalContents[index]->timeTolerance; 
442
    }  
443

444
    /// Sets \ref time
445
446
    inline double setTime(double t) 
    { 
447
      time = t; 
448
449
450
451
452
      if (time > endTime) 
	time = endTime;
      if (time < startTime) 
	time = startTime;

453
      return time;
454
    }
455

456
    /// Gets \ref time
457
458
    inline double getTime() 
    { 
459
      return time; 
460
    }  
461

462
    /// Gets \ref &time
463
464
    inline double *getTimePtr() 
    { 
465
      return &time; 
466
    }  
467

468
    /// Sets \ref timestep
469
470
    inline double setTimestep(double t) 
    { 
471
      timestep = t; 
472
      if (timestep > maxTimestep)
473
	timestep = maxTimestep;
474
      if (timestep < minTimestep)
475
	timestep = minTimestep;
476
      if (time + timestep > endTime)
477
478
	timestep = endTime - time;
      
479
      return timestep;
480
    }
481

Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
482
483
484
485
    /** \brief
     * Returns true, if the end time is reached and no more timestep
     * computations must be done.
     */
486
487
    inline bool reachedEndTime() 
    {
488
489
490
491
      if (nTimesteps > 0) 
	return !(timestepNumber < nTimesteps);

      return !(time < endTime - DBL_TOL);
Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
492
493
    }

494
    /// Gets \ref timestep
495
496
    inline double getTimestep() 
    { 
497
      return timestep; 
498
    }
499

500
    /// Sets \ref minTimestep
501
502
    inline void setMinTimestep(double t) 
    { 
503
      minTimestep = t; 
504
    }
505

506
    /// Gets \ref minTimestep
507
508
    inline double getMinTimestep() 
    { 
509
      return minTimestep; 
510
    }  
511

512
    /// Sets \ref maxTimestep
513
514
    inline void setMaxTimestep(double t) 
    { 
515
      maxTimestep = t; 
516
    }
517

518
    /// Gets \ref maxTimestep
519
520
    inline double getMaxTimestep() 
    { 
521
      return maxTimestep; 
522
    }  
523

524
    /// Gets \ref &timestep
525
526
    inline double *getTimestepPtr() 
    { 
527
      return &timestep; 
528
    }  
529

530
    /// Sets \ref startTime = time
531
532
    inline void setStartTime(double time) 
    { 
533
      startTime = time; 
534
    }
535

536
    /// Sets \ref endTime = time
537
538
    inline void setEndTime(double time) 
    { 
539
      endTime = time; 
540
    }
541

542
    /// Returns \ref startTime
543
544
    inline double getStartTime() 
    { 
545
      return startTime; 
546
    }
547

548
    /// Returns \ref endTime
549
550
    inline double getEndTime() 
    { 
551
      return endTime; 
552
    }
553

554
    /// Returns \ref timeErrLow.
555
556
    inline double getTimeErrLow(int index) 
    { 
557
      return scalContents[index]->timeErrLow; 
558
    }  
559

560
    /// Returns whether coarsening is allowed or not.
561
562
    inline bool isCoarseningAllowed(int index) 
    {
563
      return (scalContents[index]->coarsenAllowed == 1);
564
    }
565

566
    /// Returns whether coarsening is allowed or not.
567
568
    inline bool isRefinementAllowed(int index) 
    {
569
      return (scalContents[index]->refinementAllowed == 1);
570
    }
571

572
    ///
573
574
    inline void allowRefinement(bool allow, int index) 
    {
575
      scalContents[index]->refinementAllowed = allow;
576
    }
577

578
    ///
579
580
    inline void allowCoarsening(bool allow, int index) 
    {
581
      scalContents[index]->coarsenAllowed = allow;
582
    }
583

584
    /// Returns \ref refineBisections
585
586
    inline const int getRefineBisections(int index) const 
    {
587
      return scalContents[index]->refineBisections;
588
    }
589

590
    /// Returns \ref coarseBisections
591
592
    inline const int getCoarseBisections(int index) const 
    {
593
      return scalContents[index]->coarseBisections;
594
    }
595

596
597
    inline int getSize() 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
598
      return scalContents.size(); 
599
    }
600

601
602
    inline void setSolverIterations(int it) 
    {
603
      solverIterations = it;
604
    }
605
  
606
607
    inline int getSolverIterations() 
    {
608
      return solverIterations;
609
    }
610
  
611
612
    inline void setMaxSolverIterations(int it) 
    {
613
      maxSolverIterations = it;
614
    }
615
  
616
617
    inline int getMaxSolverIterations() 
    {
618
      return maxSolverIterations;
619
    }
620
  
621
622
    inline void setSolverTolerance(double tol) 
    {
623
      solverTolerance = tol;
624
    }
625
  
626
627
    inline double getSolverTolerance() 
    {
628
      return solverTolerance;
629
    }
630
  
631
632
    inline void setSolverResidual(double res) 
    {
633
      solverResidual = res;
634
    }
635
  
636
637
    inline double getSolverResidual() 
    {
638
      return solverResidual;
639
    }
640

641
    /// Returns true, if the adaptive procedure was deserialized from a file.
642
643
    const bool isDeserialized() const 
    {
644
645
646
      return isDeserialized_;
    }

647
648
    inline void setIsDeserialized(bool b) 
    {
649
650
651
      isDeserialized_ = b;
    }

652
    /// Creates new scalContents with the given size.
653
654
    void setScalContents(int newSize);

655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
    /** \brief
     * Resets timestep, current time and time boundaries without
     * any check. Is used by the parareal algorithm.
     */
    void resetTimeValues(double newTimeStep,
			 double newStartTime,
			 double newEndTime)
    {
      time = newStartTime;
      startTime = newStartTime;
      endTime = newEndTime;
      timestep = newTimeStep;
      timestepNumber = 0;
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
670
    void serialize(std::ostream& out);
671

Thomas Witkowski's avatar
Thomas Witkowski committed
672
    void deserialize(std::istream& in);
673
674

  protected:
675
    /// Name.
Thomas Witkowski's avatar
Thomas Witkowski committed
676
    std::string name;
677

678
    /// Current space iteration
679
680
681
682
683
684
685
686
    int spaceIteration;

    /** \brief
     * maximal allowed number of iterations of the adaptive procedure; if 
     * maxIteration <= 0, no iteration bound is used
     */
    int maxSpaceIteration;

687
    /// Current timestep iteration
688
689
    int timestepIteration;

690
    /// Maximal number of iterations for choosing a timestep
691
692
    int maxTimestepIteration;

693
    /// Current time iteration
694
695
    int timeIteration;

696
    /// Maximal number of time iterations
697
698
    int maxTimeIteration;

699
    /// Actual time, end of time interval for current time step
700
701
    double time;

702
    /// Initial time
703
704
    double startTime;

705
    /// Final time
706
707
    double endTime;

708
    /// Current time step size
709
710
    double timestep;

711
    /// Minimal step size
712
713
    double minTimestep;

714
    /// Maximal step size
715
716
    double maxTimestep;

717
    /// Number of current time step
718
    int timestepNumber;
Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
719
720
721
722
723
724
725

    /** \brief
     * Per default this value is 0 and not used. If it is set to a non-zero value,
     * the computation of the stationary problem is done nTimesteps times with a
     * fixed timestep.
     */
    int nTimesteps;
726
  
727
    /// number of iterations needed of linear or nonlinear solver
728
729
    int solverIterations;

730
    /// maximal number of iterations needed of linear or nonlinear solver
731
732
    int maxSolverIterations;

733
    ///
734
735
    double solverTolerance;

736
    ///
737
738
    double solverResidual;

739
    /// Scalar adapt infos.
Thomas Witkowski's avatar
Thomas Witkowski committed
740
    std::vector<ScalContent*> scalContents;
741

742
    /// Is true, if the adaptive procedure was deserialized from a file.
743
    bool isDeserialized_;
744
745
746
747
748
  };

}

#endif //  AMDIS_ADAPTINFO_H