Teko  Version of the Day
 All Classes Files Functions Variables Pages
Public Member Functions | List of all members
Teko::NS::InvLSCStrategy Class Reference

A strategy that takes a single inverse factory and uses that for all inverses. If no mass matrix is passed in the diagonal of the 1,1 block is used. More...

#include <Teko_InvLSCStrategy.hpp>

Inheritance diagram for Teko::NS::InvLSCStrategy:
Inheritance graph
[legend]

Public Member Functions

virtual void setPressureStabMatrix (const Teko::LinearOp &psm)
 
virtual void initializeState (const BlockedLinearOp &A, LSCPrecondState *state) const
 Initialize the state object using this blocked linear operator. More...
 
void computeInverses (const BlockedLinearOp &A, LSCPrecondState *state) const
 
virtual void setEigSolveParam (int sz)
 Set the number of power series iterations to use when finding the spectral radius. More...
 
virtual int getEigSolveParam ()
 Return the number of power series iterations to use when finding the spectral radius. More...
 
virtual void setUseFullLDU (bool val)
 Set to true to use the Full LDU decomposition, false otherwise. More...
 
virtual void setRowZeroing (bool val)
 Set to true to zero the rows of F when computing the spectral radius. More...
 
virtual void setMassMatrix (const LinearOp &mass)
 set the mass matrix to use in computing the scaling More...
 
virtual void setHScaling (const LinearOp &hScaling)
 
virtual void setHScaling (const MultiVector &hScaling)
 
virtual void setWScaling (const MultiVector &wScaling)
 
Constructors
 InvLSCStrategy ()
 
 InvLSCStrategy (const Teuchos::RCP< InverseFactory > &factory, bool rzn=false)
 
 InvLSCStrategy (const Teuchos::RCP< InverseFactory > &factory, LinearOp &mass, bool rzn=false)
 
 InvLSCStrategy (const Teuchos::RCP< InverseFactory > &invFactF, const Teuchos::RCP< InverseFactory > &invFactS, bool rzn=false)
 
 InvLSCStrategy (const Teuchos::RCP< InverseFactory > &invFactF, const Teuchos::RCP< InverseFactory > &invFactS, LinearOp &mass, bool rzn=false)
 
virtual void buildState (BlockedLinearOp &A, BlockPreconditionerState &state) const
 Functions inherited from LSCStrategy. More...
 
virtual LinearOp getInvBQBt (const BlockedLinearOp &A, BlockPreconditionerState &state) const
 
virtual LinearOp getInvBHBt (const BlockedLinearOp &A, BlockPreconditionerState &state) const
 
virtual LinearOp getInvF (const BlockedLinearOp &A, BlockPreconditionerState &state) const
 
virtual LinearOp getOuterStabilization (const BlockedLinearOp &A, BlockPreconditionerState &state) const
 
virtual LinearOp getInnerStabilization (const BlockedLinearOp &, BlockPreconditionerState &) const
 
virtual LinearOp getInvMass (const BlockedLinearOp &A, BlockPreconditionerState &state) const
 
virtual LinearOp getHScaling (const BlockedLinearOp &A, BlockPreconditionerState &state) const
 
virtual bool useFullLDU () const
 
virtual void setSymmetric (bool isSymmetric)
 
virtual void initializeFromParameterList (const Teuchos::ParameterList &pl, const InverseLibrary &invLib)
 Initialize from a parameter list. More...
 
virtual Teuchos::RCP
< Teuchos::ParameterList > 
getRequestedParameters () const
 For assiting in construction of the preconditioner. More...
 
virtual bool updateRequestedParameters (const Teuchos::ParameterList &pl)
 For assiting in construction of the preconditioner. More...
 
- Public Member Functions inherited from Teko::NS::LSCStrategy
void setRequestHandler (const Teuchos::RCP< RequestHandler > &rh)
 This method sets the request handler for this object. More...
 
Teuchos::RCP< RequestHandlergetRequestHandler () const
 This method gets the request handler uses by this object. More...
 

Detailed Description

A strategy that takes a single inverse factory and uses that for all inverses. If no mass matrix is passed in the diagonal of the 1,1 block is used.

A strategy that takes a single inverse factory and uses that for all inverses. Optionally the mass matrix can be passed in, if it is the diagonal is extracted and that is used to form the inverse approximation.

Definition at line 66 of file Teko_InvLSCStrategy.hpp.

Member Function Documentation

