Anasazi  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
List of all members
Anasazi::LOBPCG< ScalarType, MV, OP > Class Template Reference

This class provides the Locally Optimal Block Preconditioned Conjugate Gradient (LOBPCG) iteration, a preconditioned iteration for solving linear Hermitian eigenproblems. More...

#include <AnasaziLOBPCG.hpp>

Inheritance diagram for Anasazi::LOBPCG< ScalarType, MV, OP >:
Anasazi::Eigensolver< ScalarType, MV, OP >

Public Member Functions

Constructor/Destructor
 LOBPCG (const Teuchos::RCP< Eigenproblem< ScalarType, MV, OP > > &problem, const Teuchos::RCP< SortManager< typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > > &sorter, const Teuchos::RCP< OutputManager< ScalarType > > &printer, const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &tester, const Teuchos::RCP< MatOrthoManager< ScalarType, MV, OP > > &ortho, Teuchos::ParameterList &params)
 LOBPCG constructor with eigenproblem, solver utilities, and parameter list of solver options. More...
 
virtual ~LOBPCG ()
 LOBPCG destructor More...
 
Solver methods
void iterate ()
 This method performs LOBPCG iterations until the status test indicates the need to stop or an error occurs (in which case, an exception is thrown). More...
 
void initialize (LOBPCGState< ScalarType, MV > &newstate)
 Initialize the solver to an iterate, optionally providing the Ritz values, residual, and search direction. More...
 
void initialize ()
 Initialize the solver with the initial vectors from the eigenproblem or random data. More...
 
bool isInitialized () const
 Indicates whether the solver has been initialized or not. More...
 
LOBPCGState< ScalarType, MV > getState () const
 Get the current state of the eigensolver. More...
 
Status methods
int getNumIters () const
 Get the current iteration count. More...
 
void resetNumIters ()
 Reset the iteration count. More...
 
Teuchos::RCP< const MV > getRitzVectors ()
 Get the Ritz vectors from the previous iteration. More...
 
std::vector< Value< ScalarType > > getRitzValues ()
 Get the Ritz values from the previous iteration. More...
 
std::vector< int > getRitzIndex ()
 Get the index used for extracting Ritz vectors from getRitzVectors(). More...
 
std::vector< typename
Teuchos::ScalarTraits
< ScalarType >::magnitudeType > 
getResNorms ()
 Get the current residual norms. More...
 
std::vector< typename
Teuchos::ScalarTraits
< ScalarType >::magnitudeType > 
getRes2Norms ()
 Get the current residual 2-norms. More...
 
std::vector< typename
Teuchos::ScalarTraits
< ScalarType >::magnitudeType > 
getRitzRes2Norms ()
 Get the 2-norms of the residuals. More...
 
int getCurSubspaceDim () const
 Get the dimension of the search subspace used to generate the current eigenvectors and eigenvalues. More...
 
int getMaxSubspaceDim () const
 Get the maximum dimension allocated for the search subspace. For LOBPCG, this always returns 3*getBlockSize(), the dimension of the subspace colspan([X H P]). More...
 
Accessor routines from Eigensolver
void setStatusTest (Teuchos::RCP< StatusTest< ScalarType, MV, OP > > test)
 Set a new StatusTest for the solver. More...
 
Teuchos::RCP< StatusTest
< ScalarType, MV, OP > > 
getStatusTest () const
 Get the current StatusTest used by the solver. More...
 
const Eigenproblem< ScalarType,
MV, OP > & 
getProblem () const
 Get a constant reference to the eigenvalue problem. More...
 
void setBlockSize (int blockSize)
 Set the blocksize to be used by the iterative solver in solving this eigenproblem. More...
 
int getBlockSize () const
 Get the blocksize to be used by the iterative solver in solving this eigenproblem. More...
 
void setAuxVecs (const Teuchos::Array< Teuchos::RCP< const MV > > &auxvecs)
 Set the auxiliary vectors for the solver. More...
 
Teuchos::Array< Teuchos::RCP
< const MV > > 
getAuxVecs () const
 Get the current auxiliary vectors. More...
 
LOBPCG-specific accessor routines
void setFullOrtho (bool fullOrtho)
 Instruct the LOBPCG iteration to use full orthogonality. More...
 
bool getFullOrtho () const
 Determine if the LOBPCG iteration is using full orthogonality. More...
 
bool hasP ()
 Indicates whether the search direction given by getState() is valid. More...
 
