AdaptInfo.h 18.8 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
24
25
26
27
28

/** \file AdaptInfo.h */

#ifndef AMDIS_ADAPTINFO_H
#define AMDIS_ADAPTINFO_H

#include "MatrixVector.h"
29
#include "Initfile.h"
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
#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:
50
      /// Constructor.
51
      ScalContent(std::string prefix)
52
	: est_sum(0.0),
53
54
55
56
57
58
59
60
61
62
63
64
	  est_t_sum(0.0),
	  est_max(0.0),
	  est_t_max(0.0),
	  fac_max(0.0),
	  fac_sum(1.0),
	  spaceTolerance(0.0),
	  timeTolerance(0.0),
	  timeErrLow(0.0),
	  coarsenAllowed(0),
	  refinementAllowed(1),
	  refineBisections(1),
	  coarseBisections(1)	  	
65
      {
66
67
68
69
70
71
72
73
	Parameters::get(prefix + "->tolerance", spaceTolerance);
	Parameters::get(prefix + "->time tolerance", timeTolerance);
	Parameters::get(prefix + "->coarsen allowed", coarsenAllowed);
	Parameters::get(prefix + "->refinement allowed", refinementAllowed);
	Parameters::get(prefix + "->refine bisections", refineBisections);
	Parameters::get(prefix + "->coarsen bisections", coarseBisections);
	Parameters::get(prefix + "->sum factor", fac_sum);
	Parameters::get(prefix + "->max factor", fac_max);
74

75
	timeErrLow = timeTolerance * 0.3;
76
      }
77

78
      /// Sum of all error estimates
79
80
      double est_sum;

81
      /// Sum of all time error estimates
82
83
      double est_t_sum;

84
      /// maximal local error estimate.
85
86
      double est_max;

87
      /// Maximum of all time error estimates
88
      double est_t_max;
89
90
91
      
      /// factors to combine max and integral time estimate
      double fac_max, fac_sum;
92

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

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

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

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

105
      /// true if refinement is allowed, false otherwise.
106
107
108
109
110
111
112
      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
       */
113
      int refineBisections;
114
115
116
117
118
119

      /** \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
       */                          
120
      int coarseBisections;    
121
122
123
    };

  public:
124
    /// Constructor.
Thomas Witkowski's avatar
Thomas Witkowski committed
125
    AdaptInfo(std::string name_, int size = 1) 
126
127
128
129
130
131
132
133
134
135
136
      : 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),
137
	lastProcessedTimestep(0.0),
138
139
140
	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
        scalContents(size),
147
148
	deserialized(false),
	rosenbrockMode(false)
149
    {
150
      init();
151
152
      char number[5];
      for (int i = 0; i < size; i++) {
153
154
	sprintf(number, "[%d]", i);
	scalContents[i] = new ScalContent(name + std::string(number));  
155
      }
156
    }
157

158
    /// Destructor.
Thomas Witkowski's avatar
Thomas Witkowski committed
159
160
161
    virtual ~AdaptInfo() 
    {
      for (unsigned int i = 0;  i < scalContents.size(); i++)
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
  delete scalContents[i];
    }

    /// Sets initial values to time/timestep variables
    inline void init()
    {
      Parameters::get(name + "->start time", startTime);
      time = startTime;
      Parameters::get(name + "->timestep", timestep);
      Parameters::get(name + "->end time", endTime);
      Parameters::get(name + "->max iteration", maxSpaceIteration);
      Parameters::get(name + "->max timestep iteration", maxTimestepIteration);
      Parameters::get(name + "->max time iteration", maxTimeIteration);
      Parameters::get(name + "->min timestep", minTimestep);
      Parameters::get(name + "->max timestep", maxTimestep);
      Parameters::get(name + "->number of timesteps", nTimesteps);
178
    }
179

