1// This file is part of Eigen, a lightweight C++ template library
2// for linear algebra.
3//
4// Copyright (C) 2008-2009 Gael Guennebaud <gael.guennebaud@inria.fr>
5//
6// This Source Code Form is subject to the terms of the Mozilla
7// Public License v. 2.0. If a copy of the MPL was not distributed
8// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9
10/* NOTE The class IterationController has been adapted from the iteration
11 *      class of the GMM++ and ITL libraries.
12 */
13
14//=======================================================================
15// Copyright (C) 1997-2001
16// Authors: Andrew Lumsdaine <lums@osl.iu.edu>
17//          Lie-Quan Lee     <llee@osl.iu.edu>
18//
19// This file is part of the Iterative Template Library
20//
21// You should have received a copy of the License Agreement for the
22// Iterative Template Library along with the software;  see the
23// file LICENSE.
24//
25// Permission to modify the code and to distribute modified code is
26// granted, provided the text of this NOTICE is retained, a notice that
27// the code was modified is included with the above COPYRIGHT NOTICE and
28// with the COPYRIGHT NOTICE in the LICENSE file, and that the LICENSE
29// file is distributed with the modified code.
30//
31// LICENSOR MAKES NO REPRESENTATIONS OR WARRANTIES, EXPRESS OR IMPLIED.
32// By way of example, but not limitation, Licensor MAKES NO
33// REPRESENTATIONS OR WARRANTIES OF MERCHANTABILITY OR FITNESS FOR ANY
34// PARTICULAR PURPOSE OR THAT THE USE OF THE LICENSED SOFTWARE COMPONENTS
35// OR DOCUMENTATION WILL NOT INFRINGE ANY PATENTS, COPYRIGHTS, TRADEMARKS
36// OR OTHER RIGHTS.
37//=======================================================================
38
39//========================================================================
40//
41// Copyright (C) 2002-2007 Yves Renard
42//
43// This file is a part of GETFEM++
44//
45// Getfem++ is free software; you can redistribute it and/or modify
46// it under the terms of the GNU Lesser General Public License as
47// published by the Free Software Foundation; version 2.1 of the License.
48//
49// This program is distributed in the hope that it will be useful,
50// but WITHOUT ANY WARRANTY; without even the implied warranty of
51// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
52// GNU Lesser General Public License for more details.
53// You should have received a copy of the GNU Lesser General Public
54// License along with this program; if not, write to the Free Software
55// Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301,
56// USA.
57//
58//========================================================================
59
60#include "../../../../Eigen/src/Core/util/NonMPL2.h"
61
62#ifndef EIGEN_ITERATION_CONTROLLER_H
63#define EIGEN_ITERATION_CONTROLLER_H
64
65namespace Eigen {
66
67/** \ingroup IterativeSolvers_Module
68  * \class IterationController
69  *
70  * \brief Controls the iterations of the iterative solvers
71  *
72  * This class has been adapted from the iteration class of GMM++ and ITL libraries.
73  *
74  */
75class IterationController
76{
77  protected :
78    double m_rhsn;        ///< Right hand side norm
79    size_t m_maxiter;     ///< Max. number of iterations
80    int m_noise;          ///< if noise > 0 iterations are printed
81    double m_resmax;      ///< maximum residual
82    double m_resminreach, m_resadd;
83    size_t m_nit;         ///< iteration number
84    double m_res;         ///< last computed residual
85    bool m_written;
86    void (*m_callback)(const IterationController&);
87  public :
88
89    void init()
90    {
91      m_nit = 0; m_res = 0.0; m_written = false;
92      m_resminreach = 1E50; m_resadd = 0.0;
93      m_callback = 0;
94    }
95
96    IterationController(double r = 1.0E-8, int noi = 0, size_t mit = size_t(-1))
97      : m_rhsn(1.0), m_maxiter(mit), m_noise(noi), m_resmax(r) { init(); }
98
99    void operator ++(int) { m_nit++; m_written = false; m_resadd += m_res; }
100    void operator ++() { (*this)++; }
101
102    bool first() { return m_nit == 0; }
103
104    /* get/set the "noisyness" (verbosity) of the solvers */
105    int noiseLevel() const { return m_noise; }
106    void setNoiseLevel(int n) { m_noise = n; }
107    void reduceNoiseLevel() { if (m_noise > 0) m_noise--; }
108
109    double maxResidual() const { return m_resmax; }
110    void setMaxResidual(double r) { m_resmax = r; }
111
112    double residual() const { return m_res; }
113
114    /* change the user-definable callback, called after each iteration */
115    void setCallback(void (*t)(const IterationController&))
116    {
117      m_callback = t;
118    }
119
120    size_t iteration() const { return m_nit; }
121    void setIteration(size_t i) { m_nit = i; }
122
123    size_t maxIterarions() const { return m_maxiter; }
124    void setMaxIterations(size_t i) { m_maxiter = i; }
125
126    double rhsNorm() const { return m_rhsn; }
127    void setRhsNorm(double r) { m_rhsn = r; }
128
129    bool converged() const { return m_res <= m_rhsn * m_resmax; }
130    bool converged(double nr)
131    {
132      m_res = internal::abs(nr);
133      m_resminreach = (std::min)(m_resminreach, m_res);
134      return converged();
135    }
136    template<typename VectorType> bool converged(const VectorType &v)
137    { return converged(v.squaredNorm()); }
138
139    bool finished(double nr)
140    {
141      if (m_callback) m_callback(*this);
142      if (m_noise > 0 && !m_written)
143      {
144        converged(nr);
145        m_written = true;
146      }
147      return (m_nit >= m_maxiter || converged(nr));
148    }
149    template <typename VectorType>
150    bool finished(const MatrixBase<VectorType> &v)
151    { return finished(double(v.squaredNorm())); }
152
153};
154
155} // end namespace Eigen
156
157#endif // EIGEN_ITERATION_CONTROLLER_H
158