Output methods
void currentStatus (std::ostream &os)
 This method requests that the solver print out its current status to screen. More...
 
- Public Member Functions inherited from Anasazi::Eigensolver< ScalarType, MV, OP >
 Eigensolver ()
 Default Constructor. More...
 
 Eigensolver (const Teuchos::RCP< Eigenproblem< ScalarType, MV, OP > > &problem, const Teuchos::RCP< SortManager< ScalarType > > &sorter, const Teuchos::RCP< OutputManager< ScalarType > > &printer, const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &tester, const Teuchos::RCP< OrthoManager< ScalarType, MV > > &ortho, Teuchos::ParameterList &params)
 Basic Constructor. More...
 
virtual ~Eigensolver ()
 Destructor. More...
 

Detailed Description

template<class ScalarType, class MV, class OP>
class Anasazi::LOBPCG< ScalarType, MV, OP >

This class provides the Locally Optimal Block Preconditioned Conjugate Gradient (LOBPCG) iteration, a preconditioned iteration for solving linear Hermitian eigenproblems.

This implementation is a modification of the one found in A. Knyazev, "Toward the optimal preconditioned eigensolver: Locally optimal block preconditioner conjugate gradient method", SIAM J. Sci. Comput., vol 23, n 2, pp. 517-541.

The modification consists of the orthogonalization steps recommended in U. Hetmaniuk and R. Lehoucq, "Basis Selection in LOBPCG", Journal of Computational Physics.

These modifcation are referred to as full orthogonalization, and consist of also conducting the local optimization using an orthonormal basis.

Author
Chris Baker, Ulrich Hetmaniuk, Rich Lehoucq, Heidi Thornquist

Definition at line 222 of file AnasaziLOBPCG.hpp.

Constructor & Destructor Documentation

template<class ScalarType , class MV , class OP >
Anasazi::LOBPCG< ScalarType, MV, OP >::LOBPCG ( const Teuchos::RCP< Eigenproblem< ScalarType, MV, OP > > &  problem,
const Teuchos::RCP< SortManager< typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > > &  sorter,
const Teuchos::RCP< OutputManager< ScalarType > > &  printer,
const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &  tester,
const Teuchos::RCP< MatOrthoManager< ScalarType, MV, OP > > &  ortho,
Teuchos::ParameterList params 
)

LOBPCG constructor with eigenproblem, solver utilities, and parameter list of solver options.

This constructor takes pointers required by the eigensolver, in addition to a parameter list of options for the eigensolver. These options include the following:

  • "Block Size" - an int specifying the block size used by the algorithm. This can also be specified using the setBlockSize() method.
  • "Full Ortho" - a bool specifying whether the solver should employ a full orthogonalization technique. This can also be specified using the setFullOrtho() method.

Definition at line 624 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
virtual Anasazi::LOBPCG< ScalarType, MV, OP >::~LOBPCG ( )
inlinevirtual

LOBPCG destructor

Definition at line 244 of file AnasaziLOBPCG.hpp.

Member Function Documentation

template<class ScalarType , class MV , class OP >
void Anasazi::LOBPCG< ScalarType, MV, OP >::iterate ( )
virtual

This method performs LOBPCG iterations until the status test indicates the need to stop or an error occurs (in which case, an exception is thrown).

iterate() will first determine whether the solver is initialized; if not, it will call initialize() using default arguments. After initialization, the solver performs LOBPCG iterations until the status test evaluates as Passed, at which point the method returns to the caller.

The LOBPCG iteration proceeds as follows:

  1. The current residual (R) is preconditioned to form H
  2. H is orthogonalized against the auxiliary vectors and, if full orthogonalization
    is enabled, against X and P.
  3. The basis [X H P] is used to project the problem matrices.
  4. The projected eigenproblem is solved, and the desired eigenvectors and eigenvalues are selected.
  5. These are used to form the new eigenvector estimates (X) and the search directions (P).
    If full orthogonalization is enabled, these are generated to be mutually orthogonal and with orthonormal columns.
  6. The new residual (R) is formed.

The status test is queried at the beginning of the iteration.

Possible exceptions thrown include std::logic_error, std::invalid_argument or one of the LOBPCG-specific exceptions.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 1350 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
void Anasazi::LOBPCG< ScalarType, MV, OP >::initialize ( LOBPCGState< ScalarType, MV > &  newstate)

Initialize the solver to an iterate, optionally providing the Ritz values, residual, and search direction.

