basic_iteration.hpp 3.78 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
// Software License for MTL
//
// Copyright (c) 2007 The Trustees of Indiana University.
//               2008 Dresden University of Technology and the Trustees of Indiana University.
//               2010 SimuNova UG (haftungsbeschränkt), www.simunova.com.
// All rights reserved.
// Authors: Peter Gottschling and Andrew Lumsdaine
//
// This file is part of the Matrix Template Library
//
// See also license.mtl.txt in the distribution.

#ifndef ITL_BASIC_ITERATION_INCLUDE
#define ITL_BASIC_ITERATION_INCLUDE

#include <iostream>
#include <complex>
#include <string>

namespace itl {


template <class Real>
class basic_iteration
{
  public:
    typedef basic_iteration self;
    typedef Real            real;

    template <class Vector>
    basic_iteration(const Vector& r0, int max_iter_, Real t, Real a = Real(0))
      : error(0), i(0), norm_r0(std::abs(two_norm(r0))),
	max_iter(max_iter_), rtol_(t), atol_(a), is_finished(false), my_quite(false), my_suppress(false) { }

    basic_iteration(Real nb, int max_iter_, Real t, Real a = Real(0))
36
      : error(0), i(0), norm_r0(nb), max_iter(max_iter_), rtol_(t), atol_(a), is_finished(false),
37
38
39
40
41
42
	my_quite(false), my_suppress(false) {}

    virtual ~basic_iteration() {}

    bool check_max()
    {
43
	if (i >= max_iter)
44
45
46
47
48
	    error= 1, is_finished= true, err_msg= "Too many iterations.";
	return is_finished;
    }

    template <class Vector>
49
    bool finished(const Vector& r)
50
51
52
53
54
55
    {
	if (converged(two_norm(r)))
	    return is_finished= true;
	return check_max();
    }

56
    bool finished(const Real& r)
57
58
59
60
61
62
63
    {
	if (converged(r))
	    return is_finished= true;
	return check_max();
    }

    template <typename T>
64
    bool finished(const std::complex<T>& r)
65
    {
66
	if (converged(std::abs(r)))
67
68
69
70
71
72
73
74
75
76
77
	    return is_finished= true;
	return check_max();
    }

    bool finished() const { return is_finished; }

    template <class T>
    int terminate(const T& r) { finished(r); return error; }

    bool converged(const Real& r) { resid_= r; return converged(); }

78
79
    bool converged() const
    {
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
	if (norm_r0 == 0)
	    return resid_ <= atol_;  // ignore relative tolerance if |r0| is zero
	return resid_ <= rtol_ * norm_r0 || resid_ <= atol_;
    }

    self& operator++() { ++i; return *this; }

    self& operator+=(int n) { i+= n; return *this; }

    bool first() const { return i <= 1; }

    virtual operator int() const { return error; }

    virtual int error_code() const { return error; }

    bool is_converged() const { return is_finished && error == 0; }

    int iterations() const { return i; }
98

99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
    int max_iterations() const { return max_iter; }

    void set_max_iterations(int m) { max_iter= m; }

    Real resid() const { return resid_; }

    Real relresid() const { return resid_ / norm_r0; }

    Real normb() const { return norm_r0; }

    Real tol() const { return rtol_; }
    Real atol() const { return atol_; }

    int fail(int err_code) { error = err_code; return error_code(); }

    int fail(int err_code, const std::string& msg)
    { error = err_code; err_msg = msg; return error_code(); }

    void set(Real v) { norm_r0 = v; }

    void set_quite(bool q) { my_quite= q; }

    bool is_quite() const { return my_quite; }

    void suppress_resume(bool s) { my_suppress= s; }

    bool resume_suppressed() const { return my_suppress; }

    void update_progress(const basic_iteration& that)
    {
	i= that.i;
	resid_= that.resid_;
	if (that.error > 1) { // copy error except too many iterations
	    error= that.error;
	    err_msg= that.err_msg;
	    is_finished= true;
135
	} else
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
	    finished(resid_);
    }

  protected:
    int          error, i;
    Real         norm_r0;
    int          max_iter;
    Real         rtol_, atol_, resid_;
    std::string  err_msg;
    bool         is_finished, my_quite, my_suppress;
};


} // namespace itl

#endif // ITL_BASIC_ITERATION_INCLUDE