180
    /// Resets all variables to zero (or something equivalent)
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
      Parameters::get(name + "->timestep", 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
      for (unsigned int i = 0; i < scalContents.size(); i++) {
200
	if (!(scalContents[i]->est_sum < scalContents[i]->spaceTolerance))
201
	  return false;
202
      }
203

204
      return true;
205
    }
206

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

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

223
      return true;
224
    }
225

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

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

242
      return true;
243
    }
244
245
246
247
    /// Returns the time estimation as a combination 
    /// of maximal and integral time error 
    double getTimeEstCombined(unsigned i) const 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
248
249
250
      return 
	scalContents[i]->est_t_max * scalContents[i]->fac_max +
	scalContents[i]->est_t_sum * scalContents[i]->fac_sum; 
251
252
    }

253

254
    /// Print debug information about time error and its bound.
255
256
    void printTimeErrorLowInfo() 
    {
257
      for (unsigned int i = 0; i < scalContents.size(); i++){
258
	std::cout << "    Time error estimate     ["<<i<<"] = " 
259
		  << getTimeEstCombined(i) << "\n"
260
		  << "    Time error estimate sum ["<<i<<"] = " 
261
		  << scalContents[i]->est_t_sum << "\n" 
262
		  << "    Time error estimate max ["<<i<<"] = " 
263
		  << scalContents[i]->est_t_max << "\n"
264
		  << "    Time error low bound    ["<<i<<"] = " 
265
		  << scalContents[i]->timeErrLow << "\n"
266
		  << "    Time error high bound   ["<<i<<"] = " 
267
		  << scalContents[i]->timeTolerance << "\n";
268
      }
269
    }
270

271
    /// Returns \ref spaceIteration.
272
273
    inline int getSpaceIteration() 
    { 
274
      return spaceIteration; 
275
    }
276

277
    /// Sets \ref spaceIteration.
278
279
    inline void setSpaceIteration(int it) 
    { 
280
      spaceIteration = it; 
281
    }
282
  
283
    /// Returns \ref maxSpaceIteration.
284
285
    inline int getMaxSpaceIteration() 
    { 
286
      return maxSpaceIteration;
287
    }
288

289
    /// Sets \ref maxSpaceIteration.
290
291
    inline void setMaxSpaceIteration(int it) 
    { 
292
      maxSpaceIteration = it; 
293
    }
294
  
295
    /// Increments \ref spaceIteration by 1;
296
297
    inline void incSpaceIteration() 
    { 
298
      spaceIteration++; 
299
    }
300

301
    /// Sets \ref timestepIteration.
302
303
    inline void setTimestepIteration(int it) 
    { 
304
      timestepIteration = it; 
305
    }
306
  
307
    /// Returns \ref timestepIteration.
308
309
    inline int getTimestepIteration() 
    { 
310
      return timestepIteration; 
311
    }
312

313
    /// Increments \ref timestepIteration by 1;
314
315
    inline void incTimestepIteration() 
    { 
316
      timestepIteration++; 
317
    }
318

319
    /// Returns \ref maxTimestepIteration.
320
321
    inline int getMaxTimestepIteration() 
    { 
322
      return maxTimestepIteration; 
323
    }
324

325
    /// Sets \ref maxTimestepIteration.
326
327
    inline void setMaxTimestepIteration(int it) 
    { 
328
      maxTimestepIteration = it; 
329
    }
330
  
331
    /// Sets \ref timeIteration.
332
333
    inline void setTimeIteration(int it) 
    { 
334
      timeIteration = it; 
335
    }
336
  
337
    /// Returns \ref timeIteration.
338
339
    inline int getTimeIteration() 
    { 
340
      return timeIteration; 
341
    }
342

343
    /// Increments \ref timesIteration by 1;
344
345
    inline void incTimeIteration() 
    { 
346
      timeIteration++; 
347
    }
348

349
    /// Returns \ref maxTimeIteration.
350
351
    inline int getMaxTimeIteration() 
    { 
352
      return maxTimeIteration; 
353
    }
354

355
    /// Sets \ref maxTimeIteration.
