42#ifndef BELOS_GMRES_ITERATION_HPP
43#define BELOS_GMRES_ITERATION_HPP
62 template <
class ScalarType,
class MV>
71 Teuchos::RCP<const MV>
V;
74 Teuchos::RCP<const MV>
Z;
82 Teuchos::RCP<const Teuchos::SerialDenseMatrix<int,ScalarType> >
H;
85 Teuchos::RCP<const Teuchos::SerialDenseMatrix<int,ScalarType> >
R;
88 Teuchos::RCP<const Teuchos::SerialDenseMatrix<int,ScalarType> >
z;
92 H(Teuchos::null),
R(Teuchos::null),
140template<
class ScalarType,
class MV,
class OP>
198 virtual void setSize(
int blockSize,
int numBlocks) = 0;
Belos header file which uses auto-configuration information to include necessary C++ headers.
Pure virtual base class which describes the basic interface to the linear solver iteration.
Collection of types and exceptions used within the Belos solvers.
Parent class to all Belos exceptions.
GmresIterationInitFailure is thrown when the GmresIteration object is unable to generate an initial i...
GmresIterationInitFailure(const std::string &what_arg)
GmresIterationLAPACKFailure is thrown when a nonzero return value is passed back from an LAPACK routi...
GmresIterationLAPACKFailure(const std::string &what_arg)
GmresIterationOrthoFailure is thrown when the GmresIteration object is unable to compute independent ...
GmresIterationOrthoFailure(const std::string &what_arg)
virtual GmresIterationState< ScalarType, MV > getState() const =0
Get the current state of the linear solver.
virtual void initializeGmres(GmresIterationState< ScalarType, MV > &newstate)=0
Initialize the solver to an iterate, providing a complete state.
virtual int getCurSubspaceDim() const =0
Get the dimension of the search subspace used to generate the current solution to the linear problem.
virtual void updateLSQR(int dim=-1)=0
Method for updating QR factorization of upper Hessenberg matrix.
virtual int getMaxSubspaceDim() const =0
Get the maximum dimension allocated for the search subspace.
virtual void setSize(int blockSize, int numBlocks)=0
Set the blocksize and number of blocks to be used by the iterative solver in solving this linear prob...
Structure to contain pointers to GmresIteration state variables.
Teuchos::RCP< const MV > V
The current Krylov basis.
Teuchos::RCP< const Teuchos::SerialDenseMatrix< int, ScalarType > > z
The current right-hand side of the least squares system RY = Z.
Teuchos::RCP< const MV > Z
The current preconditioned Krylov basis (only used in flexible GMRES).
int curDim
The current dimension of the reduction.
Teuchos::RCP< const Teuchos::SerialDenseMatrix< int, ScalarType > > H
The current Hessenberg matrix.
Teuchos::RCP< const Teuchos::SerialDenseMatrix< int, ScalarType > > R
The current upper-triangular matrix from the QR reduction of H.