Specific times in output

Documentation/ProjectFile/prj/time_loop/output/t_fixed_output_times.md

@@ -0,0 +1,2 @@
+An option input of fixed times, at which the output of the results must be
+ conducted.

Documentation/ProjectFile/prj/time_loop/processes/process/time_stepping/EvolutionaryPIDcontroller/t_fixed_output_times.md

@@ -0,0 +1 @@
+The fixed times that must be reached in the time stepping.

NumLib/TimeStepping/Algorithms/CreateEvolutionaryPIDcontroller.cpp

@@ -42,25 +42,21 @@ std::unique_ptr<TimeStepAlgorithm> createEvolutionaryPIDcontroller(
     //! \ogs_file_param{prj__time_loop__processes__process__time_stepping__EvolutionaryPIDcontroller__rel_dt_max}
     auto const rel_h_max = config.getConfigParameter<double>("rel_dt_max");
 
-    auto specific_times =
-        //! \ogs_file_param{prj__time_loop__processes__process__time_stepping__EvolutionaryPIDcontroller__specific_times}
-        config.getConfigParameter<std::vector<double>>("specific_times",
+    auto fixed_output_times =
+        //! \ogs_file_param{prj__time_loop__processes__process__time_stepping__EvolutionaryPIDcontroller__fixed_output_times}
+        config.getConfigParameter<std::vector<double>>("fixed_output_times",
                                                        std::vector<double>{});
-    if (!specific_times.empty())
+    if (!fixed_output_times.empty())
     {
-        // Sort in descending order.
-        std::sort(specific_times.begin(),
-                  specific_times.end(),
-                  std::greater<double>());
-        // Remove possible duplicated elements.
-        BaseLib::makeVectorUnique(specific_times);
+        // Remove possible duplicated elements and sort in descending order.
+        BaseLib::makeVectorUnique(fixed_output_times, std::greater<double>());
     }
 
     //! \ogs_file_param{prj__time_loop__processes__process__time_stepping__EvolutionaryPIDcontroller__tol}
     auto const tol = config.getConfigParameter<double>("tol");
 
     return std::make_unique<EvolutionaryPIDcontroller>(
         t0, t_end, h0, h_min, h_max, rel_h_min, rel_h_max,
-        std::move(specific_times), tol);
+        std::move(fixed_output_times), tol);
 }
 }  // end of namespace NumLib

NumLib/TimeStepping/Algorithms/EvolutionaryPIDcontroller.cpp

@@ -13,10 +13,11 @@
 #include <functional>
 #include <limits>
 #include <vector>
+#include <logog/include/logog.hpp>
 
-#include "EvolutionaryPIDcontroller.h"
+#include "BaseLib/makeVectorUnique.h"
 
-#include <logog/include/logog.hpp>
+#include "EvolutionaryPIDcontroller.h"
 
 namespace NumLib
 {
@@ -142,18 +143,28 @@ double EvolutionaryPIDcontroller::limitStepSize(
 
 double EvolutionaryPIDcontroller::checkSpecificTimeReached(const double h_new)
 {
-    if (_specific_times.empty())
+    if (_fixed_output_times.empty())
         return h_new;
 
-    const double specific_time = _specific_times.back();
-    const double zero_threshold = std::numeric_limits<double>::epsilon();
+    const double specific_time = _fixed_output_times.back();
     if ((specific_time > _ts_current.current()) &&
-        (_ts_current.current() + h_new - specific_time > zero_threshold))
+        (_ts_current.current() + h_new - specific_time > 0.0))
     {
-        _specific_times.pop_back();
+        _fixed_output_times.pop_back();
         return specific_time - _ts_current.current();
     }
 
     return h_new;
 }
+
+void EvolutionaryPIDcontroller::addFixedOutputTimes(
+    std::vector<double> const& extra_fixed_output_times)
+{
+    _fixed_output_times.insert(_fixed_output_times.end(),
+                               extra_fixed_output_times.begin(),
+                               extra_fixed_output_times.end());
+
+    // Remove possible duplicated elements and sort in descending order.
+    BaseLib::makeVectorUnique(_fixed_output_times, std::greater<double>());
+}
 }  // end of namespace NumLib

NumLib/TimeStepping/Algorithms/EvolutionaryPIDcontroller.h

