Reference documentation for deal.II version Git 0096380e24 2020-08-06 21:03:09 -0400
\(\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 - 2019 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 
34 
65 class SolverControl : public Subscriptor
66 {
67 public:
72  enum State
73  {
75  iterate = 0,
80  };
81 
82 
83 
95  {
96  public:
97  NoConvergence(const unsigned int last_step, const double last_residual)
98  : last_step(last_step)
99  , last_residual(last_residual)
100  {}
101 
102  virtual ~NoConvergence() noexcept override = default;
103 
104  virtual void
105  print_info(std::ostream &out) const override
106  {
107  out
108  << "Iterative method reported convergence failure in step " << last_step
109  << ". The residual in the last step was " << last_residual << ".\n\n"
110  << "This error message can indicate that you have simply not allowed "
111  << "a sufficiently large number of iterations for your iterative solver "
112  << "to converge. This often happens when you increase the size of your "
113  << "problem. In such cases, the last residual will likely still be very "
114  << "small, and you can make the error go away by increasing the allowed "
115  << "number of iterations when setting up the SolverControl object that "
116  << "determines the maximal number of iterations you allow."
117  << "\n\n"
118  << "The other situation where this error may occur is when your matrix "
119  << "is not invertible (e.g., your matrix has a null-space), or if you "
120  << "try to apply the wrong solver to a matrix (e.g., using CG for a "
121  << "matrix that is not symmetric or not positive definite). In these "
122  << "cases, the residual in the last iteration is likely going to be large."
123  << std::endl;
124  }
125 
129  const unsigned int last_step;
130 
134  const double last_residual;
135  };
136 
137 
138 
150  SolverControl(const unsigned int n = 100,
151  const double tol = 1.e-10,
152  const bool log_history = false,
153  const bool log_result = true);
154 
159  virtual ~SolverControl() override = default;
160 
164  static void
166 
170  void
172 
190  virtual State
191  check(const unsigned int step, const double check_value);
192 
196  State
197  last_check() const;
198 
202  double
203  initial_value() const;
204 
209  double
210  last_value() const;
211 
215  unsigned int
216  last_step() const;
217 
221  unsigned int
222  max_steps() const;
223 
227  unsigned int
228  set_max_steps(const unsigned int);
229 
235  void
236  set_failure_criterion(const double rel_failure_residual);
237 
242  void
244 
248  double
249  tolerance() const;
250 
254  double
255  set_tolerance(const double);
256 
260  void
262 
266  const std::vector<double> &
267  get_history_data() const;
268 
274  double
275  average_reduction() const;
282  double
283  final_reduction() const;
284 
290  double
291  step_reduction(unsigned int step) const;
292 
296  void
297  log_history(const bool);
298 
302  bool
303  log_history() const;
304 
308  unsigned int
309  log_frequency(unsigned int);
310 
314  void
315  log_result(const bool);
316 
320  bool
321  log_result() const;
322 
329 
330 protected:
334  unsigned int maxsteps;
335 
339  double tol;
340 
345 
349  double initial_val;
350 
354  double lvalue;
355 
359  unsigned int lstep;
360 
366 
371 
379 
384 
388  unsigned int m_log_frequency;
389 
396 
401 
408  std::vector<double> history_data;
409 };
410 
411 
424 {
425 public:
431  ReductionControl(const unsigned int maxiter = 100,
432  const double tolerance = 1.e-10,
433  const double reduce = 1.e-2,
434  const bool log_history = false,
435  const bool log_result = true);
436 
442 
448  operator=(const SolverControl &c);
449 
454  virtual ~ReductionControl() override = default;
455 
459  static void
461 
465  void
467 
473  virtual State
474  check(const unsigned int step, const double check_value) override;
475 
479  double
480  reduction() const;
481 
485  double
486  set_reduction(const double);
487 
488 protected:
492  double reduce;
493 
498  double reduced_tol;
499 };
500 
512 {
513 public:
518  IterationNumberControl(const unsigned int maxiter = 100,
519  const double tolerance = 1e-12,
520  const bool log_history = false,
521  const bool log_result = true);
522 
528 
535  operator=(const SolverControl &c);
536 
541  virtual ~IterationNumberControl() override = default;
542 
548  virtual State
549  check(const unsigned int step, const double check_value) override;
550 };
551 
552 
565 {
566 public:
573  ConsecutiveControl(const unsigned int maxiter = 100,
574  const double tolerance = 1.e-10,
575  const unsigned int n_consecutive_iterations = 2,
576  const bool log_history = false,
577  const bool log_result = false);
578 
584 
591  operator=(const SolverControl &c);
592 
597  virtual ~ConsecutiveControl() override = default;
598 
603  virtual State
604  check(const unsigned int step, const double check_value) override;
605 
606 protected:
612 
617 };
618 
620 //---------------------------------------------------------------------------
621 
622 #ifndef DOXYGEN
623 
624 inline unsigned int
626 {
627  return maxsteps;
628 }
629 
630 
631 
632 inline unsigned int
633 SolverControl::set_max_steps(const unsigned int newval)
634 {
635  unsigned int old = maxsteps;
636  maxsteps = newval;
637  return old;
638 }
639 
640 
641 
642 inline void
643 SolverControl::set_failure_criterion(const double rel_failure_residual)
644 {
645  relative_failure_residual = rel_failure_residual;
646  check_failure = true;
647 }
648 
649 
650 
651 inline void
653 {
655  failure_residual = 0;
656  check_failure = false;
657 }
658 
659 
660 
661 inline double
663 {
664  return tol;
665 }
666 
667 
668 
669 inline double
670 SolverControl::set_tolerance(const double t)
671 {
672  double old = tol;
673  tol = t;
674  return old;
675 }
676 
677 
678 inline void
679 SolverControl::log_history(const bool newval)
680 {
681  m_log_history = newval;
682 }
683 
684 
685 
686 inline bool
688 {
689  return m_log_history;
690 }
691 
692 
693 inline void
694 SolverControl::log_result(const bool newval)
695 {
696  m_log_result = newval;
697 }
698 
699 
700 inline bool
702 {
703  return m_log_result;
704 }
705 
706 
707 inline double
709 {
710  return reduce;
711 }
712 
713 
714 inline double
715 ReductionControl::set_reduction(const double t)
716 {
717  double old = reduce;
718  reduce = t;
719  return old;
720 }
721 
722 #endif // DOXYGEN
723 
725 
726 #endif
NoConvergence(const unsigned int last_step, const double last_residual)
double initial_value() const
Stop iteration, goal not reached.
__global__ void reduction(Number *result, const Number *v, const size_type N)
bool log_history() const
virtual State check(const unsigned int step, const double check_value)
bool history_data_enabled
Continue iteration.
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
const std::vector< double > & get_history_data() const
static void declare_parameters(ParameterHandler &param)
unsigned int n_consecutive_iterations
double step_reduction(unsigned int step) const
double final_reduction() const
void clear_failure_criterion()
double reduction() const
double failure_residual
ExceptionBase operator=(const ExceptionBase &)=delete
const unsigned int last_step
unsigned int log_frequency(unsigned int)
Stop iteration, goal reached.
unsigned int m_log_frequency
double set_reduction(const double)
#define DeclException0(Exception0)
Definition: exceptions.h:470
static ::ExceptionBase & ExcHistoryDataRequired()
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
double relative_failure_residual
double last_value() const
double average_reduction() const
unsigned int lstep
double tolerance() const
virtual ~SolverControl() override=default
virtual ~NoConvergence() noexcept override=default
unsigned int n_converged_iterations
State last_check() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
unsigned int max_steps() const
void enable_history_data()
double set_tolerance(const double)
unsigned int maxsteps
void parse_parameters(ParameterHandler &param)
std::vector< double > history_data
virtual void print_info(std::ostream &out) const override
SolverControl(const unsigned int n=100, const double tol=1.e-10, const bool log_history=false, const bool log_result=true)
bool log_result() const
unsigned int set_max_steps(const unsigned int)
void set_failure_criterion(const double rel_failure_residual)