Fixed time step suggestion for explicit methods & give more verbose errors for...

Fixed time step suggestion for explicit methods & give more verbose errors for grid functions out of the transport time interval
Signed-off-by: Santiago Ospina <santiago.ospina@iup.uni-heidelberg.de>

dune/dorie/interface/transport_simulation.cc

@@ -139,27 +139,49 @@ void TransportSimulation<Traits>::step ()
     operator_setup();
   bool succeed = false;
   bool exception = false;
-  while (not succeed)
-  {
-    auto t = controller->getTime();
-    suggest_timestep(controller->getDT());
 
-    double courant_number = 0.5; //FIXME: add keyword in inifile
-    if (not ts_param->implicit())
-      this->suggest_timestep(courant_number*Dune::Dorie::cfl_condition(*gf_water_flux));
+  // Suggest times step for explicit methods. 
+  // This assumes that timestep is never increased on failure!
+  if (not ts_param->implicit()) {
+    double courant_number = 0.2; //FIXME: add keyword in inifile
+    const auto cfl = Dune::Dorie::cfl_condition(*gf_water_flux);
+    controller->suggest_timestep(courant_number*cfl);
+  }
 
-    // for some reason dt!=this->dt, any idea why?
+  while (not succeed)
+  {
+    // Obtain time variables
+    const auto t = controller->getTime();
+    const auto dt = controller->getDT();
     
     // Check whether time intervals of grid functions are valid
-    if (Dune::FloatCmp::gt(water_flux_interval.begin,t+this->dt))
-      DUNE_THROW(Dune::Exception," ");
-    if (Dune::FloatCmp::gt(saturation_interval.begin,t+this->dt))
-      DUNE_THROW(Dune::Exception," ");
-
-    if (Dune::FloatCmp::lt(water_flux_interval.end,t+this->dt))
-      DUNE_THROW(Dune::Exception,water_flux_interval.end<< " < "<<t+this->dt);
-    if (Dune::FloatCmp::lt(saturation_interval.end,t+this->dt))
-      DUNE_THROW(Dune::Exception," ");
+    if (Dune::FloatCmp::gt(water_flux_interval.begin,t+dt))
+      DUNE_THROW(Dune::InvalidStateException,
+        "Water Flux is in invalid state:\n" 
+        <<   "Water Flux time domain [" << water_flux_interval.begin << "," 
+                                        << water_flux_interval.end << "]\n" <<
+        "Transport time interval [" << t << "," << t+dt << "]");
+
+    if (Dune::FloatCmp::gt(saturation_interval.begin,t+dt))
+      DUNE_THROW(Dune::InvalidStateException,
+        "Saturation is in invalid state:\n" 
+        << "Saturation time domain [" << saturation_interval.begin << "," 
+                                      << saturation_interval.end << "]\n" <<
+      "Transport time interval [" << t << "," << t+dt << "]");
+
+    if (Dune::FloatCmp::lt(water_flux_interval.end,t+dt))
+      DUNE_THROW(Dune::InvalidStateException,
+        "Water Flux is in invalid state:\n" 
+        <<   "Water Flux time domain [" << water_flux_interval.begin << "," 
+                                        << water_flux_interval.end << "]\n" <<
+        "Transport time interval [" << t << "," << t+dt << "]");
+
+    if (Dune::FloatCmp::lt(saturation_interval.end,t+dt))
+      DUNE_THROW(Dune::InvalidStateException,
+        "Water Saturation is in invalid state:\n" 
+        << "Saturation time domain [" << saturation_interval.begin << "," 
+                                      << saturation_interval.end << "]\n" <<
+      "Transport time interval [" << t << "," << t+dt << "]");
 
     const bool solver_warnings = verbose > 0 && helper.rank() == 0 ?
       true : false;
@@ -170,7 +192,7 @@ void TransportSimulation<Traits>::step ()
 
       if (not solver_warnings)
         dwarn.push(false);
-      osm->apply(t, this->dt, *u, *unext);
+      osm->apply(t, dt, *u, *unext);
       if (not solver_warnings)
         dwarn.pop();
       u = unext;

dune/dorie/interface/transport_simulation.hh

@@ -240,7 +240,6 @@ protected:
   const int verbose;
 
   const Dune::VTK::OutputType output_type;
-  TimeField dt;
   bool operator_setup_flag = true;
 public:
 
@@ -296,10 +295,7 @@ public:
     return controller->getTime();
   }
 
-  void suggest_timestep(TimeField _dt) override
-  {
-    dt = std::min(_dt,controller->getDT());
-  }
+  void suggest_timestep(TimeField dt) final {controller->suggest_timestep(dt);}
 
   // TODO: Define data exchange with richards
   std::shared_ptr<GFSolute> get_solute()