NOX  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
LOCA::Epetra::AugmentedOp Class Reference

Epetra operator representing a $n+m$ bordered matrix. More...

#include <LOCA_Epetra_AugmentedOp.H>

Inheritance diagram for LOCA::Epetra::AugmentedOp:
Inheritance graph
[legend]
Collaboration diagram for LOCA::Epetra::AugmentedOp:
Collaboration graph
[legend]

Public Member Functions

 AugmentedOp (const Teuchos::RCP< LOCA::GlobalData > &global_data, const Teuchos::RCP< Epetra_Operator > &jac, const Teuchos::RCP< const Epetra_MultiVector > &a, const Teuchos::RCP< const Epetra_MultiVector > &b, const Teuchos::RCP< const NOX::Abstract::MultiVector::DenseMatrix > c)
 Constructor. More...
 
virtual ~AugmentedOp ()
 Destructor.
 
virtual int SetUseTranspose (bool UseTranspose)
 If set true, transpose of this operator will be applied. More...
 
virtual int Apply (const Epetra_MultiVector &Input, Epetra_MultiVector &Result) const
 Returns the result of a Epetra_Operator applied to a Epetra_MultiVector Input in Result. More...
 
virtual int ApplyInverse (const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
 Returns the result of a Epetra_Operator inverse applied to an Epetra_MultiVector Input in Result. More...
 
virtual double NormInf () const
 Returns the infinity norm of the bordered matrix. More...
 
virtual const char * Label () const
 Returns a character std::string describing the operator.
 
virtual bool UseTranspose () const
 Returns the current UseTranspose setting.
 
virtual bool HasNormInf () const
 Returns true if the this object can provide an approximate Inf-norm, false otherwise.
 
virtual const Epetra_CommComm () const
 Returns a reference to the Epetra_Comm communicator associated with this operator.
 
virtual const Epetra_MapOperatorDomainMap () const
 Returns the Epetra_Map object associated with the domain of this matrix operator.
 
virtual const Epetra_MapOperatorRangeMap () const
 Returns the Epetra_Map object associated with the range of this matrix operator.
 
virtual void init (const Epetra_MultiVector &x)
 Initialiazes operator for a solve.
 
virtual Teuchos::RCP
< Epetra_MultiVector
buildEpetraAugmentedMultiVec (const Epetra_MultiVector &x, const NOX::Abstract::MultiVector::DenseMatrix *y, bool doCopy) const
 Builds an extended vector from components. More...
 
virtual void setEpetraAugmentedMultiVec (Epetra_MultiVector &x, NOX::Abstract::MultiVector::DenseMatrix &y, const Epetra_MultiVector &augMultiVec) const
 Sets components from extended vector. More...
 

Protected Member Functions

void buildExtendedMap (const Epetra_BlockMap &map, Epetra_Map *&extMapPtr, bool buildImporter, bool haveParam)
 Builds extended domain, range maps.
 
int blockMap2PointMap (const Epetra_BlockMap &BlockMap, Epetra_Map *&PointMap) const
 Converts a block map to an equivalent point map.
 

Protected Attributes

Teuchos::RCP< LOCA::GlobalDataglobalData
 LOCA global data object.
 
std::string label
 Label for operator.
 
Teuchos::RCP< Epetra_OperatorjacOperator
 Stores operator representing $J$.
 
const Epetra_BlockMapunderlyingMap
 Stores underlying domain map.
 
const Epetra_CommunderlyingComm
 Stores comm.
 
Epetra_LocalMap localMap
 Local map for generating Epetra matrices.
 
Teuchos::RCP< const
Epetra_MultiVector
a
 Stores pointer to a multivector.
 
Teuchos::RCP< const
Epetra_MultiVector
b
 Stores pointer to b multivector.
 
Epetra_MultiVector c
 Stores c matrix.
 
int underlyingLength
 Stores underlying vector local length.
 
int numConstraints
 Number of constraints.
 
Epetra_MapextendedMapPtr
 Stores extended domain map.
 
Epetra_MapextendedImportMapPtr
 Stores extended turning point map for importing param component.
 
Epetra_ImportextendedImporter
 Stores importer object for importing param component.
 
Epetra_MultiVectorimportedInput
 Stores imported input multivector.
 
Epetra_MultiVectorresult_y
 Stores parameter component of result multivector.
 
Epetra_MultiVectortmp
 Stores temporary multivector used in ApplyInverse()
 
bool haveParamComponent
 Flag indicating whether we have the parameter component.
 
bool useTranspose
 Flag indicating whether to use transpose of operator.
 
Teuchos::LAPACK< int, double > dlapack
 LAPACK Wrappers.
 

Detailed Description

Epetra operator representing a $n+m$ bordered matrix.

The LOCA::Epetra::AugmentedOp is an Epetra_Operator representing the $n+m$ bordered matrix

\[ \begin{bmatrix} J & A \\ B^T & C \end{bmatrix} \]

where $J$ is an Epetra_Operator representing an $n\times n$ matrix, and $A$ and $B$ are length $n$ Epetra_MultiVector's with $m$ columns, and $C$ is an $m\times m$ dense matrix. It is assumed the Epetra_Map's for $A$, $B$, and $J$ are the same and the corresponding map for the bordered matrix is constructed from this map by storing the additional components on processor 0. The buildEpetraAugmentedMultiVec() method can be used to construct an Epetra_MultiVector using this map, a supplied length $n$ Epetra_MultiVector and an $m\times m$ matrix, while setEpetraAugmentedMultiVec() splits an extended multivector into its length $n$ and $m$ components. The Apply() method performs the $n+m\times n+m$ matrix multiplication while ApplyInverse() uses a block-elimination algorithm to compute the inverse using the ApplyInverse() method of the underlying operator $J$. In this way, linear systems of the form

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

can be solved in a matrix-free mode using the Apply() method. This operator can also represent a preconditioner of the form

\[ \begin{bmatrix} M & A \\ B^T & C \end{bmatrix} \]

using the ApplyInvese() method, where $M$ is a preconditioner for $J$. Note that if $J$ is nearly singular, the preconditioner should not be too good because otherwise the preconditining operation represented by ApplyInverse() becomes unstable.

Constructor & Destructor Documentation

LOCA::Epetra::AugmentedOp::AugmentedOp ( const Teuchos::RCP< LOCA::GlobalData > &  global_data,
const Teuchos::RCP< Epetra_Operator > &  jac,
const Teuchos::RCP< const Epetra_MultiVector > &  a,
const Teuchos::RCP< const Epetra_MultiVector > &  b,
const Teuchos::RCP< const NOX::Abstract::MultiVector::DenseMatrix c 
)

Constructor.

Builds the bordered operator using the supplied operator jac and Epetra_Vector's a and b. It is assumed a, b, and jac all have the same map.

References buildExtendedMap(), extendedImporter, extendedImportMapPtr, extendedMapPtr, haveParamComponent, Epetra_Comm::MyPID(), underlyingComm, underlyingMap, and View.

Member Function Documentation

int LOCA::Epetra::AugmentedOp::Apply ( const Epetra_MultiVector Input,
Epetra_MultiVector Result 
) const
virtual

Returns the result of a Epetra_Operator applied to a Epetra_MultiVector Input in Result.

Computes the extended matrix-vector product

\[ \begin{bmatrix} J & A \\ B^T & C \end{bmatrix} \begin{bmatrix} X \\ Y \end{bmatrix} = \begin{bmatrix} JX + AY \\ B^T X + CY \end{bmatrix} \]

or its transpose if UseTranpose() is true.

Implements Epetra_Operator.

References Insert, and View.

int LOCA::Epetra::AugmentedOp::ApplyInverse ( const Epetra_MultiVector X,
Epetra_MultiVector Y 
) const
virtual

Returns the result of a Epetra_Operator inverse applied to an Epetra_MultiVector Input in Result.

Solves the extended system

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

using the following block-elimination algorithm:

\[ \begin{split} X_1 &= J^{-1} F, \\ X_2 &= J^{-1} A, \\ Y &= (C-B^T X_2)^{-1}(G-B^T X_1), \\ X &= X_1 - X_2 Y \end{split} \]

If UseTranpose() is true, the tranpose of the system is solved.

Implements Epetra_Operator.

References Insert, and View.

Teuchos::RCP< Epetra_MultiVector > LOCA::Epetra::AugmentedOp::buildEpetraAugmentedMultiVec ( const Epetra_MultiVector x,
const NOX::Abstract::MultiVector::DenseMatrix y,
bool  doCopy 
) const
virtual

Builds an extended vector from components.

Builds an extended vector using the map representing the bordered matrix. If doCopy is true, the contents of x are copied into the extended vector, otherwise only space for the extended vector is created.

References Teuchos::rcp().

Referenced by LOCA::BorderedSolver::EpetraAugmented::applyInverse().

double LOCA::Epetra::AugmentedOp::NormInf ( ) const
virtual

Returns the infinity norm of the bordered matrix.

This is defined only if NormInf() of the underlying operator $J$ is defined and is given by $\|J\|_\infty+\|A\|_\infty+\|B\|_\infty$.

Implements Epetra_Operator.

void LOCA::Epetra::AugmentedOp::setEpetraAugmentedMultiVec ( Epetra_MultiVector x,
NOX::Abstract::MultiVector::DenseMatrix y,
const Epetra_MultiVector augMultiVec 
) const
virtual

Sets components from extended vector.

Splits the extended vector augMultiVec into components x and y by copying values out of extVec.

References Insert.

Referenced by LOCA::BorderedSolver::EpetraAugmented::applyInverse().

int LOCA::Epetra::AugmentedOp::SetUseTranspose ( bool  UseTranspose)
virtual

If set true, transpose of this operator will be applied.

Note that is only valid if the underlying operator $J$ supports a transpose.

Implements Epetra_Operator.


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