43 #include <Epetra_ConfigDefs.h> 47 #include "Epetra_MpiComm.h" 51 #include "Epetra_SerialComm.h" 52 #include "Epetra_Time.h" 53 #include "Epetra_Map.h" 54 #include "Epetra_CrsGraph.h" 55 #include "Epetra_CrsMatrix.h" 56 #include "Epetra_Vector.h" 58 #include "../epetra_test_err.h" 67 int main(
int argc,
char *argv[]) {
77 MPI_Init(&argc,&argv);
78 Epetra_MpiComm Comm(MPI_COMM_WORLD);
81 Epetra_SerialComm Comm;
86 if (argv[1][0]==
'-' && argv[1][1]==
'v') verbose =
true;
90 int verbose_int = verbose ? 1 : 0;
91 Comm.Broadcast(&verbose_int, 1, 0);
92 verbose = verbose_int==1 ? true :
false;
95 Comm.SetTracebackMode(0);
98 if (verbose && Comm.MyPID()==0)
125 int MyPID = Comm.MyPID();
126 int NumProc = Comm.NumProc();
129 bool verbose1 = verbose;
131 if (verbose) verbose = (MyPID==0);
134 cerr <<
"================check_rowpermute_crsmatrix_local_diagonal==========" 138 int NumMyElements = 5;
139 long long NumGlobalElements = ((
long long)NumMyElements)*NumProc;
141 Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
143 long long* p =
new long long[NumMyElements];
144 long long firstGlobalRow = ((
long long)MyPID)*NumMyElements;
150 cout <<
"Permutation P:"<<endl;
155 for(i=0; i<NumMyElements; ++i) {
156 p[i] = firstGlobalRow+NumMyElements-1-i;
158 cout <<
"p["<<firstGlobalRow+i<<
"]: "<<p[i]<<endl;
162 Epetra_CrsMatrix
A(
Copy, Map, 1);
170 for(i=0; i<NumMyElements; ++i) {
171 long long row = firstGlobalRow+i;
175 A.InsertGlobalValues(row, 1, &val, &col);
181 cout <<
"********** matrix A: **************"<<endl;
187 Epetra_CrsMatrix&
B = P(
A);
190 cout <<
"************ permuted matrix B: ***************"<<endl;
201 int MyPID = Comm.MyPID();
202 int NumProc = Comm.NumProc();
205 bool verbose1 = verbose;
207 if (verbose) verbose = (MyPID==0);
210 cerr <<
"================check_rowpermute_crsgraph_local_diagonal==========" 214 int NumMyElements = 5;
215 long long NumGlobalElements = ((
long long)NumMyElements)*NumProc;
217 Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
219 long long* p =
new long long[NumMyElements];
220 long long firstGlobalRow = ((
long long)MyPID)*NumMyElements;
226 cout <<
"Permutation P:"<<endl;
231 for(i=0; i<NumMyElements; ++i) {
232 p[i] = firstGlobalRow+NumMyElements-1-i;
234 cout <<
"p["<<firstGlobalRow+i<<
"]: "<<p[i]<<endl;
238 Epetra_CrsGraph Agrph(
Copy, Map, 1);
245 for(i=0; i<NumMyElements; ++i) {
246 long long row = firstGlobalRow+i;
249 Agrph.InsertGlobalIndices(row, 1, &col);
252 Agrph.FillComplete();
255 cout <<
"*************** graph Agrph: ********************"<<endl;
256 cout << Agrph << endl;
261 Epetra_CrsGraph& Bgrph = P(Agrph);
264 cout <<
"************* permuted graph Bgrph: ****************"<<endl;
265 cout << Bgrph << endl;
275 int MyPID = Comm.MyPID();
276 int NumProc = Comm.NumProc();
279 bool verbose1 = verbose;
281 if (verbose) verbose = (MyPID==0);
284 cerr <<
"================check_colpermute_crsgraph==========" 288 int NumMyElements = 5;
289 long long NumGlobalElements = ((
long long)NumMyElements)*NumProc;
291 Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
293 long long* p =
new long long[NumMyElements];
294 long long firstGlobalRow = ((
long long)MyPID)*NumMyElements;
297 cout <<
"Permutation P:"<<endl;
302 for(i=0; i<NumMyElements; ++i) {
303 long long row = firstGlobalRow+i;
304 p[i] = NumGlobalElements - row - 1;
306 cout <<
"p["<<firstGlobalRow+i<<
"]: "<<p[i]<<endl;
310 Epetra_CrsGraph Agrph(
Copy, Map, 1);
316 for(i=0; i<NumMyElements; ++i) {
317 long long row = firstGlobalRow+i;
318 col = NumGlobalElements - row - 1;
320 Agrph.InsertGlobalIndices(row, 1, &col);
323 long long colm1 = col-1;
324 Agrph.InsertGlobalIndices(row, 1, &colm1);
327 if (col < NumGlobalElements-1) {
328 long long colp1 = col+1;
329 Agrph.InsertGlobalIndices(row, 1, &colp1);
333 Agrph.FillComplete();
336 cout <<
"*************** graph Agrph: ********************"<<endl;
337 cout << Agrph << endl;
342 bool column_permutation =
true;
343 Epetra_CrsGraph& Bgrph = P(Agrph, column_permutation);
346 cout <<
"************* column-permuted graph Bgrph: ****************"<<endl;
347 cout << Bgrph << endl;
359 int MyPID = Comm.MyPID();
360 int NumProc = Comm.NumProc();
363 bool verbose1 = verbose;
365 if (verbose) verbose = (MyPID==0);
368 cerr <<
"================check_rowpermute_crsmatrix_global_diagonal==========" 372 int NumMyElements = 5;
373 long long NumGlobalElements = ((
long long)NumMyElements)*NumProc;
375 Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
377 long long* p =
new long long[NumMyElements];
378 long long firstGlobalRow = ((
long long)MyPID)*NumMyElements;
386 Epetra_CrsMatrix
A(
Copy, Map, 1);
394 for(i=0; i<NumMyElements; ++i) {
395 long long row = firstGlobalRow+i;
399 A.InsertGlobalValues(row, 1, &val, &col);
405 cout <<
"******************* matrix A: ****************************"<<endl;
410 cout <<
"Permutation P:"<<endl;
413 for(i=0; i<NumMyElements; ++i) {
414 long long globalrow = NumGlobalElements-(firstGlobalRow+i)-1;
417 cout <<
"p["<<firstGlobalRow+i<<
"]: "<<p[i]<<endl;
423 Epetra_CrsMatrix& Bglobal = Pglobal(
A);
426 cout <<
"******************* permuted matrix Bglobal: *******************" <<endl;
427 cout << Bglobal << endl;
437 int MyPID = Comm.MyPID();
438 int NumProc = Comm.NumProc();
441 bool verbose1 = verbose;
443 if (verbose) verbose = (MyPID==0);
446 cerr <<
"================check_colpermute_crsmatrix==========" 450 int NumMyElements = 5;
451 long long NumGlobalElements = ((
long long)NumMyElements)*NumProc;
453 Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
455 long long* p =
new long long[NumMyElements];
456 long long firstGlobalRow = ((
long long)MyPID)*NumMyElements;
459 cout <<
"Permutation P:"<<endl;
464 for(i=0; i<NumMyElements; ++i) {
465 long long row = firstGlobalRow+i;
466 p[i] = NumGlobalElements - row - 1;
468 cout <<
"p["<<firstGlobalRow+i<<
"]: "<<p[i]<<endl;
472 Epetra_CrsMatrix
A(
Copy, Map, 1);
479 for(i=0; i<NumMyElements; ++i) {
480 long long row = firstGlobalRow+i;
481 col = NumGlobalElements - row - 1;
484 A.InsertGlobalValues(row, 1, &val, &col);
487 long long colm1 = col-1;
489 A.InsertGlobalValues(row, 1, &val, &colm1);
492 if (col < NumGlobalElements-1) {
493 long long colp1 = col+1;
495 A.InsertGlobalValues(row, 1, &val, &colp1);
502 cout <<
"*************** matrix A: ********************"<<endl;
508 bool column_permutation =
true;
509 Epetra_CrsMatrix&
B = P(
A, column_permutation);
512 cout <<
"************* column-permuted matrix B: ****************"<<endl;
525 int MyPID = Comm.MyPID();
526 int NumProc = Comm.NumProc();
529 bool verbose1 = verbose;
531 if (verbose) verbose = (MyPID==0);
534 cerr <<
"================check_rowpermute_multivector_local==========" 538 int NumMyElements = 5;
539 long long NumGlobalElements = ((
long long)NumMyElements)*NumProc;
541 Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
543 long long* p =
new long long[NumMyElements];
544 long long firstGlobalRow = ((
long long)MyPID)*NumMyElements;
550 cout <<
"Permutation P:"<<endl;
555 for(i=0; i<NumMyElements; ++i) {
556 p[i] = firstGlobalRow+NumMyElements-1-i;
558 cout <<
"p["<<firstGlobalRow+i<<
"]: "<<p[i]<<endl;
562 Epetra_MultiVector v(Map, 3);
568 for(i=0; i<NumMyElements; ++i) {
569 v0[i] = 1.0*(firstGlobalRow+i) + 0.1;
570 v1[i] = 1.0*(firstGlobalRow+i) + 0.2;
571 v2[i] = 1.0*(firstGlobalRow+i) + 0.3;
575 cout <<
"*************** MultiVector v: ********************"<<endl;
581 Epetra_MultiVector& Pv = P(v);
584 cout <<
"************* permuted MultiVector Pv: ****************"<<endl;
int check_rowpermute_multivector_local(Epetra_Comm &Comm, bool verbose)
int check_rowpermute_crsmatrix_local_diagonal(Epetra_Comm &Comm, bool verbose)
std::string EpetraExt_Version()
int check_rowpermute_crsgraph_local_diagonal(Epetra_Comm &Comm, bool verbose)
int check_colpermute_crsgraph(Epetra_Comm &Comm, bool verbose)
int main(int argc, char *argv[])
int check_colpermute_crsmatrix(Epetra_Comm &Comm, bool verbose)
int check_rowpermute_crsmatrix_global_diagonal(Epetra_Comm &Comm, bool verbose)