NOX
Development
|
Concrete implementation for creating an Epetra_RowMatrix Jacobian via finite differencing of the residual using coloring. This method assumes the existence of a valid parallel coloring of the columns of the Jacobian (aka from Isorropia). More...
#include <NOX_Epetra_FiniteDifferenceColoringWithUpdate.H>
Public Member Functions | |
FiniteDifferenceColoringWithUpdate (Teuchos::ParameterList &printingParams, const Teuchos::RCP< Interface::Required > &i, const NOX::Epetra::Vector &initialGuess, const Teuchos::RCP< Epetra_MapColoring > &colorMap, double beta=1.0e-6, double alpha=1.0e-4) | |
Constructor with no frills. | |
FiniteDifferenceColoringWithUpdate (Teuchos::ParameterList &printingParams, const Teuchos::RCP< Interface::Required > &i, const NOX::Epetra::Vector &initialGuess, const Teuchos::RCP< Epetra_CrsGraph > &rawGraph, const Teuchos::RCP< Epetra_MapColoring > &colorMap, double beta=1.0e-6, double alpha=1.0e-4) | |
Constructor with graph. | |
FiniteDifferenceColoringWithUpdate (Teuchos::ParameterList &printingParams, const Teuchos::RCP< Interface::Required > &i, const NOX::Epetra::Vector &initialGuess, const Teuchos::RCP< Epetra_MapColoring > &colorMap, const Teuchos::RCP< Epetra_MapColoring > &updateColorMap, double beta=1.0e-6, double alpha=1.0e-4) | |
Constructor with update map. | |
FiniteDifferenceColoringWithUpdate (Teuchos::ParameterList &printingParams, const Teuchos::RCP< Interface::Required > &i, const NOX::Epetra::Vector &initialGuess, const Teuchos::RCP< Epetra_CrsGraph > &rawGraph, const Teuchos::RCP< Epetra_MapColoring > &colorMap, const Teuchos::RCP< Epetra_MapColoring > &updatecolorMap, double beta=1.0e-6, double alpha=1.0e-4) | |
Constructor with update map and graph. | |
virtual void | forceJacobianRecompute () |
virtual | ~FiniteDifferenceColoringWithUpdate () |
Pure virtual destructor. | |
virtual bool | computeJacobian (const Epetra_Vector &x, Epetra_Operator &Jac) |
Computes (or updates) the Jacobian given the specified input vector, x. Returns true if computation was successful. | |
virtual bool | computeJacobian (const Epetra_Vector &x) |
Computes (or updates) Jacobian given the specified input vector, x. Returns true if computation was successful. | |
void | SetProbingDiagnostics (bool activate) |
Disable/Enable the (low computational cost) probing diagnostics. | |
![]() | |
FiniteDifference (Teuchos::ParameterList &printingParams, const Teuchos::RCP< NOX::Epetra::Interface::Required > &i, const NOX::Epetra::Vector &initialGuess, double beta=1.0e-6, double alpha=1.0e-4) | |
Constructor with scalar beta. | |
FiniteDifference (Teuchos::ParameterList &printingParams, const Teuchos::RCP< NOX::Epetra::Interface::Required > &i, const NOX::Epetra::Vector &initialGuess, const Teuchos::RCP< const Epetra_Vector > &beta, double alpha=1.0e-4) | |
Constructor with vector beta. | |
FiniteDifference (Teuchos::ParameterList &printingParams, const Teuchos::RCP< NOX::Epetra::Interface::Required > &i, const NOX::Epetra::Vector &initialGuess, const Teuchos::RCP< Epetra_CrsGraph > &g, double beta=1.0e-6, double alpha=1.0e-4) | |
Constructor that takes a pre-constructed Epetra_CrsGraph so it does not have to determine the non-zero entries in the matrix. | |
FiniteDifference (Teuchos::ParameterList &printingParams, const Teuchos::RCP< NOX::Epetra::Interface::Required > &i, const NOX::Epetra::Vector &initialGuess, const Teuchos::RCP< Epetra_CrsGraph > &g, const Teuchos::RCP< const Epetra_Vector > &beta, double alpha=1.0e-4) | |
Constructor with output control that takes a pre-constructed Epetra_CrsGraph so it does not have to determine the non-zero entries in the matrix. | |
virtual | ~FiniteDifference () |
Pure virtual 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_Map & | OperatorDomainMap () const |
Returns the Epetra_BlockMap object associated with the domain of this matrix operator. | |
virtual const Epetra_Map & | OperatorRangeMap () 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_Comm & | Comm () const |
See Epetra_RowMatrix documentation. | |
virtual const Epetra_Map & | RowMatrixRowMap () const |
See Epetra_RowMatrix documentation. | |
virtual const Epetra_Map & | RowMatrixColMap () const |
See Epetra_RowMatrix documentation. | |
virtual const Epetra_Import * | RowMatrixImporter () const |
See Epetra_RowMatrix documentation. | |
virtual const Epetra_BlockMap & | Map () const |
See Epetra_SrcDistObj documentation. | |
virtual bool | computePreconditioner (const Epetra_Vector &x, Epetra_Operator &Prec, Teuchos::ParameterList *precParams=0) |
Compute an Epetra_RowMatrix to be used by Aztec preconditioners given the specified input vector, x. Returns true if computation was successful. | |
virtual void | setDifferenceMethod (DifferenceType type) |
Set the type of perturbation method used (default is Forward) | |
virtual Epetra_CrsMatrix & | getUnderlyingMatrix () const |
An accessor method for the underlying Epetra_CrsMatrix. | |
virtual void | Print (std::ostream &) const |
Output the underlying matrix. | |
void | setGroupForComputeF (NOX::Abstract::Group &group) |
Register a NOX::Abstract::Group derived object and use the computeF() method of that group for the perturbation instead of the NOX::Epetra::Interface::Required::computeF() method. This is required for LOCA to get the operators correct during homotopy. | |
![]() | |
Jacobian () | |
Constructor. | |
virtual | ~Jacobian () |
Destructor. | |
![]() | |
Preconditioner () | |
Constructor. | |
virtual | ~Preconditioner () |
Destructor. | |
Protected Member Functions | |
bool | differenceProbe (const Epetra_Vector &x, Epetra_CrsMatrix &Jac, const Epetra_MapColoring &colors) |
![]() | |
Teuchos::RCP< Epetra_CrsMatrix > | createGraphAndJacobian (Interface::Required &i, const Epetra_Vector &x) |
Constructs an Epetra_CrsGraph and Epetra_RowMatrix for the Jacobian. This is only called if the user does not supply an Epetra_CrsGraph. | |
bool | computeF (const Epetra_Vector &input, Epetra_Vector &result, NOX::Epetra::Interface::Required::FillType) |
Protected Attributes | |
bool | jacobianComputed |
bool | use_update |
bool | use_probing_diags |
Teuchos::RCP< Epetra_MapColoring > | colorMap_ |
Teuchos::RCP< Epetra_MapColoring > | updateColorMap_ |
![]() | |
const NOX::Utils | utils |
Printing Utilities object. | |
Teuchos::RCP< Epetra_CrsGraph > | graph |
Pointer to the Jacobian graph. | |
Teuchos::RCP< Epetra_CrsMatrix > | jacobian |
Pointer to the Jacobian. | |
Teuchos::RCP < NOX::Epetra::Interface::Required > | interface |
User provided interface function. | |
Epetra_Vector | x_perturb |
Perturbed solution vector - a work array that needs to be mutable. | |
Epetra_Vector | fo |
Function evaluation at currentX - a work array that needs to be mutable. | |
Epetra_Vector | fp |
Function evaluation at perturbX - a work array that needs to be mutable. | |
Teuchos::RCP< Epetra_Vector > | fmPtr |
Optional pointer to function evaluation at -perturbX - needed only for centered finite differencing. | |
Epetra_Vector | Jc |
Column vector of the jacobian - a work array that needs to be mutable. | |
double | alpha |
Constant for the perturbation calculation. | |
double | beta |
Constant for the perturbation calculation. | |
Teuchos::RCP< const Epetra_Vector > | betaVector |
Vector for the perturbation calculation. | |
BetaType | betaType |
Flag that sets whether ![]() | |
DifferenceType | diffType |
Define types for use of the perturbation parameter ![]() | |
std::string | label |
label for the Epetra_RowMatrix | |
bool | useGroupForComputeF |
Flag to enables the use of a group instead of the interface for the computeF() calls in the directional difference calculation. | |
Teuchos::RCP < NOX::Abstract::Group > | groupPtr |
Pointer to the group for possible use in computeF() calls. | |
Additional Inherited Members | |
![]() | |
enum | DifferenceType { Forward, Backward, Centered } |
Define types for use of the perturbation parameter ![]() | |
![]() | |
enum | BetaType { Scalar, Vector } |
Concrete implementation for creating an Epetra_RowMatrix Jacobian via finite differencing of the residual using coloring. This method assumes the existence of a valid parallel coloring of the columns of the Jacobian (aka from Isorropia).
Unlike the class NOX::FiniteDifferenceColoring, this class allows for "update" colorings, for use in situations where part of the Jacobian changes from iteration to iteration, but part does not. The first time (or any time after the forceJacobianRecompute method is called) the method uses the complete coloring. Afterwards, it uses the "update" coloring and only changes the entries that can change.
WARNING: The "update" coloring assumes that rows AND columns corresponding to uncolored (aka color 0) nodes do not change from call to call. If either the row or the column corresponding to a given node change then you must make sure it gets colored.
WARNING: Centered Difference Coloring is NOT supported as of yet.