@@ -55,15 +55,15 @@ class EvolutionaryPIDcontroller final : public TimeStepAlgorithm
                               const double h0, const double h_min,
                               const double h_max, const double rel_h_min,
                               const double rel_h_max,
-                              const std::vector<double>&& specific_times,
+                              const std::vector<double>&& fixed_output_times,
                               const double tol)
         : TimeStepAlgorithm(t0, t_end),
           _h0(h0),
           _h_min(h_min),
           _h_max(h_max),
           _rel_h_min(rel_h_min),
           _rel_h_max(rel_h_max),
-          _specific_times(std::move(specific_times)),
+          _fixed_output_times(std::move(fixed_output_times)),
           _tol(tol),
           _e_n_minus1(0.),
           _e_n_minus2(0.),
@@ -91,6 +91,11 @@ class EvolutionaryPIDcontroller final : public TimeStepAlgorithm
     /// Get a flag to indicate that this algorithm need to compute
     /// solution error.
     bool isSolutionErrorComputationNeeded() override { return true; }
+
+    /// \copydoc NumLib::TimeStepAlgorithm::addFixedOutputTimes
+    void addFixedOutputTimes(
+        std::vector<double> const& extra_fixed_output_times) override;
+
 private:
     const double _kP = 0.075;  ///< Parameter. \see EvolutionaryPIDcontroller
     const double _kI = 0.175;  ///< Parameter. \see EvolutionaryPIDcontroller
@@ -106,7 +111,7 @@ class EvolutionaryPIDcontroller final : public TimeStepAlgorithm
     const double _rel_h_max;
 
     // Given times that steps have to reach.
-    std::vector<double> _specific_times;
+    std::vector<double> _fixed_output_times;
 
     const double _tol;
 

NumLib/TimeStepping/Algorithms/TimeStepAlgorithm.h

@@ -85,6 +85,19 @@ class TimeStepAlgorithm
     /// Get a flag to indicate whether this algorithm needs to compute
     /// solution error. The default return value is false.
     virtual bool isSolutionErrorComputationNeeded() { return false; }
+
+    /**
+     * Add specified times to the existing vector of the specified times.
+     * If there are specified times, they will be used as constraints in the
+     * computing of time step size such that the time step can exactly reach at
+     * the specified times. The function is mainly used to accept the specified
+     * times from the configuration of output.
+     */
+    virtual void addFixedOutputTimes(
+        std::vector<double> const& /*fixed_output_times*/)
+    {
+    }
+
 protected:
     /// initial time
     const double _t_initial;

ProcessLib/Output/CreateOutput.cpp

@@ -15,6 +15,7 @@
 
 #include "BaseLib/ConfigTree.h"
 #include "BaseLib/FileTools.h"
+#include "BaseLib/makeVectorUnique.h"
 
 #include "Output.h"
 
@@ -43,6 +44,8 @@ std::unique_ptr<Output> createOutput(const BaseLib::ConfigTree& config,
     // Construction of output times
     std::vector<Output::PairRepeatEachSteps> repeats_each_steps;
 
+    std::vector<double> fixed_output_times;
+
     //! \ogs_file_param{prj__time_loop__output__timesteps}
     if (auto const timesteps = config.getConfigSubtreeOptional("timesteps"))
     {
@@ -71,13 +74,25 @@ std::unique_ptr<Output> createOutput(const BaseLib::ConfigTree& config,
         repeats_each_steps.emplace_back(1, 1);
     }
 
+    auto fixed_output_times_ptr =
+        //! \ogs_file_param{prj__time_loop__output__fixed_output_times}
+        config.getConfigParameterOptional<std::vector<double>>(
+            "fixed_output_times");
+    if (fixed_output_times_ptr)
+    {
+        fixed_output_times = std::move(*fixed_output_times_ptr);
+        // Remove possible duplicated elements and sort in descending order.
+        BaseLib::makeVectorUnique(fixed_output_times, std::greater<double>());
+    }
+
     bool const output_iteration_results =
         //! \ogs_file_param{prj__time_loop__output__output_iteration_results}
         config.getConfigParameter<bool>("output_iteration_results", false);
 
     return std::make_unique<Output>(output_directory, prefix, compress_output,
                                     data_mode, output_iteration_results,
-                                    std::move(repeats_each_steps));
+                                    std::move(repeats_each_steps),
+                                    std::move(fixed_output_times));
 }
 
 }  // namespace ProcessLib

