AztecOO  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
AztecOO Class Reference

AztecOO: An object-oriented wrapper for Aztec. More...

#include <AztecOO.h>

Collaboration diagram for AztecOO:
Collaboration graph
[legend]

Classes

struct  MatrixData
 
struct  OperatorData
 

Public Member Functions

int SetPreconditioner (AZ_PRECOND *Prec)
 AztecOO External Preconditioner Set (object) More...
 
int SetPreconditioner (AZ_PREC_FUN prec_function, void *prec_data)
 AztecOO External Preconditioner Set (function and data) More...
 
int SetScaling (struct AZ_SCALING *Scaling)
 AztecOO External Scaling Set. More...
 
int SetMatrixName (int label)
 AztecOO Label Matrix for Aztec. More...
 
void SetLabel (const char *const Label)
 Set Label this AztecOO object. More...
 
const char * GetLabel () const
 Get the label describing this AztecOO object. More...
 
Constructors/destructors.
 AztecOO (Epetra_Operator *A, Epetra_MultiVector *X, Epetra_MultiVector *B)
 AztecOO Constructor. More...
 
 AztecOO (Epetra_RowMatrix *A, Epetra_MultiVector *X, Epetra_MultiVector *B)
 AztecOO Constructor. More...
 
 AztecOO (const Epetra_LinearProblem &LinearProblem)
 AztecOO Constructor. More...
 
 AztecOO ()
 AztecOO Default constructor.
 
 AztecOO (const AztecOO &Solver)
 AztecOO Copy Constructor. More...
 
virtual ~AztecOO (void)
 AztecOO Destructor. More...
 
Post-construction setup methods.
int SetProblem (const Epetra_LinearProblem &prob, bool call_SetPrecMatrix=false)
 AztecOO Epetra_LinearProblem Set. More...
 
int SetUserOperator (Epetra_Operator *UserOperator)
 AztecOO User Operator Set. More...
 
int SetUserMatrix (Epetra_RowMatrix *UserMatrix, bool call_SetPrecMatrix=false)
 AztecOO User Matrix Set. More...
 
int SetLHS (Epetra_MultiVector *X)
 AztecOO LHS Set. More...
 
int SetRHS (Epetra_MultiVector *B)
 AztecOO RHS Set. More...
 
int UnsetLHSRHS ()
 AztecOO unset LHS and RHS. More...
 
int SetPrecMatrix (Epetra_RowMatrix *PrecMatrix)
 AztecOO Preconditioner Matrix Set. More...
 
int SetPrecOperator (Epetra_Operator *PrecOperator)
 AztecOO External Preconditioner Set. More...
 
int SetStatusTest (AztecOO_StatusTest *StatusTest)
 AztecOO External Convergence/Status Test Set. More...
 
void SetOutputStream (std::ostream &ostrm)
 Set std::ostream for Aztec's screen output. More...
 
void SetErrorStream (std::ostream &errstrm)
 Set std::ostream for Aztec's error output. More...
 
Explicit preconditioner construction/assessment/destruction methods.
int ConstructPreconditioner (double &condest)
 Forces explicit construction and retention of an AztecOO native preconditioner. More...
 
int DestroyPreconditioner ()
 Destroys a preconditioner computed using ConstructPreconditioner(). More...
 
double Condest () const
 Returns the condition number estimate for the current preconditioner, if one exists, returns -1.0 if no estimate.
 
Check/Attribute Access Methods.
int CheckInput () const
 Prints a summary of solver parameters, performs simple sanity checks.
 
Epetra_LinearProblemGetProblem () const
 Get a pointer to the Linear Problem used to construct this solver; returns zero if not available.
 
Epetra_OperatorGetUserOperator () const
 Get a pointer to the user operator A.
 
Epetra_RowMatrixGetUserMatrix () const
 Get a pointer to the user matrix A.
 
Epetra_OperatorGetPrecOperator () const
 Get a pointer to the preconditioner operator.
 
Epetra_RowMatrixGetPrecMatrix () const
 Get a pointer to the matrix used to construct the preconditioner.
 
Epetra_MultiVectorGetLHS () const
 Get a pointer to the left-hand-side X.
 
Epetra_MultiVectorGetRHS () const
 Get a pointer to the right-hand-side B.
 
