1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52
|
#include "AztecOO.h"
#include "AztecOO_Version.h"
#include "Epetra_SerialComm.h"
#include "Epetra_Map.h"
#include "Epetra_Vector.h"
#include "Epetra_CrsMatrix.h"
int main(int argc, char* argv[]) {
Epetra_SerialComm Comm;
int NumMyElements = 1000;
if (Comm.MyPID() == 0)
std::cout << AztecOO_Version() << std::endl;
std::cout << Comm << std::endl;
Epetra_Map Map(-1, NumMyElements, 0, Comm);
int NumGlobalElements = Map.NumGlobalElements();
Epetra_CrsMatrix A(Copy, Map, 3);
double negOne = -1.0;
double posTwo = 2.0;
for (int i = 0; i < NumMyElements; i++) {
int GlobalRow = A.GRID(i); int RowLess1 = GlobalRow - 1; int RowPlus1 = GlobalRow + 1;
if (RowLess1 != -1) A.InsertGlobalValues(GlobalRow, 1, &negOne, &RowLess1);
if (RowPlus1 != NumGlobalElements) A.InsertGlobalValues(GlobalRow, 1, &negOne, &RowPlus1);
A.InsertGlobalValues(GlobalRow, 1, &posTwo, &GlobalRow);
};
A.FillComplete();
Epetra_Vector x(Map);
Epetra_Vector b(Map);
b.Random();
Epetra_LinearProblem problem(&A, &x, &b);
AztecOO solver(problem);
solver.SetAztecOption(AZ_precond, AZ_Jacobi);
solver.Iterate(10, 1.0E-2);
std::cout << "Solver performed " << solver.NumIters() << " iterations." << std::endl;
std::cout << "Norm of true residual = " << solver.TrueResidual() << std::endl;
for (int i = 0; i < 10; i++) {
std::cout << x[i] << std::endl;
}
return 0;
}
|