Note
LOBPCGState contains fields V, KV and MV: These are ignored by initialize()

The LOBPCG eigensolver contains a certain amount of state relating to the current iterate, including the current residual, the current search direction, and the images of these spaces under the eigenproblem operators.

initialize() gives the user the opportunity to manually set these, although this must be done with caution, abiding by the rules given below. All notions of orthogonality and orthonormality are derived from the inner product specified by the orthogonalization manager.

Postcondition
  • isInitialized() == true (see post-conditions of isInitialize())
  • If newstate.P != Teuchos::null, hasP() == true.
    Otherwise, hasP() == false

The user has the option of specifying any component of the state using initialize(). However, these arguments are assumed to match the post-conditions specified under isInitialized(). Any component of the state (i.e., KX) not given to initialize() will be generated.

Definition at line 974 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
void Anasazi::LOBPCG< ScalarType, MV, OP >::initialize ( )
virtual

Initialize the solver with the initial vectors from the eigenproblem or random data.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 1310 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
bool Anasazi::LOBPCG< ScalarType, MV, OP >::isInitialized ( ) const
virtual

Indicates whether the solver has been initialized or not.

Returns
bool indicating the state of the solver.
Postcondition
If isInitialized() == true:
  • X is orthogonal to auxiliary vectors and has orthonormal columns
  • KX == Op*X
  • MX == M*X if M != Teuchos::null
    Otherwise, MX == Teuchos::null
  • getRitzValues() returns the sorted Ritz values with respect to X
  • getResNorms(), getRes2Norms(), getRitzResNorms() are correct
  • If hasP() == true,
    • P orthogonal to auxiliary vectors
    • If getFullOrtho() == true,
      • P is orthogonal to X and has orthonormal columns
    • KP == Op*P
    • MP == M*P if M != Teuchos::null
      Otherwise, MP == Teuchos::null

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2181 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
LOBPCGState< ScalarType, MV > Anasazi::LOBPCG< ScalarType, MV, OP >::getState ( ) const

Get the current state of the eigensolver.

The data is only valid if isInitialized() == true. The data for the search directions P is only meaningful if hasP() == true. Finally, the data for the preconditioned residual (H) is only meaningful in the situation where the solver throws an ::LOBPCGRitzFailure exception during iterate().

Returns
An LOBPCGState object containing const views to the current solver state.

Definition at line 2313 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
int Anasazi::LOBPCG< ScalarType, MV, OP >::getNumIters ( ) const
virtual

Get the current iteration count.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2305 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
void Anasazi::LOBPCG< ScalarType, MV, OP >::resetNumIters ( )
virtual

Reset the iteration count.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2297 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
Teuchos::RCP< const MV > Anasazi::LOBPCG< ScalarType, MV, OP >::getRitzVectors ( )
virtual

Get the Ritz vectors from the previous iteration.

Returns
A multivector with getBlockSize() vectors containing the sorted Ritz vectors corresponding to the most significant Ritz values. The i-th vector of the return corresponds to the i-th Ritz vector; there is no need to use getRitzIndex().

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2289 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
std::vector< Value< ScalarType > > Anasazi::LOBPCG< ScalarType, MV, OP >::getRitzValues ( )
virtual

Get the Ritz values from the previous iteration.

Returns
A vector of length getCurSubspaceDim() containing the Ritz values from the previous projected eigensolve.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2261 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
std::vector< int > Anasazi::LOBPCG< ScalarType, MV, OP >::getRitzIndex ( )
virtual

Get the index used for extracting Ritz vectors from getRitzVectors().

Because BlockDavidson is a Hermitian solver, all Ritz values are real and all Ritz vectors can be represented in a single column of a multivector. Therefore, getRitzIndex() is not needed when using the output from getRitzVectors().

Returns
An int vector of size getCurSubspaceDim() composed of zeros.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2252 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
std::vector< typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > Anasazi::LOBPCG< ScalarType, MV, OP >::getResNorms ( )
virtual

Get the current residual norms.

Returns
A vector of length getBlockSize() containing the norms of the residuals, with respect to the orthogonalization manager norm() method.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 1977 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
std::vector< typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > Anasazi::LOBPCG< ScalarType, MV, OP >::getRes2Norms ( )
virtual

Get the current residual 2-norms.

Returns
A vector of length getBlockSize() containing the 2-norms of the residuals.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 1991 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
std::vector< typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > Anasazi::LOBPCG< ScalarType, MV, OP >::getRitzRes2Norms ( )
virtual

