Belos::IMGSOrthoManager< ScalarType, MV, OP > Class Template Reference

An implementation of the Belos::MatOrthoManager that performs orthogonalization using multiple steps of modified Gram-Schmidt. More...

#include <BelosIMGSOrthoManager.hpp>

Inheritance diagram for Belos::IMGSOrthoManager< ScalarType, MV, OP >:

Inheritance graph
[legend]
List of all members.

Public Member Functions

Constructor/Destructor
 IMGSOrthoManager (const std::string &label="Belos", Teuchos::RCP< const OP > Op=Teuchos::null, const int max_ortho_steps=1, const MagnitudeType blk_tol=10 *MGT::squareroot(MGT::eps()), const MagnitudeType sing_tol=10 *MGT::eps())
 Constructor specifying re-orthogonalization tolerance.
 ~IMGSOrthoManager ()
 Destructor.
Accessor routines
void setBlkTol (const MagnitudeType blk_tol)
 Set parameter for block re-orthogonalization threshhold.
void setSingTol (const MagnitudeType sing_tol)
 Set parameter for singular block detection.
MagnitudeType getBlkTol () const
 Return parameter for block re-orthogonalization threshhold.
MagnitudeType getSingTol () const
 Return parameter for singular block detection.
Orthogonalization methods
void project (MV &X, Teuchos::RCP< MV > MX, Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > > C, Teuchos::Array< Teuchos::RCP< const MV > > Q) const
 Given a list of (mutually and internally) orthonormal bases Q, this method takes a multivector X and projects it onto the space orthogonal to the individual Q[i], optionally returning the coefficients of X for the individual Q[i]. All of this is done with respect to the inner product innerProd().
void project (MV &X, Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > > C, Teuchos::Array< Teuchos::RCP< const MV > > Q) const
 This method calls project(X,Teuchos::null,C,Q); see documentation for that function.
int normalize (MV &X, Teuchos::RCP< MV > MX, Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > B) const
 This method takes a multivector X and attempts to compute an orthonormal basis for $colspan(X)$, with respect to innerProd().
int normalize (MV &X, Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > B) const
 This method calls normalize(X,Teuchos::null,B); see documentation for that function.
int projectAndNormalize (MV &X, Teuchos::RCP< MV > MX, Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > > C, Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > B, Teuchos::Array< Teuchos::RCP< const MV > > Q) const
 Given a set of bases Q[i] and a multivector X, this method computes an orthonormal basis for $colspan(X) - \sum_i colspan(Q[i])$.
int projectAndNormalize (MV &X, Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > > C, Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > B, Teuchos::Array< Teuchos::RCP< const MV > > Q) const
 This method calls projectAndNormalize(X,Teuchos::null,C,B,Q); see documentation for that function.
Error methods
Teuchos::ScalarTraits< ScalarType
>::magnitudeType 
orthonormError (const MV &X) const
 This method computes the error in orthonormality of a multivector, measured as the Frobenius norm of the difference innerProd(X,Y) - I.
Teuchos::ScalarTraits< ScalarType
>::magnitudeType 
orthonormError (const MV &X, Teuchos::RCP< const MV > MX) const
 This method computes the error in orthonormality of a multivector, measured as the Frobenius norm of the difference innerProd(X,Y) - I. The method has the option of exploiting a caller-provided MX.
Teuchos::ScalarTraits< ScalarType
>::magnitudeType 
orthogError (const MV &X1, const MV &X2) const
 This method computes the error in orthogonality of two multivectors, measured as the Frobenius norm of innerProd(X,Y).
Teuchos::ScalarTraits< ScalarType
>::magnitudeType 
orthogError (const MV &X1, Teuchos::RCP< const MV > MX1, const MV &X2) const
 This method computes the error in orthogonality of two multivectors, measured as the Frobenius norm of innerProd(X,Y). The method has the option of exploiting a caller-provided MX.
Label methods
void setLabel (const std::string &label)
 This method sets the label used by the timers in the orthogonalization manager.
const std::string & getLabel () const
 This method returns the label being used by the timers in the orthogonalization manager.

Detailed Description

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

An implementation of the Belos::MatOrthoManager that performs orthogonalization using multiple steps of modified Gram-Schmidt.

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

Definition at line 54 of file BelosIMGSOrthoManager.hpp.


Constructor & Destructor Documentation

template<class ScalarType, class MV, class OP>
Belos::IMGSOrthoManager< ScalarType, MV, OP >::IMGSOrthoManager ( const std::string &  label = "Belos",
Teuchos::RCP< const OP >  Op = Teuchos::null,
const int  max_ortho_steps = 1,
const MagnitudeType  blk_tol = 10*MGT::squarerootMGT::eps() ),
const MagnitudeType  sing_tol = 10*MGT::eps() 
) [inline]

Constructor specifying re-orthogonalization tolerance.

Definition at line 68 of file BelosIMGSOrthoManager.hpp.

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

Destructor.

Definition at line 93 of file BelosIMGSOrthoManager.hpp.


