AdaptInfo.h 18.5 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
      : 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),
138
	lastProcessedTimestep(0.0),
139
140
141
	minTimestep(0.0),
	maxTimestep(1.0),
	timestepNumber(0),
Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
142
	nTimesteps(0),
143
144
145
146
	solverIterations(0),
	maxSolverIterations(0),
	solverTolerance(1e-8),
	solverResidual(0.0),
147
148
        scalContents(size),
	isDeserialized_(false)
149
150
151
152
153
154
155
156
157
158
159
160
    {
      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
161
162
      GET_PARAMETER(0, name_ + "->number of timesteps", "%d", &nTimesteps);

163
164
165
166
167
168
      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
169
	  scalContents[i] = new ScalContent(name + std::string(number)); 
170
171
	}
      }
172
    }
173

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

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

      GET_PARAMETER(0, name + "->timestep", "%f", &timestep);
193
      lastProcessedTimestep=timestep;
194
    }
195

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

207
      return true;
208
    }
209

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

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

226
      return true;
227
    }
228

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

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

245
      return true;
246
247
    }

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

257
    /// Returns \ref spaceIteration.
258
259
    inline int getSpaceIteration() 
    { 
260
      return spaceIteration; 
261
    }
262

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

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

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

299
    /// Increments \ref timestepIteration by 1;
300
301
    inline void incTimestepIteration() 
    { 
302
      timestepIteration++; 
303
    }
304

305
    /// Returns \ref maxTimestepIteration.
306
307
    inline int getMaxTimestepIteration() 
    { 
308
      return maxTimestepIteration; 
309
    }
310

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

329
    /// Increments \ref timesIteration by 1;
330
331
    inline void incTimeIteration() 
    { 
332
      timeIteration++; 
333
    }
334

335
    /// Returns \ref maxTimeIteration.
336
337
    inline int getMaxTimeIteration() 
    { 
338
      return maxTimeIteration; 
339
    }
340

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

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

359
    /// Increments \ref timestepNumber by 1;
360
361
    inline void incTimestepNumber() 
    { 
362
      timestepNumber++; 
363
    }
364

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

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

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

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

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

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

397
      return scalContents[index]->est_sum; 
398
    }
399

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

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

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

414
      return scalContents[index]->est_max; 
415
    }
416

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

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

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

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

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

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

456
      return time;
457
    }
458

459
    /// Gets \ref time
460
461
    inline double getTime() 
    { 
462
      return time; 
463
    }  
464

465
    /// Gets \ref &time
466
467
    inline double *getTimePtr() 
    { 
468
      return &time; 
469
    }  
470

471
    /// Sets \ref timestep
472
473
    inline double setTimestep(double t) 
    { 
474
      timestep = t; 
475
      if (timestep > maxTimestep)
476
	timestep = maxTimestep;
477
      if (timestep < minTimestep)
478
	timestep = minTimestep;
479
      if (time + timestep > endTime)
480
481
	timestep = endTime - time;
      
482
      return timestep;
483
    }
484
485
486
487
488
489
490
491
492
493
494
495
496
    /// Gets \ref timestep
    inline double getTimestep() 
    { 
      return timestep; 
    }

    inline void setLastProcessedTimestep(double t){
	lastProcessedTimestep=t;
    } 

    inline double getLastProcessedTimestep(){
	return lastProcessedTimestep;
    } 
497

Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
498
499
500
501
    /** \brief
     * Returns true, if the end time is reached and no more timestep
     * computations must be done.
     */
