LOCA::BorderedSystem::EpetraHouseholder Class Reference

Bordered system solver strategy based on Householder transformations. More...

#include <LOCA_BorderedSystem_EpetraHouseholder.H>

Inheritance diagram for LOCA::BorderedSystem::EpetraHouseholder:

[legend]
Collaboration diagram for LOCA::BorderedSystem::EpetraHouseholder:
[legend]
List of all members.

Public Member Functions

 EpetraHouseholder (const Teuchos::RefCountPtr< LOCA::GlobalData > &global_data, const Teuchos::RefCountPtr< LOCA::Parameter::SublistParser > &topParams, const Teuchos::RefCountPtr< NOX::Parameter::List > &solverParams)
 Constructor.
virtual ~EpetraHouseholder ()
 Destructor.
virtual void setIsContiguous (bool flag)
 Set flag indicating whether F and A are stored in a continguous multivector.
virtual void setMatrixBlocks (const Teuchos::RefCountPtr< const NOX::Abstract::Group > &group, const Teuchos::RefCountPtr< const NOX::Abstract::MultiVector > &blockA, const Teuchos::RefCountPtr< const LOCA::MultiContinuation::ConstraintInterface > &blockB, const Teuchos::RefCountPtr< const NOX::Abstract::MultiVector::DenseMatrix > &blockC)
 Set blocks.
virtual NOX::Abstract::Group::ReturnType apply (const NOX::Abstract::MultiVector &X, const NOX::Abstract::MultiVector::DenseMatrix &Y, NOX::Abstract::MultiVector &U, NOX::Abstract::MultiVector::DenseMatrix &V) const
 Computed extended matrix-multivector product.
virtual NOX::Abstract::Group::ReturnType applyTranspose (const NOX::Abstract::MultiVector &X, const NOX::Abstract::MultiVector::DenseMatrix &Y, NOX::Abstract::MultiVector &U, NOX::Abstract::MultiVector::DenseMatrix &V) const
 Computed extended matrix transpose-multivector product.
virtual NOX::Abstract::Group::ReturnType applyInverse (NOX::Parameter::List &params, const NOX::Abstract::MultiVector *F, const NOX::Abstract::MultiVector::DenseMatrix *G, NOX::Abstract::MultiVector &X, NOX::Abstract::MultiVector::DenseMatrix &Y) const
 Solves the extended system using the technique described above.

Protected Member Functions

void factorConstraints ()
 Factors the constraint derivatives using the QR-factorization and generates the Y, T, and R matrices as described above.
void computeHouseholderVector (int col, const NOX::Abstract::MultiVector &A_x, const NOX::Abstract::MultiVector::DenseMatrix &A_p, NOX::Abstract::MultiVector &V_x, NOX::Abstract::MultiVector::DenseMatrix &V_p, double &beta)
 Computes the Householder vector V_x, V_p for column col of A_x and A_p, starting at row col of A_p.
void applyHouseholderVector (const NOX::Abstract::MultiVector &V_x, const NOX::Abstract::MultiVector::DenseMatrix &V_p, double beta, NOX::Abstract::MultiVector &A_x, NOX::Abstract::MultiVector::DenseMatrix &A_p)
 Applies the Householder vector V_x, V_p to the matrix sub-block represented by A_x and A_p.
void applyCompactWY (NOX::Abstract::MultiVector &x, NOX::Abstract::MultiVector::DenseMatrix &y, bool isZeroX, bool isZeroY, bool useTranspose) const
 Applies the operator Q as described above overwriting x and y. If either of x or y are zero on input, set the corresponding isZeroX or isZeroY flags. Set\ useTranspose to true to instead apply the transpose of Q.
void applyCompactWY (const NOX::Abstract::MultiVector *input_x, const NOX::Abstract::MultiVector::DenseMatrix *input_y, NOX::Abstract::MultiVector &result_x, NOX::Abstract::MultiVector::DenseMatrix &result_y, bool useTranspose) const
 Another version of applyCompactWY() that does not overwrite its inputs. If either input is zero, set the corresponding pointer to NULL.
NOX::Abstract::Group::ReturnType solveAZero (NOX::Parameter::List &params, const LOCA::MultiContinuation::ConstraintInterface *BB, const NOX::Abstract::MultiVector::DenseMatrix *CC, const NOX::Abstract::MultiVector *F, const NOX::Abstract::MultiVector::DenseMatrix *G, NOX::Abstract::MultiVector &X, NOX::Abstract::MultiVector::DenseMatrix &Y) const
 Solves the extended system of equations when the A block is zero using a simple block elimination.