Member Function Documentation

template<class ScalarType, class MV, class OP>
void Belos::IMGSOrthoManager< ScalarType, MV, OP >::setBlkTol ( const MagnitudeType  blk_tol  )  [inline]

Set parameter for block re-orthogonalization threshhold.

Definition at line 101 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
void Belos::IMGSOrthoManager< ScalarType, MV, OP >::setSingTol ( const MagnitudeType  sing_tol  )  [inline]

Set parameter for singular block detection.

Definition at line 104 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
MagnitudeType Belos::IMGSOrthoManager< ScalarType, MV, OP >::getBlkTol (  )  const [inline]

Return parameter for block re-orthogonalization threshhold.

Definition at line 107 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
MagnitudeType Belos::IMGSOrthoManager< ScalarType, MV, OP >::getSingTol (  )  const [inline]

Return parameter for singular block detection.

Definition at line 110 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
void Belos::IMGSOrthoManager< ScalarType, MV, OP >::project ( MV &  X,
Teuchos::RCP< MV >  MX,
Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > >  C,
Teuchos::Array< Teuchos::RCP< const MV > >  Q 
) const [virtual]

Given a list of (mutually and internally) orthonormal bases Q, this method takes a multivector X and projects it onto the space orthogonal to the individual Q[i], optionally returning the coefficients of X for the individual Q[i]. All of this is done with respect to the inner product innerProd().

After calling this routine, X will be orthogonal to each of the Q[i].

The method uses either one or two steps of modified Gram-Schmidt. The algebraically equivalent projection matrix is $P_Q = I - Q Q^H Op$, if Op is the matrix specified for use in the inner product. Note, this is not an orthogonal projector.

Parameters:
X [in/out] The multivector to be modified. On output, X will be orthogonal to Q[i] with respect to innerProd().
MX [in/out] The image of X under the operator Op. If $ MX != 0$: On input, this is expected to be consistent with X. On output, this is updated consistent with updates to X. If $ MX == 0$ or $ Op == 0$: MX is not referenced.
C [out] The coefficients of X in the *Q[i], with respect to innerProd(). If C[i] is a non-null pointer and *C[i] matches the dimensions of X and *Q[i], then the coefficients computed during the orthogonalization routine will be stored in the matrix *C[i]. If C[i] is a non-null pointer whose size does not match the dimensions of X and *Q[i], then a std::invalid_argument std::exception will be thrown. Otherwise, if C.size() < i or C[i] is a null pointer, then the orthogonalization manager will declare storage for the coefficients and the user will not have access to them.
Q [in] A list of multivector bases specifying the subspaces to be orthogonalized against. Each Q[i] is assumed to have orthonormal columns, and the Q[i] are assumed to be mutually orthogonal.

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

Definition at line 529 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
void Belos::IMGSOrthoManager< ScalarType, MV, OP >::project ( MV &  X,
Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > >  C,
Teuchos::Array< Teuchos::RCP< const MV > >  Q 
) const [inline, virtual]

This method calls project(X,Teuchos::null,C,Q); see documentation for that function.

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

Definition at line 152 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
int Belos::IMGSOrthoManager< ScalarType, MV, OP >::normalize ( MV &  X,
Teuchos::RCP< MV >  MX,
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > >  B 
) const [virtual]

This method takes a multivector X and attempts to compute an orthonormal basis for $colspan(X)$, with respect to innerProd().

The method uses modified Gram-Schmidt, so that the coefficient matrix B is upper triangular.

This routine returns an integer rank stating the rank of the computed basis. If X does not have full rank and the normalize() routine does not attempt to augment the subspace, then rank may be smaller than the number of columns in X. In this case, only the first rank columns of output X and first rank rows of B will be valid.

The method attempts to find a basis with dimension the same as the number of columns in X. It does this by augmenting linearly dependant vectors in X with random directions. A finite number of these attempts will be made; therefore, it is possible that the dimension of the computed basis is less than the number of vectors in X.

Parameters:
X [in/out] The multivector to the modified. On output, X will have some number of orthonormal columns (with respect to innerProd()).
MX [in/out] The image of X under the operator Op. If $ MX != 0$: On input, this is expected to be consistent with X. On output, this is updated consistent with updates to X. If $ MX == 0$ or $ Op == 0$: MX is not referenced.
B [out] The coefficients of the original X with respect to the computed basis. The first rows in B corresponding to the valid columns in X will be upper triangular.
Returns:
Rank of the basis computed by this method.

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

Definition at line 515 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
int Belos::IMGSOrthoManager< ScalarType, MV, OP >::normalize ( MV &  X,
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > >  B 
) const [inline, virtual]

This method calls normalize(X,Teuchos::null,B); see documentation for that function.

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

Definition at line 190 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
int Belos::IMGSOrthoManager< ScalarType, MV, OP >::projectAndNormalize ( MV &  X,
Teuchos::RCP< MV >  MX,
Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > >  C,
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > >  B,
Teuchos::Array< Teuchos::RCP< const MV > >  Q 
) const [virtual]