ProcessLib/Output/Output.cpp

@@ -22,29 +22,6 @@
 
 namespace
 {
-//! Determines if there should be output at the given \c timestep.
-template <typename CountsSteps>
-bool shallDoOutput(unsigned timestep, CountsSteps const& repeats_each_steps)
-{
-    unsigned each_steps = 1;
-
-    for (auto const& pair : repeats_each_steps)
-    {
-        each_steps = pair.each_steps;
-
-        if (timestep > pair.repeat * each_steps)
-        {
-            timestep -= pair.repeat * each_steps;
-        }
-        else
-        {
-            break;
-        }
-    }
-
-    return timestep % each_steps == 0;
-}
-
 //! Converts a vtkXMLWriter's data mode string to an int. See
 /// Output::_output_file_data_mode.
 int convertVtkDataMode(std::string const& data_mode)
@@ -70,16 +47,52 @@ int convertVtkDataMode(std::string const& data_mode)
 
 namespace ProcessLib
 {
+bool Output::shallDoOutput(unsigned timestep, double const t)
+{
+    unsigned each_steps = 1;
+
+    for (auto const& pair : _repeats_each_steps)
+    {
+        each_steps = pair.each_steps;
+
+        if (timestep > pair.repeat * each_steps)
+        {
+            timestep -= pair.repeat * each_steps;
+        }
+        else
+        {
+            break;
+        }
+    }
+
+    bool make_output = timestep % each_steps == 0;
+
+    if (_fixed_output_times.empty())
+        return make_output;
+
+    const double specific_time = _fixed_output_times.back();
+    const double zero_threshold = std::numeric_limits<double>::min();
+    if (std::fabs(specific_time - t) < zero_threshold)
+    {
+        _fixed_output_times.pop_back();
+        make_output = true;
+    }
+
+    return make_output;
+}
+
 Output::Output(std::string output_directory, std::string prefix,
                bool const compress_output, std::string const& data_mode,
                bool const output_nonlinear_iteration_results,
-               std::vector<PairRepeatEachSteps> repeats_each_steps)
+               std::vector<PairRepeatEachSteps> repeats_each_steps,
+               std::vector<double>&& fixed_output_times)
     : _output_directory(std::move(output_directory)),
       _output_file_prefix(std::move(prefix)),
       _output_file_compression(compress_output),
       _output_file_data_mode(convertVtkDataMode(data_mode)),
       _output_nonlinear_iteration_results(output_nonlinear_iteration_results),
-      _repeats_each_steps(std::move(repeats_each_steps))
+      _repeats_each_steps(std::move(repeats_each_steps)),
+      _fixed_output_times(std::move(fixed_output_times))
 {
 }
 
@@ -168,7 +181,7 @@ void Output::doOutput(Process const& process,
                       const double t,
                       GlobalVector const& x)
 {
-    if (shallDoOutput(timestep, _repeats_each_steps))
+    if (shallDoOutput(timestep, t))
     {
         doOutputAlways(process, process_id, process_output, timestep, t, x);
     }
@@ -186,7 +199,7 @@ void Output::doOutputLastTimestep(Process const& process,
                                   const double t,
                                   GlobalVector const& x)
 {
-    if (!shallDoOutput(timestep, _repeats_each_steps))
+    if (!shallDoOutput(timestep, t))
     {
         doOutputAlways(process, process_id, process_output, timestep, t, x);
     }

ProcessLib/Output/Output.h

@@ -42,7 +42,8 @@ class Output
     Output(std::string output_directory, std::string prefix,
            bool const compress_output, std::string const& data_mode,
            bool const output_nonlinear_iteration_results,
-           std::vector<PairRepeatEachSteps> repeats_each_steps);
+           std::vector<PairRepeatEachSteps> repeats_each_steps,
+           std::vector<double>&& fixed_output_times);
 
     //! TODO doc. Opens a PVD file for each process.
     void addProcess(ProcessLib::Process const& process, const int process_id);
@@ -77,6 +78,8 @@ class Output
                                     GlobalVector const& x,
                                     const unsigned iteration);
 
+    std::vector<double> getFixedOutputTimes() {return _fixed_output_times;}
+
 private:
     struct ProcessData
     {
@@ -100,6 +103,9 @@ class Output
     //! Describes after which timesteps to write output.
     std::vector<PairRepeatEachSteps> _repeats_each_steps;
 
+    //! Given times that steps have to reach.
+    std::vector<double> _fixed_output_times;
+
     std::multimap<Process const*, ProcessData> _process_to_process_data;
 
     /**
@@ -109,6 +115,11 @@ class Output
      * @return Address of a ProcessData.
      */
     ProcessData* findProcessData(Process const& process, const int process_id);
+
+    //! Determines if there should be output at the given \c timestep or \c t.
+    bool shallDoOutput(unsigned timestep, double const t);
 };
 
+
+
 }  // namespace ProcessLib

