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

A concrete implementation of a Broyden-type operator for NOX. More...

#include <NOX_Epetra_BroydenOperator.H>

Inheritance diagram for NOX::Epetra::BroydenOperator:
Inheritance graph
[legend]
Collaboration diagram for NOX::Epetra::BroydenOperator:
Collaboration graph
[legend]

Classes

class  ReplacementInterface
 

Public Member Functions

 BroydenOperator (Teuchos::ParameterList &nlParams, const Teuchos::RCP< NOX::Utils > &utils, Epetra_Vector &solnVec, const Teuchos::RCP< Epetra_CrsMatrix > &broydMat0, bool verbose=false)
 Constructor taking an initial matrix to be updated.
 
 BroydenOperator (const BroydenOperator &)
 Copy Constructor.
 
virtual ~BroydenOperator ()
 Destructor.
 
virtual const char * Label () const
 Returns a character std::string describing the name of the operator.
 
virtual int SetUseTranspose (bool UseTranspose)
 If set true, the transpose of this operator will be applied.
 
virtual int Apply (const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
 Return the result on an Epetra_Operator applied to an Epetra_MultiVector X in Y.
 
virtual int ApplyInverse (const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
 Return the result on an Epetra_Operator inverse applied to an Epetra_MultiVector X in Y.
 
virtual bool UseTranspose () const
 Returns the current use transpose setting.
 
virtual bool HasNormInf () const
 Returns true if the this object can provide an approximate Inf-norm, false otherwise.
 
virtual const Epetra_MapOperatorDomainMap () const
 Returns the Epetra_BlockMap object associated with the domain of this matrix operator.
 
virtual const Epetra_MapOperatorRangeMap () const
 Returns the Epetra_BlockMap object associated with the range of this matrix operator.
 
virtual bool Filled () const
 See Epetra_RowMatrix documentation.
 
virtual int NumMyRowEntries (int MyRow, int &NumEntries) const
 See Epetra_RowMatrix documentation.
 
virtual int MaxNumEntries () const
 See Epetra_RowMatrix documentation.
 
virtual int ExtractMyRowCopy (int MyRow, int Length, int &NumEntries, double *Values, int *Indices) const
 See Epetra_RowMatrix documentation.
 
virtual int ExtractDiagonalCopy (Epetra_Vector &Diagonal) const
 See Epetra_RowMatrix documentation.
 
virtual int Multiply (bool TransA, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
 See Epetra_RowMatrix documentation.
 
virtual int Solve (bool Upper, bool Trans, bool UnitDiagonal, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
 See Epetra_RowMatrix documentation.
 
virtual int InvRowSums (Epetra_Vector &x) const
 See Epetra_RowMatrix documentation.
 
virtual int LeftScale (const Epetra_Vector &x)
 See Epetra_RowMatrix documentation.
 
virtual int InvColSums (Epetra_Vector &x) const
 See Epetra_RowMatrix documentation.
 
virtual int RightScale (const Epetra_Vector &x)
 See Epetra_RowMatrix documentation.
 
virtual double NormInf () const
 See Epetra_RowMatrix documentation.
 
virtual double NormOne () const
 See Epetra_RowMatrix documentation.
 
virtual int NumGlobalNonzeros () const
 See Epetra_RowMatrix documentation.
 
virtual long long NumGlobalNonzeros64 () const
 
virtual int NumGlobalRows () const
 See Epetra_RowMatrix documentation.
 
virtual long long NumGlobalRows64 () const
 
virtual int NumGlobalCols () const
 See Epetra_RowMatrix documentation.
 
virtual long long NumGlobalCols64 () const
 
virtual int NumGlobalDiagonals () const
 See Epetra_RowMatrix documentation.
 
virtual long long NumGlobalDiagonals64 () const
 
virtual int NumMyNonzeros () const
 See Epetra_RowMatrix documentation.
 
virtual int NumMyRows () const
 See Epetra_RowMatrix documentation.
 
virtual int NumMyCols () const
 See Epetra_RowMatrix documentation.
 
virtual int NumMyDiagonals () const
 See Epetra_RowMatrix documentation.
 
virtual bool LowerTriangular () const
 See Epetra_RowMatrix documentation.
 
virtual bool UpperTriangular () const
 See Epetra_RowMatrix documentation.
 
virtual const Epetra_CommComm () const
 See Epetra_RowMatrix documentation.
 
virtual const Epetra_MapRowMatrixRowMap () const
 See Epetra_RowMatrix documentation.
 
virtual const Epetra_MapRowMatrixColMap () const
 See Epetra_RowMatrix documentation.
 
virtual const Epetra_ImportRowMatrixImporter () const
 See Epetra_RowMatrix documentation.
 
virtual const Epetra_BlockMapMap () const
 See Epetra_SrcDistObj documentation.
 
void setStepVector (Epetra_Vector &vec)
 Set the current step vector,

\[ y_k = x_{k+1} - x_k \]

.

 
void setStepVector (NOX::Epetra::Vector &vec)
 Set the current step vector,

\[ y_k = x_{k+1} - x_k \]

.

 
void setYieldVector (Epetra_Vector &vec)
 Set the current yield vector,

\[ y_k = F_{k+1} - F_k \]

.

 
void setYieldVector (NOX::Epetra::Vector &vec)
 Set the current yield vector,

\[ y_k = F_{k+1} - F_k \]

.

 
bool computeSparseBroydenUpdate ()
 Compute the sparse Broyden update.
 
void removeEntriesFromBroydenUpdate (const Epetra_CrsGraph &graph)
 Remove entries from being involved in Broyden updates.
 
const Epetra_CrsMatrixgetBroydenMatrix ()
 Return a const reference to the Broyden matrix. The matrix is not owned but is obtained from the client at construction.
 
void resetBroydenMatrix (const Epetra_CrsMatrix &mat)
 Reset the values of our matrix.
 
void addReplacementInterface (ReplacementInterface *i)
 Register replacement interface.
 
"Is" functions

Checks to see if various objects have been computed. Returns true if the corresponding "compute" function has been called since the last update to the solution vector (via instantiation or computeX).

virtual bool isStep () const
 
virtual bool isYield () const
 
virtual bool isBroyden () const
 
- Public Member Functions inherited from NOX::Observer
 Observer ()
 Constructor.
 
virtual ~Observer ()
 Destructor.
 
virtual void runPreSolutionUpdate (const NOX::Abstract::Vector &, const NOX::Solver::Generic &)
 User defined method that will be executed prior to the update of the solution vector during a call to NOX::Solver::Generic::step(). This is intended to allow users to adjust the direction before the solution update, typically based on knowledge of the problem formulation. The direction is const as we can't guarantee that changes to the direction won't violate assumptions of the solution algorithm. Users can change the update/direciton after a const cast, but NOX may not function as expected. Use at your own risk! More...
 
virtual void runPostSolutionUpdate (const NOX::Solver::Generic &)
 User defined method that will be executed after the update of the solution vector during a call to NOX::Solver::Generic::step(). This is intended to allow users to adjust the direction after the solution update, typically based on knowledge of the problem formulation (e.g. clipping negative mass fractions). The direction is const as we can't guarantee that changes to the direction won't violate assumptions of the solution algorithm. Users can change the update/direciton after a const cast, but NOX may not function as expected. Use at your own risk! More...
 
virtual void runPreLineSearch (const NOX::Solver::Generic &)
 User defined method that will be executed before a call to NOX::LineSearch::Generic::compute(). Only to be used in NOX::Solver::LineSearchBased!
 
virtual void runPostLineSearch (const NOX::Solver::Generic &)
 User defined method that will be executed after a call to NOX::LineSearch::Generic::compute(). Only to be used in NOX::Solver::LineSearchBased!
 
- Public Member Functions inherited from NOX::Epetra::Interface::Jacobian
 Jacobian ()
 Constructor.
 
virtual ~Jacobian ()
 Destructor.
 
- Public Member Functions inherited from NOX::Epetra::Interface::Preconditioner
 Preconditioner ()
 Constructor.
 
virtual ~Preconditioner ()
 Destructor.
 

Protected Member Functions

virtual bool initialize (Teuchos::ParameterList &nlParams, const Epetra_Vector &x)
 Initialize operator and data.
 
virtual bool computeJacobian (const Epetra_Vector &x, Epetra_Operator &Jac)
 Allow for fill of Jacobian matrix whose values will replace corresponding entries in the Broyden matrix.
 
virtual bool computePreconditioner (const Epetra_Vector &, Epetra_Operator &Prec, Teuchos::ParameterList *params=0)
 Allow for fill of preconditioning matrix whose values will replace corresponding entries in the Broyden matrix.
 
virtual void runPreSolve (const NOX::Solver::Generic &solver)
 Set a flag to skip the first call to computeSparseBroydenUpdate as no valid step vector or residual vector is available.
 
virtual void runPreIterate (const NOX::Solver::Generic &solver)
 Needed to preserve any existing call through to a user-defined Pre/Preost Operator.
 
virtual void runPostIterate (const NOX::Solver::Generic &solver)
 Update the Broyden matrix using changes in residuals the solution vector from the most recent nnlinear iteration.
 
virtual void runPostSolve (const NOX::Solver::Generic &solver)
 Needed to preserve any existing call through to a user-defined Pre/Preost Operator.
 
void replaceBroydenMatrixValues (const Epetra_CrsMatrix &mat)
 Replace values in Broyden matrix with either Jacobian or preconditioning matrix entries.
 

Protected Attributes

bool verbose
 
Teuchos::RCP< NOX::Epetra::VectorstepVec
 
Teuchos::RCP< NOX::Epetra::VectoryieldVec
 
Teuchos::RCP< NOX::Epetra::VectorworkVec
 
Teuchos::RCP< NOX::Epetra::VectoroldX
 
Teuchos::RCP< NOX::Epetra::VectoroldF
 
Teuchos::RCP< Epetra_CrsMatrixcrsMatrix
 
Teuchos::RCP
< NOX::Epetra::Interface::Jacobian
jacIntPtr
 
Teuchos::RCP< Epetra_CrsMatrixjacMatrixPtr
 
Teuchos::RCP
< NOX::Epetra::Interface::Preconditioner
precIntPtr
 
Teuchos::RCP< Epetra_CrsMatrixprecMatrixPtr
 
Teuchos::ParameterListnlParams
 Reference to top-level nonlinear solver parameters list.
 
const Teuchos::RCP< NOX::Utils > & utils
 Reference to NOX::Utils object.
 
Teuchos::RCP< NOX::Observerobserver
 Pointer to a user defined NOX::Observer object.
 
std::string label
 label for the Epetra_RowMatrix
 
std::vector< bool > entriesRemoved
 Flag to signal removal of some entries.
 
std::map< int, std::list< int > > retainedEntries
 Container of entries to omit from Broyden updates.
 
std::vector
< ReplacementInterface * > 
replacementInterfaces
 Container of entries to omit from Broyden updates.
 
IsValid flags

True if objects are current with respect to the currect stepVec.

bool isValidStep
 
bool isValidYield
 
bool isValidBroyden
 

Detailed Description

A concrete implementation of a Broyden-type operator for NOX.

This operator is intended to allow cheap updates to an existing Jacobian or preconditioning matrix that would otherwise be difficult or impossible to obtain by other means. It computes updates using secant approximations emobdied in truncated Broyden updates that preserve matrix sparsity.

This class derives from NOX::Abstract::PrePostOperator in order to perform a Broyden-type update on an existing matrix that it holds but does not own. This update is performed after each nonlinear iteration within method runPostIterate(...) according to the recursive formula:

\[ \tilde{B}_{k+1} = \tilde{B}_k + \frac{({y_k - \tilde{B}_k s_k})s_k^T}{s^T s} \]

where

\[ y_k = F_{k+1} - F_k \]

and

\[ s_k = x_{k+1} - x_k \]

The tilde on the matrices $ B $ indicates that the updates are constrained so that the nonzero structure of the original matrix passed into the constructor is preserved. Inasmuch as unconstrained Broyden updates produce dense matrices, these constrained updates lead to a loss of Broyden-matrix properties, e.g.

\[ \tilde{B}_{k+1} s_k \ne \tilde{B}_k + s_k \]

\[ \tilde{B}_{k+1} q \ne \tilde{B}_k q \quad \forall q : s_k^T q = 0 \]

One could recover these properties by passing into the constructor a dense Epetra_CrsMatrix, though the cost of typical use of this matrix, e.g. applying ILU to it, would be significant. Additionally, "better" values obtained from another Jacobian or preconditioning matrix can be used to replace corresponding values in the updated Broyden matrix by passing the Jacobian or preconditioning matrix and its associated interface to the constructor. The structure of the Jacobain or preconditioning matrix typically represents a subset of the Broyden matrix, e.g. a block diagonal matrix.


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