From c69b63f07899e308c8e5d90e85a4b33eaf0b260b Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Wed, 20 Nov 2013 09:54:38 +0000 Subject: [PATCH] Introduce IterationNumberControl to perform a fixed number of solver steps without throwing an exception. git-svn-id: https://svn.dealii.org/trunk@31722 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/lac/solver_control.h | 52 ++++++++++++ deal.II/source/lac/solver_control.cc | 37 +++++++++ tests/lac/iteration_number_control.cc | 85 ++++++++++++++++++++ tests/lac/iteration_number_control.output | 6 ++ 4 files changed, 180 insertions(+) create mode 100644 tests/lac/iteration_number_control.cc create mode 100644 tests/lac/iteration_number_control.output diff --git a/deal.II/include/deal.II/lac/solver_control.h b/deal.II/include/deal.II/lac/solver_control.h index f303cddcf8..80d998fd86 100644 --- a/deal.II/include/deal.II/lac/solver_control.h +++ b/deal.II/include/deal.II/lac/solver_control.h @@ -462,6 +462,58 @@ protected: double reduced_tol; }; +/** + * Specialization of @p SolverControl which returns @p success if a given + * number of iteration was performed, irrespective of the actual + * residual. This is useful in cases where you don't want to solve exactly, + * but rather want to perform a fixed number of iterations, e.g. in an inner + * solver. The arguments given to this class are exactly the same as for the + * SolverControl class and the solver terminates similarly when one of the + * given tolerance or the maximum iteration count were reached. The only + * difference to SolverControl is that the solver returns success in the + * latter case. + * + * @author Martin Kronbichler + */ +class IterationNumberControl : public SolverControl +{ +public: + /** + * Constructor. Provide exactly the same arguments as the constructor of + * the SolverControl class. + */ + IterationNumberControl (const unsigned int maxiter = 100, + const double tolerance = 1e-12, + const bool log_history = false, + const bool log_result = true); + + /** + * Initialize with a SolverControl object. The result will emulate + * SolverControl by setting #reduce to zero. + */ + IterationNumberControl (const SolverControl &c); + + /** + * Assign a SolverControl object to ReductionControl. The result of the + * assignment will emulate SolverControl by setting #reduce to zero. + */ + IterationNumberControl &operator= (const SolverControl &c); + + /** + * Virtual destructor is needed as there are virtual functions in this + * class. + */ + virtual ~IterationNumberControl(); + + /** + * Decide about success or failure of an iteration. This function bases + * success solely on the fact if a given number of iterations was reached or + * the check value reached exactly zero. + */ + virtual State check (const unsigned int step, + const double check_value); +}; + /*@}*/ //--------------------------------------------------------------------------- diff --git a/deal.II/source/lac/solver_control.cc b/deal.II/source/lac/solver_control.cc index 8f40d154a4..6eab7e4d4f 100644 --- a/deal.II/source/lac/solver_control.cc +++ b/deal.II/source/lac/solver_control.cc @@ -303,4 +303,41 @@ ReductionControl::parse_parameters (ParameterHandler ¶m) set_reduction (param.get_double("Reduction")); } + +/*---------------------- IterationNumberControl -----------------------------*/ + + +IterationNumberControl::IterationNumberControl(const unsigned int n, + const double tolerance, + const bool m_log_history, + const bool m_log_result) + : + SolverControl (n, tolerance, m_log_history, m_log_result){} + + +IterationNumberControl::~IterationNumberControl() +{} + + +SolverControl::State +IterationNumberControl::check (const unsigned int step, + const double check_value) +{ + // check whether the given number of iterations was reached, and return + // success in that case. Otherwise, go on to the check of the base class. + if (step >= this->maxsteps) + { + if (m_log_result) + deallog << "Convergence step " << step + << " value " << check_value << std::endl; + lstep = step; + lvalue = check_value; + + lcheck = success; + return success; + } + else + return SolverControl::check(step, check_value); +} + DEAL_II_NAMESPACE_CLOSE diff --git a/tests/lac/iteration_number_control.cc b/tests/lac/iteration_number_control.cc new file mode 100644 index 0000000000..26ff41d21a --- /dev/null +++ b/tests/lac/iteration_number_control.cc @@ -0,0 +1,85 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 1998 - 2013 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +// check that IterationNumber control works as expected, i.e., it terminates +// with success when reaching the maximum number of iterations + + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include + + +void check() +{ + // Check 1: should still terminate after one iteration when using a matrix + // with only one eigenvalue + const unsigned int size = 20; + FullMatrix A(size, size); + for (unsigned int i=0; i in(size), out(size); + in = 1.; + + { + IterationNumberControl control(5); + SolverCG<> solver(control); + solver.solve(A, out, in, PreconditionIdentity()); + Assert (control.last_step() == 1, ExcInternalError()); + } + for (unsigned int i=0; i solver(control); + solver.solve(A, out, in, PreconditionIdentity()); + Assert (control.last_step() == 5, ExcInternalError()); + } + bool solved_exactly = true; + for (unsigned int i=0; i 1e-8) + solved_exactly = false; + Assert(solved_exactly == false, ExcInternalError()); + + deallog << "OK" << std::endl; +} + +int main() +{ + std::ofstream logfile("output"); +// logfile.setf(std::ios::fixed); + deallog << std::setprecision(4); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + check(); +} + diff --git a/tests/lac/iteration_number_control.output b/tests/lac/iteration_number_control.output new file mode 100644 index 0000000000..e2c05a877e --- /dev/null +++ b/tests/lac/iteration_number_control.output @@ -0,0 +1,6 @@ + +DEAL:cg::Starting value 4.472 +DEAL:cg::Convergence step 1 value 0 +DEAL:cg::Starting value 4.472 +DEAL:cg::Convergence step 5 value 0.6308 +DEAL::OK -- 2.39.5