356
357
    inline void setMaxTimeIteration(int it) 
    { 
358
      maxTimeIteration = it; 
359
    }
360
  
361
    /// Returns \ref timestepNumber.
362
363
    inline int getTimestepNumber() 
    { 
364
      return timestepNumber; 
365
    }
366

367
368
369
    /// Sets \ref timestepNumber.
    inline void setTimestepNumber(int num) 
    {
370
      timestepNumber = std::min(nTimesteps, num);
371
372
    }
    
Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
373
    /// Returns \ref nTimesteps.
374
375
    inline int getNumberOfTimesteps() 
    {
Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
376
377
378
      return nTimesteps;
    }

379
380
381
    /// Sets \ref nTimesteps.
    inline void setNumberOfTimesteps(int num) 
    {
382
      nTimesteps = std::max(0, num);
383
384
    }

385
    /// Increments \ref timestepNumber by 1;
386
387
    inline void incTimestepNumber() 
    { 
388
      timestepNumber++; 
389
    }
390

391
    /// Sets \ref est_sum.
392
393
    inline void setEstSum(double e, int index) 
    {
394
      scalContents[index]->est_sum = e;
395
    }
396

397
    /// Sets \ref est_max.
398
399
    inline void setEstMax(double e, int index) 
    {
400
      scalContents[index]->est_max = e;
401
    }
402

403
    /// Sets \ref est_max.
404
405
    inline void setTimeEstMax(double e, int index) 
    {
406
      scalContents[index]->est_t_max = e;
407
    }
408

409
    /// Sets \ref est_t_sum.
410
411
    inline void setTimeEstSum(double e, int index) 
    {
412
      scalContents[index]->est_t_sum = e;
413
    }
414

415
    /// Returns \ref est_sum.
416
417
    inline double getEstSum(int index) 
    { 
418
      FUNCNAME_DBG("AdaptInfo::getEstSum()");
Thomas Witkowski's avatar
Thomas Witkowski committed
419
420
421
422

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

423
      return scalContents[index]->est_sum; 
424
    }
425

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

432
    /// Returns \ref est_max.
433
434
    inline double getEstMax(int index) 
    { 
435
      FUNCNAME_DBG("AdaptInfo::getEstSum()");
Thomas Witkowski's avatar
Thomas Witkowski committed
436
437
438
439

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

440
      return scalContents[index]->est_max; 
441
    }
442

443
    /// Returns \ref est_max.
444
445
    inline double getTimeEstMax(int index) 
    { 
446
      return scalContents[index]->est_t_max; 
447
    }
448

449
    /// Returns \ref est_t_sum.
450
451
    inline double getTimeEstSum(int index) 
    { 
452
      return scalContents[index]->est_t_sum; 
453
    }
Praetorius, Simon's avatar
Praetorius, Simon committed
454
455
456
457
458
459
460
461
462
463
464
    
    /// Returns \ref est_t the estimated overall time error
    inline double getTimeEst()
    {
      return est_t;
    }
    
    inline void setTimeEst(double value)
    {
      est_t = value;
    }
465

466
    /// Returns \ref spaceTolerance.
467
468
    inline double getSpaceTolerance(int index) 
    { 
469
      return scalContents[index]->spaceTolerance; 
470
    }  
471

472
    /// Sets \ref spaceTolerance.
473
474
    inline void setSpaceTolerance(int index, double tol) 
    { 
475
      scalContents[index]->spaceTolerance = tol; 
476
    }  
477

478
    /// Returns \ref timeTolerance.
479
480
    inline double getTimeTolerance(int index) 
    { 
481
      return scalContents[index]->timeTolerance; 
482
    }  
483

484
    /// Sets \ref time
485
486
    inline double setTime(double t) 
    { 
487
      time = t; 
488
489
490
491
492
      if (time > endTime) 
	time = endTime;
      if (time < startTime) 
	time = startTime;

493
      return time;
494
    }
495

496
    /// Gets \ref time