Get the 2-norms of the residuals.

The Ritz residuals are not defined for the LOBPCG iteration. Hence, this method returns the 2-norms of the direct residuals, and is equivalent to calling getRes2Norms().

Returns
A vector of length getBlockSize() containing the 2-norms of the direct residuals.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2243 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
int Anasazi::LOBPCG< ScalarType, MV, OP >::getCurSubspaceDim ( ) const
virtual

Get the dimension of the search subspace used to generate the current eigenvectors and eigenvalues.

LOBPCG employs a sequential subspace iteration, maintaining a fixed-rank basis, as opposed to an expanding subspace mechanism employed by Krylov-subspace solvers like BlockKrylovSchur and BlockDavidson.

Returns
An integer specifying the rank of the subspace generated by the eigensolver. If isInitialized() == false, the return is 0. Otherwise, the return will be 2*getBlockSize() or 3*getBlockSize().

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2233 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
int Anasazi::LOBPCG< ScalarType, MV, OP >::getMaxSubspaceDim ( ) const
virtual

Get the maximum dimension allocated for the search subspace. For LOBPCG, this always returns 3*getBlockSize(), the dimension of the subspace colspan([X H P]).

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2226 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
void Anasazi::LOBPCG< ScalarType, MV, OP >::setStatusTest ( Teuchos::RCP< StatusTest< ScalarType, MV, OP > >  test)
virtual

Set a new StatusTest for the solver.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2273 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
Teuchos::RCP< StatusTest< ScalarType, MV, OP > > Anasazi::LOBPCG< ScalarType, MV, OP >::getStatusTest ( ) const
virtual

Get the current StatusTest used by the solver.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2282 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
const Eigenproblem< ScalarType, MV, OP > & Anasazi::LOBPCG< ScalarType, MV, OP >::getProblem ( ) const
virtual

Get a constant reference to the eigenvalue problem.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2218 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
void Anasazi::LOBPCG< ScalarType, MV, OP >::setBlockSize ( int  blockSize)
virtual

Set the blocksize to be used by the iterative solver in solving this eigenproblem.

If the block size is reduced, then the new iterate (and residual and search direction) are chosen as the subset of the current iterate preferred by the sort manager. Otherwise, the solver state is set to uninitialized.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 701 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
int Anasazi::LOBPCG< ScalarType, MV, OP >::getBlockSize ( ) const
virtual

Get the blocksize to be used by the iterative solver in solving this eigenproblem.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2211 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
void Anasazi::LOBPCG< ScalarType, MV, OP >::setAuxVecs ( const Teuchos::Array< Teuchos::RCP< const MV > > &  auxvecs)
virtual

Set the auxiliary vectors for the solver.

Because the current iterate X and search direction P cannot be assumed orthogonal to the new auxiliary vectors, a call to setAuxVecs() with a non-empty argument will reset the solver to the uninitialized state.

In order to preserve the current state, the user will need to extract it from the solver using getState(), orthogonalize it against the new auxiliary vectors, and manually reinitialize the solver using initialize().

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 929 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
Teuchos::Array< Teuchos::RCP< const MV > > Anasazi::LOBPCG< ScalarType, MV, OP >::getAuxVecs ( ) const
virtual

Get the current auxiliary vectors.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2204 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
void Anasazi::LOBPCG< ScalarType, MV, OP >::setFullOrtho ( bool  fullOrtho)

Instruct the LOBPCG iteration to use full orthogonality.

If the getFullOrtho() == false and isInitialized() == true and hasP() == true, then P will be invalidated by setting full orthogonalization to true.

Definition at line 1320 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
bool Anasazi::LOBPCG< ScalarType, MV, OP >::getFullOrtho ( ) const

Determine if the LOBPCG iteration is using full orthogonality.

Definition at line 2196 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
bool Anasazi::LOBPCG< ScalarType, MV, OP >::hasP ( )

Indicates whether the search direction given by getState() is valid.

Definition at line 2189 of file AnasaziLOBPCG.hpp.

template<class ScalarType , class MV , class OP >
void Anasazi::LOBPCG< ScalarType, MV, OP >::currentStatus ( std::ostream &  os)
virtual

This method requests that the solver print out its current status to screen.

Implements Anasazi::Eigensolver< ScalarType, MV, OP >.

Definition at line 2136 of file AnasaziLOBPCG.hpp.


The documentation for this class was generated from the following file: