diff --git a/Documentation/ProjectFile/prj/time_loop/time_stepping/EvolutionaryPIDcontroller/t_norm_type.md b/Documentation/ProjectFile/prj/time_loop/time_stepping/EvolutionaryPIDcontroller/t_norm_type.md deleted file mode 100644 index de46a911c18183e2aa51bd519b3a5bfe54c72d67..0000000000000000000000000000000000000000 --- a/Documentation/ProjectFile/prj/time_loop/time_stepping/EvolutionaryPIDcontroller/t_norm_type.md +++ /dev/null @@ -1 +0,0 @@ -The norm type of the global vector for error calculation. diff --git a/NumLib/ODESolver/ConvergenceCriterion.h b/NumLib/ODESolver/ConvergenceCriterion.h index 154b23b4de4c4d7514a4507d1ad2670d45495528..59babdd8450e4a88623e1a2c5d183b7372f56b4d 100644 --- a/NumLib/ODESolver/ConvergenceCriterion.h +++ b/NumLib/ODESolver/ConvergenceCriterion.h @@ -10,6 +10,7 @@ #pragma once #include <memory> +#include "MathLib/LinAlg/LinAlg.h" // For MathLib::VecNormType #include "NumLib/NumericsConfig.h" namespace BaseLib { @@ -27,6 +28,11 @@ namespace NumLib class ConvergenceCriterion { public: + ConvergenceCriterion(const MathLib::VecNormType norm_type) + : _norm_type(norm_type) + { + } + //! Tells if the change of the solution between iterations is checked. //! //! \remark @@ -76,11 +82,14 @@ public: //! Tell if the convergence criterion is satisfied. virtual bool isSatisfied() const { return _satisfied; }; + MathLib::VecNormType getVectorNormType () const { return _norm_type; } + virtual ~ConvergenceCriterion() = default; protected: bool _satisfied = true; bool _is_first_iteration = true; + const MathLib::VecNormType _norm_type; }; //! Creates a convergence criterion from the given configuration. diff --git a/NumLib/ODESolver/ConvergenceCriterionDeltaX.cpp b/NumLib/ODESolver/ConvergenceCriterionDeltaX.cpp index 88de56434264079d10a16d5fb528cb703bf0a6ee..1bbece408622b75a3d9f16963e343b707dfc350c 100644 --- a/NumLib/ODESolver/ConvergenceCriterionDeltaX.cpp +++ b/NumLib/ODESolver/ConvergenceCriterionDeltaX.cpp @@ -18,10 +18,9 @@ namespace NumLib ConvergenceCriterionDeltaX::ConvergenceCriterionDeltaX( boost::optional<double>&& absolute_tolerance, boost::optional<double>&& relative_tolerance, - MathLib::VecNormType norm_type) - : _abstol(std::move(absolute_tolerance)), - _reltol(std::move(relative_tolerance)), - _norm_type(norm_type) + const MathLib::VecNormType norm_type) + : ConvergenceCriterion(norm_type), _abstol(std::move(absolute_tolerance)), + _reltol(std::move(relative_tolerance)) { if ((!_abstol) && (!_reltol)) OGS_FATAL( diff --git a/NumLib/ODESolver/ConvergenceCriterionDeltaX.h b/NumLib/ODESolver/ConvergenceCriterionDeltaX.h index 75ede7ec830da063c8cb22d4b466c8e17e12a87f..71715ca1950c2e25469e5ddace7d37d9a6ca1efd 100644 --- a/NumLib/ODESolver/ConvergenceCriterionDeltaX.h +++ b/NumLib/ODESolver/ConvergenceCriterionDeltaX.h @@ -27,7 +27,7 @@ public: ConvergenceCriterionDeltaX( boost::optional<double>&& absolute_tolerance, boost::optional<double>&& relative_tolerance, - MathLib::VecNormType norm_type); + const MathLib::VecNormType norm_type); bool hasDeltaXCheck() const override { return true; } bool hasResidualCheck() const override { return false; } @@ -40,7 +40,6 @@ public: private: const boost::optional<double> _abstol; const boost::optional<double> _reltol; - const MathLib::VecNormType _norm_type; }; std::unique_ptr<ConvergenceCriterionDeltaX> createConvergenceCriterionDeltaX( diff --git a/NumLib/ODESolver/ConvergenceCriterionPerComponent.h b/NumLib/ODESolver/ConvergenceCriterionPerComponent.h index 0ff42fdebf1c5c6f302110ce0df9158242e64f90..8da246ad60d1856e59dec0f516a5061665b79194 100644 --- a/NumLib/ODESolver/ConvergenceCriterionPerComponent.h +++ b/NumLib/ODESolver/ConvergenceCriterionPerComponent.h @@ -11,6 +11,8 @@ #include "ConvergenceCriterion.h" +#include "MathLib/LinAlg/LinAlg.h" // For MathLib::VecNormType + namespace MeshLib { class Mesh; @@ -26,6 +28,9 @@ class LocalToGlobalIndexMap; class ConvergenceCriterionPerComponent : public ConvergenceCriterion { public: + ConvergenceCriterionPerComponent(const MathLib::VecNormType norm_type) + :ConvergenceCriterion(norm_type) {} + //! Sets the d.o.f. table used to extract data for a specific component. virtual void setDOFTable(NumLib::LocalToGlobalIndexMap const& dof_table, MeshLib::Mesh const& mesh) = 0; diff --git a/NumLib/ODESolver/ConvergenceCriterionPerComponentDeltaX.cpp b/NumLib/ODESolver/ConvergenceCriterionPerComponentDeltaX.cpp index b9faaeed41ed0f4ae6e37633a7185bee5e30b930..f5ef433565e803e7d0f9a1ee6141cdbe1ee3f51e 100644 --- a/NumLib/ODESolver/ConvergenceCriterionPerComponentDeltaX.cpp +++ b/NumLib/ODESolver/ConvergenceCriterionPerComponentDeltaX.cpp @@ -19,10 +19,10 @@ namespace NumLib ConvergenceCriterionPerComponentDeltaX::ConvergenceCriterionPerComponentDeltaX( std::vector<double>&& absolute_tolerances, std::vector<double>&& relative_tolerances, - MathLib::VecNormType norm_type) - : _abstols(std::move(absolute_tolerances)), - _reltols(std::move(relative_tolerances)), - _norm_type(norm_type) + const MathLib::VecNormType norm_type) + : ConvergenceCriterionPerComponent(norm_type), + _abstols(std::move(absolute_tolerances)), + _reltols(std::move(relative_tolerances)) { if (_abstols.size() != _reltols.size()) OGS_FATAL( diff --git a/NumLib/ODESolver/ConvergenceCriterionPerComponentDeltaX.h b/NumLib/ODESolver/ConvergenceCriterionPerComponentDeltaX.h index a35ea708527e8c83dc54775637b05e33692b2180..3622a32e286d9b979749181d94a014e4b2eebdd6 100644 --- a/NumLib/ODESolver/ConvergenceCriterionPerComponentDeltaX.h +++ b/NumLib/ODESolver/ConvergenceCriterionPerComponentDeltaX.h @@ -29,7 +29,7 @@ public: ConvergenceCriterionPerComponentDeltaX( std::vector<double>&& absolute_tolerances, std::vector<double>&& relative_tolerances, - MathLib::VecNormType norm_type); + const MathLib::VecNormType norm_type); bool hasDeltaXCheck() const override { return true; } bool hasResidualCheck() const override { return false; } @@ -46,7 +46,6 @@ public: private: const std::vector<double> _abstols; const std::vector<double> _reltols; - const MathLib::VecNormType _norm_type; LocalToGlobalIndexMap const* _dof_table = nullptr; MeshLib::Mesh const* _mesh = nullptr; }; diff --git a/NumLib/ODESolver/ConvergenceCriterionPerComponentResidual.cpp b/NumLib/ODESolver/ConvergenceCriterionPerComponentResidual.cpp index 29600a819952605a9669b730bef056398eaa357e..acd515f7ccab2dc91f921428ac29f554d13ff568 100644 --- a/NumLib/ODESolver/ConvergenceCriterionPerComponentResidual.cpp +++ b/NumLib/ODESolver/ConvergenceCriterionPerComponentResidual.cpp @@ -20,10 +20,10 @@ ConvergenceCriterionPerComponentResidual:: ConvergenceCriterionPerComponentResidual( std::vector<double>&& absolute_tolerances, std::vector<double>&& relative_tolerances, - MathLib::VecNormType norm_type) - : _abstols(std::move(absolute_tolerances)), + const MathLib::VecNormType norm_type) + : ConvergenceCriterionPerComponent(norm_type), + _abstols(std::move(absolute_tolerances)), _reltols(std::move(relative_tolerances)), - _norm_type(norm_type), _residual_norms_0(_abstols.size()) { if (_abstols.size() != _reltols.size()) diff --git a/NumLib/ODESolver/ConvergenceCriterionPerComponentResidual.h b/NumLib/ODESolver/ConvergenceCriterionPerComponentResidual.h index f54c1a6a932e04919a0d4f4dcc8f5d5f7dd755a2..743a6fc5da5b63e9428de7d94cd63d312e721257 100644 --- a/NumLib/ODESolver/ConvergenceCriterionPerComponentResidual.h +++ b/NumLib/ODESolver/ConvergenceCriterionPerComponentResidual.h @@ -30,7 +30,7 @@ public: ConvergenceCriterionPerComponentResidual( std::vector<double>&& absolute_tolerances, std::vector<double>&& relative_tolerances, - MathLib::VecNormType norm_type); + const MathLib::VecNormType norm_type); bool hasDeltaXCheck() const override { return true; } bool hasResidualCheck() const override { return true; } @@ -47,7 +47,6 @@ public: private: const std::vector<double> _abstols; const std::vector<double> _reltols; - const MathLib::VecNormType _norm_type; LocalToGlobalIndexMap const* _dof_table = nullptr; MeshLib::Mesh const* _mesh = nullptr; std::vector<double> _residual_norms_0; diff --git a/NumLib/ODESolver/ConvergenceCriterionResidual.cpp b/NumLib/ODESolver/ConvergenceCriterionResidual.cpp index 0424d80288737849f5b4c1b538786cb9cb4d5323..e160ac6a4bf7a4385d711f0056fca15d41ddf2bb 100644 --- a/NumLib/ODESolver/ConvergenceCriterionResidual.cpp +++ b/NumLib/ODESolver/ConvergenceCriterionResidual.cpp @@ -18,10 +18,9 @@ namespace NumLib ConvergenceCriterionResidual::ConvergenceCriterionResidual( boost::optional<double>&& absolute_tolerance, boost::optional<double>&& relative_tolerance, - MathLib::VecNormType norm_type) - : _abstol(std::move(absolute_tolerance)), - _reltol(std::move(relative_tolerance)), - _norm_type(norm_type) + const MathLib::VecNormType norm_type) + : ConvergenceCriterion(norm_type), _abstol(std::move(absolute_tolerance)), + _reltol(std::move(relative_tolerance)) { if ((!_abstol) && (!_reltol)) OGS_FATAL( diff --git a/NumLib/ODESolver/ConvergenceCriterionResidual.h b/NumLib/ODESolver/ConvergenceCriterionResidual.h index b598cfbc83a3920dcbcdf4df61ba4473b5d57b41..b63af27e9b9af5010a7e7f297341b567b4db28d0 100644 --- a/NumLib/ODESolver/ConvergenceCriterionResidual.h +++ b/NumLib/ODESolver/ConvergenceCriterionResidual.h @@ -27,7 +27,7 @@ public: ConvergenceCriterionResidual( boost::optional<double>&& absolute_tolerance, boost::optional<double>&& relative_tolerance, - MathLib::VecNormType norm_type); + const MathLib::VecNormType norm_type); bool hasDeltaXCheck() const override { return true; } bool hasResidualCheck() const override { return true; } @@ -41,7 +41,6 @@ public: private: const boost::optional<double> _abstol; const boost::optional<double> _reltol; - const MathLib::VecNormType _norm_type; double _residual_norm_0; }; diff --git a/NumLib/ODESolver/TimeDiscretization.cpp b/NumLib/ODESolver/TimeDiscretization.cpp index c6619d7641893a94518f08d3936b0880b2adb8f7..e7cf70c71188bf60c952ccd3af9c1d07cf47ae3b 100644 --- a/NumLib/ODESolver/TimeDiscretization.cpp +++ b/NumLib/ODESolver/TimeDiscretization.cpp @@ -135,4 +135,4 @@ void BackwardDifferentiationFormula::getWeightedOldX(GlobalVector& y) const LinAlg::scale(y, 1.0 / _delta_t); } -} // end of namespace NumLib \ No newline at end of file +} // end of namespace NumLib diff --git a/NumLib/TimeStepping/Algorithms/CreateEvolutionaryPIDcontroller.cpp b/NumLib/TimeStepping/Algorithms/CreateEvolutionaryPIDcontroller.cpp index b677ea5af52ec6658c2bc8caafe26d682e92a56b..84b63208c16ec97563017fcaeec4d48e414b0abc 100644 --- a/NumLib/TimeStepping/Algorithms/CreateEvolutionaryPIDcontroller.cpp +++ b/NumLib/TimeStepping/Algorithms/CreateEvolutionaryPIDcontroller.cpp @@ -58,15 +58,9 @@ std::unique_ptr<TimeStepAlgorithm> createEvolutionaryPIDcontroller( //! \ogs_file_param{prj__time_loop__time_stepping__EvolutionaryPIDcontroller__tol} auto const tol = config.getConfigParameter<double>("tol"); - auto const norm_type_opt = - //! \ogs_file_param{prj__time_loop__time_stepping__EvolutionaryPIDcontroller__norm_type} - config.getConfigParameterOptional<std::string>("norm_type"); - const MathLib::VecNormType norm_type = - (norm_type_opt) ? MathLib::convertStringToVecNormType(*norm_type_opt) - : MathLib::VecNormType::NORM2; return std::make_unique<EvolutionaryPIDcontroller>( t0, t_end, h0, h_min, h_max, rel_h_min, rel_h_max, - std::move(specific_times), tol, norm_type); + std::move(specific_times), tol); } } // end of namespace NumLib diff --git a/NumLib/TimeStepping/Algorithms/EvolutionaryPIDcontroller.h b/NumLib/TimeStepping/Algorithms/EvolutionaryPIDcontroller.h index 73eb690ffc607a85fea82201d36b7587317f2503..5e5226667f64bc5624bbbd2968be39d5247cd14a 100644 --- a/NumLib/TimeStepping/Algorithms/EvolutionaryPIDcontroller.h +++ b/NumLib/TimeStepping/Algorithms/EvolutionaryPIDcontroller.h @@ -56,8 +56,7 @@ public: const double h_max, const double rel_h_min, const double rel_h_max, const std::vector<double>&& specific_times, - const double tol, - const MathLib::VecNormType norm_type) + const double tol) : TimeStepAlgorithm(t0, t_end), _h0(h0), _h_min(h_min), @@ -68,8 +67,7 @@ public: _tol(tol), _e_n_minus1(0.), _e_n_minus2(0.), - _is_accepted(true), - _norm_type(norm_type) + _is_accepted(true) { } @@ -86,10 +84,6 @@ public: /// Get a flag to indicate that this algorithm need to compute /// solution error. bool isSolutionErrorComputationNeeded() override { return true; } - MathLib::VecNormType getSolutionNormType() const override - { - return _norm_type; - } private: const double _kP = 0.075; ///< Parameter. \see EvolutionaryPIDcontroller @@ -115,9 +109,6 @@ private: bool _is_accepted; - /// Type of the norm of the solution vector. - const MathLib::VecNormType _norm_type; - double limitStepSize(const double h_new, const double h_n) const { const double h_in_range = std::max(_h_min, std::min(h_new, _h_max)); diff --git a/NumLib/TimeStepping/Algorithms/TimeStepAlgorithm.h b/NumLib/TimeStepping/Algorithms/TimeStepAlgorithm.h index 4c00eed7a434ed7bad6ea68a43758a9752e24ba0..1c459e5b7972963f64f992988d2e963cd53a32c1 100644 --- a/NumLib/TimeStepping/Algorithms/TimeStepAlgorithm.h +++ b/NumLib/TimeStepping/Algorithms/TimeStepAlgorithm.h @@ -16,8 +16,6 @@ #include "NumLib/TimeStepping/TimeStep.h" -#include "MathLib/LinAlg/LinAlg.h" // For MathLib::VecNormType - namespace NumLib { /** @@ -80,11 +78,6 @@ public: return _dt_vector; } - virtual MathLib::VecNormType getSolutionNormType() const - { - return MathLib::VecNormType::NORM2; - } - /// Get a flag to indicate whether this algorithm needs to compute /// solution error. The default return value is false. virtual bool isSolutionErrorComputationNeeded() { return false; } diff --git a/ProcessLib/UncoupledProcessesTimeLoop.cpp b/ProcessLib/UncoupledProcessesTimeLoop.cpp index 7c58ee812fd5095fe8a77eb361b98293f965a7a4..e8b122ca6a79706d77a373a1f1d0657cff094ff9 100644 --- a/ProcessLib/UncoupledProcessesTimeLoop.cpp +++ b/ProcessLib/UncoupledProcessesTimeLoop.cpp @@ -559,12 +559,17 @@ double UncoupledProcessesTimeLoop::computeTimeStepping( auto& time_disc = ppd.time_disc; auto const& x = *_process_solutions[i]; + auto const& conv_crit = ppd.conv_crit; + const MathLib::VecNormType norm_type = + (conv_crit) ? conv_crit->getVectorNormType() + : MathLib::VecNormType::NORM2; + const double solution_error = (timestepper->isSolutionErrorComputationNeeded()) ? ((t == timestepper->begin()) - ? 0. // Always accepts the first step + ? 0. // Always accepts the zeroth step : time_disc->getRelativeChangeFromPreviousTimestep( - x, timestepper->getSolutionNormType())) + x, norm_type)) : 0.; if (!timestepper->next(solution_error) && // In case of FixedTimeStepping, which makes timestepper->next(...) diff --git a/Tests/Data b/Tests/Data index cf593d3828f045ce0341c3e1c2ed5c8046ba9126..ca1db9c602b485f66cead0e82924181c26b86539 160000 --- a/Tests/Data +++ b/Tests/Data @@ -1 +1 @@ -Subproject commit cf593d3828f045ce0341c3e1c2ed5c8046ba9126 +Subproject commit ca1db9c602b485f66cead0e82924181c26b86539 diff --git a/Tests/NumLib/TestTimeSteppingEvolutionaryPIDcontroller.cpp b/Tests/NumLib/TestTimeSteppingEvolutionaryPIDcontroller.cpp index 2fe8060513c0c4dbdd030850ce0fafe5c1577a01..1e509148fcf6b0f8545a63f5141ba51edc79e4e8 100644 --- a/Tests/NumLib/TestTimeSteppingEvolutionaryPIDcontroller.cpp +++ b/Tests/NumLib/TestTimeSteppingEvolutionaryPIDcontroller.cpp @@ -46,7 +46,6 @@ TEST(NumLibTimeStepping, testEvolutionaryPIDcontroller) " <rel_dt_min> 0.01 </rel_dt_min>" " <rel_dt_max> 5 </rel_dt_max>" " <tol> 1.e-3 </tol>" - " <norm_type> NORM2 </norm_type>" "</time_stepping>"; auto const PIDStepper = createTestTimeStepper(xml);