55 #include "Epetra_Map.h" 57 #include "Epetra_CrsGraph.h" 59 #include "Teuchos_RefCountPtr.hpp" 61 #include "Epetra_MpiComm.h" 63 #include "Epetra_SerialComm.h" 69 int NumGlobalRows1,
int NumMyRows1,
int LevelFill1,
bool verbose);
71 int main(
int argc,
char *argv[])
83 MPI_Init(&argc,&argv);
102 if (argc>1)
if (argv[1][0]==
'-' && argv[1][1]==
'v') {
112 int MyPID = Comm.
MyPID();
118 if (
verbose) cout << Comm <<endl;
125 if (verbose1 && argc != 4) {
128 cout <<
"Setting nx = " << nx <<
", ny = " << ny << endl;
130 else if (!verbose1 && argc != 3) {
135 nx = atoi(argv[nextarg++]);
136 if (nx<3) {cout <<
"nx = " << nx <<
": Must be greater than 2 for meaningful graph." << endl; exit(1);}
137 ny = atoi(argv[nextarg++]);
138 if (ny<3) {cout <<
"ny = " << ny <<
": Must be greater than 2 for meaningful graph." << endl; exit(1);}
141 int NumGlobalPoints = nx*ny;
145 cout <<
"\n\n*****Building 5 point matrix, Level 1 and 2 filled matrices for" << endl
146 <<
" nx = " << nx <<
", ny = " << ny << endl<< endl;
151 Epetra_Map Map(NumGlobalPoints, IndexBase, Comm);
165 std::vector<int> Indices(4);
167 for (j=0; j<ny; j++) {
168 for (i=0; i<nx; i++) {
170 if (Map.
MyGID(Row)) {
177 if (i>0) Indices[k++] = i-1 + j *nx;
178 if (j>0) Indices[k++] = i +(j-1)*nx;
185 if ((i<nx-1) && (j>0 )) Indices[k++] = i+1 +(j-1)*nx;
193 if ((i<nx-2) && (j>0 )) Indices[k++] = i+2 +(j-1)*nx;
207 if (i<nx-1) Indices[k++] = i+1 + j *nx;
208 if (j<ny-1) Indices[k++] = i +(j+1)*nx;
216 if ((i>0 ) && (j<ny-1)) Indices[k++] = i-1 +(j+1)*nx;
223 if ((i>1 ) && (j<ny-1)) Indices[k++] = i-2 +(j+1)*nx;
232 if (
verbose) cout <<
"\n\nCompleting A" << endl<< endl;
234 if (
verbose) cout <<
"\n\nCompleting L0" << endl<< endl;
236 if (
verbose) cout <<
"\n\nCompleting U0" << endl<< endl;
238 if (
verbose) cout <<
"\n\nCompleting L1" << endl<< endl;
240 if (
verbose) cout <<
"\n\nCompleting U1" << endl<< endl;
242 if (
verbose) cout <<
"\n\nCompleting L2" << endl<< endl;
244 if (
verbose) cout <<
"\n\nCompleting U2" << endl<< endl;
247 if (
verbose) cout <<
"\n\n*****Testing ILU(0) constructor on A" << endl<< endl;
252 assert(
check(L0, U0, ILU0, NumGlobalPoints, NumMyPoints, 0,
verbose)==0);
254 if (
verbose) cout <<
"\n\n*****Testing ILU(1) constructor on A" << endl<< endl;
259 assert(
check(L1, U1, ILU1, NumGlobalPoints, NumMyPoints, 1,
verbose)==0);
261 if (
verbose) cout <<
"\n\n*****Testing ILU(2) constructor on A" << endl<< endl;
266 assert(
check(L2, U2, ILU2, NumGlobalPoints, NumMyPoints, 2,
verbose)==0);
268 if (
verbose) cout <<
"\n\n*****Testing copy constructor" << endl<< endl;
272 assert(
check(L2, U2, ILUC, NumGlobalPoints, NumMyPoints, 2,
verbose)==0);
274 if (
verbose) cout <<
"\n\n*****Testing copy constructor" << endl<< endl;
276 Teuchos::RefCountPtr<Ifpack_IlukGraph> OverlapGraph;
277 for (
int overlap = 1; overlap < 4; overlap++) {
278 if (
verbose) cout <<
"\n\n*********************************************" << endl;
279 if (
verbose) cout <<
"\n\nConstruct Level 1 fill with Overlap = " << overlap <<
".\n\n" << endl;
282 assert(OverlapGraph->ConstructFilledGraph()==0);
285 cout <<
"Number of Global Rows = " << OverlapGraph->NumGlobalRows() << endl;
286 cout <<
"Number of Global Nonzeros = " << OverlapGraph->NumGlobalNonzeros() << endl;
287 cout <<
"Number of Local Rows = " << OverlapGraph->NumMyRows() << endl;
288 cout <<
"Number of Local Nonzeros = " << OverlapGraph->NumMyNonzeros() << endl;
296 int NumElements1 = 6;
297 int NumPoints1 = NumElements1;
302 std::vector<int> NumNz1(NumPoints1);
307 for (i=0; i<NumPoints1; i++)
308 if (i==0 || i == NumPoints1-1)
315 Epetra_Map Map1(NumPoints1, NumPoints1, 1, Comm);
323 std::vector<int> Indices1(2);
326 for (i=0; i<NumPoints1; i++)
333 else if (i == NumPoints1-1)
335 Indices1[0] = NumPoints1-1;
352 if (
verbose) cout <<
"\n\nPrint out tridiagonal matrix with IndexBase = 1.\n\n" << endl;
355 if (
verbose) cout <<
"\n\nConstruct Level 1 fill with IndexBase = 1.\n\n" << endl;
360 if (
verbose) cout <<
"\n\nPrint out Level 1 ILU graph of tridiagonal matrix with IndexBase = 1.\n\n" << endl;
361 if (verbose1) cout << ILU11 << endl;
376 int NumGlobalRows1,
int NumMyRows1,
int LevelFill1,
bool verbose) {
381 int NumIndices, * Indices;
382 int NumIndices1, * Indices1;
397 assert(NumIndices==NumIndices1);
398 for (j=0; j<NumIndices1; j++) {
399 if (debug &&(Indices[j]!=Indices1[j])) {
401 cout <<
"Proc " << MyPID
402 <<
" Local Row = " << i
403 <<
" L.Indices["<< j <<
"] = " << Indices[j]
404 <<
" L1.Indices["<< j <<
"] = " << Indices1[j] << endl;
406 assert(Indices[j]==Indices1[j]);
408 Nout += (NumIndices-NumIndices1);
412 assert(NumIndices==NumIndices1);
413 for (j=0; j<NumIndices1; j++) {
414 if (debug &&(Indices[j]!=Indices1[j])) {
416 cout <<
"Proc " << MyPID
417 <<
" Local Row = " << i
418 <<
" U.Indices["<< j <<
"] = " << Indices[j]
419 <<
" U1.Indices["<< j <<
"] = " << Indices1[j] << endl;
421 assert(Indices[j]==Indices1[j]);
423 Nout += (NumIndices-NumIndices1);
429 if (
verbose) cout <<
"\n\nNumber of Global Rows = " << NumGlobalRows << endl<< endl;
431 assert(NumGlobalRows==NumGlobalRows1);
434 if (
verbose) cout <<
"\n\nNumber of Global Nonzero entries = " 435 << NumGlobalNonzeros << endl<< endl;
444 if (
verbose) cout <<
"\n\nNumber of Rows = " << NumMyRows << endl<< endl;
446 assert(NumMyRows==NumMyRows1);
449 if (
verbose) cout <<
"\n\nNumber of Nonzero entries = " << NumMyNonzeros << endl<< endl;
453 if (
verbose) cout <<
"\n\nLU check OK" << endl<< endl;
int NumGlobalRows() const
Returns the number of global matrix rows.
int NumMyRows() const
Returns the number of local matrix rows.
int NumGlobalNonzeros() const
Returns the number of nonzero entries in the global graph.
virtual Epetra_CrsGraph & L_Graph()
Returns the graph of lower triangle of the ILU(k) graph as a Epetra_CrsGraph.
int NumMyNonzeros() const
Ifpack_IlukGraph: A class for constructing level filled graphs for use with ILU(k) class precondition...
const Epetra_BlockMap & RowMap() const
int InsertGlobalIndices(int_type GlobalRow, int NumIndices, int_type *Indices)
virtual int MyPID() const=0
int main(int argc, char *argv[])
virtual int SumAll(double *PartialSums, double *GlobalSums, int Count) const=0
virtual Epetra_CrsGraph & U_Graph()
Returns the graph of upper triangle of the ILU(k) graph as a Epetra_CrsGraph.
bool MyGID(int GID_in) const
int check(Epetra_CrsGraph &L, Epetra_CrsGraph &U, Ifpack_IlukGraph &LU, int NumGlobalRows1, int NumMyRows1, int LevelFill1, bool verbose)
std::string Ifpack_Version()
const Epetra_Comm & Comm() const
int NumGlobalNonzeros() const
int ExtractMyRowView(int LocalRow, int &NumIndices, int *&Indices) const
int NumMyNonzeros() const
Returns the number of nonzero entries in the local graph.
virtual int ConstructFilledGraph()
Does the actual construction of the graph.