Reference documentation for deal.II version GIT relicensing-233-g802318d791 2024-03-28 20:20:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
solver_control.cc
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 1998 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
18
20
21#include <cmath>
22#include <sstream>
23
25
26/*----------------------- SolverControl ---------------------------------*/
27
28
29SolverControl::SolverControl(const unsigned int maxiter,
30 const double tolerance,
31 const bool m_log_history,
32 const bool m_log_result)
33 : maxsteps(maxiter)
34 , tol(tolerance)
35 , lcheck(failure)
36 , initial_val(numbers::signaling_nan<double>())
37 , lvalue(numbers::signaling_nan<double>())
38 , lstep(numbers::invalid_unsigned_int)
39 , check_failure(false)
40 , relative_failure_residual(0)
41 , failure_residual(0)
42 , m_log_history(m_log_history)
43 , m_log_frequency(1)
44 , m_log_result(m_log_result)
45 , history_data_enabled(false)
46{}
47
48
49
51SolverControl::check(const unsigned int step, const double check_value)
52{
53 // if this is the first time we
54 // come here, then store the
55 // residual for later comparisons
56 if (step == 0)
57 {
58 initial_val = check_value;
59
61 history_data.clear();
62 }
63
64 if (m_log_history && ((step % m_log_frequency) == 0))
65 deallog << "Check " << step << "\t" << check_value << std::endl;
66
67 lstep = step;
68 lvalue = check_value;
69
70 if (step == 0)
71 {
72 if (check_failure)
74
75 if (m_log_result)
76 deallog << "Starting value " << check_value << std::endl;
77 }
78
80 history_data.push_back(check_value);
81
82 if (check_value <= tol)
83 {
84 if (m_log_result)
85 deallog << "Convergence step " << step << " value " << check_value
86 << std::endl;
88 return success;
89 }
90
91 if ((step >= maxsteps) || numbers::is_nan(check_value) ||
92 (check_failure && (check_value > failure_residual)))
93 {
94 if (m_log_result)
95 deallog << "Failure step " << step << " value " << check_value
96 << std::endl;
98 return failure;
99 }
100
101 lcheck = iterate;
102 return iterate;
103}
104
105
106
109{
110 return lcheck;
111}
112
113
114double
116{
117 return initial_val;
118}
119
120
121double
123{
124 return lvalue;
125}
126
127
128unsigned int
130{
131 return lstep;
132}
133
134
135unsigned int
137{
138 if (f == 0)
139 f = 1;
140 unsigned int old = m_log_frequency;
141 m_log_frequency = f;
142 return old;
143}
144
145
146void
151
152
153
154const std::vector<double> &
156{
158 Assert(
159 history_data.size() > 0,
161 "The SolverControl object was asked for the solver history "
162 "data, but there is no data. Possibly you requested the data before the "
163 "solver was run."));
164
165 return history_data;
166}
167
168
169
170double
183
184
185
186double
187SolverControl::step_reduction(unsigned int step) const
188{
191 Assert(step <= lstep, ExcIndexRange(step, 1, lstep + 1));
192 Assert(step > 0, ExcIndexRange(step, 1, lstep + 1));
193
194 return history_data[step] / history_data[step - 1];
195}
196
197
198double
203
204
205void
207{
208 param.declare_entry("Max steps", "100", Patterns::Integer());
209 param.declare_entry("Tolerance", "1.e-10", Patterns::Double());
210 param.declare_entry("Log history", "false", Patterns::Bool());
211 param.declare_entry("Log frequency", "1", Patterns::Integer());
212 param.declare_entry("Log result", "true", Patterns::Bool());
213}
214
215
216void
218{
219 set_max_steps(param.get_integer("Max steps"));
220 set_tolerance(param.get_double("Tolerance"));
221 log_history(param.get_bool("Log history"));
222 log_result(param.get_bool("Log result"));
223 log_frequency(param.get_integer("Log frequency"));
224}
225
226/*----------------------- ReductionControl ---------------------------------*/
227
228
230 const double tol,
231 const double red,
232 const bool m_log_history,
233 const bool m_log_result)
234 : SolverControl(n, tol, m_log_history, m_log_result)
235 , reduce(red)
236 , reduced_tol(numbers::signaling_nan<double>())
237{}
238
239
241 : SolverControl(c)
242 , reduce(numbers::signaling_nan<double>())
243 , reduced_tol(numbers::signaling_nan<double>())
244{
245 set_reduction(0.);
246}
247
248
251{
253 set_reduction(0.);
254 return *this;
255}
256
257
258
260ReductionControl::check(const unsigned int step, const double check_value)
261{
262 // if this is the first time we
263 // come here, then store the
264 // residual for later comparisons
265 if (step == 0)
266 {
267 initial_val = check_value;
268 reduced_tol = check_value * reduce;
269
271 history_data.clear();
272 }
273
274 // check whether desired reduction
275 // has been achieved. also check
276 // for equality in case initial
277 // residual already was zero
278 if (check_value <= reduced_tol)
279 {
280 if (m_log_result)
281 deallog << "Convergence step " << step << " value " << check_value
282 << std::endl;
283 lstep = step;
284 lvalue = check_value;
285 lcheck = success;
286
288 history_data.push_back(check_value);
289
290 return success;
291 }
292 else
293 return SolverControl::check(step, check_value);
294}
295
296
297
298void
304
305
306void
312
313
314/*---------------------- IterationNumberControl -----------------------------*/
315
316
318 const double tolerance,
319 const bool m_log_history,
320 const bool m_log_result)
321 : SolverControl(n, tolerance, m_log_history, m_log_result)
322{}
323
324
325
327IterationNumberControl::check(const unsigned int step, const double check_value)
328{
329 // check whether the given number of iterations was reached, and return
330 // success in that case. Otherwise, go on to the check of the base class.
331 if (step >= this->maxsteps)
332 {
333 if (m_log_result)
334 deallog << "Convergence step " << step << " value " << check_value
335 << std::endl;
336 lstep = step;
337 lvalue = check_value;
338
339 lcheck = success;
340 return success;
341 }
342 else
343 return SolverControl::check(step, check_value);
344}
345
346/*------------------------ ConsecutiveControl -------------------------------*/
347
348
350 const unsigned int n,
351 const double tolerance,
352 const unsigned int n_consecutive_iterations,
353 const bool m_log_history,
354 const bool m_log_result)
355 : SolverControl(n, tolerance, m_log_history, m_log_result)
356 , n_consecutive_iterations(n_consecutive_iterations)
357 , n_converged_iterations(0)
358{
360 ExcMessage("n_consecutive_iterations should be positive"));
361}
362
363
364
366 : SolverControl(c)
367 , n_consecutive_iterations(1)
368 , n_converged_iterations(0)
369{}
370
371
372
381
382
383
385ConsecutiveControl::check(const unsigned int step, const double check_value)
386{
387 // reset the counter if ConsecutiveControl is being reused
388 if (step == 0)
390 else
391 {
392 // check two things:
393 // (i) steps are ascending without repetitions
394 // (ii) user started from zero even when solver is being reused.
395 Assert(step - 1 == lstep,
396 ExcMessage("steps should be ascending integers."));
397 }
398
399 SolverControl::State state = SolverControl::check(step, check_value);
400 // check if we need to override the success:
401 if (state == success)
402 {
405 {
406 return success;
407 }
408 else
409 {
410 lcheck = iterate;
411 return iterate;
412 }
413 }
414 else
415 {
417 return state;
418 }
419}
420
ConsecutiveControl & operator=(const SolverControl &c)
ConsecutiveControl(const unsigned int maxiter=100, const double tolerance=1.e-10, const unsigned int n_consecutive_iterations=2, const bool log_history=false, const bool log_result=false)
unsigned int n_consecutive_iterations
unsigned int n_converged_iterations
virtual State check(const unsigned int step, const double check_value) override
virtual State check(const unsigned int step, const double check_value) override
IterationNumberControl(const unsigned int maxiter=100, const double tolerance=1e-12, const bool log_history=false, const bool log_result=true)
long int get_integer(const std::string &entry_string) const
bool get_bool(const std::string &entry_name) const
void declare_entry(const std::string &entry, const std::string &default_value, const Patterns::PatternBase &pattern=Patterns::Anything(), const std::string &documentation="", const bool has_to_be_set=false)
double get_double(const std::string &entry_name) const
double set_reduction(const double)
void parse_parameters(ParameterHandler &param)
virtual State check(const unsigned int step, const double check_value) override
static void declare_parameters(ParameterHandler &param)
ReductionControl & operator=(const SolverControl &c)
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)
static void declare_parameters(ParameterHandler &param)
State last_check() const
double average_reduction() const
const std::vector< double > & get_history_data() const
unsigned int lstep
unsigned int last_step() const
double last_value() const
unsigned int m_log_frequency
unsigned int log_frequency(unsigned int)
virtual State check(const unsigned int step, const double check_value)
unsigned int maxsteps
bool log_history() const
void enable_history_data()
double step_reduction(unsigned int step) const
std::vector< double > history_data
double failure_residual
double relative_failure_residual
unsigned int set_max_steps(const unsigned int)
double set_tolerance(const double)
SolverControl(const unsigned int n=100, const double tol=1.e-10, const bool log_history=false, const bool log_result=true)
double initial_value() const
bool log_result() const
void parse_parameters(ParameterHandler &param)
double final_reduction() const
@ iterate
Continue iteration.
@ success
Stop iteration, goal reached.
@ failure
Stop iteration, goal not reached.
Subscriptor & operator=(const Subscriptor &)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
#define Assert(cond, exc)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcHistoryDataRequired()
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
LogStream deallog
Definition logstream.cc:36
bool is_nan(const double x)
Definition numbers.h:530
::VectorizedArray< Number, width > pow(const ::VectorizedArray< Number, width > &, const Number p)