497
498
    inline double getTime() 
    { 
499
      return time; 
500
    }  
501

502
    /// Gets \ref &time
503
504
    inline double *getTimePtr() 
    { 
505
      return &time; 
506
    }  
507

508
    /// Sets \ref timestep
509
510
    inline double setTimestep(double t) 
    { 
511
      timestep = t; 
512
      if (timestep > maxTimestep)
513
	timestep = maxTimestep;
514
      if (timestep < minTimestep)
515
	timestep = minTimestep;
516
      if (time + timestep > endTime)
517
518
	timestep = endTime - time;
      
519
      return timestep;
520
    }
521
522
523
524
525
526
    /// Gets \ref timestep
    inline double getTimestep() 
    { 
      return timestep; 
    }

527
528
529
    inline void setLastProcessedTimestep(double t)
    {
      lastProcessedTimestep = t;
530
531
    } 

532
533
    inline double getLastProcessedTimestep()
    {
534
535
	return lastProcessedTimestep;
    } 
536

537
538
    /// Returns true, if the end time is reached and no more timestep
    /// computations must be done.
539
540
    inline bool reachedEndTime() 
    {
541
542
543
      if (nTimesteps > 0) 
	return !(timestepNumber < nTimesteps);

544
      return !(fabs(time - endTime) > DBL_TOL);
Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
545
546
    }

547

548
    /// Sets \ref minTimestep
549
550
    inline void setMinTimestep(double t) 
    { 
551
      minTimestep = t; 
552
    }
553

554
    /// Gets \ref minTimestep
555
556
    inline double getMinTimestep() 
    { 
557
      return minTimestep; 
558
    }  
559

560
    /// Sets \ref maxTimestep
561
562
    inline void setMaxTimestep(double t) 
    { 
563
      maxTimestep = t; 
564
    }
565

566
    /// Gets \ref maxTimestep
567
568
    inline double getMaxTimestep() 
    { 
569
      return maxTimestep; 
570
    }  
571
    
572
    /// Gets \ref &timestep
573
574
    inline double *getTimestepPtr() 
    { 
575
      return &timestep; 
576
    }  
577

578
    /// Sets \ref startTime = time
579
580
    inline void setStartTime(double time) 
    { 
581
      startTime = time; 
582
    }
583

584
    /// Sets \ref endTime = time
585
586
    inline void setEndTime(double time) 
    { 
587
      endTime = time; 
588
    }
589

590
    /// Returns \ref startTime
591
592
    inline double getStartTime() 
    { 
593
      return startTime; 
594
    }
595

596
    /// Returns \ref endTime
597
598
    inline double getEndTime() 
    { 
599
      return endTime; 
600
    }
601

602
    /// Returns \ref timeErrLow.
603
604
    inline double getTimeErrLow(int index) 
    { 
605
      return scalContents[index]->timeErrLow; 
606
    }  
607

608
    /// Returns whether coarsening is allowed or not.
609
610
    inline bool isCoarseningAllowed(int index) 
    {
611
      return (scalContents[index]->coarsenAllowed == 1);
612
    }
613

614
    /// Returns whether coarsening is allowed or not.
615
616
    inline bool isRefinementAllowed(int index) 
    {
617
      return (scalContents[index]->refinementAllowed == 1);
618
    }
619

620
    ///
621
622
    inline void allowRefinement(bool allow, int index) 
    {
623
      scalContents[index]->refinementAllowed = allow;
624
    }
625

626
    ///
627
628
    inline void allowCoarsening(bool allow, int index) 
    {
629
      scalContents[index]->coarsenAllowed = allow;
630
    }
631

632
    /// Returns \ref refineBisections
633
    inline int getRefineBisections(int index) const 
634
    {
635
      return scalContents[index]->refineBisections;
636
    }
637

638
    /// Returns \ref coarseBisections
639
    inline int getCoarseBisections(int index) const 
640
    {
641
      return scalContents[index]->coarseBisections;
642
    }    
643

