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

Concrete implementation for creating an Epetra_RowMatrix Jacobian via finite differencing of the residual using coloring. More...

#include <NOX_Epetra_FiniteDifferenceColoring.H>

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

Public Member Functions

 FiniteDifferenceColoring (Teuchos::ParameterList &printingParams, const Teuchos::RCP< Interface::Required > &i, const NOX::Epetra::Vector &initialGuess, const Teuchos::RCP< Epetra_MapColoring > &colorMap, const Teuchos::RCP< std::vector< Epetra_IntVector > > &columns, bool parallelColoring=false, bool distance1=false, double beta=1.0e-6, double alpha=1.0e-4)
 Constructor with output control.
 
 FiniteDifferenceColoring (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< std::vector< Epetra_IntVector > > &columns, bool parallelColoring=false, bool distance1=false, double beta=1.0e-6, double alpha=1.0e-4)
 Constructor with output control.
 
virtual ~FiniteDifferenceColoring ()
 Pure virtual destructor.
 
virtual bool computeJacobian (const Epetra_Vector &x, Epetra_Operator &Jac)
 Compute Jacobian given the specified input vector, x. Returns true if computation was successful.
 
virtual bool computeJacobian (const Epetra_Vector &x)
 Compute Jacobian given the specified input vector, x. Returns true if computation was successful.
 
virtual void createColorContainers ()
 Output the coloring map, index map and underlying matrix. More...
 
- Public Member Functions inherited from NOX::Epetra::FiniteDifference
 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_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.
 
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_CrsMatrixgetUnderlyingMatrix () 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.
 
- 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 Types

enum  ColoringType { NOX_SERIAL, NOX_PARALLEL }
 
- Protected Types inherited from NOX::Epetra::FiniteDifference
enum  BetaType { Scalar, Vector }
 

Protected Attributes

ColoringType coloringType
 Enum flag for type of coloring being used.
 
bool distance1
 bool flag for specifying special case of distance-1 coloring
 
Teuchos::RCP< const
Epetra_MapColoring
colorMap
 Color Map created by external algorithm.
 
Teuchos::RCP< std::vector
< Epetra_IntVector > > 
columns
 vector of Epetra_IntVectors containing columns corresponding to a given row and color
 
int numColors
 Number of colors in Color Map.
 
int maxNumColors
 Max Number of colors on all procs in Color Map.
 
int * colorList
 List of colors in Color Map.
 
std::list< int > listOfAllColors
 List of colors in Color Map.
 
std::map< int, int > colorToNumMap
 Inverse mapping from color id to its position in colorList.
 
Epetra_MapcMap
 Coloring Map created by external algorithm.
 
Epetra_ImportImporter
 Importer needed for mapping Color Map to unColored Map.
 
Epetra_VectorcolorVect
 Color vector based on Color Map containing perturbations.
 
Epetra_VectorbetaColorVect
 Color vector based on Color Map containing beta value(s)
 
Epetra_VectormappedColorVect
 Color vector based on unColorred Map containing perturbations.
 
Epetra_VectorxCol_perturb
 Perturbed solution vector based on column map.
 
const Epetra_BlockMapcolumnMap
 Overlap Map (Column Map of Matrix Graph) needed for parallel.
 
Epetra_ImportrowColImporter
 An Import object needed in parallel to map from row-space to column-space.
 
- Protected Attributes inherited from NOX::Epetra::FiniteDifference
const NOX::Utils utils
 Printing Utilities object.
 
Teuchos::RCP< Epetra_CrsGraphgraph
 Pointer to the Jacobian graph.
 
Teuchos::RCP< Epetra_CrsMatrixjacobian
 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_VectorfmPtr
 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_VectorbetaVector
 Vector for the perturbation calculation.
 
BetaType betaType
 Flag that sets whether $ \beta $ is a scalar or a vector.
 
DifferenceType diffType
 Define types for use of the perturbation parameter $ \delta$.
 
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

- Public Types inherited from NOX::Epetra::FiniteDifference
enum  DifferenceType { Forward, Backward, Centered }
 Define types for use of the perturbation parameter $ \delta$.
 
