Belos Version of the Day
Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP > Class Template Reference

The Belos::PseudoBlockGmresSolMgr provides a powerful and fully-featured solver manager over the pseudo-block GMRES iteration. More...

#include <BelosPseudoBlockGmresSolMgr.hpp>

Inheritance diagram for Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >:
Inheritance graph
[legend]

List of all members.

Public Member Functions

Constructors/Destructor
 PseudoBlockGmresSolMgr ()
 Empty constructor for BlockGmresSolMgr. This constructor takes no arguments and sets the default values for the solver. The linear problem must be passed in using setProblem() before solve() is called on this object. The solver values can be changed using setParameters().
 PseudoBlockGmresSolMgr (const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem, const Teuchos::RCP< Teuchos::ParameterList > &pl)
 Basic constructor for PseudoBlockGmresSolMgr.
virtual ~PseudoBlockGmresSolMgr ()
 Destructor.
Accessor methods
const LinearProblem
< ScalarType, MV, OP > & 
getProblem () const
 Return a reference to the linear problem being solved by this solver manager.
Teuchos::RCP< const
Teuchos::ParameterList
getValidParameters () const
 Get a parameter list containing the valid parameters for this object.
Teuchos::RCP< const
Teuchos::ParameterList
getCurrentParameters () const
 Get a parameter list containing the current parameters for this object.
Teuchos::Array< Teuchos::RCP
< Teuchos::Time > > 
getTimers () const
 Return the timers for this object.
int getNumIters () const
 Get the iteration count for the most recent call to solve().
bool isLOADetected () const
 Whether a "loss of accuracy" was detected during the last solve().
Set methods
void setProblem (const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem)
 Set the linear problem to solve.
void setParameters (const Teuchos::RCP< Teuchos::ParameterList > &params)
 Set the parameters the solver manager should use to solve the linear problem.
virtual void setUserConvStatusTest (const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &userConvStatusTest)
 Set a custom status test.
Reset methods
void reset (const ResetType type)
 Performs a reset of the solver manager specified by the ResetType. This informs the solver manager that the solver should prepare for the next call to solve by resetting certain elements of the iterative solver strategy.
Solver application methods
ReturnType solve ()
 This method performs possibly repeated calls to the underlying linear solver's iterate() routine until the problem has been solved (as decided by the solver manager) or the solver manager decides to quit.
Overridden from Teuchos::Describable
std::string description () const
 Return a description of the pseudoblock GMRES solver manager.

Detailed Description

template<class ScalarType, class MV, class OP>
class Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >

The Belos::PseudoBlockGmresSolMgr provides a powerful and fully-featured solver manager over the pseudo-block GMRES iteration.

Author:
Heidi Thornquist, Chris Baker, and Teri Barth
Examples:

BlockGmres/PseudoBlockGmresEpetraExFile.cpp, and BlockGmres/PseudoBlockPrecGmresEpetraExFile.cpp.

Definition at line 112 of file BelosPseudoBlockGmresSolMgr.hpp.


Constructor & Destructor Documentation

template<class ScalarType , class MV , class OP >
Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::PseudoBlockGmresSolMgr ( )

Empty constructor for BlockGmresSolMgr. This constructor takes no arguments and sets the default values for the solver. The linear problem must be passed in using setProblem() before solve() is called on this object. The solver values can be changed using setParameters().

Definition at line 447 of file BelosPseudoBlockGmresSolMgr.hpp.

template<class ScalarType , class MV , class OP >
Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::PseudoBlockGmresSolMgr ( const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &  problem,
const Teuchos::RCP< Teuchos::ParameterList > &  pl 
)

Basic constructor for PseudoBlockGmresSolMgr.

