Teko  Version of the Day
 All Classes Files Functions Variables Pages
Public Member Functions | Protected Attributes | Related Functions | List of all members
Teko::LU2x2InverseOp Class Reference

This linear operator approximates the inverse of a block $ 2\times 2 $ operator using a block $ LDU $ decomposition. More...

#include <Teko_LU2x2InverseOp.hpp>

Inheritance diagram for Teko::LU2x2InverseOp:
Inheritance graph
[legend]

Public Member Functions

 LU2x2InverseOp (const BlockedLinearOp &A, const LinearOp &invA00, const LinearOp &invS)
 This constructor explicitly takes the parts of $ A $ required to build the inverse operator. More...
 
 LU2x2InverseOp (const BlockedLinearOp &A, const LinearOp &hatInvA00, const LinearOp &tildeInvA00, const LinearOp &invS)
 This constructor explicitly takes the parts of $ A $ required to build the inverse operator. More...
 
Inherited methods from Thyra::LinearOpBase
virtual VectorSpace range () const
 Range space of this operator. More...
 
virtual VectorSpace domain () const
 Domain space of this operator. More...
 
virtual void implicitApply (const BlockedMultiVector &x, BlockedMultiVector &y, const double alpha=1.0, const double beta=0.0) const
 Perform a matrix vector multiply with this operator. More...
 
- Public Member Functions inherited from Teko::BlockImplicitLinearOp
virtual void implicitApply (const Thyra::EOpTransp M_trans, const BlockedMultiVector &x, BlockedMultiVector &y, const double alpha=1.0, const double beta=0.0) const
 Perform a matrix vector multiply with this implicitly defined blocked operator. More...
 

Protected Attributes

const BlockedLinearOp A_
 operator $ A $ More...
 
const LinearOp hatInvA00_
 inverse of $ A_{00} $ More...
 
const LinearOp tildeInvA00_
 inverse of $ A_{00} $ More...
 
const LinearOp invS_
 inverse of $ S $ More...
 
const LinearOp A10_
 operator $ A_{10} $ More...
 
const LinearOp A01_
 operator $ A_{01} $ More...
 
Teuchos::RCP< const
Thyra::ProductVectorSpaceBase
< double > > 
productRange_
 Range vector space. More...
 
Teuchos::RCP< const
Thyra::ProductVectorSpaceBase
< double > > 
productDomain_
 Domain vector space. More...
 

Related Functions

(Note that these are not member functions.)

ModifiableLinearOp createDiagnosticLinearOp (const Teuchos::RCP< std::ostream > &os, const ModifiableLinearOp &A, const std::string &label)
 Constructor method for building DiagnosticLinearOp. More...
 
ModifiableLinearOp createDiagnosticLinearOp (const Teuchos::RCP< std::ostream > &os, const LinearOp &A, const std::string &label)
 Constructor method for building DiagnosticLinearOp. More...
 
ModifiableLinearOp createDiagnosticLinearOp (const Teuchos::RCP< std::ostream > &os, const Teko::LinearOp &fwdOp, const ModifiableLinearOp &A, const std::string &label)
 Constructor method for building DiagnosticLinearOp. More...
 
LinearOp createLU2x2InverseOp (BlockedLinearOp &A, LinearOp &invA00, LinearOp &invS)
 Constructor method for building LU2x2InverseOp. More...
 
LinearOp createLU2x2InverseOp (BlockedLinearOp &A, LinearOp &invA00, LinearOp &invS, const std::string &str)
 Constructor method for building LU2x2InverseOp. More...
 
LinearOp createLU2x2InverseOp (BlockedLinearOp &A, LinearOp &hatInvA00, LinearOp &tildeInvA00, LinearOp &invS)
 Constructor method for building LU2x2InverseOp. More...
 
LinearOp createLU2x2InverseOp (BlockedLinearOp &A, LinearOp &hatInvA00, LinearOp &tildeInvA00, LinearOp &invS, const std::string &str)
 Constructor method for building LU2x2InverseOp. More...
 

Additional Inherited Members

- Protected Member Functions inherited from Teko::BlockImplicitLinearOp
virtual bool opSupportedImpl (const Thyra::EOpTransp M_trans) const
 Functions required by Thyra::LinearOpBase. More...
 

Detailed Description

This linear operator approximates the inverse of a block $ 2\times 2 $ operator using a block $ LDU $ decomposition.

For a matrix that is blocked like

$ A = \left[\begin{array}{cc} A_{00} & A_{01} \\ A_{10} & A_{11} \end{array}\right] $

this class evaluates the $A^{-1}$ given $A_{00}^{-1}$ and the inverse of the Schur complement. The $ LDU $ factorization is defined as

$ A = \left[ \begin{array}{cc} I & 0 \\ A_{10} A_{00}^{-1} & I \end{array} \right] \left[ \begin{array}{cc} A_{00} & 0 \\ 0 & -S \end{array} \right] \left[ \begin{array}{cc} I & A_{00}^{-1} A_{01} \\ 0 & I \end{array} \right]$