502
503
    inline bool reachedEndTime() 
    {
504
505
506
507
      if (nTimesteps > 0) 
	return !(timestepNumber < nTimesteps);

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

510

511
    /// Sets \ref minTimestep
512
513
    inline void setMinTimestep(double t) 
    { 
514
      minTimestep = t; 
515
    }
516

517
    /// Gets \ref minTimestep
518
519
    inline double getMinTimestep() 
    { 
520
      return minTimestep; 
521
    }  
522

523
    /// Sets \ref maxTimestep
524
525
    inline void setMaxTimestep(double t) 
    { 
526
      maxTimestep = t; 
527
    }
528

529
    /// Gets \ref maxTimestep
530
531
    inline double getMaxTimestep() 
    { 
532
      return maxTimestep; 
533
    }  
534

535
    /// Gets \ref &timestep
536
537
    inline double *getTimestepPtr() 
    { 
538
      return &timestep; 
539
    }  
540

541
    /// Sets \ref startTime = time
542
543
    inline void setStartTime(double time) 
    { 
544
      startTime = time; 
545
    }
546

547
    /// Sets \ref endTime = time
548
549
    inline void setEndTime(double time) 
    { 
550
      endTime = time; 
551
    }
552

553
    /// Returns \ref startTime
554
555
    inline double getStartTime() 
    { 
556
      return startTime; 
557
    }
558

559
    /// Returns \ref endTime
560
561
    inline double getEndTime() 
    { 
562
      return endTime; 
563
    }
564

565
    /// Returns \ref timeErrLow.
566
567
    inline double getTimeErrLow(int index) 
    { 
568
      return scalContents[index]->timeErrLow; 
569
    }  
570

571
    /// Returns whether coarsening is allowed or not.
572
573
    inline bool isCoarseningAllowed(int index) 
    {
574
      return (scalContents[index]->coarsenAllowed == 1);
575
    }
576

577
    /// Returns whether coarsening is allowed or not.
578
579
    inline bool isRefinementAllowed(int index) 
    {
580
      return (scalContents[index]->refinementAllowed == 1);
581
    }
582

583
    ///
584
585
    inline void allowRefinement(bool allow, int index) 
    {
586
      scalContents[index]->refinementAllowed = allow;
587
    }
588

589
    ///
590
591
    inline void allowCoarsening(bool allow, int index) 
    {
592
      scalContents[index]->coarsenAllowed = allow;
593
    }
594

595
    /// Returns \ref refineBisections
596
597
    inline const int getRefineBisections(int index) const 
    {
598
      return scalContents[index]->refineBisections;
599
    }
600

601
    /// Returns \ref coarseBisections
602
603
    inline const int getCoarseBisections(int index) const 
    {
604
      return scalContents[index]->coarseBisections;
605
    }
606

607
608
    inline int getSize() 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
609
      return scalContents.size(); 
610
    }
611

612
613
    inline void setSolverIterations(int it) 
    {
614
      solverIterations = it;
615
    }
616
  
617
618
    inline int getSolverIterations() 
    {
619
      return solverIterations;
620
    }
621
  
622
623
    inline void setMaxSolverIterations(int it) 
    {
624
      maxSolverIterations = it;
625
    }
626
  
627
628
    inline int getMaxSolverIterations() 
    {
629
      return maxSolverIterations;
630
    }
631
  
632
633
    inline void setSolverTolerance(double tol) 
    {
634
      solverTolerance = tol;
635
    }
636
  
637
638
    inline double getSolverTolerance() 
    {
639
      return solverTolerance;
640
    }
641
  
642
643
    inline void setSolverResidual(double res) 
    {
644
      solverResidual = res;
645
    }
646
  
647
648
    inline double getSolverResidual() 
    {
649
      return solverResidual;
650
    }
651

652
    /// Returns true, if the adaptive procedure was deserialized from a file.
653
654
    const bool isDeserialized() const 
    {
655
656
657
      return isDeserialized_;
    }

658
659
    inline void setIsDeserialized(bool b) 
    {
660
661
662
      isDeserialized_ = b;
    }

663
    /// Creates new scalContents with the given size.
664
665
    void setScalContents(int newSize);

666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
    /** \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
681
    void serialize(std::ostream& out);
682

Thomas Witkowski's avatar
Thomas Witkowski committed
683
    void deserialize(std::istream& in);
684
685

  protected:
686
    /// Name.
Thomas Witkowski's avatar
Thomas Witkowski committed
687
    std::string name;
688

689
    /// Current space iteration
690
691
692
693
694
695
696
697
    int spaceIteration;

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

698
    /// Current timestep iteration
699
700
    int timestepIteration;

701
    /// Maximal number of iterations for choosing a timestep
702
703
    int maxTimestepIteration;

704
    /// Current time iteration
705
706
    int timeIteration;

707
    /// Maximal number of time iterations
708
709
    int maxTimeIteration;

710
    /// Actual time, end of time interval for current time step
711
712
    double time;

713
    /// Initial time
714
715
    double startTime;

716
    /// Final time
717
718
    double endTime;

719
    ///Time step size to be used
720
721
    double timestep;

722
723
724
    /// Last processed time step size of finished iteration
    double lastProcessedTimestep;

725
    /// Minimal step size
726
727
    double minTimestep;

728
    /// Maximal step size
729
730
    double maxTimestep;

731
    /// Number of current time step
732
    int timestepNumber;
Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
733
734
735
736
737
738
739

    /** \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;
740
  
741
    /// number of iterations needed of linear or nonlinear solver
742
743
    int solverIterations;

744
    /// maximal number of iterations needed of linear or nonlinear solver
745
746
    int maxSolverIterations;

747
    ///
748
749
    double solverTolerance;

750
    ///
751
752
    double solverResidual;

753
    /// Scalar adapt infos.
Thomas Witkowski's avatar
Thomas Witkowski committed
754
    std::vector<ScalContent*> scalContents;
755

756
    /// Is true, if the adaptive procedure was deserialized from a file.
757
    bool isDeserialized_;
758
759
760
761
762
  };

}

#endif //  AMDIS_ADAPTINFO_H