diff --git a/NumLib/ODESolver/TimeDiscretization.h b/NumLib/ODESolver/TimeDiscretization.h
index db3fb84d8b7dcf4df55d617bb1df5306f867235e..49ce64cb86c93218e40ddcbe50d918ccaaa6341a 100644
--- a/NumLib/ODESolver/TimeDiscretization.h
+++ b/NumLib/ODESolver/TimeDiscretization.h
@@ -142,6 +142,18 @@ public:
virtual void pushState(const double t, GlobalVector const& x,
InternalMatrixStorage const& strg) = 0;
+ /*!
+ * Restores the given vector x to its old value.
+ * Used only for repeating of the time step with new time step size when
+ * the current time step is rejected. The restored x is only used as
+ * an initial guess for linear solver in the first Picard nonlinear
+ * iteration.
+ *
+ * \param x The solution at the current time step, which is going to be
+ * reset to its previous value.
+ */
+ virtual void popState(GlobalVector& x) = 0;
+
/*! Indicate that the computation of a new timestep is being started now.
*
* \warning Currently changing timestep sizes are not supported. Thus,
@@ -262,6 +274,11 @@ public:
MathLib::LinAlg::copy(x, _x_old);
}
+ void popState(GlobalVector& x) override
+ {
+ MathLib::LinAlg::copy(_x_old, x);
+ }
+
void nextTimestep(const double t, const double delta_t) override
{
_t = t;
@@ -314,6 +331,11 @@ public:
MathLib::LinAlg::copy(x, _x_old);
}
+ void popState(GlobalVector& x) override
+ {
+ MathLib::LinAlg::copy(_x_old, x);
+ }
+
void nextTimestep(const double t, const double delta_t) override
{
_t_old = _t;
@@ -392,6 +414,11 @@ public:
strg.pushMatrices();
}
+ void popState(GlobalVector& x) override
+ {
+ MathLib::LinAlg::copy(_x_old, x);
+ }
+
void nextTimestep(const double t, const double delta_t) override
{
_t = t;
@@ -463,6 +490,11 @@ public:
void pushState(const double, GlobalVector const& x,
InternalMatrixStorage const&) override;
+ void popState(GlobalVector& x) override
+ {
+ MathLib::LinAlg::copy(*_xs_old[_xs_old.size()-1], x);
+ }
+
void nextTimestep(const double t, const double delta_t) override
{
_t = t;
diff --git a/NumLib/TimeStepping/Algorithms/EvolutionaryPIDcontroller.cpp b/NumLib/TimeStepping/Algorithms/EvolutionaryPIDcontroller.cpp
index e1571221440478f3ffafb9ceebfde2de01626fe6..63ef785620750e543ce7246e945ac57edd38091f 100644
--- a/NumLib/TimeStepping/Algorithms/EvolutionaryPIDcontroller.cpp
+++ b/NumLib/TimeStepping/Algorithms/EvolutionaryPIDcontroller.cpp
@@ -16,6 +16,8 @@
#include "EvolutionaryPIDcontroller.h"
+#include <logog/include/logog.hpp>
+
namespace NumLib
{
bool EvolutionaryPIDcontroller::next(const double solution_error)
@@ -35,6 +37,13 @@ bool EvolutionaryPIDcontroller::next(const double solution_error)
_ts_current = _ts_prev;
_ts_current += h_new;
+
+ WARN("\tThis step is rejected due to the relative change from the\n"
+ "\tsolution of the previous time step to the current solution\n"
+ "\t exceeds the given tolerance of %g.\n"
+ "\tThis time step will be repeated with a new time step size of"
+ " %g or the simulation will be halted." , _tol, h_new);
+
return false;
}
diff --git a/ProcessLib/LiquidFlow/Tests.cmake b/ProcessLib/LiquidFlow/Tests.cmake
index d2933143d52615305e4acdd59e0ac982addc226f..e55c59e92002c1ade861139cde6bbba0b01257d6 100644
--- a/ProcessLib/LiquidFlow/Tests.cmake
+++ b/ProcessLib/LiquidFlow/Tests.cmake
@@ -105,10 +105,10 @@ AddTest(
REQUIREMENTS NOT OGS_USE_MPI
ABSTOL 1.1 RELTOL 1e-6
DIFF_DATA
- mesh2D.vtu gravity_driven_adaptive_dt_pcs_1_ts_11_t_300.000000.vtu OGS5_PRESSURE1 pressure
- mesh2D.vtu gravity_driven_adaptive_dt_pcs_1_ts_11_t_300.000000.vtu OGS5_TEMPERATURE1 temperature
- mesh2D.vtu gravity_driven_adaptive_dt_pcs_1_ts_11_t_300.000000.vtu v_x_ref v_x
- mesh2D.vtu gravity_driven_adaptive_dt_pcs_1_ts_11_t_300.000000.vtu v_y_ref v_y
+ mesh2D.vtu gravity_driven_adaptive_dt_pcs_1_ts_10_t_300.000000.vtu OGS5_PRESSURE1 pressure
+ mesh2D.vtu gravity_driven_adaptive_dt_pcs_1_ts_10_t_300.000000.vtu OGS5_TEMPERATURE1 temperature
+ mesh2D.vtu gravity_driven_adaptive_dt_pcs_1_ts_10_t_300.000000.vtu v_x_ref v_x
+ mesh2D.vtu gravity_driven_adaptive_dt_pcs_1_ts_10_t_300.000000.vtu v_y_ref v_y
)
# PETSc/MPI
@@ -218,8 +218,8 @@ AddTest(
REQUIREMENTS OGS_USE_MPI
ABSTOL 1.1 RELTOL 1e-7
DIFF_DATA
- mesh2D.vtu gravity_driven_adaptive_dt_pcs_1_ts_11_t_300_000000_0.vtu OGS5_PRESSURE1 pressure
- mesh2D.vtu gravity_driven_adaptive_dt_pcs_1_ts_11_t_300_000000_0.vtu OGS5_TEMPERATURE1 temperature
+ mesh2D.vtu gravity_driven_adaptive_dt_pcs_1_ts_10_t_300_000000_0.vtu OGS5_PRESSURE1 pressure
+ mesh2D.vtu gravity_driven_adaptive_dt_pcs_1_ts_10_t_300_000000_0.vtu OGS5_TEMPERATURE1 temperature
# To be activated when the output of velocity is correct under PETSc version.
# mesh2D.vtu gravity_driven_pcs_1_ts_5_t_300_000000_0.vtu v_x_ref v_x
# mesh2D.vtu gravity_driven_pcs_1_ts_5_t_300_000000_0.vtu v_y_ref v_y
diff --git a/ProcessLib/RichardsFlow/Tests.cmake b/ProcessLib/RichardsFlow/Tests.cmake
index 35b569a8d578bc4b292ee4d98cadba46c1c6acb0..51b95a26429a61fedc27d2d3e99cbd465ea6f859 100644
--- a/ProcessLib/RichardsFlow/Tests.cmake
+++ b/ProcessLib/RichardsFlow/Tests.cmake
@@ -31,3 +31,29 @@ AddTest(
ref_t_20000.000000.vtu richards_pcs_0_ts_18200_t_20000.000000.vtu pressure pressure
REQUIREMENTS NOT OGS_USE_MPI
)
+
+AddTest(
+ NAME 2D_RichardsFlow_h_us_quad_small_Adpative_dt
+ PATH Parabolic/Richards
+ EXECUTABLE ogs
+ EXECUTABLE_ARGS RichardsFlow_2d_small_adaptive_dt.prj
+ TESTER vtkdiff
+ ABSTOL 1e-8 RELTOL 1e-3
+ DIFF_DATA
+ ref_t_1600.000000.vtu richards_pcs_0_ts_805_t_1600.000000.vtu pressure pressure
+ REQUIREMENTS NOT OGS_USE_MPI
+)
+
+#PETSc/MPI
+AddTest(
+ NAME 2D_RichardsFlow_h_us_quad_small_Adpative_dt
+ PATH Parabolic/Richards
+ EXECUTABLE_ARGS RichardsFlow_2d_small_adaptive_dt.prj
+ WRAPPER mpirun
+ WRAPPER_ARGS -np 1
+ TESTER vtkdiff
+ REQUIREMENTS OGS_USE_MPI
+ ABSTOL 1e-8 RELTOL 1e-3
+ DIFF_DATA
+ ref_t_1600.000000.vtu richards_pcs_0_ts_805_t_1600_000000_0.vtu pressure pressure
+)
diff --git a/ProcessLib/UncoupledProcessesTimeLoop.cpp b/ProcessLib/UncoupledProcessesTimeLoop.cpp
index 60a7afc6c0ec3775c31f8b17bb725d46cc40625a..d39069ec0cc494f8cd3b82f86778d3b80543bf2d 100644
--- a/ProcessLib/UncoupledProcessesTimeLoop.cpp
+++ b/ProcessLib/UncoupledProcessesTimeLoop.cpp
@@ -620,21 +620,23 @@ double UncoupledProcessesTimeLoop::computeTimeStepping(
continue;
}
+ auto& time_disc = ppd.time_disc;
+ auto& mat_strg = *ppd.mat_strg;
+ auto& x = *_process_solutions[i];
if (all_process_steps_accepted)
{
- auto& time_disc = ppd.time_disc;
- auto& mat_strg = *ppd.mat_strg;
- auto const& x = *_process_solutions[i];
time_disc->pushState(t, x, mat_strg);
+ _repeating_times_of_rejected_step = 0;
}
else
{
if (t < _end_time)
{
WARN(
- "Time step %d is rejected. "
- "The computation is back to the previous time.",
- accepted_steps + rejected_steps);
+ "Time step %d was rejected %d times "
+ "and it will be repeated with a reduced step size.",
+ accepted_steps + 1, _repeating_times_of_rejected_step++);
+ time_disc->popState(x);
}
}
}
@@ -642,13 +644,28 @@ double UncoupledProcessesTimeLoop::computeTimeStepping(
if (!is_initial_step)
{
if (all_process_steps_accepted)
+ {
accepted_steps++;
+ _last_step_rejected = false;
+ }
else
{
+ if (std::abs(dt -prev_dt) < std::numeric_limits<double>::min()
+ && _last_step_rejected)
+ {
+ OGS_FATAL("\tThis time step is rejected and the new computed"
+ " step size is the same as\n"
+ "\tthat was just used.\n"
+ "\tSuggest to adjust the parameters of the time"
+ " stepper or try other time stepper.\n"
+ "\tThe program stops");
+ }
+
if (t < _end_time)
{
t -= prev_dt;
rejected_steps++;
+ _last_step_rejected = true;
}
}
}
@@ -729,7 +746,7 @@ bool UncoupledProcessesTimeLoop::loop()
t += dt;
const double prev_dt = dt;
- const std::size_t timesteps = accepted_steps + rejected_steps + 1;
+ const std::size_t timesteps = accepted_steps + 1;
// TODO, input option for time unit.
INFO("=== Time stepping at step #%u and time %g with step size %g",
timesteps, t, dt);
diff --git a/ProcessLib/UncoupledProcessesTimeLoop.h b/ProcessLib/UncoupledProcessesTimeLoop.h
index df0c9e058b3fc1ba66587f458fbe733ee990053c..fe62579f90c291e370c095b19973a1b10d6ce969 100644
--- a/ProcessLib/UncoupledProcessesTimeLoop.h
+++ b/ProcessLib/UncoupledProcessesTimeLoop.h
@@ -66,6 +66,8 @@ private:
std::unique_ptr<Output> _output;
std::vector<std::unique_ptr<SingleProcessData>> _per_process_data;
+ bool _last_step_rejected = false;
+ bool _repeating_times_of_rejected_step = 0;
const double _start_time;
const double _end_time;
diff --git a/Tests/Data b/Tests/Data
index d482ae358213a45b0887765ba6105e57194dba4d..e9754efad6c6a2467a66f6d0138c94e999c3972e 160000
--- a/Tests/Data
+++ b/Tests/Data
@@ -1 +1 @@
-Subproject commit d482ae358213a45b0887765ba6105e57194dba4d
+Subproject commit e9754efad6c6a2467a66f6d0138c94e999c3972e