where the Schur complement is $ S=-A_{11}+A_{10} A_{00}^{-1} A_{01} $ . In order to do this 2 evaluations of $ A_{00}^{-1} $ and a single evalution of $ S^{-1} $ are needed. For increased flexibility both evaluations of $A_{00}^{-1}$ can be specified independently. For righthand side vector $[f, g]^T$ and solution vector $[u,v]^T$ the two inverses ( $A$-hat and $A$-tilde) are needed to evaluate

$\hat{A}_{00} u^* = f$,

$\tilde{A}_{00} v = A_{01} v$

where $u^*$ is an intermediate step.

Definition at line 101 of file Teko_LU2x2InverseOp.hpp.

Constructor & Destructor Documentation

Teko::LU2x2InverseOp::LU2x2InverseOp ( const BlockedLinearOp &  A,
const LinearOp &  invA00,
const LinearOp &  invS 
)

This constructor explicitly takes the parts of $ A $ required to build the inverse operator.

This constructor explicitly takes the parts of $ A $ required to build the inverse operator.

Parameters
[in]AThe block $ 2 \times 2 $ $A$ operator.
[in]invA00An approximate inverse of $ A_{00} $, used for both $\hat{A}_{00}$ and $\tilde{A}_{00}$
[in]invSAn approximate inverse of $ S = -A_{11} + A_{10} A_{00}^{-1} A_{01} $.

Definition at line 63 of file Teko_LU2x2InverseOp.cpp.

Teko::LU2x2InverseOp::LU2x2InverseOp ( const BlockedLinearOp &  A,
const LinearOp &  hatInvA00,
const LinearOp &  tildeInvA00,
const LinearOp &  invS 
)

This constructor explicitly takes the parts of $ A $ required to build the inverse operator.

This constructor explicitly takes the parts of $ A $ required to build the inverse operator.

Parameters
[in]AThe block $ 2 \times 2 $ $A$ operator.
[in]hatInvA00An approximate inverse of $ \hat{A}_{00} $
[in]tildeInvA00An approximate inverse of $ \tilde{A}_{00} $
[in]invSAn approximate inverse of $ S = -A_{11} + A_{10} A_{00}^{-1} A_{01} $.

Definition at line 90 of file Teko_LU2x2InverseOp.cpp.

Member Function Documentation

virtual VectorSpace Teko::LU2x2InverseOp::range ( ) const
inlinevirtual

Range space of this operator.

Implements Teko::BlockImplicitLinearOp.

Definition at line 137 of file Teko_LU2x2InverseOp.hpp.

virtual VectorSpace Teko::LU2x2InverseOp::domain ( ) const
inlinevirtual

Domain space of this operator.

Implements Teko::BlockImplicitLinearOp.

Definition at line 140 of file Teko_LU2x2InverseOp.hpp.

void Teko::LU2x2InverseOp::implicitApply ( const BlockedMultiVector &  x,
BlockedMultiVector &  y,
const double  alpha = 1.0,
const double  beta = 0.0 
) const
virtual

Perform a matrix vector multiply with this operator.

The apply function takes one vector as input and applies the inverse $ LDU $ decomposition. The result is returned in $y$. If this operator is reprsented as $M$ then $ y = \alpha M x + \beta y $ (ignoring conjugation!).

Parameters
[in]x
[in,out]y
[in]alpha(default=1)
[in]beta(default=0)

Implements Teko::BlockImplicitLinearOp.

Definition at line 107 of file Teko_LU2x2InverseOp.cpp.

Friends And Related Function Documentation

ModifiableLinearOp createDiagnosticLinearOp ( const Teuchos::RCP< std::ostream > &  os,
const ModifiableLinearOp &  A,
const std::string &  label 
)
related

Constructor method for building DiagnosticLinearOp.

Constructor method for building DiagnosticLinearOp.

Parameters
[in]osOutput stream to print diagnostics to
[in]AOperator to be wrapped
[in]labelString for outputing with diagnostics
Returns
A linear operator that wrapping A that will print diagnostics on descruction.

Definition at line 172 of file Teko_DiagnosticLinearOp.hpp.

ModifiableLinearOp createDiagnosticLinearOp ( const Teuchos::RCP< std::ostream > &  os,
const LinearOp &  A,
const std::string &  label 
)
related

Constructor method for building DiagnosticLinearOp.

Constructor method for building DiagnosticLinearOp.

Parameters
[in]osOutput stream to print diagnostics to
[in]AOperator to be wrapped
[in]labelString for outputing with diagnostics
Returns
A linear operator that wrapping A that will print diagnostics on descruction.

Definition at line 190 of file Teko_DiagnosticLinearOp.hpp.