void PrintLinearSystem (const char *name)
 Print linear-system to files. More...
 
Standard AztecOO option and parameter setting methods.
int SetParameters (Teuchos::ParameterList &parameterlist, bool cerr_warning_if_unused=false)
 Method to set options/parameters using a ParameterList object. More...
 
int SetAztecDefaults ()
 AztecOO function to restore default options/parameter settings. More...
 
int SetAztecOption (int option, int value)
 AztecOO option setting function. More...
 
int GetAztecOption (int option)
 AztecOO option getting function. More...
 
int SetAztecParam (int param, double value)
 AztecOO param setting function. More...
 
const int * GetAllAztecOptions () const
 AztecOO option setting function. More...
 
const double * GetAllAztecParams () const
 AztecOO param setting function. More...
 
int SetAllAztecOptions (const int *options)
 AztecOO option setting function. More...
 
int SetAllAztecParams (const double *params)
 AztecOO param setting function. More...
 
Standard AztecOO solve methods.
int Iterate (long long MaxIters, double Tolerance)
 AztecOO iteration function. More...
 
int Iterate (Epetra_RowMatrix *A, Epetra_MultiVector *X, Epetra_MultiVector *B, int MaxIters, double Tolerance)
 AztecOO iteration function. More...
 
Specialist AztecOO solve method.
int recursiveIterate (int MaxIters, double Tolerance)
 AztecOO iteration functions. More...
 
const double * GetAztecStatus () const
 Return the Aztec status after iterating. More...
 
Adaptive Solve methods.
int SetUseAdaptiveDefaultsTrue ()
 Force the AdaptiveIterate() method to use default adaptive strategy.
 
int SetAdaptiveParams (int NumTrials, double *athresholds, double *rthresholds, double condestThreshold, double maxFill, int maxKspace)
 Set the parameter that control the AdaptiveIterate() method. More...
 
int AdaptiveIterate (int MaxIters, int MaxSolveAttempts, double Tolerance)
 Attempts to solve the given linear problem using an adaptive strategy.
 
Post-solve access functions
int NumIters () const
 Returns the total number of iterations performed on this problem.
 
double TrueResidual () const
 Returns the true unscaled residual for this problem.
 
double ScaledResidual () const
 Returns the true scaled residual for this problem.
 
double RecursiveResidual () const
 Returns the recursive residual for this problem.
 
double SolveTime () const
 Returns the solve time.
 
int GetAllAztecStatus (double *status)
 AztecOO status extraction function. More...
 

Protected Member Functions

int AllocAzArrays ()
 
void DeleteAzArrays ()
 
int SetAztecVariables ()
 
int SetProblemOptions (ProblemDifficultyLevel PDL, bool ProblemSymmetric)
 
int SetProcConfig (const Epetra_Comm &Comm)
 
void DeleteMemory ()
 

Protected Attributes

char * Label_
 
Epetra_LinearProblemProblem_
 
Epetra_MultiVectorX_
 
Epetra_MultiVectorB_
 
Epetra_VectorResidualVector_
 
int N_local_
 
int x_LDA_
 
double * x_
 
int b_LDA_
 
double * b_
 
int * proc_config_
 
int * options_
 
double * params_
 
double * status_
 
AZ_MATRIXAmat_
 
AZ_MATRIXPmat_
 
AZ_PRECONDPrec_
 
struct AZ_SCALINGScaling_
 
bool Scaling_created_
 
AztecOO_StatusTestStatusTest_
 
struct AZ_CONVERGE_STRUCTconv_info_
 
double condest_
 
bool useAdaptiveDefaults_
 
int NumTrials_
 
double maxFill_
 
int maxKspace_
 
double * athresholds_
 
double * rthresholds_
 
double condestThreshold_
 
bool inConstructor_
 
bool procConfigSet_
 
MatrixDataUserMatrixData_
 
MatrixDataPrecMatrixData_
 
OperatorDataUserOperatorData_
 
OperatorDataPrecOperatorData_
 
std::ostream * out_stream_
 
std::ostream * err_stream_
 

Detailed Description

AztecOO: An object-oriented wrapper for Aztec.

