Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
solver_control.h
1 // ---------------------------------------------------------------------
2 // @f$Id: solver_control.h 31722 2013-11-20 09:54:38Z kronbichler @f$
3 //
4 // Copyright (C) 1998 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__solver_control_h
18 #define __deal2__solver_control_h
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/base/subscriptor.h>
23 
24 #include <vector>
25 
26 DEAL_II_NAMESPACE_OPEN
27 
28 class ParameterHandler;
29 
32 
66 class SolverControl : public Subscriptor
67 {
68 public:
69 
74  enum State
75  {
77  iterate = 0,
82  };
83 
84 
85 
97  {
98  public:
99  NoConvergence (const unsigned int last_step,
100  const double last_residual)
101  : last_step (last_step), last_residual(last_residual)
102  {}
103 
104  virtual ~NoConvergence () throw () {}
105 
106  virtual void print_info (std::ostream &out) const
107  {
108  out << "Iterative method reported convergence failure in step "
109  << last_step << " with residual " << last_residual << std::endl;
110  }
111 
115  const unsigned int last_step;
116 
120  const double last_residual;
121  };
122 
123 
124 
136  SolverControl (const unsigned int n = 100,
137  const double tol = 1.e-10,
138  const bool log_history = false,
139  const bool log_result = true);
140 
145  virtual ~SolverControl();
146 
150  static void declare_parameters (ParameterHandler &param);
151 
155  void parse_parameters (ParameterHandler &param);
156 
179  virtual State check (const unsigned int step,
180  const double check_value);
181 
185  State last_check() const;
186 
190  double initial_value() const;
191 
196  double last_value() const;
197 
201  unsigned int last_step() const;
202 
206  unsigned int max_steps () const;
207 
211  unsigned int set_max_steps (const unsigned int);
212 
218  void set_failure_criterion (const double rel_failure_residual);
219 
224  void clear_failure_criterion ();
225 
229  double tolerance () const;
230 
234  double set_tolerance (const double);
235 
240  void enable_history_data();
241 
247  double average_reduction() const;
254  double final_reduction() const;
255 
261  double step_reduction(unsigned int step) const;
262 
266  void log_history (const bool);
267 
271  bool log_history () const;
272 
276  unsigned int log_frequency (unsigned int);
277 
281  void log_result (const bool);
282 
286  bool log_result () const;
287 
293  DeclException0(ExcHistoryDataRequired);
294 
295 protected:
299  unsigned int maxsteps;
300 
304  double tol;
305 
310 
314  double initial_val;
315 
319  double lvalue;
320 
324  unsigned int lstep;
325 
331 
336 
344 
349 
353  unsigned int m_log_frequency;
354 
361 
367 
374  std::vector<double> history_data;
375 };
376 
377 
393 {
394 public:
400  ReductionControl (const unsigned int maxiter = 100,
401  const double tolerance = 1.e-10,
402  const double reduce = 1.e-2,
403  const bool log_history = false,
404  const bool log_result = true);
405 
410  ReductionControl (const SolverControl &c);
411 
417 
422  virtual ~ReductionControl();
423 
427  static void declare_parameters (ParameterHandler &param);
428 
432  void parse_parameters (ParameterHandler &param);
433 
439  virtual State check (const unsigned int step,
440  const double check_value);
441 
445  double reduction () const;
446 
450  double set_reduction (const double);
451 
452 protected:
456  double reduce;
457 
462  double reduced_tol;
463 };
464 
479 {
480 public:
485  IterationNumberControl (const unsigned int maxiter = 100,
486  const double tolerance = 1e-12,
487  const bool log_history = false,
488  const bool log_result = true);
489 
495 
501 
506  virtual ~IterationNumberControl();
507 
513  virtual State check (const unsigned int step,
514  const double check_value);
515 };
516 
518 //---------------------------------------------------------------------------
519 
520 #ifndef DOXYGEN
521 
522 inline unsigned int
524 {
525  return maxsteps;
526 }
527 
528 
529 
530 inline unsigned int
531 SolverControl::set_max_steps (const unsigned int newval)
532 {
533  unsigned int old = maxsteps;
534  maxsteps = newval;
535  return old;
536 }
537 
538 
539 
540 inline void
541 SolverControl::set_failure_criterion (const double rel_failure_residual)
542 {
543  relative_failure_residual=rel_failure_residual;
544  check_failure=true;
545 }
546 
547 
548 
549 inline void
551 {
554  check_failure=false;
555 }
556 
557 
558 
559 inline double
561 {
562  return tol;
563 }
564 
565 
566 
567 inline double
568 SolverControl::set_tolerance (const double t)
569 {
570  double old = tol;
571  tol = t;
572  return old;
573 }
574 
575 
576 inline void
577 SolverControl::log_history (const bool newval)
578 {
579  m_log_history = newval;
580 }
581 
582 
583 
584 inline bool
586 {
587  return m_log_history;
588 }
589 
590 
591 inline void
592 SolverControl::log_result (const bool newval)
593 {
594  m_log_result = newval;
595 }
596 
597 
598 inline bool
600 {
601  return m_log_result;
602 }
603 
604 
605 inline double
607 {
608  return reduce;
609 }
610 
611 
612 inline double
613 ReductionControl::set_reduction (const double t)
614 {
615  double old = reduce;
616  reduce = t;
617  return old;
618 }
619 
620 #endif // DOXYGEN
621 
622 DEAL_II_NAMESPACE_CLOSE
623 
624 #endif
Stop iteration, goal not reached.
double step_reduction(unsigned int step) const
virtual State check(const unsigned int step, const double check_value)
virtual void print_info(std::ostream &out) const
static void declare_parameters(ParameterHandler &param)
bool history_data_enabled
bool log_history() const
Continue iteration.
void parse_parameters(ParameterHandler &param)
double tolerance() const
void clear_failure_criterion()
double average_reduction() const
double initial_value() const
IterationNumberControl & operator=(const SolverControl &c)
virtual State check(const unsigned int step, const double check_value)
double failure_residual
virtual ~IterationNumberControl()
double last_value() const
const unsigned int last_step
unsigned int log_frequency(unsigned int)
Stop iteration, goal reached.
unsigned int max_steps() const
unsigned int m_log_frequency
double set_reduction(const double)
ReductionControl & operator=(const SolverControl &c)
double relative_failure_residual
unsigned int last_step() const
unsigned int lstep
virtual ~ReductionControl()
static void declare_parameters(ParameterHandler &param)
double reduction() const
IterationNumberControl(const unsigned int maxiter=100, const double tolerance=1e-12, const bool log_history=false, const bool log_result=true)
virtual ~SolverControl()
State last_check() const
double final_reduction() const
void enable_history_data()
double set_tolerance(const double)
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)
unsigned int maxsteps
virtual State check(const unsigned int step, const double check_value)
void parse_parameters(ParameterHandler &param)
std::vector< double > history_data
bool log_result() const
SolverControl(const unsigned int n=100, const double tol=1.e-10, const bool log_history=false, const bool log_result=true)
unsigned int set_max_steps(const unsigned int)
DeclException0(ExcHistoryDataRequired)
void set_failure_criterion(const double rel_failure_residual)