Given a set of bases Q[i] and a multivector X, this method computes an orthonormal basis for $colspan(X) - \sum_i colspan(Q[i])$.

This routine returns an integer rank stating the rank of the computed basis. If the subspace $colspan(X) - \sum_i colspan(Q[i])$ does not have dimension as large as the number of columns of X and the orthogonalization manager doe not attempt to augment the subspace, then rank may be smaller than the number of columns of X. In this case, only the first rank columns of output X and first rank rows of B will be valid.

The method attempts to find a basis with dimension the same as the number of columns in X. It does this by augmenting linearly dependant vectors with random directions. A finite number of these attempts will be made; therefore, it is possible that the dimension of the computed basis is less than the number of vectors in X.

Parameters:
X [in/out] The multivector to the modified. On output, the relevant rows of X will be orthogonal to the Q[i] and will have orthonormal columns (with respect to innerProd()).
MX [in/out] The image of X under the operator Op. If $ MX != 0$: On input, this is expected to be consistent with X. On output, this is updated consistent with updates to X. If $ MX == 0$ or $ Op == 0$: MX is not referenced.
C [out] The coefficients of the original X in the *Q[i], with respect to innerProd(). If C[i] is a non-null pointer and *C[i] matches the dimensions of X and *Q[i], then the coefficients computed during the orthogonalization routine will be stored in the matrix *C[i]. If C[i] is a non-null pointer whose size does not match the dimensions of X and *Q[i], then a std::invalid_argument std::exception will be thrown. Otherwise, if C.size() < i<> or C[i] is a null pointer, then the orthogonalization manager will declare storage for the coefficients and the user will not have access to them.
B [out] The coefficients of the original X with respect to the computed basis. The first rows in B corresponding to the valid columns in X will be upper triangular.
Q [in] A list of multivector bases specifying the subspaces to be orthogonalized against. Each Q[i] is assumed to have orthonormal columns, and the Q[i] are assumed to be mutually orthogonal.
Returns:
Rank of the basis computed by this method.

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

Definition at line 375 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
int Belos::IMGSOrthoManager< ScalarType, MV, OP >::projectAndNormalize ( MV &  X,
Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > >  C,
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > >  B,
Teuchos::Array< Teuchos::RCP< const MV > >  Q 
) const [inline, virtual]

This method calls projectAndNormalize(X,Teuchos::null,C,B,Q); see documentation for that function.

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

Definition at line 234 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
Teuchos::ScalarTraits<ScalarType>::magnitudeType Belos::IMGSOrthoManager< ScalarType, MV, OP >::orthonormError ( const MV &  X  )  const [inline, virtual]

This method computes the error in orthonormality of a multivector, measured as the Frobenius norm of the difference innerProd(X,Y) - I.

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

Definition at line 250 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
Teuchos::ScalarTraits< ScalarType >::magnitudeType Belos::IMGSOrthoManager< ScalarType, MV, OP >::orthonormError ( const MV &  X,
Teuchos::RCP< const MV >  MX 
) const [virtual]

This method computes the error in orthonormality of a multivector, measured as the Frobenius norm of the difference innerProd(X,Y) - I. The method has the option of exploiting a caller-provided MX.

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

Definition at line 349 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
Teuchos::ScalarTraits<ScalarType>::magnitudeType Belos::IMGSOrthoManager< ScalarType, MV, OP >::orthogError ( const MV &  X1,
const MV &  X2 
) const [inline, virtual]

This method computes the error in orthogonality of two multivectors, measured as the Frobenius norm of innerProd(X,Y).

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

Definition at line 265 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
Teuchos::ScalarTraits< ScalarType >::magnitudeType Belos::IMGSOrthoManager< ScalarType, MV, OP >::orthogError ( const MV &  X1,
Teuchos::RCP< const MV >  MX1,
const MV &  X2 
) const [virtual]

This method computes the error in orthogonality of two multivectors, measured as the Frobenius norm of innerProd(X,Y). The method has the option of exploiting a caller-provided MX.

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

Definition at line 364 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
void Belos::IMGSOrthoManager< ScalarType, MV, OP >::setLabel ( const std::string &  label  )  [virtual]

This method sets the label used by the timers in the orthogonalization manager.

Implements Belos::OrthoManager< ScalarType, MV >.

Definition at line 327 of file BelosIMGSOrthoManager.hpp.

template<class ScalarType, class MV, class OP>
const std::string& Belos::IMGSOrthoManager< ScalarType, MV, OP >::getLabel (  )  const [inline, virtual]

This method returns the label being used by the timers in the orthogonalization manager.

Implements Belos::OrthoManager< ScalarType, MV >.

Definition at line 287 of file BelosIMGSOrthoManager.hpp.


The documentation for this class was generated from the following file:
Generated on Wed May 12 21:45:52 2010 for Belos by  doxygen 1.4.7