Currently it accepts a Petra matrix, initial guess and RHS as separate arguments, or alternatively, accepts a Epetra_LinearProblem. If constructed using a Epetra_LinearProblem, AztecOO will infer some solver/preconditioner, etc., options and parameters. Users may override these choices and manually choose from among the full set of Aztec options using the SetAztecOption() and SetAztecParam() functions.

AztecOO will solve a linear systems of equations: $ AX=B $, using Epetra objects and the Aztec solver library, where $ A$ is an Epetra_Operator or Epetra_RowMatrix (note that the Epetra_Operator class is a base class for Epetra_RowMatrix so that Epetra_RowMatrix isa Epetra_Operator.) $ X$ and $ B$ are Epetra_MultiVector objects.

Warning
AztecOO does not presently support solution of more than one simultaneous right-hand-side.

Constructor & Destructor Documentation

AztecOO::AztecOO ( Epetra_Operator A,
Epetra_MultiVector X,
Epetra_MultiVector B 
)

AztecOO Constructor.

Creates a AztecOO instance, passing in already-defined objects for the linear operator (as an Epetra_Operator), left-hand-side and right-hand-side.

Note: Use of this constructor may prohibit use of native AztecOO preconditioners, since an Epetra_Operator is not necessarily an Epetra_RowMatrix and all AztecOO incomplete factorization preconditioners are based on having explicit access to matrix coefficients. Polynomial preconditioners are available if the Epetra_Operator passed in here has a non-trivial definition of the NormInf() method and HasNormInf() returns true.

References SetAztecDefaults(), SetLHS(), SetRHS(), SetUserMatrix(), and SetUserOperator().

AztecOO::AztecOO ( Epetra_RowMatrix A,
Epetra_MultiVector X,
Epetra_MultiVector B 
)

AztecOO Constructor.

Creates a AztecOO instance, passing in already-defined objects for the linear operator (as an Epetra_RowMatrix), left-hand-side and right-hand-side.

Note: Use of this constructor allows full access to native AztecOO preconditioners, using the Epetra_RowMatrix A passed in here as the basis for computing the preconditioner. All AztecOO incomplete factorization preconditioners are based on having explicit access to matrix coefficients. Polynomial preconditioners are also available. It is possible to change the matrix used for computing incomplete factorization by calling the SetPrecMatrix() method. It is also possible to provide a user-supplied preconditioner via SetPrecOperator().

References SetAztecDefaults(), SetLHS(), SetRHS(), and SetUserMatrix().

AztecOO::AztecOO ( const Epetra_LinearProblem LinearProblem)

AztecOO Constructor.

Creates a AztecOO instance, using a Epetra_LinearProblem, passing in an already-defined Epetra_LinearProblem object. The Epetra_LinearProblem class is the preferred method for passing in the linear problem to AztecOO because this class provides scaling capabilities and self-consistency checks that are not available when using other constructors.

Note: If the Epetra_LinearProblem passed in here has a non-trivial pointer to an Epetra_Matrix then use of this constructor allows full access to native AztecOO preconditioners, using the Epetra_RowMatrix A passed in here as the basis for computing the preconditioner. All AztecOO incomplete factorization preconditioners are based on having explicit access to matrix coefficients. Polynomial preconditioners are also available. It is possible to change the matrix used for computing incomplete factorization by calling the SetPrecMatrix() method. It is also possible to provide a user-supplied preconditioner by call SetPrecOperator().

If the Epetra_LinearProblems passed in here has only an Epetra_Operator, then use of this constructor may prohibit use of native AztecOO preconditioners, since an Epetra_Operator is not necessarily an Epetra_RowMatrix and all AztecOO incomplete factorization preconditioners are based on having explicit access to matrix coefficients. Polynomial preconditioners are available if the Epetra_Operator passed in here has a non-trivial definition of the NormInf() method and HasNormInf() returns true.

References SetAztecDefaults(), and SetProblem().

AztecOO::AztecOO ( const AztecOO Solver)
AztecOO::~AztecOO ( void  )
virtual

AztecOO Destructor.

Completely deletes a AztecOO object.

Member Function Documentation

int AztecOO::ConstructPreconditioner ( double &  condest)

Forces explicit construction and retention of an AztecOO native preconditioner.