NOX::Abstract::Group::ReturnType solveBZero (NOX::Parameter::List &params, const NOX::Abstract::MultiVector *AA, const NOX::Abstract::MultiVector::DenseMatrix *CC, const NOX::Abstract::MultiVector *F, const NOX::Abstract::MultiVector::DenseMatrix *G, NOX::Abstract::MultiVector &X, NOX::Abstract::MultiVector::DenseMatrix &Y) const
 Solves the extended system of equations when the B block is zero using a simple block elimination.

Protected Attributes

Teuchos::RefCountPtr< LOCA::GlobalDataglobalData
 Global data object.
Teuchos::RefCountPtr< NOX::Parameter::ListsolverParams
 Solver parameters.
Teuchos::RefCountPtr< const
LOCA::EpetraNew::Group
grp
 Pointer to group storing J.
Teuchos::RefCountPtr< const
NOX::Abstract::MultiVector
A
 Pointer to A block.
Teuchos::RefCountPtr< const
NOX::Abstract::MultiVector
B
 Pointer to B block.
Teuchos::RefCountPtr< const
NOX::Abstract::MultiVector::DenseMatrix
C
 Pointer to C block.
Teuchos::RefCountPtr< const
LOCA::MultiContinuation::ConstraintInterfaceMVDX
constraints
 Pointer to constraint interface.
Teuchos::RefCountPtr< NOX::Abstract::MultiVectorhouse_x
 Solution component of Householder multivec.
NOX::Abstract::MultiVector::DenseMatrix house_p
 Parameter component of Householder multivec.
NOX::Abstract::MultiVector::DenseMatrix T
 "T" matrix in compact WY representation
NOX::Abstract::MultiVector::DenseMatrix R
 R matrix in QR factorization.
Teuchos::RefCountPtr< NOX::Abstract::MultiVectorv_x
 Vector store 1 single Householder vector.
int numConstraints
 Number of constraint equations.
bool isZeroA
 flag indicating whether A block is zero
bool isZeroB
 flag indicating whether B block is zero
bool isZeroC
 flag indicating whether C block is zero
bool isContiguous
 flag indicating whether F and A are contiguous
Teuchos::BLAS< int, double > dblas
 BLAS Wrappers.

Detailed Description

Bordered system solver strategy based on Householder transformations.

This class solves the extended system of equations

\[ \begin{bmatrix} J & A \\ B^T & C \end{bmatrix} \begin{bmatrix} X \\ Y \end{bmatrix} = \begin{bmatrix} F \\ G \end{bmatrix} \]

using Householder tranformations. The algorithm works as follows: First consider a slightly rearranged version of the extended system of equations:

\[ \begin{bmatrix} C & B^T \\ A & J \end{bmatrix} \begin{bmatrix} Y \\ X \end{bmatrix} = \begin{bmatrix} G \\ F \end{bmatrix}. \]

Let

\[ Q^T \begin{bmatrix} C^T \\ B \end{bmatrix} = \begin{bmatrix} R \\ 0 \end{bmatrix} \]

be the QR decomposition of the constraints matrix where $Q\in\Re^{n+m\times n+m}$ and $R\in\Re^{m\times m}$. Define

\[ \begin{bmatrix} Z_Y \\ Z_X \end{bmatrix} = Q^T \begin{bmatrix} Y \\ X \end{bmatrix}, \]

then the extended system of equations is equivalent to

\[ \begin{bmatrix} R^T & 0 \\ [A & J] Q \end{bmatrix} \begin{bmatrix} Z_Y \\ Z_X \end{bmatrix} = \begin{bmatrix} G \\ F \end{bmatrix} \]

and hence

\[ \begin{split} Z_Y &= R^{-T} G \\ [A \;\; J] Q \begin{bmatrix} 0 \\ Z_X \end{bmatrix} &= F - [A \;\; J] Q \begin{bmatrix} Z_Y \\ 0 \end{bmatrix}. \end{split} \]

This last equation equation can be written

\[ P Z_X = \tilde{F} \]

where $P\in\Re^{n\times n}$ is given by

\[ P Z_X = [A \;\; J] Q \begin{bmatrix} 0 \\ Z_X \end{bmatrix} \]

and

\[ \tilde{F} = F - [A \;\; J] Q \begin{bmatrix} Z_Y \\ 0 \end{bmatrix}. \]

We then recover $X$ and $Y$ by

\[ \begin{bmatrix} Y \\ X \end{bmatrix} = Q \begin{bmatrix} Z_Y \\ Z_X \end{bmatrix}. \]