void Teko::NS::InvLSCStrategy::buildState ( BlockedLinearOp &  A,
BlockPreconditionerState state 
) const
virtual

Functions inherited from LSCStrategy.

This informs the strategy object to build the state associated with this operator.

Parameters
[in]AThe linear operator to be preconditioned by LSC.
[in]stateState object for storying reusable information about the operator A.

Implements Teko::NS::LSCStrategy.

Definition at line 167 of file Teko_InvLSCStrategy.cpp.

LinearOp Teko::NS::InvLSCStrategy::getInvBQBt ( const BlockedLinearOp &  A,
BlockPreconditionerState state 
) const
virtual

Get the inverse of $B Q_u^{-1} B^T$.

Parameters
[in]AThe linear operator to be preconditioned by LSC.
[in]stateState object for storying reusable information about the operator A.
Returns
An (approximate) inverse of $B Q_u^{-1} B^T$.

Implements Teko::NS::LSCStrategy.

Definition at line 202 of file Teko_InvLSCStrategy.cpp.

LinearOp Teko::NS::InvLSCStrategy::getInvBHBt ( const BlockedLinearOp &  A,
BlockPreconditionerState state 
) const
virtual

Get the inverse of $B H B^T - \gamma C$.

Parameters
[in]AThe linear operator to be preconditioned by LSC.
[in]stateState object for storying reusable information about the operator A.
Returns
An (approximate) inverse of $B H B^T - \gamma C$.

Implements Teko::NS::LSCStrategy.

Definition at line 207 of file Teko_InvLSCStrategy.cpp.

LinearOp Teko::NS::InvLSCStrategy::getInvF ( const BlockedLinearOp &  A,
BlockPreconditionerState state 
) const
virtual

Get the inverse of the $F$ block.

Parameters
[in]AThe linear operator to be preconditioned by LSC.
[in]stateState object for storying reusable information about the operator A.
Returns
An (approximate) inverse of $F$.

Implements Teko::NS::LSCStrategy.

Definition at line 212 of file Teko_InvLSCStrategy.cpp.

LinearOp Teko::NS::InvLSCStrategy::getOuterStabilization ( const BlockedLinearOp &  A,
BlockPreconditionerState state 
) const
virtual

Get the inverse for stabilizing the whole schur complement approximation.

Parameters
[in]AThe linear operator to be preconditioned by LSC.
[in]stateState object for storying reusable information about the operator A.
Returns
The operator to stabilize the whole Schur complement.

Implements Teko::NS::LSCStrategy.

Definition at line 217 of file Teko_InvLSCStrategy.cpp.

virtual LinearOp Teko::NS::InvLSCStrategy::getInnerStabilization ( const BlockedLinearOp &  A,
BlockPreconditionerState state 
) const
inlinevirtual

Get the inverse to stablized stabilizing the Schur complement approximation using a placement on the ``inside''. That is what is the value for $C_I$. This quantity may be null.

Parameters
[in]AThe linear operator to be preconditioned by LSC.
[in]stateState object for storying reusable information about the operator A.
Returns
The operator to stabilize the whole Schur complement.

Implements Teko::NS::LSCStrategy.

Definition at line 135 of file Teko_InvLSCStrategy.hpp.

LinearOp Teko::NS::InvLSCStrategy::getInvMass ( const BlockedLinearOp &  A,
BlockPreconditionerState state 
) const
virtual

Get the inverse mass matrix.

Parameters
[in]AThe linear operator to be preconditioned by LSC.
[in]stateState object for storying reusable information about the operator A.
Returns
The inverse of the mass matrix $Q_u$.

Implements Teko::NS::LSCStrategy.

Definition at line 226 of file Teko_InvLSCStrategy.cpp.

LinearOp Teko::NS::InvLSCStrategy::getHScaling ( const BlockedLinearOp &  A,
BlockPreconditionerState state 
) const
virtual

Get the $H$ scaling matrix.

Parameters
[in]AThe linear operator to be preconditioned by LSC.
[in]stateState object for storying reusable information about the operator A.
Returns
The $H$ scaling matrix.

Implements Teko::NS::LSCStrategy.

Definition at line 235 of file Teko_InvLSCStrategy.cpp.

virtual bool Teko::NS::InvLSCStrategy::useFullLDU ( ) const
inlinevirtual

Should the approximation of the inverse use a full LDU decomposition, or is a upper triangular approximation sufficient.

Returns
True if the full LDU decomposition should be used, otherwise only an upper triangular version is used.

Implements Teko::NS::LSCStrategy.

