Belos Version of the Day
Loading...
Searching...
No Matches
BelosGmresIteration.hpp
Go to the documentation of this file.
1//@HEADER
2// ************************************************************************
3//
4// Belos: Block Linear Solvers Package
5// Copyright 2004 Sandia Corporation
6//
7// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8// the U.S. Government retains certain rights in this software.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38//
39// ************************************************************************
40//@HEADER
41
42#ifndef BELOS_GMRES_ITERATION_HPP
43#define BELOS_GMRES_ITERATION_HPP
44
49#include "BelosConfigDefs.hpp"
50#include "BelosTypes.hpp"
51#include "BelosIteration.hpp"
52
53namespace Belos {
54
56
57
62 template <class ScalarType, class MV>
68 int curDim;
69
71 Teuchos::RCP<const MV> V;
72
74 Teuchos::RCP<const MV> Z;
75
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;
89
90
91 GmresIterationState() : curDim(0), V(Teuchos::null), Z(Teuchos::null),
92 H(Teuchos::null), R(Teuchos::null),
93 z(Teuchos::null)
94 {}
95 };
96
98
100
101
113 class GmresIterationInitFailure : public BelosError {public:
114 GmresIterationInitFailure(const std::string& what_arg) : BelosError(what_arg)
115 {}};
116
124 GmresIterationOrthoFailure(const std::string& what_arg) : BelosError(what_arg)
125 {}};
126
134 GmresIterationLAPACKFailure(const std::string& what_arg) : BelosError(what_arg)
135 {}};
136
138
139
140template<class ScalarType, class MV, class OP>
141class GmresIteration : virtual public Iteration<ScalarType,MV,OP> {
142
143 public:
144
146
147
162
171
173
174
176
179 virtual void updateLSQR( int dim = -1 ) = 0;
180
182 virtual int getCurSubspaceDim() const = 0;
183
185 virtual int getMaxSubspaceDim() const = 0;
186
188
190
191
198 virtual void setSize(int blockSize, int numBlocks) = 0;
199
201
202};
203
204} // end Belos namespace
205
206#endif /* BELOS_GMRES_ITERATION_HPP */
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.

Generated for Belos by doxygen 1.10.0