ProcessLib/RichardsFlow/Tests.cmake

@@ -37,6 +37,12 @@ AddTest(
     TESTER vtkdiff
     DIFF_DATA
     ref_t_1600.000000.vtu richards_pcs_0_ts_803_t_1600.000000.vtu pressure pressure 1e-8 1e-3
+# The following three comparisons are used just to check whether the output is
+# made at the fixed times of 50, 100 and 500, which are given in the project
+# file of RichardsFlow_2d_small_adaptive_dt.prj
+    richards_pcs_0_ts_28_spec_t_50.000000.vtu richards_pcs_0_ts_28_t_50.000000.vtu pressure pressure 1e-10 1e-10
+    richards_pcs_0_ts_53_spec_t_100.000000.vtu richards_pcs_0_ts_53_t_100.000000.vtu pressure pressure 1e-10 1e-10
+    richards_pcs_0_ts_253_spec_t_500.000000.vtu richards_pcs_0_ts_253_t_500.000000.vtu pressure pressure 1e-10 1e-10
     REQUIREMENTS NOT OGS_USE_MPI
 )
 

ProcessLib/UncoupledProcessesTimeLoop.cpp

@@ -515,7 +515,7 @@ bool UncoupledProcessesTimeLoop::loop()
 {
     // initialize output, convergence criterion, etc.
     {
-        unsigned process_id = 0;
+        int process_id = 0;
         for (auto& process_data : _per_process_data)
         {
             auto& pcs = process_data->process;
@@ -532,6 +532,13 @@ bool UncoupledProcessesTimeLoop::loop()
                                        pcs.getMesh());
             }
 
+            // Add the fixed times of output to time stepper in order that
+            // the time stepping is performed and the results are output at
+            // these times. Note: only the adaptive time steppers can have the
+            // the fixed times.
+            auto& timestepper = process_data->timestepper;
+            timestepper->addFixedOutputTimes(_output->getFixedOutputTimes());
+
             ++process_id;
         }
     }

Tests/Data/Parabolic/HT/StaggeredCoupling/ConstViscosity/square_5500x5500_staggered_scheme_adaptive_dt.prj

@@ -100,7 +100,7 @@
                     <rel_dt_min> 0.1 </rel_dt_min>
                     <rel_dt_max> 10 </rel_dt_max>
                     <tol> 10.0 </tol>
-                    <specific_times>5.0e9 2.e10</specific_times>
+                    <fixed_output_times>5.0e9 2.e10</fixed_output_times>
                 </time_stepping>
             </process>
             <process ref="ConstViscosityThermalConvection">
@@ -128,7 +128,7 @@
                     <rel_dt_min> 0.1 </rel_dt_min>
                     <rel_dt_max> 10 </rel_dt_max>
                     <tol> 10.0 </tol>
-                    <specific_times>5.0e9 2.e10</specific_times>
+                    <fixed_output_times>5.0e9 2.e10</fixed_output_times>
                 </time_stepping>
             </process>
         </processes>

Tests/Data/Parabolic/Richards/RichardsFlow_2d_small_adaptive_dt.prj

@@ -109,6 +109,7 @@
                     <each_steps>100000000</each_steps>
                 </pair>
             </timesteps>
+            <fixed_output_times> 50.0 100.0 500.</fixed_output_times>
             <output_iteration_results>false</output_iteration_results>
         </output>
     </time_loop>