43 #include "Ifpack_ConfigDefs.h" 44 #include "Ifpack_Condest.h" 45 #include "Ifpack_CondestType.h" 46 #include "Ifpack_Preconditioner.h" 47 #include "Epetra_Vector.h" 48 #include "Epetra_LinearProblem.h" 49 #include "Epetra_Map.h" 50 #include "Epetra_RowMatrix.h" 51 #ifdef HAVE_IFPACK_AZTECOO 56 const Ifpack_CondestType CT,
59 Epetra_RowMatrix* Matrix)
61 double ConditionNumberEstimate = -1.0;
63 if (CT == Ifpack_Cheap) {
66 Epetra_Vector Ones(IFP.OperatorDomainMap());
69 Epetra_Vector OnesResult(IFP.OperatorRangeMap());
73 IFPACK_CHK_ERR(OnesResult.Abs(OnesResult));
75 IFPACK_CHK_ERR(OnesResult.MaxValue(&ConditionNumberEstimate));
78 else if (CT == Ifpack_CG) {
80 #ifdef HAVE_IFPACK_AZTECOO 82 Matrix = (Epetra_RowMatrix*)&(IFP.
Matrix());
84 Epetra_Vector LHS(IFP.OperatorDomainMap());
86 Epetra_Vector RHS(IFP.OperatorRangeMap());
88 Epetra_LinearProblem Problem;
89 Problem.SetOperator(Matrix);
93 AztecOO Solver(Problem);
94 Solver.SetAztecOption(AZ_output,AZ_none);
95 Solver.SetAztecOption(AZ_solver,AZ_cg_condnum);
96 Solver.Iterate(MaxIters,Tol);
98 const double* status = Solver.GetAztecStatus();
99 ConditionNumberEstimate = status[AZ_condnum];
102 }
else if (CT == Ifpack_GMRES) {
104 #ifdef HAVE_IFPACK_AZTECOO 106 Matrix = (Epetra_RowMatrix*)&(IFP.
Matrix());
108 Epetra_Vector LHS(IFP.OperatorDomainMap());
110 Epetra_Vector RHS(IFP.OperatorRangeMap());
112 Epetra_LinearProblem Problem;
113 Problem.SetOperator(Matrix);
114 Problem.SetLHS(&LHS);
115 Problem.SetRHS(&RHS);
117 AztecOO Solver(Problem);
118 Solver.SetAztecOption(AZ_solver,AZ_gmres_condnum);
119 Solver.SetAztecOption(AZ_output,AZ_none);
123 Solver.SetAztecOption(AZ_kspace,MaxIters);
124 Solver.Iterate(MaxIters,Tol);
126 const double* status = Solver.GetAztecStatus();
127 ConditionNumberEstimate = status[AZ_condnum];
131 return(ConditionNumberEstimate);
virtual const Epetra_RowMatrix & Matrix() const =0
Returns a pointer to the matrix to be preconditioned.
Ifpack_Preconditioner: basic class for preconditioning in Ifpack.
virtual int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const =0
Applies the preconditioner to vector X, returns the result in Y.