Anasazi Version of the Day
Public Member Functions | Static Public Member Functions
Anasazi::TsqrOrthoManagerImpl< ScalarType, MV > Class Template Reference

TSQR-based OrthoManager subclass implementation. More...

#include <AnasaziTsqrOrthoManagerImpl.hpp>

List of all members.

Public Member Functions

 TsqrOrthoManagerImpl (const Teuchos::RCP< const Teuchos::ParameterList > &params, const std::string &label)
void innerProd (const MV &X, const MV &Y, serial_matrix_type &Z) const
 Euclidean inner product.
void norm (const MV &X, std::vector< MagnitudeType > &normvec) const
void project (MV &X, prev_coeffs_type C, const_prev_mvs_type Q)
 Compute $C := Q^* X$ and $X := X - Q C$.
int normalize (MV &X, serial_matrix_ptr B)
 Orthogonalize the columns of X in place.
int normalizeNoCopy (MV &X, MV &Q, serial_matrix_ptr B)
 Normalize X into Q*B, overwriting X with invalid values.
int projectAndNormalize (MV &X, prev_coeffs_type C, serial_matrix_ptr B, const_prev_mvs_type Q)
int projectAndNormalizeNoCopy (MV &X_in, MV &X_out, prev_coeffs_type C, serial_matrix_ptr B, const_prev_mvs_type Q)
 Project and normalize X_in into X_out; overwrite X_in.
MagnitudeType orthonormError (const MV &X) const
 Return $ \| I - X^* \cdot X \|_F $.
MagnitudeType orthogError (const MV &X1, const MV &X2) const
 Return the Frobenius norm of the inner product of X1 and X1.
MagnitudeType blockReorthogThreshold () const
MagnitudeType relativeRankTolerance () const

Static Public Member Functions

static Teuchos::RCP< const
Teuchos::ParameterList
getDefaultParameters ()
 Get default parameters for TsqrOrthoManagerImpl.

Detailed Description

template<class ScalarType, class MV>
class Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >

TSQR-based OrthoManager subclass implementation.

TsqrOrthoManagerImpl implements the interface defined by OrthoManager. It doesn't actually inherit from OrthoManager, which gives us a bit more freedom when defining the actual subclass of OrthoManager (TsqrOrthoManager).

This class uses a combination of Tall Skinny QR (TSQR) and Block Gram-Schmidt (BGS) to orthogonalize multivectors.