644
645
    inline int getSize() 
    { 
Thomas Witkowski's avatar
Thomas Witkowski committed
646
      return scalContents.size(); 
647
    }
648

649
650
651
652
653
    inline bool getRosenbrockMode()
    {
      return rosenbrockMode;
    }

654
655
    inline void setSolverIterations(int it) 
    {
656
      solverIterations = it;
657
    }
658
  
659
660
    inline int getSolverIterations() 
    {
661
      return solverIterations;
662
    }
663
  
664
665
    inline void setMaxSolverIterations(int it) 
    {
666
      maxSolverIterations = it;
667
    }
668
  
669
670
    inline int getMaxSolverIterations() 
    {
671
      return maxSolverIterations;
672
    }
673
  
674
675
    inline void setSolverTolerance(double tol) 
    {
676
      solverTolerance = tol;
677
    }
678
  
679
680
    inline double getSolverTolerance() 
    {
681
      return solverTolerance;
682
    }
683
  
684
685
    inline void setSolverResidual(double res) 
    {
686
      solverResidual = res;
687
    }
688
  
689
690
    inline double getSolverResidual() 
    {
691
      return solverResidual;
692
    }
693

694
    /// Returns true, if the adaptive procedure was deserialized from a file.
695
    bool isDeserialized() const 
696
    {
697
      return deserialized;
698
699
    }

700
701
    inline void setIsDeserialized(bool b) 
    {
702
703
704
705
706
707
      deserialized = b;
    }

    inline void setRosenbrockMode(bool b)
    {
      rosenbrockMode = b;
708
709
    }

710
    /// Creates new scalContents with the given size.
711
712
    void setScalContents(int newSize);

713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
    /** \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
728
    void serialize(std::ostream& out);
729

Thomas Witkowski's avatar
Thomas Witkowski committed
730
    void deserialize(std::istream& in);
731
732

  protected:
733
    /// Name.
Thomas Witkowski's avatar
Thomas Witkowski committed
734
    std::string name;
735

736
    /// Current space iteration
737
738
739
740
741
742
743
744
    int spaceIteration;

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

745
    /// Current timestep iteration
746
747
    int timestepIteration;

748
    /// Maximal number of iterations for choosing a timestep
749
750
    int maxTimestepIteration;

751
    /// Current time iteration
752
753
    int timeIteration;

754
    /// Maximal number of time iterations
755
756
    int maxTimeIteration;

757
    /// Actual time, end of time interval for current time step
758
759
    double time;

760
    /// Initial time
761
762
    double startTime;

763
    /// Final time
764
765
    double endTime;

766
    ///Time step size to be used
767
768
    double timestep;

769
770
771
    /// Last processed time step size of finished iteration
    double lastProcessedTimestep;

772
    /// Minimal step size
773
774
    double minTimestep;

775
    /// Maximal step size
776
    double maxTimestep;
777
    
778
    /// Number of current time step
779
    int timestepNumber;
Thomas Witkowski's avatar
* Bla    
Thomas Witkowski committed
780
781
782
783
784
785
786

    /** \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;
787
  
788
    /// number of iterations needed of linear or nonlinear solver
789
790
    int solverIterations;

791
    /// maximal number of iterations needed of linear or nonlinear solver
792
793
    int maxSolverIterations;

794
    ///
795
796
    double solverTolerance;

797
    ///
798
799
    double solverResidual;

800
    /// Scalar adapt infos.
Thomas Witkowski's avatar
Thomas Witkowski committed
801
    std::vector<ScalContent*> scalContents;
802

803
    /// Is true, if the adaptive procedure was deserialized from a file.
804
805
806
807
    bool deserialized;

    /// Is true, if the time adaption is controlled by a Rosenbrock Method.
    bool rosenbrockMode;
Praetorius, Simon's avatar
Praetorius, Simon committed
808
809
810
    
    /// overall time error estimate
    double est_t;
811
812
813
814
815
  };

}

#endif //  AMDIS_ADAPTINFO_H