Reference documentation for deal.II version GIT c14369f203 2023-10-01 07:40: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.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2022 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 #ifndef dealii_solver_control_h
17 #define dealii_solver_control_h
18 
19 
20 #include <deal.II/base/config.h>
21 
23 
24 #include <vector>
25 
27 
28 #ifndef DOXYGEN
29 class ParameterHandler;
30 #endif
31 
67 class SolverControl : public Subscriptor
68 {
69 public:
74  enum State
75  {
77  iterate = 0,
81  failure
82  };
83 
84 
85 
96  {
97  public:
98  NoConvergence(const unsigned int last_step, const double last_residual)
101  {}
102 
103  virtual ~NoConvergence() noexcept override = default;
104 
105  virtual void
106  print_info(std::ostream &out) const override
107  {
108  out
109  << "Iterative method reported convergence failure in step " << last_step
110  << ". The residual in the last step was " << last_residual << ".\n\n"
111  << "This error message can indicate that you have simply not allowed "
112  << "a sufficiently large number of iterations for your iterative solver "
113  << "to converge. This often happens when you increase the size of your "
114  << "problem. In such cases, the last residual will likely still be very "
115  << "small, and you can make the error go away by increasing the allowed "
116  << "number of iterations when setting up the SolverControl object that "
117  << "determines the maximal number of iterations you allow."
118  << "\n\n"
119  << "The other situation where this error may occur is when your matrix "
120  << "is not invertible (e.g., your matrix has a null-space), or if you "
121  << "try to apply the wrong solver to a matrix (e.g., using CG for a "
122  << "matrix that is not symmetric or not positive definite). In these "
123  << "cases, the residual in the last iteration is likely going to be large."
124  << std::endl;
125  }
126 
130  const unsigned int last_step;
131 
135  const double last_residual;
136  };
137 
138 
139 
151  explicit SolverControl(const unsigned int n = 100,
152  const double tol = 1.e-10,
153  const bool log_history = false,
154  const bool log_result = true);
155 
160  virtual ~SolverControl() override = default;
161 
165  static void
167 
171  void
173 
191  virtual State
192  check(const unsigned int step, const double check_value);
193 
197  State
198  last_check() const;
199 
203  double
204  initial_value() const;
205 
210  double
211  last_value() const;
212 
216  unsigned int
217  last_step() const;
218 
222  unsigned int
223  max_steps() const;
224 
228  unsigned int
229  set_max_steps(const unsigned int);
230 
236  void
237  set_failure_criterion(const double rel_failure_residual);
238 
243  void
245 
249  double
250  tolerance() const;
251 
255  double
256  set_tolerance(const double);
257 
261  void
263 
267  const std::vector<double> &
268  get_history_data() const;
269 
275  double
276  average_reduction() const;
283  double
284  final_reduction() const;
285 
291  double
292  step_reduction(unsigned int step) const;
293 
297  void
298  log_history(const bool);
299 
303  bool
304  log_history() const;
305 
309  unsigned int
310  log_frequency(unsigned int);
311 
315  void
316  log_result(const bool);
317 
321  bool
322  log_result() const;
323 
330 
331 protected:
335  unsigned int maxsteps;
336 
340  double tol;
341 
346 
350  double initial_val;
351 
355  double lvalue;
356 
360  unsigned int lstep;
361 
367 
372 
380 
385 
389  unsigned int m_log_frequency;
390 
397 
402 
409  std::vector<double> history_data;
410 };
411 
412 
425 {
426 public:
432  explicit ReductionControl(const unsigned int maxiter = 100,
433  const double tolerance = 1.e-10,
434  const double reduce = 1.e-2,
435  const bool log_history = false,
436  const bool log_result = true);
437 
442  explicit ReductionControl(const SolverControl &c);
443 
449  operator=(const SolverControl &c);
450 
455  virtual ~ReductionControl() override = default;
456 
460  static void
462 
466  void
468 
474  virtual State
475  check(const unsigned int step, const double check_value) override;
476 
480  double
481  reduction() const;
482 
486  double
487  set_reduction(const double);
488 
489 protected:
493  double reduce;
494 
499  double reduced_tol;
500 };
501 
513 {
514 public:
519  explicit IterationNumberControl(const unsigned int maxiter = 100,
520  const double tolerance = 1e-12,
521  const bool log_history = false,
522  const bool log_result = true);
523 
529 
537 
542  virtual ~IterationNumberControl() override = default;
543 
549  virtual State
550  check(const unsigned int step, const double check_value) override;
551 };
552 
553 
566 {
567 public:
574  explicit ConsecutiveControl(const unsigned int maxiter = 100,
575  const double tolerance = 1.e-10,
576  const unsigned int n_consecutive_iterations = 2,
577  const bool log_history = false,
578  const bool log_result = false);
579 
584  explicit ConsecutiveControl(const SolverControl &c);
585 
592  operator=(const SolverControl &c);
593 
598  virtual ~ConsecutiveControl() override = default;
599 
604  virtual State
605  check(const unsigned int step, const double check_value) override;
606 
607 protected:
613 
618 };
619 
621 //---------------------------------------------------------------------------
622 
623 #ifndef DOXYGEN
624 
625 inline unsigned int
627 {
628  return maxsteps;
629 }
630 
631 
632 
633 inline unsigned int
634 SolverControl::set_max_steps(const unsigned int newval)
635 {
636  unsigned int old = maxsteps;
637  maxsteps = newval;
638  return old;
639 }
640 
641 
642 
643 inline void
644 SolverControl::set_failure_criterion(const double rel_failure_residual)
645 {
646  relative_failure_residual = rel_failure_residual;
647  check_failure = true;
648 }
649 
650 
651 
652 inline void
654 {
656  failure_residual = 0;
657  check_failure = false;
658 }
659 
660 
661 
662 inline double
664 {
665  return tol;
666 }
667 
668 
669 
670 inline double
671 SolverControl::set_tolerance(const double t)
672 {
673  double old = tol;
674  tol = t;
675  return old;
676 }
677 
678 
679 inline void
680 SolverControl::log_history(const bool newval)
681 {
682  m_log_history = newval;
683 }
684 
685 
686 
687 inline bool
689 {
690  return m_log_history;
691 }
692 
693 
694 inline void
695 SolverControl::log_result(const bool newval)
696 {
697  m_log_result = newval;
698 }
699 
700 
701 inline bool
703 {
704  return m_log_result;
705 }
706 
707 
708 inline double
710 {
711  return reduce;
712 }
713 
714 
715 inline double
716 ReductionControl::set_reduction(const double t)
717 {
718  double old = reduce;
719  reduce = t;
720  return old;
721 }
722 
723 #endif // DOXYGEN
724 
726 
727 #endif
ConsecutiveControl & operator=(const SolverControl &c)
virtual ~ConsecutiveControl() override=default
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)
IterationNumberControl(const SolverControl &c)
IterationNumberControl & operator=(const SolverControl &c)
virtual ~IterationNumberControl() override=default
double set_reduction(const double)
void parse_parameters(ParameterHandler &param)
virtual ~ReductionControl() override=default
double reduction() const
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)
const unsigned int last_step
NoConvergence(const unsigned int last_step, const double last_residual)
virtual void print_info(std::ostream &out) const override
virtual ~NoConvergence() noexcept override=default
void log_history(const bool)
static void declare_parameters(ParameterHandler &param)
State last_check() const
double average_reduction() const
virtual ~SolverControl() override=default
const std::vector< double > & get_history_data() const
void clear_failure_criterion()
unsigned int lstep
unsigned int last_step() const
double last_value() const
void log_result(const bool)
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
void set_failure_criterion(const double rel_failure_residual)
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)
unsigned int max_steps() const
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
double tolerance() 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
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
#define DeclException0(Exception0)
Definition: exceptions.h:467
static ::ExceptionBase & ExcHistoryDataRequired()
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)