The Block Gram-Schmidt procedure used here is inspired by that of G. W. Stewart ("Block Gram-Schmidt Orthogonalization", SISC vol 31 #1 pp. 761--775, 2008), except that we use TSQR+SVD instead of standard Gram-Schmidt with orthogonalization to handle the current block. "Orthogonalization faults" may still happen, but we do not handle them by default. Rather, we make one BGS pass, do TSQR+SVD, check the resulting column norms, and make a second BGS pass (+ TSQR+SVD) if necessary. If we then detect an orthogonalization fault, we throw TsqrOrthoFault.

Note:
Despite the "Impl" part of the name of this class, we don't actually use it for the "pImpl" idiom.

Definition at line 107 of file AnasaziTsqrOrthoManagerImpl.hpp.


Constructor & Destructor Documentation

template<class ScalarType , class MV >
Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >::TsqrOrthoManagerImpl ( const Teuchos::RCP< const Teuchos::ParameterList > &  params,
const std::string &  label 
)

Constructor

Parameters:
params[in] Configuration parameters, both for this orthogonalization manager, and for TSQR itself (as an RCP<const ParameterList> under "TsqrImpl"). Call the getDefaultParameters() class method for default parameters and their documentation, including TSQR implementation parameters.
label[in] Label for timers. This only matters if the compile-time option for enabling timers is set.

Definition at line 831 of file AnasaziTsqrOrthoManagerImpl.hpp.


Member Function Documentation

template<class ScalarType , class MV >
Teuchos::RCP< const Teuchos::ParameterList > Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >::getDefaultParameters ( ) [static]

Get default parameters for TsqrOrthoManagerImpl.

Get a (pointer to a) default list of parameters for configuring a TsqrOrthoManager or TsqrMatOrthoManager instance. The same parameters work for both.

Note:
To get nondefault behavior, a good thing to do is to make a deep copy of the returned parameter list, and then modify individual entries as desired.
TSQR implementation configuration options are stored under "TsqrImpl" as an RCP<const ParameterList>. (Don't call sublist() to get them; call get().)
Warning:
This method is not reentrant. It should only be called by one thread at a time.

Definition at line 1739 of file AnasaziTsqrOrthoManagerImpl.hpp.

template<class ScalarType , class MV >
void Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >::innerProd ( const MV &  X,
const MV &  Y,
serial_matrix_type Z 
) const [inline]

Euclidean inner product.

Compute the Euclidean block inner product X^* Y, and store the result in Z.

Parameters:
X[in]
Y[in]
Z[out] On output, $X^* Y$

Definition at line 170 of file AnasaziTsqrOrthoManagerImpl.hpp.

template<class ScalarType , class MV >
void Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >::norm ( const MV &  X,
std::vector< MagnitudeType > &  normvec 
) const [inline]

Compute the 2-norm of each column j of X

Parameters:
X[in] Multivector for which to compute column norms
normvec[out] On output: normvec[j] is the 2-norm of column j of X. normvec is resized if necessary so that it has at least as many entries as there are columns of X.
Note:
Performance of this method depends on how MultiVecTraits implements column norm computation for the given multivector type MV. It may or may not be the case that a reduction is performed for every column of X. Furthermore, whether or not the columns of X are contiguous (as opposed to a view of noncontiguous columns) may also affect performance. The computed results should be the same regardless, except perhaps for small rounding differences due to a different order of operations.

Definition at line 192 of file AnasaziTsqrOrthoManagerImpl.hpp.

template<class ScalarType , class MV >
void Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >::project ( MV &  X,
prev_coeffs_type  C,
const_prev_mvs_type  Q 
)

Compute $C := Q^* X$ and $X := X - Q C$.

Project X against the span of the (Euclidean) orthogonal vectors Q, and store the resulting coefficients in C.

Parameters:
X[in/out] On input: the vectors to project. On output: $X := X - Q C$ where $C := Q^* X$.
C[out] The projection coefficients $C := Q^* X$
Q[in] The orthogonal basis against which to project

FIXME (mfh 12 Jan 2011): Original: project(X, Q, C = tuple(serial_matrix_ptr(null))) const

Definition at line 859 of file AnasaziTsqrOrthoManagerImpl.hpp.

template<class ScalarType , class MV >
int Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >::normalize ( MV &  X,
serial_matrix_ptr  B 
)

Orthogonalize the columns of X in place.

Orthogonalize the columns of X in place, storing the resulting coefficients in B. Return the rank of X. If X is full rank, then X*B on output is a QR factorization of X on input. If X is not full rank, then the first rank columns of X on output form a basis for the column space of X (on input). Additional options control randomization of the null space basis.

Parameters:
X[in/out]
B[out]
Returns:
Rank of X

FIXME (mfh 12 Jan 2011): Original: normalize(X, B = null)

Definition at line 962 of file AnasaziTsqrOrthoManagerImpl.hpp.

template<class ScalarType , class MV >
int Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >::normalizeNoCopy ( MV &  X,
MV &  Q,
serial_matrix_ptr  B 
)

Normalize X into Q*B, overwriting X with invalid values.

Normalize X into Q*B, overwriting X with invalid values.

Note:
We expose this interface to applications because TSQR is not able to compute an orthogonal basis in place; it needs scratch space. Applications can exploit this interface to avoid excessive copying of vectors when using TSQR for orthogonalization.
Parameters:
X[in/out] Vector(s) to orthogonalize
B[out] Orthogonalization coefficients
Returns:
Rank of X
Note:
Q must have at least as many columns as X. It may have more columns than X; those columns are ignored.

Definition at line 1316 of file AnasaziTsqrOrthoManagerImpl.hpp.

template<class ScalarType , class MV >
int Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >::projectAndNormalize ( MV &  X,
prev_coeffs_type  C,
serial_matrix_ptr  B,
const_prev_mvs_type  Q 
)

Equivalent (in exact arithmetic) to project(X,C,Q) followed by normalize(X,B). However, this method performs reorthogonalization more efficiently and accurately.

Parameters:
X[in/out] The vectors to project against Q and normalize.
C[out] The projection coefficients
B[out] The normalization coefficients
Q[in] The orthogonal basis against which to project
Returns:
Rank of X after projection

FIXME (mfh 12 Jan 2011): Original: projectAndNormalize(X, Q, C = tuple(serial_matrix_ptr(null)), B = null)

Definition at line 1033 of file AnasaziTsqrOrthoManagerImpl.hpp.

template<class ScalarType , class MV >
int Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >::projectAndNormalizeNoCopy ( MV &  X_in,
MV &  X_out,
prev_coeffs_type  C,
serial_matrix_ptr  B,
const_prev_mvs_type  Q 
)

Project and normalize X_in into X_out; overwrite X_in.

Project X_in against Q, storing projection coefficients in C, and normalize X_in into X_out, storing normalization coefficients in B. On output, X_out has the resulting orthogonal vectors and X_in is overwritten with invalid values.

Parameters:
X_in[in/out] On input: The vectors to project against Q and normalize. Overwritten with invalid values on output.
X_out[out] On output: the normalized input vectors after projection against Q.
C[out] The projection coefficients
B[out] The normalization coefficients
Q[in] The orthogonal basis against which to project
Returns:
Rank of X_in after projection
Note:
We expose this interface to applications for the same reason that we expose normalizeNoCopy().

Definition at line 1451 of file AnasaziTsqrOrthoManagerImpl.hpp.

template<class ScalarType , class MV >
MagnitudeType Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >::orthonormError ( const MV &  X) const [inline]

Return $ \| I - X^* \cdot X \|_F $.

Return the Frobenius norm of I - X^* X, which is an absolute measure of the orthogonality of the columns of X.

Definition at line 302 of file AnasaziTsqrOrthoManagerImpl.hpp.

template<class ScalarType , class MV >
MagnitudeType Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >::orthogError ( const MV &  X1,
const MV &  X2 
) const [inline]

Return the Frobenius norm of the inner product of X1 and X1.

Definition at line 315 of file AnasaziTsqrOrthoManagerImpl.hpp.

template<class ScalarType , class MV >
MagnitudeType Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >::blockReorthogThreshold ( ) const [inline]

Relative tolerance for triggering a block reorthogonalization. If any column norm in a block decreases by this amount, then we reorthogonalize.

Definition at line 328 of file AnasaziTsqrOrthoManagerImpl.hpp.

template<class ScalarType , class MV >
MagnitudeType Anasazi::TsqrOrthoManagerImpl< ScalarType, MV >::relativeRankTolerance ( ) const [inline]

Relative tolerance for determining (via the SVD) whether a block is of full numerical rank.

Definition at line 332 of file AnasaziTsqrOrthoManagerImpl.hpp.


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