Reference documentation for deal.II version GIT d6cf33b98c 2023-09-22 19:50: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\}}\)
solver_control.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2023 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #include <deal.II/base/logstream.h>
19 
21 
22 #include <cmath>
23 #include <sstream>
24 
26 
27 /*----------------------- SolverControl ---------------------------------*/
28 
29 
30 SolverControl::SolverControl(const unsigned int maxiter,
31  const double tolerance,
32  const bool m_log_history,
33  const bool m_log_result)
34  : maxsteps(maxiter)
35  , tol(tolerance)
36  , lcheck(failure)
37  , initial_val(numbers::signaling_nan<double>())
38  , lvalue(numbers::signaling_nan<double>())
40  , check_failure(false)
41  , relative_failure_residual(0)
42  , failure_residual(0)
43  , m_log_history(m_log_history)
44  , m_log_frequency(1)
45  , m_log_result(m_log_result)
46  , history_data_enabled(false)
47 {}
48 
49 
50 
52 SolverControl::check(const unsigned int step, const double check_value)
53 {
54  // if this is the first time we
55  // come here, then store the
56  // residual for later comparisons
57  if (step == 0)
58  {
59  initial_val = check_value;
60 
62  history_data.clear();
63  }
64 
65  if (m_log_history && ((step % m_log_frequency) == 0))
66  deallog << "Check " << step << "\t" << check_value << std::endl;
67 
68  lstep = step;
69  lvalue = check_value;
70 
71  if (step == 0)
72  {
73  if (check_failure)
75 
76  if (m_log_result)
77  deallog << "Starting value " << check_value << std::endl;
78  }
79 
81  history_data.push_back(check_value);
82 
83  if (check_value <= tol)
84  {
85  if (m_log_result)
86  deallog << "Convergence step " << step << " value " << check_value
87  << std::endl;
88  lcheck = success;
89  return success;
90  }
91 
92  if ((step >= maxsteps) || numbers::is_nan(check_value) ||
93  (check_failure && (check_value > failure_residual)))
94  {
95  if (m_log_result)
96  deallog << "Failure step " << step << " value " << check_value
97  << std::endl;
98  lcheck = failure;
99  return failure;
100  }
101 
102  lcheck = iterate;
103  return iterate;
104 }
105 
106 
107 
110 {
111  return lcheck;
112 }
113 
114 
115 double
117 {
118  return initial_val;
119 }
120 
121 
122 double
124 {
125  return lvalue;
126 }
127 
128 
129 unsigned int
131 {
132  return lstep;
133 }
134 
135 
136 unsigned int
138 {
139  if (f == 0)
140  f = 1;
141  unsigned int old = m_log_frequency;
142  m_log_frequency = f;
143  return old;
144 }
145 
146 
147 void
149 {
150  history_data_enabled = true;
151 }
152 
153 
154 
155 const std::vector<double> &
157 {
159  Assert(
160  history_data.size() > 0,
161  ExcMessage(
162  "The SolverControl object was asked for the solver history "
163  "data, but there is no data. Possibly you requested the data before the "
164  "solver was run."));
165 
166  return history_data;
167 }
168 
169 
170 
171 double
173 {
174  if (lstep == 0)
175  return 0.;
176 
181 
182  return std::pow(history_data[lstep] / history_data[0], 1. / lstep);
183 }
184 
185 
186 
187 double
188 SolverControl::step_reduction(unsigned int step) const
189 {
192  Assert(step <= lstep, ExcIndexRange(step, 1, lstep + 1));
193  Assert(step > 0, ExcIndexRange(step, 1, lstep + 1));
194 
195  return history_data[step] / history_data[step - 1];
196 }
197 
198 
199 double
201 {
202  return step_reduction(lstep);
203 }
204 
205 
206 void
208 {
209  param.declare_entry("Max steps", "100", Patterns::Integer());
210  param.declare_entry("Tolerance", "1.e-10", Patterns::Double());
211  param.declare_entry("Log history", "false", Patterns::Bool());
212  param.declare_entry("Log frequency", "1", Patterns::Integer());
213  param.declare_entry("Log result", "true", Patterns::Bool());
214 }
215 
216 
217 void
219 {
220  set_max_steps(param.get_integer("Max steps"));
221  set_tolerance(param.get_double("Tolerance"));
222  log_history(param.get_bool("Log history"));
223  log_result(param.get_bool("Log result"));
224  log_frequency(param.get_integer("Log frequency"));
225 }
226 
227 /*----------------------- ReductionControl ---------------------------------*/
228 
229 
231  const double tol,
232  const double red,
233  const bool m_log_history,
234  const bool m_log_result)
235  : SolverControl(n, tol, m_log_history, m_log_result)
236  , reduce(red)
237  , reduced_tol(numbers::signaling_nan<double>())
238 {}
239 
240 
242  : SolverControl(c)
244  , reduced_tol(numbers::signaling_nan<double>())
245 {
246  set_reduction(0.);
247 }
248 
249 
252 {
254  set_reduction(0.);
255  return *this;
256 }
257 
258 
259 
261 ReductionControl::check(const unsigned int step, const double check_value)
262 {
263  // if this is the first time we
264  // come here, then store the
265  // residual for later comparisons
266  if (step == 0)
267  {
268  initial_val = check_value;
269  reduced_tol = check_value * reduce;
270 
272  history_data.clear();
273  }
274 
275  // check whether desired reduction
276  // has been achieved. also check
277  // for equality in case initial
278  // residual already was zero
279  if (check_value <= reduced_tol)
280  {
281  if (m_log_result)
282  deallog << "Convergence step " << step << " value " << check_value
283  << std::endl;
284  lstep = step;
285  lvalue = check_value;
286  lcheck = success;
287 
289  history_data.push_back(check_value);
290 
291  return success;
292  }
293  else
294  return SolverControl::check(step, check_value);
295 }
296 
297 
298 
299 void
301 {
303  param.declare_entry("Reduction", "1.e-2", Patterns::Double());
304 }
305 
306 
307 void
309 {
311  set_reduction(param.get_double("Reduction"));
312 }
313 
314 
315 /*---------------------- IterationNumberControl -----------------------------*/
316 
317 
319  const double tolerance,
320  const bool m_log_history,
321  const bool m_log_result)
322  : SolverControl(n, tolerance, m_log_history, m_log_result)
323 {}
324 
325 
326 
328 IterationNumberControl::check(const unsigned int step, const double check_value)
329 {
330  // check whether the given number of iterations was reached, and return
331  // success in that case. Otherwise, go on to the check of the base class.
332  if (step >= this->maxsteps)
333  {
334  if (m_log_result)
335  deallog << "Convergence step " << step << " value " << check_value
336  << std::endl;
337  lstep = step;
338  lvalue = check_value;
339 
340  lcheck = success;
341  return success;
342  }
343  else
344  return SolverControl::check(step, check_value);
345 }
346 
347 /*------------------------ ConsecutiveControl -------------------------------*/
348 
349 
351  const unsigned int n,
352  const double tolerance,
353  const unsigned int n_consecutive_iterations,
354  const bool m_log_history,
355  const bool m_log_result)
356  : SolverControl(n, tolerance, m_log_history, m_log_result)
357  , n_consecutive_iterations(n_consecutive_iterations)
358  , n_converged_iterations(0)
359 {
361  ExcMessage("n_consecutive_iterations should be positive"));
362 }
363 
364 
365 
367  : SolverControl(c)
368  , n_consecutive_iterations(1)
369  , n_converged_iterations(0)
370 {}
371 
372 
373 
376 {
380  return *this;
381 }
382 
383 
384 
386 ConsecutiveControl::check(const unsigned int step, const double check_value)
387 {
388  // reset the counter if ConsecutiveControl is being reused
389  if (step == 0)
391  else
392  {
393  // check two things:
394  // (i) steps are ascending without repetitions
395  // (ii) user started from zero even when solver is being reused.
396  Assert(step - 1 == lstep,
397  ExcMessage("steps should be ascending integers."));
398  }
399 
400  SolverControl::State state = SolverControl::check(step, check_value);
401  // check if we need to override the success:
402  if (state == success)
403  {
406  {
407  return success;
408  }
409  else
410  {
411  lcheck = iterate;
412  return iterate;
413  }
414  }
415  else
416  {
418  return state;
419  }
420 }
421 
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.
bool history_data_enabled
Subscriptor & operator=(const Subscriptor &)
Definition: subscriptor.h:291
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
Definition: exceptions.h:1616
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)
Definition: exceptions.h:1705
LogStream deallog
Definition: logstream.cc:37
T reduce(const T &local_value, const MPI_Comm comm, const std::function< T(const T &, const T &)> &combiner, const unsigned int root_process=0)
static const unsigned int invalid_unsigned_int
Definition: types.h:213
T signaling_nan()
bool is_nan(const double x)
Definition: numbers.h:531