#include <base/subscriptor.h>
+class ParameterHandler;
/**
* 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);
*/
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
*/
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.
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;
};
* 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);
*/
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
*/
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.
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
+
+
#include <base/logstream.h>
+#include <base/parameter_handler.h>
#include <lac/solver_control.h>
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)
{};
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;
if (check_value <= tol)
{
- if (log_result)
+ if (_log_result)
deallog << "Convergence step " << step
<< " value " << check_value << endl;
return success;
}
return iterate;
-};
+}
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)
{};
if (check_value < reduced_tol)
{
- if (log_result)
+ if (_log_result)
deallog << "Convergence step " << step
<< " value " << check_value << endl;
return success;
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"));
+}
+