- Protected Member Functions inherited from NOX::Epetra::FiniteDifference
Teuchos::RCP< Epetra_CrsMatrixcreateGraphAndJacobian (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)
 

Detailed Description

Concrete implementation for creating an Epetra_RowMatrix Jacobian via finite differencing of the residual using coloring.

The Jacobian entries are calculated via 1st or 2nd order finite differencing. This requires $ N + 1 $ or $ 2N + 1 $ calls to computeF(), respectively, where $ N $ is the number of colors.

\[ J_{ij} = \frac{\partial F_i}{\partial x_j} = \frac{F_i(x+\delta\mathbf{e}_j) - F_i(x)}{\delta} \]

where $ J$ is the Jacobian, $ F$ is the function evaluation, $ x$ is the solution vector, and $\delta$ is a small perturbation to the $ x_j$ entry.

Instead of perturbing each $ N_{dof} $ problem degrees of freedom sequentially and then evaluating all $ N_{dof} $ functions for each perturbation, coloring allows several degrees of freedom (all belonging to the same color) to be perturbed at the same time. This reduces the total number of function evaluations needed to compute $\mathbf{J}$ from $ N_{dof}^2 $ as is required using FiniteDifference to $ N\cdot N_{dof} $, often representing substantial computational savings.

Coloring is based on a user-supplied color map generated using an appropriate algorithm, eg greedy-algorithm - Y. Saad, "Iterative Methods for Sparse %Linear Systems, 2<sup>nd</sup> ed.," chp. 3, SIAM, 2003.. Use can be made of the coloring algorithm provided by the EpetraExt package in Trilinos. The 1Dfem_nonlinearColoring and Brusselator example problems located in the nox/epetra-examples subdirectory demonstrate use of the EpetraExt package, and the 1Dfem_nonlinearColoring directory also contains a stand-alone coloring algorithm very similar to that in EpetraExt.

The perturbation, $ \delta $, is calculated using the following equation:

\[ \delta = \alpha * | x_j | + \beta \]

where $ \alpha $ is a scalar value (defaults to 1.0e-4) and $ \beta $ is another scalar (defaults to 1.0e-6).

Since both FiniteDifferenceColoring and FiniteDifference inherit from the Epetra_RowMatrix class, they can be used as preconditioning matrices for AztecOO preconditioners.

As for FiniteDifference, 1st order accurate Forward and Backward differences as well as 2nd order accurate Centered difference can be specified using setDifferenceMethod with the appropriate enumerated type passed as the argument.

Using FiniteDifferenceColoring in Parallel

Two ways of using this class in a distributed parallel environment are currently supported. From an application standpoint, the two approaches differ only in the status of the solution iterate used in the residual fill. If an object of this class is contructed with parallelColoring = true the solution iterate will be passe back in a non-ghosted form. On the contrary, setting this parameter to false in the constructor will cause the solution iterate to be in a ghosted form when calling back for a residual fill. When using the second approach, the user should be aware that the perturbed vector used to compute residuals has already been scattered to a form consistent with the column space of the Epetra_CrsGraph. In practice, this means that the perturbed vector used by computeF() has already been scattered to a ghosted or overlapped state. The application should then not perform this step but rather simply use the vector provided with the possible exception of requiring a local index reordering to bring the column-space based vector in sync with a potentially different ghosted index ordering. See the Brusselator and %1Dfem_nonlinearColoring example problems for details.

Special Case for Approximate Jacobian Construction

Provision is made for a simplified and cheaper use of coloring that currently provides only for the diagonal of the Jacobian to be computed. This is based on using a first-neighbors coloring of the original Jacobian graph using the Epetra_Ext MapColoring class with the distance1 argument set to true. This same argument should also be set to true in the constructor to this class. The result will be a diagonal Jacobian filled in a much more efficient manner.

Member Function Documentation

void FiniteDifferenceColoring::createColorContainers ( )
virtual

Output the coloring map, index map and underlying matrix.

Create containers for using color and index maps in parallel coloring

References colorList, colorMap, colorToNumMap, listOfAllColors, maxNumColors, and numColors.

Referenced by FiniteDifferenceColoring().


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