The equation $P Z_X = \tilde{F}$ is solved using an iterative solver using the definition of $P Z_X$ above, in this case AztecOO. The system is preconditioned using the preconditioner for $J$. The operator $Q$ is generated using the standard Householder QR algorithm (Algorithm 5.2.1, G. Golub and C. Van Loan, "Matrix Computations," 3rd Edition, Johns Hopkins, Baltimore, 1996) and is stored using the compact WY representation: $Q = I + Y^T T Y$ (see R. Schreiver and C. Van Loan, "A Storage-Efficient WY Represntation for Products of Householder Transformations," SIAM J. Sci. Stat. Comput., Vol. 10, No. 1, pp. 53-57, January 1989). The operator representing $P$ is encapsulated in the class LOCA::Epetra::CompactWYOp.


Constructor & Destructor Documentation

LOCA::BorderedSystem::EpetraHouseholder::EpetraHouseholder const Teuchos::RefCountPtr< LOCA::GlobalData > &  global_data,
const Teuchos::RefCountPtr< LOCA::Parameter::SublistParser > &  topParams,
const Teuchos::RefCountPtr< NOX::Parameter::List > &  solverParams
 

Constructor.

Parameters:
global_data [in] Global data object
topParams [in] Parsed top-level parameter list
solverParams [in] Bordered solver parameters. Currently none are referenced.


Member Function Documentation

NOX::Abstract::Group::ReturnType LOCA::BorderedSystem::EpetraHouseholder::apply const NOX::Abstract::MultiVector X,
const NOX::Abstract::MultiVector::DenseMatrix Y,
NOX::Abstract::MultiVector U,
NOX::Abstract::MultiVector::DenseMatrix V
const [virtual]
 

Computed extended matrix-multivector product.

Computes

\[ \begin{bmatrix} U \\ V \end{bmatrix} = \begin{bmatrix} J & A \\ B^T & C \end{bmatrix} \begin{bmatrix} X \\ Y \end{bmatrix} = \begin{bmatrix} J*X + A*Y \\ B^T*X + C*Y \end{bmatrix}. \]

Implements LOCA::BorderedSystem::AbstractStrategy.

NOX::Abstract::Group::ReturnType LOCA::BorderedSystem::EpetraHouseholder::applyInverse NOX::Parameter::List params,
const NOX::Abstract::MultiVector F,
const NOX::Abstract::MultiVector::DenseMatrix G,
NOX::Abstract::MultiVector X,
NOX::Abstract::MultiVector::DenseMatrix Y
const [virtual]
 

Solves the extended system using the technique described above.

The params argument is the linear solver parameters. If isZeroF or isZeroG is true, than the corresponding F or G pointers may be NULL.

Note that if either the A or B blocks are zero, the system is solved using a simple block elimination scheme instead of the Householder scheme.

Implements LOCA::BorderedSystem::AbstractStrategy.

NOX::Abstract::Group::ReturnType LOCA::BorderedSystem::EpetraHouseholder::applyTranspose const NOX::Abstract::MultiVector X,
const NOX::Abstract::MultiVector::DenseMatrix Y,
NOX::Abstract::MultiVector U,
NOX::Abstract::MultiVector::DenseMatrix V
const [virtual]
 

Computed extended matrix transpose-multivector product.

Computes

\[ \begin{bmatrix} U \\ V \end{bmatrix} = \begin{bmatrix} J^T & B \\ A^T & C \end{bmatrix} \begin{bmatrix} X \\ Y \end{bmatrix} = \begin{bmatrix} J^T*X + B*Y \\ A^T*X + C^T*Y \end{bmatrix}. \]

Implements LOCA::BorderedSystem::AbstractStrategy.

void LOCA::BorderedSystem::EpetraHouseholder::computeHouseholderVector int  col,
const NOX::Abstract::MultiVector A_x,
const NOX::Abstract::MultiVector::DenseMatrix A_p,
NOX::Abstract::MultiVector V_x,
NOX::Abstract::MultiVector::DenseMatrix V_p,
double &  beta
[protected]
 

Computes the Householder vector V_x, V_p for column col of A_x and A_p, starting at row col of A_p.

The algorithm implemented here is essentially algorithm 5.1.1 of Golub and Van Loan.

void LOCA::BorderedSystem::EpetraHouseholder::setMatrixBlocks const Teuchos::RefCountPtr< const NOX::Abstract::Group > &  group,
const Teuchos::RefCountPtr< const NOX::Abstract::MultiVector > &  blockA,
const Teuchos::RefCountPtr< const LOCA::MultiContinuation::ConstraintInterface > &  blockB,
const Teuchos::RefCountPtr< const NOX::Abstract::MultiVector::DenseMatrix > &  blockC
[virtual]
 

Set blocks.

The blockA or blockC pointer may be null if either is zero. Whether block B is zero will be determined by querying blockB via ConstraintInterface::isConstraintDerivativesXZero.

Implements LOCA::BorderedSystem::AbstractStrategy.


The documentation for this class was generated from the following files:
Generated on Thu Sep 18 12:40:43 2008 for NOX by doxygen 1.3.9.1