Definition at line 166 of file Teko_InvLSCStrategy.hpp.

virtual void Teko::NS::InvLSCStrategy::setSymmetric ( bool  isSymmetric)
inlinevirtual

Tell strategy that this operator is supposed to be symmetric. Behavior of LSC is slightly different for non-symmetric case.

Parameters
[in]isSymmetricIs this operator symmetric?

Implements Teko::NS::LSCStrategy.

Definition at line 173 of file Teko_InvLSCStrategy.hpp.

void Teko::NS::InvLSCStrategy::initializeFromParameterList ( const Teuchos::ParameterList &  pl,
const InverseLibrary &  invLib 
)
virtual

Initialize from a parameter list.

Reimplemented from Teko::NS::LSCStrategy.

Definition at line 483 of file Teko_InvLSCStrategy.cpp.

Teuchos::RCP< Teuchos::ParameterList > Teko::NS::InvLSCStrategy::getRequestedParameters ( ) const
virtual

For assiting in construction of the preconditioner.

Reimplemented from Teko::NS::LSCStrategy.

Definition at line 550 of file Teko_InvLSCStrategy.cpp.

bool Teko::NS::InvLSCStrategy::updateRequestedParameters ( const Teuchos::ParameterList &  pl)
virtual

For assiting in construction of the preconditioner.

Reimplemented from Teko::NS::LSCStrategy.

Definition at line 580 of file Teko_InvLSCStrategy.cpp.

virtual void Teko::NS::InvLSCStrategy::setPressureStabMatrix ( const Teko::LinearOp &  psm)
inlinevirtual

When computing the Schur complement, use the passed in matrix instead of $C$ to stabilize the gradient operator.

Definition at line 188 of file Teko_InvLSCStrategy.hpp.

void Teko::NS::InvLSCStrategy::initializeState ( const BlockedLinearOp &  A,
LSCPrecondState state 
) const
virtual

Initialize the state object using this blocked linear operator.

Definition at line 242 of file Teko_InvLSCStrategy.cpp.

void Teko::NS::InvLSCStrategy::computeInverses ( const BlockedLinearOp &  A,
LSCPrecondState state 
) const

Compute the inverses required for the LSC Schur complement

Note
This method assumes that the BQBt and BHBt operators have been constructed.

Definition at line 422 of file Teko_InvLSCStrategy.cpp.

virtual void Teko::NS::InvLSCStrategy::setEigSolveParam ( int  sz)
inlinevirtual

Set the number of power series iterations to use when finding the spectral radius.

Definition at line 204 of file Teko_InvLSCStrategy.hpp.

virtual int Teko::NS::InvLSCStrategy::getEigSolveParam ( )
inlinevirtual

Return the number of power series iterations to use when finding the spectral radius.

Definition at line 207 of file Teko_InvLSCStrategy.hpp.

virtual void Teko::NS::InvLSCStrategy::setUseFullLDU ( bool  val)
inlinevirtual

Set to true to use the Full LDU decomposition, false otherwise.

Definition at line 210 of file Teko_InvLSCStrategy.hpp.

virtual void Teko::NS::InvLSCStrategy::setRowZeroing ( bool  val)
inlinevirtual

Set to true to zero the rows of F when computing the spectral radius.

Definition at line 213 of file Teko_InvLSCStrategy.hpp.

virtual void Teko::NS::InvLSCStrategy::setMassMatrix ( const LinearOp &  mass)
inlinevirtual

set the mass matrix to use in computing the scaling

Definition at line 216 of file Teko_InvLSCStrategy.hpp.

virtual void Teko::NS::InvLSCStrategy::setHScaling ( const LinearOp &  hScaling)
inlinevirtual

Set the $H$-Scaling operator used in $B H B^T$. It is expected that this will be a diagonal matrix.

Definition at line 221 of file Teko_InvLSCStrategy.hpp.

virtual void Teko::NS::InvLSCStrategy::setHScaling ( const MultiVector &  hScaling)
inlinevirtual

Set the $H$-Scaling operator used in $B H B^T$. This method takes a vector and constructs the diagonal matrix.

Definition at line 226 of file Teko_InvLSCStrategy.hpp.

virtual void Teko::NS::InvLSCStrategy::setWScaling ( const MultiVector &  wScaling)
inlinevirtual

Set the $W$-Scaling vector used in $B H B^T$. This method takes a vector.

Definition at line 233 of file Teko_InvLSCStrategy.hpp.


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