Flow123d  master-f44eb46
time_governor.cc
Go to the documentation of this file.
1 /*!
2  *
3  * Copyright (C) 2015 Technical University of Liberec. All rights reserved.
4  *
5  * This program is free software; you can redistribute it and/or modify it under
6  * the terms of the GNU General Public License version 3 as published by the
7  * Free Software Foundation. (http://www.gnu.org/licenses/gpl-3.0.en.html)
8  *
9  * This program is distributed in the hope that it will be useful, but WITHOUT
10  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
11  * FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
12  *
13  *
14  * @file time_governor.cc
15  * @ingroup application
16  * @brief Basic time management class.
17  */
18 
19 #include <limits>
20 #include "system/system.hh"
21 #include "input/accessors.hh"
22 #include "time_governor.hh"
23 #include "time_marks.hh"
24 #include "unit_si.hh"
25 #include "unit_converter.hh"
26 
27 /*******************************************************************
28  * implementation of TimeGovernor static values and methods
29  */
30 
31 //initialize constant pointer to TimeMarks object
33 
34 //const double TimeGovernor::time_step_lower_bound = numeric_limits<double>::epsilon();
35 
36 #define MAX_END_TIME 5.0e+17
37 #define MAX_END_TIME_STR "5.0e+17"
38 
39 const double TimeGovernor::inf_time = numeric_limits<double>::infinity();
40 const double TimeGovernor::max_end_time = MAX_END_TIME; // more then age of universe in seconds.
41 const double TimeGovernor::time_step_precision = 16*numeric_limits<double>::epsilon();
42 
43 
44 using namespace Input::Type;
45 
46 
47 const Tuple & TimeGovernor::get_input_time_type(double lower_bound, double upper_bound)
48 {
49  return Tuple("TimeValue", "A time with optional unit specification.")
50  .declare_key("time", Double(lower_bound, upper_bound), Default::obligatory(),
51  "The time value." )
53  "Predefined units include: `s` seconds, `min` minutes, `h` hours, `d` days, `y` years.\n"
54  "The default time unit is set from the equation's time governor, see the key `common_time_unit`"
55  "in the equation's time record.\n\n"
56  "User can benefit from the Unit Convertor funcionality and create different time units.\n"
57  "Year length example considering leap years (Gregorian calendar): `year; year = 365.2425*d`.\n"
58  "Miliseconds example : `milisec; milisec = 0.001*s`.")
59  .close();
60 }
61 
62 
64  static const Tuple &dt_step =
65  Tuple("DtLimits", "Time dependent changes in min_dt and max_dt limits.")
67  "The start time of dt step set.")
68  .declare_key("min_dt", TimeGovernor::get_input_time_type(), Default::read_time("'min_dt' value of TimeGovernor."),
69  "Soft lower limit for the time step.")
70  .declare_key("max_dt", TimeGovernor::get_input_time_type(), Default::read_time("'max_dt' value of TimeGovernor."),
71  "Whole time of the simulation if specified, infinity else.")
72  .close();
73 
74  return Record("TimeGovernor",
75  "Time axis settings of the simulation.\n"
76  "The settings is specific to a particular equation.\n"
77  "TimeGovernor allows to:\n"
78  " - define start time and end time of simulation\n"
79  " - define lower and upper limits of time steps\n"
80  " - direct fixed time marks of whole simulation\n"
81  " - set global time unit of equation (see 'common_time_unit' key)\n"
82  "Limits of time steps are defined by keys 'min_dt', 'max_dt', 'init_dt' and 'dt_limits'. Key "
83  "'init_dt' has the highest priority and allows set fix size of time steps. Pair of keys 'min_dt' "
84  "and 'max_dt' define interval of time steps. Both previous cases ('init_dt' or pair 'min_dt' "
85  "and 'max_dt') set global limits of whole simulation. In contrasts, 'dt_limits' allow set "
86  "time-dependent function of min_dt/max_dt. Used time steps of simulation can be printed to YAML "
87  "output file (see 'write_used_timesteps'.\n"
88  "Fixed time marks define exact values of time steps. They are defined in:\n"
89  " - start time and end time of simulation\n"
90  " - output times printed to output mesh file\n"
91  " - times defined in 'dt_limits' table (optional, see 'add_dt_limits_time_marks' key)")
92  .allow_auto_conversion("max_dt")
93  .declare_key("start_time", TimeGovernor::get_input_time_type(), Default("0.0"),
94  "Start time of the simulation.")
96  "End time of the simulation.\n"
97  "The default value is higher than the age of the Universe (given in seconds).")
98  .declare_key("init_dt", TimeGovernor::get_input_time_type(0.0), Default("0.0"),
99  "Initial guess for the time step.\n"
100  "It applies to equations that use an adaptive time stepping. "
101  "If set to 0.0, the time step is determined in fully autonomous "
102  "way, assuming the equation supports it.")
104  Default::read_time("Machine precision."),
105  "Soft lower limit for the time step.\n"
106  "Equation using an adaptive time stepping cannot suggest smaller time step. "
107  "The actual time step can only decrease below the limit in order to match "
108  "the prescribed input or output times.")
110  Default::read_time("Whole time of the simulation if specified, infinity else."),
111  "Hard upper limit for the time step.\n"
112  "The actual time step can only increase above the limit in order to match "
113  "the prescribed input or output times.")
114  .declare_key("dt_limits", Array(dt_step), Default::optional(),
115  "Allow to set a time dependent changes in ``min_dt`` and ``max_dt`` limits. This list is processed "
116  "at individual times overwriting previous values of ``min_dt``/``max_dt``. Limits equal to 0 are "
117  "ignored and replaced with ``min_dt``/``max_dt`` values.")
118  .declare_key("add_dt_limits_time_marks", Bool(), Default("false"), "Add all times defined in ``dt_limits`` "
119  "table to the list of fixed TimeMarks.")
120  .declare_key("write_used_timesteps", FileName::output(), Default::optional(),
121  "Write used time steps to the given file in YAML format corresponding with the format of ``dt_limits``.")
122  .declare_key("common_time_unit", UnitConverter::get_input_type(), Default("\"s\""),
123  "Common time unit of the equation.\nThis unit will be used for all time inputs and outputs "
124  "within the equation. Individually, the common time unit can be overwritten for every declared time.\n"
125  "Time units are used in the following cases:\n"
126  "1) Time units of time value keys in: TimeGovernor, FieldDescriptors.\n"
127  " The common time unit can be overwritten for every declared time.\n"
128  "2) Time units in: \n"
129  " a) input fields: FieldFE and FieldTimeFunction\n"
130  " b) time steps definition of OutputTimeSet\n"
131  " Common time unit can be overwritten by one unit value for every whole mesh data file or time function.\n"
132  "3) Time units in output files: observation times, balance times, frame times of VTK and GMSH\n"
133  " Common time unit cannot be overwritten in these cases."
134  )
135  .close();
136 }
137 
138 
139 
140 /*******************************************************************
141  * implementation of TimeUnitConversion
142  */
143 
145 : coef_(1.0), unit_string_("s") {}
146 
148 {
149  unit_string_ = input.val<std::string>("unit_formula");
151 }
152 
154 {
155  std::string unit_str = input.val<std::string>("unit_formula");
156  try {
157  return UnitSI().s().convert_unit_from(unit_str);
158  } catch (ExcInvalidUnit &e) {
159  e << input.ei_address();
160  throw;
161  } catch (ExcNoncorrespondingUnit &e) {
162  e << input.ei_address();
163  throw;
164  }
165 }
166 
167 
168 
169 double TimeUnitConversion::read_time(Input::Iterator<Input::Tuple> time_it, double default_time) const {
170  if (time_it) {
171  double time = time_it->val<double>("time");
172 
173  Input::Record unit_record;
174  if ( time_it->opt_val("unit", unit_record) ) {
175  double unit_coef = read_unit_coef_from_input(unit_record);
176  return ( time * unit_coef );
177  }
178  else {
179  return ( time * coef_ );
180  }
181  } else {
182  ASSERT(default_time!=std::numeric_limits<double>::quiet_NaN()).error("Undefined default time!");
183  return default_time;
184  }
185 }
186 
187 
188 
190  if (unit_it) {
191  return read_unit_coef_from_input(*unit_it);
192  } else {
193  return coef_;
194  }
195 }
196 
197 
198 
199 /*******************************************************************
200  * implementation of TimeStep
201  */
202 
203 TimeStep::TimeStep(double init_time, std::shared_ptr<TimeUnitConversion> time_unit_conversion) :
204 index_(0),
205 length_(1.0),
206 end_(init_time),
207 time_unit_conversion_(time_unit_conversion)
208 {}
209 
210 
211 
213 index_(0),
214 length_(TimeGovernor::inf_time),
215 end_(-TimeGovernor::inf_time)
216 {
217  time_unit_conversion_ = std::make_shared<TimeUnitConversion>();
218 }
219 
220 
221 
222 
224 index_(other.index_),
225 length_(other.length_),
226 end_(other.end_),
227 time_unit_conversion_(other.time_unit_conversion_)
228 {}
229 
230 
231 
232 TimeStep TimeStep::make_next(double new_length) const
233 {
234  return make_next(new_length, this->end_+new_length);
235 }
236 
237 
238 
239 TimeStep TimeStep::make_next(double new_lenght, double end_time) const
240 {
241  TimeStep ts;
242  ts.index_=this->index_ +1;
243  ts.length_=new_lenght;
244  ts.end_=end_time;
246  return ts;
247 }
248 
249 
250 
251 bool TimeStep::safe_compare(double t1, double t0) const
252 {
253  return t1 >= t0
254  - 16*numeric_limits<double>::epsilon()*(1.0+max(abs(t1),abs(t0)));
255 }
256 
257 
258 
259 double TimeStep::read_time(Input::Iterator<Input::Tuple> time_it, double default_time) const {
260  return time_unit_conversion_->read_time(time_it, default_time);
261 }
262 
263 
264 
266  return time_unit_conversion_->read_coef(unit_it);
267 }
268 
269 
270 
271 double TimeStep::get_coef() const {
272  return time_unit_conversion_->get_coef();
273 }
274 
275 
276 
277 std::shared_ptr<TimeUnitConversion> TimeStep::get_unit_conversion() const
278 {
279  return time_unit_conversion_;
280 }
281 
282 ostream& operator<<(ostream& out, const TimeStep& t_step) {
283  out << "time: " << t_step.end() << " dt: " << t_step.length() << endl;
284  return out;
285 }
286 
287 
288 
289 /*******************************************************************
290  * implementation of TimeGovernor
291  */
292 
293 TimeGovernor::TimeGovernor(const Input::Record &input, TimeMark::Type eq_mark_type, bool timestep_output)
294 : timestep_output_(timestep_output)
295 {
296  // use new mark type as default
297  if (eq_mark_type == TimeMark::none_type) eq_mark_type = marks().new_mark_type();
298 
299  try {
300 
301  Input::Record common_unit_record = input.val<Input::Record>("common_time_unit");
302  time_unit_conversion_ = std::make_shared<TimeUnitConversion>(common_unit_record);
303 
304  limits_time_marks_ = input.val<bool>("add_dt_limits_time_marks");
305 
306  // Get rid of rounding errors.
307  double end_time = read_time( input.find<Input::Tuple>("end_time") );
309 
310  // set permanent limits
311  init_common(read_time( input.find<Input::Tuple>("start_time") ),
312  end_time,
313  eq_mark_type);
314  Input::Array limits_array = Input::Array();
315  input.opt_val("dt_limits", limits_array);
317  read_time( input.find<Input::Tuple>("min_dt"), min_time_step_),
318  read_time( input.find<Input::Tuple>("max_dt"), max_time_step_),
319  limits_array
320  );
321 
322  // check key write_used_timesteps, open YAML file, print first time step
323  if (timestep_output_)
324  if (input.opt_val("write_used_timesteps", timesteps_output_file_) ) {
325  try {
327  } INPUT_CATCH(FilePath::ExcFileOpen, FilePath::EI_Address_String, input)
328  timesteps_output_ << "- [ " << t() << ", " << dt_limits_table_[0].min_dt << ", " << dt_limits_table_[0].max_dt << " ]\n";
330  }
331 
332  double init_dt=read_time( input.find<Input::Tuple>("init_dt") );
333  if (init_dt > 0.0) {
334  // set first time step suggested by user
335  //time_step_=min(init_dt, time_step_);
336  lower_constraint_=init_dt;
337  lower_constraint_message_ = "Initial time step set by user.";
338  upper_constraint_=init_dt;
339  upper_constraint_message_ = "Initial time step set by user.";
340  }
341 
342 
343  } catch(ExcTimeGovernorMessage &exc) {
344  exc << input.ei_address();
345  throw;
346  }
347 
348 }
349 
350 TimeGovernor::TimeGovernor(double init_time, double dt)
351 : dt_limits_pos_(0), timestep_output_(false)
352 {
353  time_unit_conversion_ = std::make_shared<TimeUnitConversion>();
355  // fixed time step
356  if (dt < time_step_precision)
357  THROW(ExcTimeGovernorMessage() << EI_Message("Fixed time step smaller then machine precision. \n") );
358 
360  is_time_step_fixed_=true;
361  time_step_changed_=true;
363 
364  // fill table limits with two records (start time, end time)
365  dt_limits_table_.push_back( DtLimitRow(init_time, dt, dt) );
366  dt_limits_table_.push_back( DtLimitRow(inf_time, dt, dt) );
367 
369  lower_constraint_message_ = "Initial time step set by user.";
371  upper_constraint_message_ = "Initial time step set by user.";
372 
373  //time_step_=dt;
374 }
375 
376 
377 // steady time governor constructor
378 TimeGovernor::TimeGovernor(double init_time, TimeMark::Type eq_mark_type)
379 : dt_limits_pos_(0), timestep_output_(false)
380 {
381  // use new mark type as default
382  if (eq_mark_type == TimeMark::none_type) eq_mark_type = marks().new_mark_type();
383 
384  time_unit_conversion_ = std::make_shared<TimeUnitConversion>();
385  init_common(init_time, inf_time, eq_mark_type);
386 
387  // fill table limits with two records (start time, end time)
390 
391  steady_ = true;
392 }
393 
394 
396 {
398  timesteps_output_.close();
399  }
400 }
401 
402 
403 // common part of constructors
404 void TimeGovernor::init_common(double init_time, double end_time, TimeMark::Type type)
405 {
406 
407  if (init_time < 0.0) {
408  THROW(ExcTimeGovernorMessage()
409  << EI_Message("Start time has to be greater or equal to 0.0\n")
410 
411  );
412  }
413 
414  recent_steps_.set_capacity(size_of_recent_steps_);
417 
418  if (end_time < init_time) {
419  THROW(ExcTimeGovernorMessage() << EI_Message("End time must be greater than start time.\n") );
420  } else {
422  }
423 
424  //if (dt == 0.0) {
425  // variable time step
426  fixed_time_step_=0.0;
427  is_time_step_fixed_=false;
428  time_step_changed_=true;
430 
432  lower_constraint_message_ = "Permanent minimal constraing, default, time_step_precision.";
434  upper_constraint_message_ = "Permanent maximal constraint, default, total simulation time.";
435  // choose maximum possible time step
436  //time_step_=max_time_step_;
437  /*} else {
438  // fixed time step
439  if (dt < time_step_lower_bound)
440  THROW(ExcTimeGovernorMessage() << EI_Message("Fixed time step smaller then machine precision. \n") );
441 
442  fixed_time_step_=dt;
443  is_time_step_fixed_=true;
444  time_step_changed_=true;
445  end_of_fixed_dt_interval_ = inf_time;
446 
447  upper_constraint_=max_time_step_=dt;
448  lower_constraint_=min_time_step_=dt;
449  time_step_=dt;
450 
451  }*/
452 
453  eq_mark_type_=type;
454  steady_=false;
456  //if (end_time_ != inf_time)
458 }
459 
460 
461 
462 void TimeGovernor::set_dt_limits( double min_dt, double max_dt, Input::Array dt_limits_list)
463 {
464  dt_limits_table_.clear();
465 
466  // check min_dt and max_dt set by user
467  if (min_dt < time_step_precision) {
468  THROW(ExcTimeGovernorMessage() << EI_Message("'min_dt' smaller than machine precision.\n") );
469  }
470  if (max_dt < min_dt) {
471  THROW(ExcTimeGovernorMessage() << EI_Message("'max_dt' smaller than 'min_dt'.\n") );
472  }
473 
474  bool first_step = true;
475  if (dt_limits_list.size())
476  for(auto it = dt_limits_list.begin<Input::Tuple>(); it != dt_limits_list.end(); ++it) {
477  double time = read_time( it->find<Input::Tuple>("time"));
478 
479  if (first_step) { // we need special check before setting first time step to the table
480  if (time > init_time_) { // table starts later than simulation, we need to add start time to the table
481  dt_limits_table_.push_back( DtLimitRow(init_time_, min_dt, max_dt) );
482  }
483  first_step = false;
484  }
485 
486  // next cases will be skipped
487  if (time < init_time_) {
488  WarningOut().fmt("Time {} define in 'dt_limits' table at address {} is lesser than start time of simulation "
489  "and can be skipped.\n", time, dt_limits_list.address_string());
490  }
491  if (dt_limits_table_.size() && (time <= dt_limits_table_[dt_limits_table_.size()-1].time) ) {
492  WarningOut().fmt("Time {} define in 'dt_limits' table at address {} is in incorrect order "
493  "and will be skipped.\n", time, dt_limits_list.address_string());
494  continue;
495  }
496  if ((time > end_time_) ) {
497  WarningOut().fmt("Time {} define in 'dt_limits' table at address {} is greater than end time of simulation "
498  "and will be skipped.\n", time, dt_limits_list.address_string());
499  continue;
500  }
501 
502  double min = read_time( it->find<Input::Tuple>("min_dt"), 0.0);
503  if (min == 0.0) min = min_dt;
504  double max = read_time( it->find<Input::Tuple>("max_dt"), 0.0);
505  if (max == 0.0) max = max_dt;
506 
507  if (min < time_step_precision) {
508  THROW(ExcTimeGovernorMessage() << EI_Message("'min_dt' in 'dt_limits' smaller than machine precision.\n") );
509  }
510  if (max < min) {
511  THROW(ExcTimeGovernorMessage() << EI_Message("'max_dt' in 'dt_limits' smaller than 'min_dt'.\n") );
512  }
513 
514  dt_limits_table_.push_back( DtLimitRow(time, min, max) );
515  if (limits_time_marks_) this->marks().add(TimeMark(time, this->equation_fixed_mark_type()));
516  }
517 
518  if (dt_limits_table_.size() == 0) {
519  // add start time to limit table if it is empty
520  dt_limits_table_.push_back( DtLimitRow(init_time_, min_dt, max_dt) );
521  }
522  if (dt_limits_table_[dt_limits_table_.size()-1].time < end_time_) {
523  // add time == end_time_ to limits table, we need only for check time, not for limits
524  dt_limits_table_.push_back( DtLimitRow(end_time_, min_dt, max_dt) );
525  }
526 
527  dt_limits_pos_ = 0;
528  while (dt_limits_table_[dt_limits_pos_+1].time <= init_time_) {
529  ++dt_limits_pos_;
530  }
531 
533 }
534 
535 
537 {
539  lower_constraint_message_ = "Permanent minimal constraint, custom.";
541  upper_constraint_message_ = "Permanent maximal constraint, custom.";
542  ++dt_limits_pos_;
543 }
544 
545 
546 // int set_constrain - dle nastaveni constraint
547 // interval - constraint - jako v cmp u Qsortu
548 // -1 vetsi nez interval (min_dt,max_dt)
549 // +1 mensi
550 // 0 OK
551 
552 int TimeGovernor::set_upper_constraint (double upper, std::string message)
553 {
554  if (upper_constraint_ < upper) {
555  //do not change upper_constraint_
556  return -1;
557  } else
558  if (lower_constraint_ > upper) {
559  // set upper constraint to the lower constraint
561  upper_constraint_message_ = "Forced lower constraint. " + message;
562  return 1;
563  } else {
564  //change upper_constraint_ to upper
565  upper_constraint_ = upper;
566  upper_constraint_message_ = message;
567  return 0;
568  }
569 }
570 
571 
572 
573 int TimeGovernor::set_lower_constraint (double lower, std::string message)
574 {
575  if (upper_constraint_ < lower) {
576  // set lower constraint to the upper constraint
578  return -1;
579  } else
580  if (lower_constraint_ > lower) {
581  //do not change lower_constraint_
582  return 1;
583  } else {
584  //change lower_constraint_ to lower
585  lower_constraint_ = lower;
586  lower_constraint_message_ = message;
587  return 0;
588  }
589 }
590 
591 
592 
593 
595  if (steady_) return 0.0;
596  end_of_fixed_dt_interval_=-inf_time; // release previous fixed interval
598  is_time_step_fixed_ = true; //flag means fixed step has been set since now
600 }
601 
602 
603 
604 
605 void TimeGovernor::add_time_marks_grid(double step, TimeMark::Type mark_type) const
606 {
607  if (end_time() == inf_time) {
608  THROW(ExcTimeGovernorMessage()
609  << EI_Message("Missing end time for making output grid required by key 'time_step' of the output stream.\n")
610  );
611  }
612 
614  // always add start time and end time
615  marks().add(TimeMark(init_time_, mark_type | eq_mark_type_));
616  marks().add(TimeMark(end_time(), mark_type | eq_mark_type_));
617 }
618 
619 
621 {
622  TimeMark::Type type = equation_mark_type() | mask;
623  return time_marks_.current(step(), type) != time_marks_.end(type);
624 }
625 
626 
627 
629  if (is_end()) return 0.0;
630 
631  if (this->step().lt(end_of_fixed_dt_interval_)) return fixed_time_step_;
632 
633  // jump to the first future fix time
635  // compute step to next fix time and apply constraints
636  double full_step = fix_time_it->time() - t();
637 
638  double step_estimate = min(full_step, upper_constraint_);
639  step_estimate = max(step_estimate, lower_constraint_); //these two must be in this order
640 
641  if (step_estimate == inf_time) return step_estimate;
642 
643  // round the time step to have integer number of steps till next fix time
644  // this always selects shorter time step,
645  // but allows time step larger then constraint by a number close to machine epsilon
646  //
647  double n_steps = ceil( full_step / step_estimate - time_step_precision);
648  step_estimate = full_step / n_steps;
649 
650  // try to avoid time_step changes
651  if (n_steps > 1 && abs(step_estimate - step().length()) < time_step_precision) {
652  step_estimate=step().length();
653  }
654 
655  // if the step estimate gets by rounding lower than lower constraint program will not crash
656  // will just output a user warning.
657  if (step_estimate < lower_constraint_) {
658  DebugOut().fmt("Time step estimate is below the lower constraint of time step. The difference is: {:.16f}.\n",
659  lower_constraint_ - step_estimate);
660  }
661 
662  return step_estimate;
663 }
664 
665 
666 
668 {
669  ASSERT_LE(0.0, t());
670  if (is_end()) return;
671 
672 
673  if (this->step().lt(end_of_fixed_dt_interval_)) {
674  // this is done for fixed step
675  // make tiny correction of time step in order to avoid big rounding errors
676  // tiny correction means that dt_changed 'is NOT changed'
679  }
680 
681  recent_steps_.push_front(recent_steps_.front().make_next(fixed_time_step_));
682 
683 
684  //checking whether fixed time step has been changed (by fix_dt_until_mark() method) since last time
686  {
687  is_time_step_fixed_ = false;
688 
689  //is true unless new fixed_dt is not equal previous time_step
690  time_step_changed_ = (step(-2).length() != step().length());
691  }
692  else
693  time_step_changed_ = false;
694  }
695  else
696  {
697  // this is done if end_of_fixed_dt_interval is not set (means it is equal to -infinity)
698  double dt=estimate_dt();
699  TimeStep step_ = recent_steps_.front().make_next(dt);
700  //DebugOut().fmt("step: {}, end: {}\n", step_.length(), step_.end());
701  recent_steps_.push_front(step_);
702  //DebugOut().fmt("last: {}, new: {}\n",step(-1).length(),step().length());
703  time_step_changed_= (step(-2).length() != step().length());
704  }
705 
708  // refreshing the upper_constraint_
711  lower_constraint_message_ = "Permanent minimal constraint, in next time.";
712  upper_constraint_message_ = "Permanent maximal constraint, in next time.";
713 
715 
716  // write time step to YAML file
718  double time = t();
719  if (time > last_printed_timestep_) {
720  if (is_end()) timesteps_output_ << "- [ " << time << ", 0, 0 ]\n";
721  else timesteps_output_ << "- [ " << time << ", " << lower_constraint_ << ", " << upper_constraint_ << " ]\n";
722  last_printed_timestep_ = time;
723  }
724  }
725 }
726 
727 
728 double TimeGovernor::reduce_timestep(double factor) {
729  double prior_dt = dt();
730  double new_upper_constraint = factor * dt();
731 
732  // Revert time.
733 // DebugOut().fmt("tg idx: {}\n", recent_steps_.front().index());
734  recent_steps_.pop_front();
735 // DebugOut().fmt("tg idx: {}\n", recent_steps_.back().index());
738 
739  // Set constraint.
740  int current_minus_new = set_upper_constraint(new_upper_constraint, "Reduce time step.");
741  if (current_minus_new < 0)
742  // current constraint < reduced dt, should not happen
743  THROW(ExcMessage() << EI_Message("Internal error."));
744 
745  next_time();
746 
747  // Return false if we hit lower time step constraint.
748  return dt() / prior_dt;
749 }
750 
751 
752 
753 const TimeStep &TimeGovernor::step(int index) const {
754  unsigned int back_idx;
755  if (index < 0) {
756  back_idx = static_cast<unsigned int>(-index-1);
757  } else {
758  back_idx = static_cast<unsigned int>(recent_steps_[0].index() - index);
759  }
760  if ( back_idx >= recent_steps_.size())
761  THROW(ExcMissingTimeStep() << EI_Index(index) << EI_BackIndex(back_idx) << EI_HistorySize(recent_steps_.size()));
762 
763  return recent_steps_[back_idx];
764 }
765 
766 
767 
768 
769 void TimeGovernor::view(const char *name) const
770 {
771 #ifdef FLOW123D_DEBUG_MESSAGES
772  MessageOut().fmt(
773  "TG[{}]:{:06d} t:{:10.4f} dt:{:10.6f} dt_int<{:10.6f},{:10.6f}> "
774  "end_time: {} end_fixed_time: {} type: {:#x}\n",
777 
778  MessageOut().fmt("Lower time step constraint [{}]: {} \nUpper time step constraint [{}]: {} \n",
781 #else
782  MessageOut().fmt(
783  "TG[{}]:{:06d} t:{:10.4f} dt:{:10.6f} dt_int<{:10.6f},{:10.6f}>\n",
784  name, tlevel(), t(), dt(), lower_constraint_, upper_constraint_);
785 #endif
786 }
787 
788 
789 
790 double TimeGovernor::read_time(Input::Iterator<Input::Tuple> time_it, double default_time) const {
791  return time_unit_conversion_->read_time(time_it, default_time);
792 }
793 
794 
795 
797  return time_unit_conversion_->read_coef(unit_it);
798 }
799 
800 
801 
802 double TimeGovernor::get_coef() const {
803  return time_unit_conversion_->get_coef();
804 }
805 
806 
807 std::shared_ptr<TimeUnitConversion> TimeGovernor::get_unit_conversion() const
808 {
809  return time_unit_conversion_;
810 }
811 
812 
813 ostream& operator<<(ostream& out, const TimeGovernor& tg)
814 {
815  static char buffer[1024];
816  sprintf(buffer, "\n%06d t:%10.4f dt:%10.6f dt_int<%10.6f,%10.6f>\n",
817  tg.tlevel(), tg.t(), tg.dt(), tg.lower_constraint(), tg.upper_constraint());
818  return (out << buffer);
819 }
820 
821 
#define ASSERT(expr)
Definition: asserts.hh:351
#define ASSERT_LE(a, b)
Definition of comparative assert macro (Less or Equal) only for debug mode.
Definition: asserts.hh:309
Dedicated class for storing path to input and output files.
Definition: file_path.hh:54
void open_stream(Stream &stream) const
Definition: file_path.cc:211
Accessor to input data conforming to declared Array.
Definition: accessors.hh:566
string address_string() const
Definition: accessors.cc:321
Iterator< ValueType > begin() const
unsigned int size() const
IteratorBase end() const
Accessor to the data with type Type::Record.
Definition: accessors.hh:291
EI_Address ei_address() const
Definition: accessors.cc:178
bool opt_val(const string &key, Ret &value) const
const Ret val(const string &key) const
Iterator< Ret > find(const string &key) const
Accessor to the data with type Type::Tuple.
Definition: accessors.hh:411
Class for declaration of inputs sequences.
Definition: type_base.hh:339
Class for declaration of the input of type Bool.
Definition: type_base.hh:452
Class Input::Type::Default specifies default value of keys of a Input::Type::Record.
Definition: type_record.hh:61
static Default obligatory()
The factory function to make an empty default value which is obligatory.
Definition: type_record.hh:110
static Default read_time(const std::string &description)
The factory function to make an default value that will be specified at the time when a key will be r...
Definition: type_record.hh:97
static Default optional()
The factory function to make an empty default value which is optional.
Definition: type_record.hh:124
Class for declaration of the input data that are floating point numbers.
Definition: type_base.hh:534
static FileName output()
The factory function for declaring type FileName for input files.
Definition: type_base.cc:531
Record type proxy class.
Definition: type_record.hh:182
virtual Record & allow_auto_conversion(const string &from_key)
Allows shorter input of the Record providing only value of the from_key given as the parameter.
Definition: type_record.cc:133
Record & close() const
Close the Record for further declarations of keys.
Definition: type_record.cc:304
Record & declare_key(const string &key, std::shared_ptr< TypeBase > type, const Default &default_value, const string &description, TypeBase::attribute_map key_attributes=TypeBase::attribute_map())
Declares a new key of the Record.
Definition: type_record.cc:503
Tuple type proxy class.
Definition: type_tuple.hh:45
Tuple & declare_key(const string &key, std::shared_ptr< TypeBase > type, const Default &default_value, const string &description, TypeBase::attribute_map key_attributes=TypeBase::attribute_map())
Overrides Record::declare_key(const string &, std::shared_ptr<TypeBase>, const Default &,...
Definition: type_tuple.cc:173
const Tuple & close() const
Close the Tuple for further declarations of keys.
Definition: type_tuple.cc:70
Basic time management functionality for unsteady (and steady) solvers (class Equation).
TimeMark::Type equation_fixed_mark_type() const
std::shared_ptr< TimeUnitConversion > time_unit_conversion_
Conversion unit of all time values within the equation.
bool time_step_changed_
Flag is set if the time step has been changed (lasts only one time step).
std::shared_ptr< TimeUnitConversion > get_unit_conversion() const
Getter for time unit conversion object.
double reduce_timestep(double factor)
Force timestep reduction in particular in the case of failure of the non-linear solver.
double lower_constraint() const
std::ofstream timesteps_output_
Handle for file for output time steps to YAML format.
double read_time(Input::Iterator< Input::Tuple > time_it, double default_time=std::numeric_limits< double >::quiet_NaN()) const
void add_time_marks_grid(double step, TimeMark::Type mark_type=TimeMark::none_type) const
void set_dt_limits(double min_dt, double max_dt, Input::Array dt_limits_list)
Sets dt limits for time dependent DT limits in simulation.
double min_time_step_
Permanent lower limit for the time step.
double last_upper_constraint_
Upper constraint used for choice of current time.
double dt() const
static const Input::Type::Tuple & get_input_time_type(double lower_bound=-std::numeric_limits< double >::max(), double upper_bound=std::numeric_limits< double >::max())
double end_time_
End time of the simulation.
bool steady_
True if the time governor is used for steady problem.
FilePath timesteps_output_file_
File path for timesteps_output_ stream.
bool limits_time_marks_
Allows add all times defined in dt_limits_table_ to list of TimeMarks.
double lower_constraint_
Lower constraint for the choice of the next time step.
double get_coef() const
static const unsigned int size_of_recent_steps_
double end_of_fixed_dt_interval_
End of interval if fixed time step.
double t() const
static const Input::Type::Record & get_input_type()
std::string upper_constraint_message_
Description of the upper constraint.
int set_lower_constraint(double lower, std::string message)
Sets lower constraint for the next time step estimating.
std::string lower_constraint_message_
Description of the upper constraint.
TimeMark::Type eq_mark_type_
TimeMark type of the equation.
std::vector< DtLimitRow > dt_limits_table_
Table of DT limits.
double upper_constraint_
Upper constraint for the choice of the next time step.
bool is_end() const
Returns true if the actual time is greater than or equal to the end time.
int set_upper_constraint(double upper, std::string message)
Sets upper constraint for the next time step estimating.
double fixed_time_step_
Next fixed time step.
double upper_constraint() const
double read_coef(Input::Iterator< Input::Record > unit_it) const
static TimeMarks time_marks_
void view(const char *name="") const
boost::circular_buffer< TimeStep > recent_steps_
Circular buffer of recent time steps. Implicit size is 3.
double end_time() const
End time.
bool is_time_step_fixed_
Flag that is set when the fixed step is set (lasts only one time step).
double last_lower_constraint_
Lower constraint used for choice of current time.
double estimate_dt() const
Estimate choice of next time step according to actual setting of constraints.
static const double time_step_precision
TimeGovernor(const Input::Record &input, TimeMark::Type fixed_time_mask=TimeMark::none_type, bool timestep_output=true)
Constructor for unsteady solvers.
static const double max_end_time
void set_permanent_constraint()
Sets permanent constraints for actual time step.
TimeMark::Type equation_mark_type() const
double last_printed_timestep_
Store last printed time to YAML output, try multiplicity output of one time.
double init_time() const
double fix_dt_until_mark()
Fixing time step until fixed time mark.
double init_time_
Initial time.
bool timestep_output_
Special flag allows forbid output time steps during multiple initialization of TimeGovernor.
bool is_current(const TimeMark::Type &mask) const
const TimeStep & step(int index=-1) const
static const double inf_time
Infinity time used for steady case.
static TimeMarks & marks()
unsigned int dt_limits_pos_
Index to actual position of DT limits.
double max_time_step_
Permanent upper limit for the time step.
int tlevel() const
void init_common(double init_time, double end_time, TimeMark::Type type)
Common part of the constructors. Set most important parameters, check they are valid and set default ...
void next_time()
Proceed to the next time according to current estimated time step.
Class used for marking specified times at which some events occur.
Definition: time_marks.hh:45
static const Type every_type
Mark Type with all bits set.
Definition: time_marks.hh:93
double time() const
Getter for the time of the TimeMark.
Definition: time_marks.hh:115
static const Type none_type
Mark Type with all bits unset.
Definition: time_marks.hh:95
Iterator over TimeMark objects in TimeMarks object (database of TimeMark objects).
Definition: time_marks.hh:351
This class is a collection of time marks to manage various events occurring during simulation time.
Definition: time_marks.hh:206
void add_time_marks(double time, double dt, double end_time, TimeMark::Type type)
Definition: time_marks.cc:106
TimeMarks::iterator end(TimeMark::Type mask) const
Iterator for the end mimics container-like of TimeMarks.
Definition: time_marks.cc:206
TimeMarks::iterator current(const TimeStep &time_step, const TimeMark::Type &mask) const
Definition: time_marks.cc:141
TimeMark::Type new_mark_type()
Definition: time_marks.cc:68
TimeMark add(const TimeMark &mark)
Definition: time_marks.cc:81
TimeMarks::iterator next(const TimeGovernor &tg, const TimeMark::Type &mask) const
Definition: time_marks.cc:149
Representation of one time step..
double read_time(Input::Iterator< Input::Tuple > time_it, double default_time=std::numeric_limits< double >::quiet_NaN()) const
double end() const
std::shared_ptr< TimeUnitConversion > time_unit_conversion_
Conversion unit of all time values within the equation.
TimeStep make_next(double new_length) const
double end_
End time point of the time step.
double read_coef(Input::Iterator< Input::Record > unit_it) const
double length() const
std::shared_ptr< TimeUnitConversion > get_unit_conversion() const
Getter for time unit conversion object.
bool safe_compare(double t1, double t0) const
unsigned int index_
Index of the step is index if the end time. Zero time step is artificial.
double get_coef() const
double length_
double read_unit_coef_from_input(const Input::Record &input) const
Reads the Unit record and computes the coef.
double coef_
Conversion coefficient of all time values within the equation.
double read_time(Input::Iterator< Input::Tuple > time_it, double default_time) const
double read_coef(Input::Iterator< Input::Record > unit_it) const
std::string unit_string_
String representation of global unit of all time values within the equation.
static Input::Type::Default get_input_default()
static const Input::Type::Record & get_input_type()
Class for representation SI units of Fields.
Definition: unit_si.hh:40
UnitSI & s(int exp=1)
Definition: unit_si.cc:76
double convert_unit_from(std::string actual_unit) const
Convert and check user-defined unit.
Definition: unit_si.cc:217
#define THROW(whole_exception_expr)
Wrapper for throw. Saves the throwing point.
Definition: exceptions.hh:53
#define INPUT_CATCH(ExceptionType, AddressEITag, input_accessor)
Definition: accessors.hh:63
#define WarningOut()
Macro defining 'warning' record of log.
Definition: logger.hh:278
#define MessageOut()
Macro defining 'message' record of log.
Definition: logger.hh:275
#define DebugOut()
Macro defining 'debug' record of log.
Definition: logger.hh:284
std::string sprintf(CStringRef format, ArgList args)
Definition: printf.h:457
Structure that stores one record of DT limit.
unsigned long int bitmap_
Definition: time_marks.hh:88
ostream & operator<<(ostream &out, const TimeStep &t_step)
#define MAX_END_TIME_STR
#define MAX_END_TIME
Basic time management class.