From e71fbd080ef2b160e7d9bde2bdabb32aa6c36a75 Mon Sep 17 00:00:00 2001 From: guido Date: Sun, 21 May 2000 13:40:02 +0000 Subject: [PATCH] parameter for solver git-svn-id: https://svn.dealii.org/trunk@2912 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/include/lac/solver_control.h | 169 +++++++++++++++++++++-- deal.II/lac/source/solver_control.cc | 80 +++++++++-- 2 files changed, 220 insertions(+), 29 deletions(-) diff --git a/deal.II/lac/include/lac/solver_control.h b/deal.II/lac/include/lac/solver_control.h index 811dd969e7..7662cc2fac 100644 --- a/deal.II/lac/include/lac/solver_control.h +++ b/deal.II/lac/include/lac/solver_control.h @@ -15,6 +15,7 @@ #include +class ParameterHandler; /** @@ -82,7 +83,8 @@ class SolverControl : public Subscriptor * specifies the whether the final result is logged * to #deallog#. Default is yes. */ - SolverControl (const unsigned int n, const double tol, + SolverControl (const unsigned int n = 100, + const double tol = 1.e-10, const bool log_history = false, const bool log_result = true); @@ -93,6 +95,16 @@ class SolverControl : public Subscriptor */ virtual ~SolverControl(); + /** + * Interface to parameter file. + */ + static void declare_parameters (ParameterHandler& param); + + /** + * Read parameters from file. + */ + void parse_parameters (ParameterHandler& param); + /** * Decide about success or failure * of an iteration. This function @@ -135,16 +147,53 @@ class SolverControl : public Subscriptor */ unsigned int last_step() const; + /** + * Maximum number of steps. + */ + unsigned int max_steps () const; + + /** + * Change maximum number of steps. + */ + unsigned int set_max_steps (const unsigned int); + + /** + * Tolerance. + */ + double tolerance () const; + + /** + * Change tolerance. + */ + double set_tolerance (const double); + + /** + * Log each iteration step. Use + * @p{log_frequency} for skipping + * steps. + */ + void log_history (const bool); + + /** + * Set logging frequency. + */ + unsigned int log_frequency (unsigned int); + + /** + * Log start and end step. + */ + void log_result (const bool); + protected: /** * Maximum number of steps. */ - const unsigned int maxsteps; + unsigned int maxsteps; /** * Prescribed tolerance to be achieved. */ - const double tol; + double tol; /** * Last value of the convergence criterion. @@ -157,16 +206,23 @@ class SolverControl : public Subscriptor unsigned int lstep; /** - * Log convergence history to #deallog#? + * Log convergence history to #deallog#. + */ + bool _log_history; + /** + * Log only every nth step. */ - const bool log_history; + unsigned int _log_frequency; + /** - * Log iteration result to #deallog#? - * If true, after finishing the iteration, a - * statement about failure or success together with - * #lstep# and #lvalue# are logged. + * Log iteration result to + * #deallog#. If true, after + * finishing the iteration, a + * statement about failure or + * success together with #lstep# + * and #lvalue# are logged. */ - const bool log_result; + bool _log_result; }; @@ -191,9 +247,9 @@ class ReductionControl : public SolverControl * the arguments of the Control * constructor. */ - ReductionControl (const unsigned int maxiter, - const double tolerance, - const double reduce, + ReductionControl (const unsigned int maxiter = 100, + const double tolerance = 1.e-10, + const double reduce = 1.e-2, const bool log_history = false, const bool log_result = true); @@ -204,6 +260,16 @@ class ReductionControl : public SolverControl */ virtual ~ReductionControl(); + /** + * Interface to parameter file. + */ + static void declare_parameters (ParameterHandler& param); + + /** + * Read parameters from file. + */ + void parse_parameters (ParameterHandler& param); + /** * Decide about success or failure * of an iteration. This function @@ -221,12 +287,21 @@ class ReductionControl : public SolverControl */ double initial_value() const; + /** + * Reduction factor. + */ + double reduction () const; + + /** + * Change reduction factor. + */ + double set_reduction (const double); protected: /** * Desired reduction factor. */ - const double reduce; + double reduce; /** * Initial value. @@ -242,5 +317,71 @@ protected: double reduced_tol; }; +//------------------------------------------------------------// + + +inline unsigned int +SolverControl::max_steps () const +{ + return maxsteps; +} + + +inline unsigned int +SolverControl::set_max_steps (unsigned int newval) +{ + unsigned int old = maxsteps; + maxsteps = newval; + return old; +} + + +inline double +SolverControl::tolerance () const +{ + return tol; +} + + +inline double +SolverControl::set_tolerance (double t) +{ + double old = tol; + tol = t; + return old; +} + + +inline void +SolverControl::log_history (bool newval) +{ + _log_history = newval; +} + + + +inline void +SolverControl::log_result (bool newval) +{ + _log_result = newval; +} + + +inline double +ReductionControl::reduction () const +{ + return reduce; +} + + +inline double +ReductionControl::set_reduction (double t) +{ + double old = reduce; + reduce = t; + return old; +} #endif + + diff --git a/deal.II/lac/source/solver_control.cc b/deal.II/lac/source/solver_control.cc index 81d27b2a13..9cb22dfe4f 100644 --- a/deal.II/lac/source/solver_control.cc +++ b/deal.II/lac/source/solver_control.cc @@ -13,6 +13,7 @@ #include +#include #include @@ -21,15 +22,16 @@ SolverControl::SolverControl (const unsigned int maxiter, const double tolerance, - const bool log_history, - const bool log_result) + const bool _log_history, + const bool _log_result) : maxsteps(maxiter), tol(tolerance), lvalue(1.e300), lstep(0), - log_history(log_history), - log_result(log_result) + _log_history(_log_history), + _log_frequency(1), + _log_result(_log_result) {}; @@ -43,18 +45,18 @@ SolverControl::State SolverControl::check (const unsigned int step, const double check_value) { - if (log_history) + if (_log_history && ((step % _log_frequency) == 0)) deallog << "Check " << step << "\t" << check_value << endl; lstep = step; lvalue = check_value; - if ((step==0) && log_result) + if ((step==0) && _log_result) deallog << "Starting value " << check_value << endl; if (step >= maxsteps) { - if (log_result) + if (_log_result) deallog << "Failure step " << step << " value " << check_value << endl; return failure; @@ -62,13 +64,13 @@ SolverControl::check (const unsigned int step, if (check_value <= tol) { - if (log_result) + if (_log_result) deallog << "Convergence step " << step << " value " << check_value << endl; return success; } return iterate; -}; +} @@ -76,26 +78,57 @@ double SolverControl::last_value() const { return lvalue; -}; +} unsigned int SolverControl::last_step() const { return lstep; -}; +} + + +unsigned int +SolverControl::log_frequency (unsigned int f) +{ + if (f==0) + f = 1; + unsigned int old = _log_frequency; + _log_frequency = f; + return old; +} + + +void +SolverControl::declare_parameters (ParameterHandler& param) +{ + param.declare_entry ("Max steps", "100", Patterns::Integer()); + param.declare_entry ("Tolerance", "1.e-10", Patterns::Double()); + param.declare_entry ("Log history", "false", Patterns::Bool()); + param.declare_entry ("Log frequency", "1", Patterns::Integer()); + param.declare_entry ("Log result", "true", Patterns::Bool()); +} +void SolverControl::parse_parameters (ParameterHandler& param) +{ + set_max_steps (param.get_integer("Max steps")); + set_tolerance (param.get_double("Tolerance")); + log_history (param.get_bool("Log history")); + log_result (param.get_bool("Log result")); + log_frequency (param.get_integer("Log frequency")); +} + /*----------------------- ReductionControl ---------------------------------*/ ReductionControl::ReductionControl(const unsigned int n, const double tol, const double red, - const bool log_history, - const bool log_result) + const bool _log_history, + const bool _log_result) : - SolverControl (n, tol, log_history, log_result), + SolverControl (n, tol, _log_history, _log_result), reduce(red) {}; @@ -122,7 +155,7 @@ ReductionControl::check (const unsigned int step, if (check_value < reduced_tol) { - if (log_result) + if (_log_result) deallog << "Convergence step " << step << " value " << check_value << endl; return success; @@ -131,3 +164,20 @@ ReductionControl::check (const unsigned int step, return SolverControl::check(step, check_value); }; + + +void +ReductionControl::declare_parameters (ParameterHandler& param) +{ + SolverControl::declare_parameters (param); + param.declare_entry("Reduction", "1.e-2", Patterns::Double()); +} + + +void +ReductionControl::parse_parameters (ParameterHandler& param) +{ + SolverControl::parse_parameters (param); + set_reduction (param.get_double("Reduction")); +} + -- 2.39.5