This constructor accepts the LinearProblem to be solved in addition to a parameter list of options for the solver manager. These options include the following:

  • "Num Blocks" - a int specifying the number of blocks allocated for the Krylov basis.
  • "Maximum Iterations" - a int specifying the maximum number of iterations the underlying solver is allowed to perform.
  • "Maximum Restarts" - a int specifying the maximum number of restarts the underlying solver is allowed to perform.
  • "Orthogonalization" - a std::string specifying the desired orthogonalization: DGKS, ICGS, and IMGS. Default: "DGKS"
  • "Verbosity" - a sum of MsgType specifying the verbosity. Default: Belos::Errors
  • "Output Style" - a OutputType specifying the style of output. Default: Belos::General
  • "Convergence Tolerance" - a MagnitudeType specifying the level that residual norms must reach to decide convergence.

Definition at line 472 of file BelosPseudoBlockGmresSolMgr.hpp.

template<class ScalarType , class MV , class OP >
virtual Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::~PseudoBlockGmresSolMgr ( ) [inline, virtual]

Destructor.

Definition at line 149 of file BelosPseudoBlockGmresSolMgr.hpp.


Member Function Documentation

template<class ScalarType , class MV , class OP >
const LinearProblem<ScalarType,MV,OP>& Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::getProblem ( ) const [inline, virtual]

Return a reference to the linear problem being solved by this solver manager.

Implements Belos::SolverManager< ScalarType, MV, OP >.

Definition at line 155 of file BelosPseudoBlockGmresSolMgr.hpp.

template<class ScalarType , class MV , class OP >
Teuchos::RCP< const Teuchos::ParameterList > Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::getValidParameters ( ) const [virtual]

Get a parameter list containing the valid parameters for this object.

Implements Belos::SolverManager< ScalarType, MV, OP >.

Definition at line 803 of file BelosPseudoBlockGmresSolMgr.hpp.

template<class ScalarType , class MV , class OP >
Teuchos::RCP<const Teuchos::ParameterList> Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::getCurrentParameters ( ) const [inline, virtual]

Get a parameter list containing the current parameters for this object.

Implements Belos::SolverManager< ScalarType, MV, OP >.

Definition at line 165 of file BelosPseudoBlockGmresSolMgr.hpp.

template<class ScalarType , class MV , class OP >
Teuchos::Array<Teuchos::RCP<Teuchos::Time> > Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::getTimers ( ) const [inline]

Return the timers for this object.

The timers are ordered as follows:

Definition at line 172 of file BelosPseudoBlockGmresSolMgr.hpp.

template<class ScalarType , class MV , class OP >
int Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::getNumIters ( ) const [inline, virtual]

Get the iteration count for the most recent call to solve().

Implements Belos::SolverManager< ScalarType, MV, OP >.

Definition at line 177 of file BelosPseudoBlockGmresSolMgr.hpp.

template<class ScalarType , class MV , class OP >
bool Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::isLOADetected ( ) const [inline, virtual]

Whether a "loss of accuracy" was detected during the last solve().

The GMRES algorithm (which Pseudoblock GMRES implements) uses two different residual norms to predict convergence: "implicit" (also called "native") and "explicit" (also called "exact," not to be confused with "exact arithmetic"). The "implicit" residuals are computed by the solver via a recurrence relation (the Arnoldi relation, in the case of GMRES). The "explicit" residuals are computed directly as $B

  • A X_k$. Implicit residuals are much cheaper to compute, since they are available almost "for free" from the recurrence relation. In contrast, computing exact residuals requires computing the current approximate solution $X_k$, applying the global operator $A$ to $X_k$, and then computing the norm of the resulting vector(s) via a global reduction. Thus, GMRES favors using the cheaper implicit residuals to predict convergence. Users typically want convergence with respect to explicit residuals, though.

Implicit and explicit residuals may differ due to rounding error. However, the difference between implicit and explicit residuals matters most when using a left (or split) preconditioner. In that case, the implicit residuals are those of the left-preconditioned problem $M_L^{-1} A X = M_L^{-1} B$ instead of the original problem $A X = B$. The implicit residual norms may thus differ significantly from the explicit residual norms, even if one could compute without rounding error.