ModifiableLinearOp createDiagnosticLinearOp ( const Teuchos::RCP< std::ostream > &  os,
const Teko::LinearOp &  fwdOp,
const ModifiableLinearOp &  A,
const std::string &  label 
)
related

Constructor method for building DiagnosticLinearOp.

Constructor method for building DiagnosticLinearOp.

Parameters
[in]osOutput stream to print diagnostics to
[in]fwdOpForward operator to compute residual with
[in]AOperator to be wrapped
[in]labelString for outputing with diagnostics
Returns
A linear operator that wrapping A that will print diagnostics on descruction.

Definition at line 209 of file Teko_DiagnosticLinearOp.hpp.

LinearOp createLU2x2InverseOp ( BlockedLinearOp &  A,
LinearOp &  invA00,
LinearOp &  invS 
)
related

Constructor method for building LU2x2InverseOp.

Constructor method for building LU2x2InverseOp.

Parameters
[in]A2x2 Operator to be decomposed
[in]invA00Approximate inverse of the operators $(0,0)$ block.
[in]invSApproximate inverse of the Schur complement
Returns
A linear operator that behaves like the inverse of the LU decomposition.

Definition at line 196 of file Teko_LU2x2InverseOp.hpp.

LinearOp createLU2x2InverseOp ( BlockedLinearOp &  A,
LinearOp &  invA00,
LinearOp &  invS,
const std::string &  str 
)
related

Constructor method for building LU2x2InverseOp.

Constructor method for building LU2x2InverseOp.

Parameters
[in]A2x2 Operator to be decomposed
[in]invA00Approximate inverse of the operators $(0,0)$ block.
[in]invSApproximate inverse of the Schur complement
[in]strString to label the operator
Returns
A linear operator that behaves like the inverse of the LU decomposition.

Definition at line 215 of file Teko_LU2x2InverseOp.hpp.

LinearOp createLU2x2InverseOp ( BlockedLinearOp &  A,
LinearOp &  hatInvA00,
LinearOp &  tildeInvA00,
LinearOp &  invS 
)
related

Constructor method for building LU2x2InverseOp.

Constructor method for building LU2x2InverseOp.

Parameters
[in]A2x2 Operator to be decomposed
[in]hatInvA00First approximate inverse of the operators $(0,0)$ block.
[in]tildeInvA00Second approximate inverse of the operators $(0,0)$ block.
[in]invSApproximate inverse of the Schur complement
Returns
A linear operator that behaves like the inverse of the LU decomposition.

Definition at line 237 of file Teko_LU2x2InverseOp.hpp.

LinearOp createLU2x2InverseOp ( BlockedLinearOp &  A,
LinearOp &  hatInvA00,
LinearOp &  tildeInvA00,
LinearOp &  invS,
const std::string &  str 
)
related

Constructor method for building LU2x2InverseOp.

Constructor method for building LU2x2InverseOp.

Parameters
[in]A2x2 Operator to be decomposed
[in]hatInvA00First approximate inverse of the operators $(0,0)$ block.
[in]tildeInvA00Second approximate inverse of the operators $(0,0)$ block.
[in]invSApproximate inverse of the Schur complement
[in]strString to label the operator
Returns
A linear operator that behaves like the inverse of the LU decomposition.

Definition at line 257 of file Teko_LU2x2InverseOp.hpp.

Member Data Documentation

const BlockedLinearOp Teko::LU2x2InverseOp::A_
protected

operator $ A $

Definition at line 165 of file Teko_LU2x2InverseOp.hpp.

const LinearOp Teko::LU2x2InverseOp::hatInvA00_
protected

inverse of $ A_{00} $

Definition at line 166 of file Teko_LU2x2InverseOp.hpp.

const LinearOp Teko::LU2x2InverseOp::tildeInvA00_
protected

inverse of $ A_{00} $

Definition at line 167 of file Teko_LU2x2InverseOp.hpp.

const LinearOp Teko::LU2x2InverseOp::invS_
protected

inverse of $ S $

Definition at line 168 of file Teko_LU2x2InverseOp.hpp.

const LinearOp Teko::LU2x2InverseOp::A10_
protected

operator $ A_{10} $

Definition at line 171 of file Teko_LU2x2InverseOp.hpp.

const LinearOp Teko::LU2x2InverseOp::A01_
protected

operator $ A_{01} $

Definition at line 172 of file Teko_LU2x2InverseOp.hpp.

Teuchos::RCP<const Thyra::ProductVectorSpaceBase<double> > Teko::LU2x2InverseOp::productRange_
protected

Range vector space.

Definition at line 174 of file Teko_LU2x2InverseOp.hpp.

Teuchos::RCP<const Thyra::ProductVectorSpaceBase<double> > Teko::LU2x2InverseOp::productDomain_
protected

Domain vector space.

Definition at line 175 of file Teko_LU2x2InverseOp.hpp.


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