AztecOO typically constructs the preconditioner on the first call to the solve function. However, there are situations where we would like to compute the preconditioner ahead of time. One particular case is when we want to confirm that the preconditioner well-conditioned. This method allows us to precompute the preconditioner. It also provides a estimate of the condition number of the preconditioner. If condest is large, e.g., > 1.0e+14, it is likely the preconditioner will fail. In this case, using threshold values (available in the incomplete factorizations) can be used to reduce the condition number.

Note: This method does not work for user-defined preconditioners (defined via calls to SetPrecOperator(). It will return with an error code of -1 for this case.

References Epetra_Operator::Comm(), Epetra_Comm::MaxAll(), Epetra_Comm::MinAll(), and Epetra_RowMatrix::NumMyCols().

Referenced by AdaptiveIterate().

int AztecOO::DestroyPreconditioner ( )

Destroys a preconditioner computed using ConstructPreconditioner().

The ConstructPreconditioner() method creates a persistent preconditioner. In other words the preconditioner will be used by all calls to the Iterate() method. DestroyPreconditioner() deletes the current preconditioner and restores AztecOO to a state where the preconditioner will computed on first use of the preconditioner solve.

Referenced by AdaptiveIterate(), and SetPreconditioner().

const int* AztecOO::GetAllAztecOptions ( ) const
inline

AztecOO option setting function.

Return a pointer to an array (size AZ_OPTIONS_SIZE) of all of the currently set aztec options.

const double* AztecOO::GetAllAztecParams ( ) const
inline

AztecOO param setting function.

Return a pointer to an array (size AZ_PARAMS_SIZE) of all of the currently set aztec parameters.

int AztecOO::GetAllAztecStatus ( double *  status)
inline

AztecOO status extraction function.

Extract Aztec status array into user-provided array. The array must be of length AZ_STATUS_SIZE as defined in the az_aztec.h header file.

Referenced by AZOO_iterate().

int AztecOO::GetAztecOption ( int  option)
inline

AztecOO option getting function.

Get a specific Aztec optioin value. Example: problem.GetAztecOption(AZ_precond)

See the Aztec 2.1 User Guide for a complete list of these options.

const double* AztecOO::GetAztecStatus ( ) const
inline

Return the Aztec status after iterating.

Returns pointer to the underlying Aztec Status array (of length AZ_STATUS_SIZE). See the Aztec documenation.

Referenced by AztecOOConditionNumber::computeConditionNumber().

const char * AztecOO::GetLabel ( ) const

Get the label describing this AztecOO object.

Returns the string used to define this object.

int AztecOO::Iterate ( long long  MaxIters,
double  Tolerance 
)

AztecOO iteration function.

Iterates on the current problem until MaxIters or Tolerance is reached.

References GetUserMatrix(), SetAztecOption(), and SetAztecParam().

Referenced by AdaptiveIterate(), AZOO_iterate(), AztecOOConditionNumber::computeConditionNumber(), and Iterate().

int AztecOO::Iterate ( Epetra_RowMatrix A,
Epetra_MultiVector X,
Epetra_MultiVector B,
int  MaxIters,
double  Tolerance 
)

AztecOO iteration function.

Iterates on the specified matrix and vectors until MaxIters or Tolerance is reached..

References Iterate(), SetLHS(), SetRHS(), and SetUserMatrix().

void AztecOO::PrintLinearSystem ( const char *  name)

Print linear-system to files.

Parameters
namePrint the matrix to the file A_'name', and print the solution and rhs vectors to files X_'name' and B_'name', respectively. Will only produce a matrix file if the run-time-type of the matrix is either Epetra_CrsMatrix or Epetra_VbrMatrix.

References Epetra_Comm::Barrier(), Epetra_BlockMap::Comm(), GetUserMatrix(), Epetra_Comm::MyPID(), Epetra_Comm::NumProc(), Epetra_DistObject::Print(), Epetra_CrsMatrix::Print(), and Epetra_RowMatrix::RowMatrixRowMap().

int AztecOO::recursiveIterate ( int  MaxIters,
double  Tolerance 
)

AztecOO iteration functions.

Iterates on the current problem until MaxIters or Tolerance is reached.. This one should be suitable for recursive invocations of Aztec.

References SetAztecOption(), and SetAztecParam().

Referenced by AztecOO_Operator::ApplyInverse().

int AztecOO::SetAdaptiveParams ( int  NumTrials,
double *  athresholds,
double *  rthresholds,
double  condestThreshold,
double  maxFill,
int  maxKspace 
)

Set the parameter that control the AdaptiveIterate() method.

The AdaptiveIterate() method attempts to solve a given problem using multiple preconditioner and iterative method tuning parameters. There are defaults that are coded into AdaptiveIterate() method, but the defaults can be over-ridden by the use of the SetAdaptiveParams() method. Details of condition number management follow:

Parameters
NumTrialsIn The number of Athresh and Rthresh pairs that should be tried when attempting to stabilize the preconditioner.
athresholdsIn The list of absolute threshold values that should be tried when attempting to stabilize the preconditioner.
rthresholdsIn The list of relative threshold values that should be tried when attempting to stabilize the preconditioner.
condestThresholdIn If the condition number estimate of the preconditioner is above this number, no attempt will be made to try iterations. Instead a new preconditioner will be computed using the next threshold pair.
maxFillIn In addition to managing the condest, the AdaptiveIterate() method will also try to increase the preconditioner fill if it is determined that this might help. maxFill specifies the maximum fill allowed.
maxKspaceIn In addition to managing the condest, the AdaptiveIterate() method will also try to increase the Krylov subspace size if GMRES is being used and it is determined that this might help. maxKspace specifies the maximum Krylov subspace allowed.

Referenced by AdaptiveIterate().

int AztecOO::SetAllAztecOptions ( const int *  options)
inline

AztecOO option setting function.

Set all Aztec option values using an existing Aztec options array.

Referenced by AZOO_iterate(), and AztecOO().

int AztecOO::SetAllAztecParams ( const double *  params)
inline

AztecOO param setting function.

Set all Aztec parameter values using an existing Aztec params array.

Referenced by AZOO_iterate(), and AztecOO().

int AztecOO::SetAztecDefaults ( )

AztecOO function to restore default options/parameter settings.

This function is called automatically within AztecOO's constructor, but if constructed using a Epetra_LinearProblem object, some options are reset based on the ProblemDifficultyLevel associated with the Epetra_LinearProblem.

See the Aztec 2.1 User Guide for a complete list of these options.

Warning
In AztecOO, the default value of options[AZ_poly_ord] is set to 1. This is different than Aztec 2.1, but the preferred value since Jacobi preconditioning is used much more often than polynomial preconditioning and one step of Jacobi is far more effective than 3 steps.

References SetLabel().

Referenced by AztecOO().

int AztecOO::SetAztecOption ( int  option,
int  value 
)
inline

AztecOO option setting function.

Set a specific Aztec option value. Example: problem.SetAztecOption(AZ_precond, AZ_Jacobi)

See the Aztec 2.1 User Guide for a complete list of these options.

Referenced by AdaptiveIterate(), AZOO_iterate(), AztecOOConditionNumber::initialize(), Iterate(), and recursiveIterate().

int AztecOO::SetAztecParam ( int  param,
double  value 
)
inline

AztecOO param setting function.

Set a specific Aztec parameter value. Example: problem.SetAztecParam(AZ_drop, 1.0E-6)

See the Aztec 2.1 User Guide for a complete list of these parameters.

Referenced by AdaptiveIterate(), Iterate(), and recursiveIterate().

void AztecOO::SetErrorStream ( std::ostream &  errstrm)

Set std::ostream for Aztec's error output.

This sets the destination for output that Aztec would normally send to stderr.

void AztecOO::SetLabel ( const char *const  Label)

Set Label this AztecOO object.

Defines the label used to describe the this object.

Referenced by SetAztecDefaults().

int AztecOO::SetLHS ( Epetra_MultiVector X)

AztecOO LHS Set.

Associates an already defined Epetra_MultiVector (or Epetra_Vector) as the initial guess and location where the solution will be return.

Referenced by AztecOO_Operator::ApplyInverse(), AztecOO(), Iterate(), and SetProblem().

int AztecOO::SetMatrixName ( int  label)

AztecOO Label Matrix for Aztec.

This is used to label individual matrices within Aztec. This might be useful if several Aztec invocations are involved corresponding to different matrices.

void AztecOO::SetOutputStream ( std::ostream &  ostrm)

Set std::ostream for Aztec's screen output.

This sets the destination for output that Aztec would normally send to stdout.

int AztecOO::SetParameters ( Teuchos::ParameterList parameterlist,
bool  cerr_warning_if_unused = false 
)

Method to set options/parameters using a ParameterList object.

This method extracts any mixture of options and parameters from a ParameterList object and uses them to set values in AztecOO's internal options and params arrays. This method may be called repeatedly. This method does not reset default values or previously-set values unless those values are contained in the current ParameterList argument. Note that if the method SetAztecDefaults() is called after this method has been called, any parameters set by this method will be lost.

A ParameterList is a collection of named ParameterEntry objects. AztecOO recognizes names which mirror the macros defined in az_aztec_defs.h. In addition, it recognizes case insensitive versions of those names, with or without the prepended 'AZ_'. So the following are equivalent and valid: "AZ_solver", "SOLVER", "Solver". To set an entry in the Aztec options array, the type of the ParameterEntry value may be either a string or an int in some cases. E.g., if selecting the solver, the following are equivalent and valid: AZ_gmres (which is an int), "AZ_gmres" (which is a string) or "GMRES" (case-insensitive, 'AZ_' is optional).

To set an entry in the Aztec params array, the type of the ParameterEntry value must be double.

By default, this method will silently ignore parameters which have unrecognized names or invalid types. Users may set the optional argument specifying that warnings be printed for unused parameters. Alternatively, users may iterate the ParameterList afterwards and check the isUsed attribute on the ParameterEntry objects.

Parameters
parameterlistObject containing parameters to be parsed by AztecOO.
cerr_warning_if_unusedOptional argument, default value is false. If true, a warning is printed to cerr stating which parameters are not used due to having unrecognized names or values of the wrong type. Default behavior is to silently ignore unused parameters.

References Teuchos::ParameterList::begin(), and Teuchos::ParameterList::end().

int AztecOO::SetPrecMatrix ( Epetra_RowMatrix PrecMatrix)

AztecOO Preconditioner Matrix Set.

Associates an already defined Epetra_Matrix as the matrix that will be used by AztecOO when constructing native AztecOO preconditioners. By default, if AztecOO native preconditioners are used, the original operator matrix will be used as the source for deriving the preconditioner. However, there are instances where a user would like to have the preconditioner be defined using a different matrix than the original operator matrix. Another common situation is where the user may not have the operator in matrix form but has a matrix that approximates the operator and can be used as the basis for an incomplete factorization. This set method allows the user to pass any Epetra_RowMatrix to AztecOO for use in constructing an AztecOO native preconditioner, as long as the matrix implements the Epetra_RowMatrix pure virtual class, and has proper domain and range map dimensions. Epetra_CrsMatrix and Epetra_VbrMatrix objects can be passed in through this method.

References Epetra_Operator::Comm(), Epetra_Operator::HasNormInf(), Epetra_RowMatrix::NormInf(), Epetra_Import::NumRecv(), and Epetra_RowMatrix::RowMatrixImporter().

Referenced by AztecOO(), and SetUserMatrix().

int AztecOO::SetPreconditioner ( AZ_PRECOND Prec)
inline

AztecOO External Preconditioner Set (object)

Associates an already defined Aztec preconditioner with this solve.

int AztecOO::SetPreconditioner ( AZ_PREC_FUN  prec_function,
void *  prec_data 
)

AztecOO External Preconditioner Set (function and data)

Associates an external function and data pointer with preconditioner

References DestroyPreconditioner().

int AztecOO::SetPrecOperator ( Epetra_Operator PrecOperator)

AztecOO External Preconditioner Set.

Associates an already defined Epetra_Operator as the preconditioner that will be called during iterations. This set method allows the user to pass any type of preconditioner to AztecOO, as long as the preconditioner implements the Epetra_Operator pure virtual class, and has proper domain and range map dimensions. Ifpack preconditioners can be passed in through this method.

References Epetra_Operator::Comm(), and Epetra_Operator::Label().

Referenced by AztecOO().

int AztecOO::SetProblem ( const Epetra_LinearProblem prob,
bool  call_SetPrecMatrix = false 
)

AztecOO Epetra_LinearProblem Set.

Associates an already defined Epetra_LinearProblem as the problem that will be solved during iterations. This method allows the user to change which problem is being solved by an existing AztecOO object.

Internally calls SetUserMatrix() if the Epetra_LinearProblem's operator can be cast to Epetra_RowMatrix, otherwise calls SetUserOperator().

IMPORTANT WARNING *** This method calls SetUserMatrix(), which also sets the preconditioner matrix to the matrix passed in, by internally calling SetPrecMatrix(), but ONLY if SetPrecMatrix() hasn't previously been called. If the user wants to make sure that any pre-existing preconditioner is replaced, they must set the optional bool argument 'call_SetPrecMatrix' to true, which will force this function to call SetPrecMatrix().

Warning
The first time this method is called, the default options and parameters are set. Therefore, this method should be called before setting any individual options or parameters.
If a preconditioner has been pre-built and associated with this AztecOO object, the Epetra_LinearProblem being passed in to this method must have compatible domain and range maps.

References Epetra_LinearProblem::GetLHS(), Epetra_LinearProblem::GetOperator(), Epetra_LinearProblem::GetPDL(), Epetra_LinearProblem::GetRHS(), Epetra_LinearProblem::IsOperatorSymmetric(), SetLHS(), SetRHS(), SetUserMatrix(), and SetUserOperator().

Referenced by AztecOO().

int AztecOO::SetRHS ( Epetra_MultiVector B)

AztecOO RHS Set.

Associates an already defined Epetra_MultiVector (or Epetra_Vector) as the right-hand-side of the linear system.

Referenced by AztecOO_Operator::ApplyInverse(), AztecOO(), Iterate(), and SetProblem().

int AztecOO::SetScaling ( struct AZ_SCALING Scaling)
inline

AztecOO External Scaling Set.

Associates an already defined Aztec scaling object with this solve.

int AztecOO::SetStatusTest ( AztecOO_StatusTest StatusTest)

AztecOO External Convergence/Status Test Set.

Assigns an already defined AztecOO_StatusTest object as the class that will determine when iterations should stop, either because convergence was reached or the iteration failed. This method allows a large variety of convergence tests to be used with AztecOO. The AztecOO_StatusTest class is a pure virtual class, so any class that implements its interface can be passed in to this set method. A number of pre-defined AztecOO_StatusTest derived classes are already available, including AztecOO_StatusTestCombo, a class that allows logical combinations of other status test objects for sophisticated convergence testing.

References Epetra_Operator::OperatorRangeMap(), and View.

int AztecOO::SetUserMatrix ( Epetra_RowMatrix UserMatrix,
bool  call_SetPrecMatrix = false 
)

AztecOO User Matrix Set.

Associates an already defined Epetra_Matrix as the matrix that will be used by AztecOO as the linear operator when solving the linear system. Epetra_CrsMatrix and Epetra_VbrMatrix objects can be passed in through this method.

IMPORTANT WARNING ***

This method sets the preconditioner matrix to the matrix passed in here, by internally calling SetPrecMatrix(), but ONLY if SetPrecMatrix() hasn't previously been called. If the user wants to make sure that any pre-existing preconditioner is replaced, they must set the optional bool argument 'call_SetPrecMatrix' to true, which will force this function to call SetPrecMatrix().

References Epetra_Operator::Comm(), Epetra_Operator::Label(), Epetra_Import::NumRecv(), Epetra_RowMatrix::RowMatrixImporter(), SetPrecMatrix(), and SetUserOperator().

Referenced by AztecOO(), Iterate(), and SetProblem().

int AztecOO::SetUserOperator ( Epetra_Operator UserOperator)

AztecOO User Operator Set.

Associates an already defined Epetra_Operator as the linear operator for the linear system system that will be solved during iterations. This set method allows the user to pass any type of linear operator to AztecOO, as long as the operator implements the Epetra_Operator pure virtual class, and has proper domain and range map dimensions. Epetra_CrsMatrix and Epetra_VbrMatrix objects can be passed in through this method.

References Epetra_Operator::Comm(), Epetra_Operator::HasNormInf(), Epetra_Operator::Label(), Epetra_Operator::NormInf(), Epetra_BlockMap::NumMyPoints(), and Epetra_Operator::OperatorRangeMap().

Referenced by AztecOO(), SetProblem(), and SetUserMatrix().

int AztecOO::UnsetLHSRHS ( )

AztecOO unset LHS and RHS.

Sets to null the previously objects..


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