When using a left preconditioner, Pseudoblock GMRES tries to detect if the implicit residuals have converged but the explicit residuals have not. In that case, it will reduce the convergence tolerance and iterate a little while longer to attempt to reduce the explicit residual norm. However, if that doesn't work, it declares a "loss of accuracy" for the affected right-hand side(s), and stops iterating on them. (Not all right-hand sides may have experienced a loss of accuracy.) Thus, the affected right-hand sides may or may not have converged to the desired residual norm tolerance. Calling this method tells you whether a "loss of accuracy" (LOA) occurred during the last solve() invocation.

When not using a left preconditioner, Pseudoblock GMRES will iterate until both the implicit and explicit residuals converge. (It does not start testing the explicit residuals until the implicit residuals have converged. This avoids whenever possible the cost of computing explicit residuals.) Implicit and explicit residuals may differ due to rounding error, even though they are identical when no rounding error occurs. In this case, the algorithm does not report a "loss of accuracy," since it continues iterating until the explicit residuals converge.

Note:
Calling solve() again resets the flag that reports whether a loss of accuracy was detected. Thus, you should call this method immediately after calling solve().

Implements Belos::SolverManager< ScalarType, MV, OP >.

Definition at line 237 of file BelosPseudoBlockGmresSolMgr.hpp.

template<class ScalarType , class MV , class OP >
void Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::setProblem ( const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &  problem) [inline, virtual]

Set the linear problem to solve.

Implements Belos::SolverManager< ScalarType, MV, OP >.

Definition at line 245 of file BelosPseudoBlockGmresSolMgr.hpp.

template<class ScalarType , class MV , class OP >
void Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::setParameters ( const Teuchos::RCP< Teuchos::ParameterList > &  params) [virtual]

Set the parameters the solver manager should use to solve the linear problem.

Implements Belos::SolverManager< ScalarType, MV, OP >.

Definition at line 507 of file BelosPseudoBlockGmresSolMgr.hpp.

template<class ScalarType , class MV , class OP >
void Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::setUserConvStatusTest ( const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &  userConvStatusTest) [virtual]

Set a custom status test.

A custom status test is not required. If you decide to set one, the current implementation will apply it sequentially (short-circuiting OR, like the || operator in C++) after Pseudoblock GMRES' standard convergence test.

Reimplemented from Belos::SolverManager< ScalarType, MV, OP >.

Definition at line 793 of file BelosPseudoBlockGmresSolMgr.hpp.

template<class ScalarType , class MV , class OP >
void Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::reset ( const ResetType  type) [inline, virtual]

Performs a reset of the solver manager specified by the ResetType. This informs the solver manager that the solver should prepare for the next call to solve by resetting certain elements of the iterative solver strategy.

Implements Belos::SolverManager< ScalarType, MV, OP >.

Definition at line 270 of file BelosPseudoBlockGmresSolMgr.hpp.

template<class ScalarType , class MV , class OP >
ReturnType Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::solve ( ) [virtual]

This method performs possibly repeated calls to the underlying linear solver's iterate() routine until the problem has been solved (as decided by the solver manager) or the solver manager decides to quit.

This method calls PseudoBlockGmresIter::iterate(), which will return either because a specially constructed status test evaluates to Passed or an std::exception is thrown.

A return from PseudoBlockGmresIter::iterate() signifies one of the following scenarios:

  • the maximum number of restarts has been exceeded. In this scenario, the current solutions to the linear system will be placed in the linear problem and return Unconverged.
  • global convergence has been met. In this case, the current solutions to the linear system will be placed in the linear problem and the solver manager will return Converged
Returns:
ReturnType specifying:
  • Converged: the linear problem was solved to the specification required by the solver manager.
  • Unconverged: the linear problem was not solved to the specification desired by the solver manager.

Implements Belos::SolverManager< ScalarType, MV, OP >.

Definition at line 940 of file BelosPseudoBlockGmresSolMgr.hpp.

template<class ScalarType , class MV , class OP >
std::string Belos::PseudoBlockGmresSolMgr< ScalarType, MV, OP >::description ( ) const [virtual]

Return a description of the pseudoblock GMRES solver manager.

Reimplemented from Teuchos::Describable.

Definition at line 1308 of file BelosPseudoBlockGmresSolMgr.hpp.


The documentation for this class was generated from the following file:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines