Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <cstddef>
7 : #include <iosfwd>
8 : #include <limits>
9 : #include <optional>
10 : #include <string>
11 :
12 : #include "NumericalAlgorithms/Convergence/Criteria.hpp"
13 : #include "NumericalAlgorithms/Convergence/Reason.hpp"
14 :
15 : /// \cond
16 : namespace PUP {
17 : class er;
18 : } // namespace PUP
19 : /// \endcond
20 :
21 : namespace Convergence {
22 :
23 : /*!
24 : * \brief Determine whether the `criteria` are met.
25 : *
26 : * \note This function assumes the `iteration_id` is that of the latest
27 : * completed step and that it is zero-indexed, where zero indicates the initial
28 : * state of the algorithm. Therefore, the `MaxIteration` criterion will match if
29 : * the `iteration_id` is equal or higher. For example, a `MaxIteration` of 0
30 : * means the algorithm should run no iterations, so it matches if the
31 : * `iteration_id` is 0 or higher since that's the initial state before any
32 : * steps have been performed. A `MaxIteration` of 1 matches if the
33 : * `iteration_id` is 1 or higher since one iteration is complete. At this point,
34 : * also the `residual_magnitude` reflects the state of the algorithm after
35 : * completion of the first iteration. The `initial_residual_magnitude` always
36 : * refers to the state before the first iteration has begun, i.e. where the
37 : * `iteration_id` is zero.
38 : *
39 : * \returns A `Convergence::Reason` if the criteria are met, or
40 : * `std::nullopt` otherwise. The possible convergence reasons are:
41 : *
42 : * - `Convergence::Reason::AbsoluteResidual` if the `residual_magnitude`
43 : * meets the convergence criteria's `absolute_residual`.
44 : * - `Convergence::Reason::RelativeResidual` if `residual_magnitude /
45 : * initial_residual_magnitude` meets the convergence criteria's
46 : * `relative_residual`.
47 : * - `Convergence::Reason::MaxIterations` if the `iteration_id` is the
48 : * convergence_criteria's `max_iterations` or higher. This is often
49 : * interpreted as an error because the algorithm did not converge in the
50 : * alloted number of iterations.
51 : */
52 1 : std::optional<Reason> criteria_match(const Criteria& criteria,
53 : size_t iteration_id,
54 : double residual_magnitude,
55 : double initial_residual_magnitude);
56 :
57 : /*!
58 : * \brief Signals convergence or termination of the algorithm.
59 : *
60 : * \details Evaluates to `true` if the algorithm has converged or terminated and
61 : * no further iterations should be performed. In this case, the `reason()`
62 : * member function provides more information. If `false`, calling `reason()` is
63 : * an error.
64 : *
65 : * The stream operator provides a human-readable description of the convergence
66 : * status.
67 : *
68 : * This type default-constructs to a state that signals the algorithm has
69 : * not yet converged.
70 : */
71 1 : struct HasConverged {
72 : public:
73 0 : HasConverged() = default;
74 : /*!
75 : * \brief Determine whether the \p criteria are met by means of
76 : * `Convergence::criteria_match`.
77 : */
78 1 : HasConverged(const Criteria& criteria, size_t iteration_id,
79 : double residual_magnitude, double initial_residual_magnitude);
80 :
81 : /// Construct at a state where `iteration_id` iterations of a total of
82 : /// `num_iterations` have completed. Use when the algorithm is intended to run
83 : /// for a fixed number of iterations. The convergence `reason()` will be
84 : /// `Convergence::Reason::NumIterations`.
85 1 : HasConverged(size_t num_iterations, size_t iteration_id);
86 :
87 : /*!
88 : * \brief Construct a state manually.
89 : *
90 : * Currently only allows to construct a state with `Reason::Error`. Use the
91 : * other constructors for the other `Reason`s.
92 : */
93 1 : HasConverged(Reason reason, std::optional<std::string> error_message,
94 : size_t iteration_id);
95 :
96 0 : explicit operator bool() const { return static_cast<bool>(reason_); }
97 :
98 : /*!
99 : * \brief The reason the algorithm has converged.
100 : *
101 : * \warning Calling this function is an error if the algorithm has not yet
102 : * converged.
103 : */
104 1 : Reason reason() const;
105 :
106 : /*!
107 : * \brief Error message if reason is `Reason::Error`
108 : *
109 : * \warning You may only call this function if the `reason()` is
110 : * `Reason::Error`.
111 : */
112 1 : const std::string& error_message() const;
113 :
114 : /// The number of iterations the algorithm has completed
115 1 : size_t num_iterations() const;
116 :
117 : /// The residual magnitude after the last iteration. NaN if no iteration has
118 : /// completed yet.
119 1 : double residual_magnitude() const;
120 :
121 : /// The residual magnitude before the first iteration. NaN if this information
122 : /// is not available yet.
123 1 : double initial_residual_magnitude() const;
124 :
125 : /// Throw an exception if the `reason()` is `Reason::Error`.
126 1 : void check_for_error() const;
127 :
128 : // NOLINTNEXTLINE(google-runtime-references)
129 0 : void pup(PUP::er& p);
130 :
131 0 : friend bool operator==(const HasConverged& lhs, const HasConverged& rhs);
132 0 : friend bool operator!=(const HasConverged& lhs, const HasConverged& rhs);
133 :
134 0 : friend std::ostream& operator<<(std::ostream& os,
135 : const HasConverged& has_converged);
136 :
137 : private:
138 0 : std::optional<Reason> reason_{};
139 0 : std::optional<std::string> error_message_ = std::nullopt;
140 0 : Criteria criteria_{};
141 0 : size_t iteration_id_{};
142 0 : double residual_magnitude_ = std::numeric_limits<double>::quiet_NaN();
143 0 : double initial_residual_magnitude_ = std::numeric_limits<double>::quiet_NaN();
144 : };
145 :
146 : } // namespace Convergence
|