EpetraExt Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
test/Permutation_LL/cxx_main.cpp
Go to the documentation of this file.
1 //@HEADER
2 // ***********************************************************************
3 //
4 // EpetraExt: Epetra Extended - Linear Algebra Services Package
5 // Copyright (2011) Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ***********************************************************************
40 //@HEADER
41 
42 //Permutation Test routine
43 #include <Epetra_ConfigDefs.h>
44 #include "EpetraExt_Version.h"
45 
46 #ifdef EPETRA_MPI
47 #include "Epetra_MpiComm.h"
48 #include <mpi.h>
49 #endif
50 
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"
57 #include "EpetraExt_Permutation.h"
58 #include "../epetra_test_err.h"
59 
63 int check_colpermute_crsgraph(Epetra_Comm& Comm, bool verbose);
64 int check_colpermute_crsmatrix(Epetra_Comm& Comm, bool verbose);
65 int check_rowpermute_multivector_local(Epetra_Comm& Comm, bool verbose);
66 
67 int main(int argc, char *argv[]) {
68 
69  int returnierr=0;
70 
71  bool verbose = false;
72 
73 #ifdef EPETRA_MPI
74 
75  // Initialize MPI
76 
77  MPI_Init(&argc,&argv);
78  Epetra_MpiComm Comm(MPI_COMM_WORLD);
79 
80 #else
81  Epetra_SerialComm Comm;
82 #endif
83 
84  // Check if we should print results to standard out
85  if (argc>1) {
86  if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
87  }
88 
89  //Make sure the value of verbose is consistent across processors.
90  int verbose_int = verbose ? 1 : 0;
91  Comm.Broadcast(&verbose_int, 1, 0);
92  verbose = verbose_int==1 ? true : false;
93 
94  if (!verbose) {
95  Comm.SetTracebackMode(0); // This should shut down error traceback reporting
96  }
97 
98  if (verbose && Comm.MyPID()==0)
99  cout << EpetraExt::EpetraExt_Version() << endl << endl;
100 
102 
104 
106 
107  EPETRA_CHK_ERR( check_colpermute_crsgraph( Comm, verbose) );
108 
109  EPETRA_CHK_ERR( check_colpermute_crsmatrix( Comm, verbose) );
110 
112 
113 
114 #ifdef EPETRA_MPI
115  MPI_Finalize();
116 #endif
117 
118  return returnierr;
119 }
120 
121 //------------------------------------------------------------------------------
123  bool verbose)
124 {
125  int MyPID = Comm.MyPID();
126  int NumProc = Comm.NumProc();
127 
128  Comm.Barrier();
129  bool verbose1 = verbose;
130 
131  if (verbose) verbose = (MyPID==0);
132 
133  if (verbose) {
134  cerr << "================check_rowpermute_crsmatrix_local_diagonal=========="
135  <<endl;
136  }
137 
138  int NumMyElements = 5;
139  long long NumGlobalElements = ((long long)NumMyElements)*NumProc;
140 
141  Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
142 
143  long long* p = new long long[NumMyElements];
144  long long firstGlobalRow = ((long long)MyPID)*NumMyElements;
145 
146  //Set up a permutation that will reverse the order of all LOCAL rows. (i.e.,
147  //this test won't cause any inter-processor data movement.)
148 
149  if (verbose) {
150  cout << "Permutation P:"<<endl;
151  }
152 
153  int i;
154 
155  for(i=0; i<NumMyElements; ++i) {
156  p[i] = firstGlobalRow+NumMyElements-1-i;
157  if (verbose1) {
158  cout << "p["<<firstGlobalRow+i<<"]: "<<p[i]<<endl;
159  }
160  }
161 
162  Epetra_CrsMatrix A(Copy, Map, 1);
163 
164  long long col;
165  double val;
166 
167  //set up a diagonal matrix A. It's diagonal because that's the easiest
168  //to fill and to examine output before and after permutation...
169 
170  for(i=0; i<NumMyElements; ++i) {
171  long long row = firstGlobalRow+i;
172  val = 1.0*row;
173  col = row;
174 
175  A.InsertGlobalValues(row, 1, &val, &col);
176  }
177 
178  A.FillComplete();
179 
180  if (verbose1) {
181  cout << "********** matrix A: **************"<<endl;
182  cout << A << endl;
183  }
184 
186 
187  Epetra_CrsMatrix& B = P(A);
188 
189  if (verbose1) {
190  cout <<"************ permuted matrix B: ***************"<<endl;
191  cout << B << endl;
192  }
193 
194  delete [] p;
195 
196  return(0);
197 }
198 
199 //------------------------------------------------------------------------------
201  bool verbose)
202 {
203  int MyPID = Comm.MyPID();
204  int NumProc = Comm.NumProc();
205 
206  Comm.Barrier();
207  bool verbose1 = verbose;
208 
209  if (verbose) verbose = (MyPID==0);
210 
211  if (verbose) {
212  cerr << "================check_rowpermute_crsgraph_local_diagonal=========="
213  <<endl;
214  }
215 
216  int NumMyElements = 5;
217  long long NumGlobalElements = ((long long)NumMyElements)*NumProc;
218 
219  Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
220 
221  long long* p = new long long[NumMyElements];
222  long long firstGlobalRow = ((long long)MyPID)*NumMyElements;
223 
224  //Set up a permutation that will reverse the order of all LOCAL rows. (i.e.,
225  //this test won't cause any inter-processor data movement.)
226 
227  if (verbose) {
228  cout << "Permutation P:"<<endl;
229  }
230 
231  int i;
232 
233  for(i=0; i<NumMyElements; ++i) {
234  p[i] = firstGlobalRow+NumMyElements-1-i;
235  if (verbose1) {
236  cout << "p["<<firstGlobalRow+i<<"]: "<<p[i]<<endl;
237  }
238  }
239 
240  Epetra_CrsGraph Agrph(Copy, Map, 1);
241 
242  long long col;
243 
244  //set up a diagonal graph. It's diagonal because that's the easiest
245  //to fill and to examine output before and after permutation...
246 
247  for(i=0; i<NumMyElements; ++i) {
248  long long row = firstGlobalRow+i;
249  col = row;
250 
251  Agrph.InsertGlobalIndices(row, 1, &col);
252  }
253 
254  Agrph.FillComplete();
255 
256  if (verbose1) {
257  cout << "*************** graph Agrph: ********************"<<endl;
258  cout << Agrph << endl;
259  }
260 
262 
263  Epetra_CrsGraph& Bgrph = P(Agrph);
264 
265  if (verbose1) {
266  cout <<"************* permuted graph Bgrph: ****************"<<endl;
267  cout << Bgrph << endl;
268  }
269 
270  delete [] p;
271 
272  return(0);
273 }
274 
275 //------------------------------------------------------------------------------
277  bool verbose)
278 {
279  int MyPID = Comm.MyPID();
280  int NumProc = Comm.NumProc();
281 
282  Comm.Barrier();
283  bool verbose1 = verbose;
284 
285  if (verbose) verbose = (MyPID==0);
286 
287  if (verbose) {
288  cerr << "================check_colpermute_crsgraph=========="
289  <<endl;
290  }
291 
292  int NumMyElements = 5;
293  long long NumGlobalElements = ((long long)NumMyElements)*NumProc;
294 
295  Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
296 
297  long long* p = new long long[NumMyElements];
298  long long firstGlobalRow = ((long long)MyPID)*NumMyElements;
299 
300  if (verbose) {
301  cout << "Permutation P:"<<endl;
302  }
303 
304  int i;
305 
306  for(i=0; i<NumMyElements; ++i) {
307  long long row = firstGlobalRow+i;
308  p[i] = NumGlobalElements - row - 1;
309  if (verbose1) {
310  cout << "p["<<firstGlobalRow+i<<"]: "<<p[i]<<endl;
311  }
312  }
313 
314  Epetra_CrsGraph Agrph(Copy, Map, 1);
315 
316  long long col;
317 
318  //set up a tri-diagonal graph.
319 
320  for(i=0; i<NumMyElements; ++i) {
321  long long row = firstGlobalRow+i;
322  col = NumGlobalElements - row - 1;
323 
324  Agrph.InsertGlobalIndices(row, 1, &col);
325 
326  if (col > 0) {
327  long long colm1 = col-1;
328  Agrph.InsertGlobalIndices(row, 1, &colm1);
329  }
330 
331  if (col < NumGlobalElements-1) {
332  long long colp1 = col+1;
333  Agrph.InsertGlobalIndices(row, 1, &colp1);
334  }
335  }
336 
337  Agrph.FillComplete();
338 
339  if (verbose1) {
340  cout << "*************** graph Agrph: ********************"<<endl;
341  cout << Agrph << endl;
342  }
343 
345 
346  bool column_permutation = true;
347  Epetra_CrsGraph& Bgrph = P(Agrph, column_permutation);
348 
349  if (verbose1) {
350  cout <<"************* column-permuted graph Bgrph: ****************"<<endl;
351  cout << Bgrph << endl;
352  }
353 
354  delete [] p;
355 
356  return(0);
357 }
358 
359 //-------------------------------------------------------------------------------
361  bool verbose)
362 {
363  int MyPID = Comm.MyPID();
364  int NumProc = Comm.NumProc();
365 
366  Comm.Barrier();
367  bool verbose1 = verbose;
368 
369  if (verbose) verbose = (MyPID==0);
370 
371  if (verbose) {
372  cerr << "================check_rowpermute_crsmatrix_global_diagonal=========="
373  <<endl;
374  }
375 
376  int NumMyElements = 5;
377  long long NumGlobalElements = ((long long)NumMyElements)*NumProc;
378 
379  Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
380 
381  long long* p = new long long[NumMyElements];
382  long long firstGlobalRow = ((long long)MyPID)*NumMyElements;
383 
384  //Now set up a permutation that will GLOBALLY reverse the order of all rows.
385  //(i.e., if there are multiple processors, there will be inter-processor
386  //data movement as rows are migrated.)
387 
388  int i;
389 
390  Epetra_CrsMatrix A(Copy, Map, 1);
391 
392  long long col;
393  double val;
394 
395  //set up a diagonal matrix A. It's diagonal because that's the easiest
396  //to fill and to examine output before and after permutation...
397 
398  for(i=0; i<NumMyElements; ++i) {
399  long long row = firstGlobalRow+i;
400  val = 1.0*row;
401  col = row;
402 
403  A.InsertGlobalValues(row, 1, &val, &col);
404  }
405 
406  A.FillComplete();
407 
408  if (verbose1) {
409  cout << "******************* matrix A: ****************************"<<endl;
410  cout << A << endl;
411  }
412 
413  if (verbose) {
414  cout << "Permutation P:"<<endl;
415  }
416 
417  for(i=0; i<NumMyElements; ++i) {
418  long long globalrow = NumGlobalElements-(firstGlobalRow+i)-1;
419  p[i] = globalrow;
420  if (verbose1) {
421  cout << "p["<<firstGlobalRow+i<<"]: "<<p[i]<<endl;
422  }
423  }
424 
425  EpetraExt::Permutation64<Epetra_CrsMatrix> Pglobal(Copy, Map, p);
426 
427  Epetra_CrsMatrix& Bglobal = Pglobal(A);
428 
429  if (verbose1) {
430  cout << "******************* permuted matrix Bglobal: *******************" <<endl;
431  cout << Bglobal << endl;
432  }
433 
434  delete [] p;
435 
436  return(0);
437 }
438 
439 //-------------------------------------------------------------------------------
441  bool verbose)
442 {
443  int MyPID = Comm.MyPID();
444  int NumProc = Comm.NumProc();
445 
446  Comm.Barrier();
447  bool verbose1 = verbose;
448 
449  if (verbose) verbose = (MyPID==0);
450 
451  if (verbose) {
452  cerr << "================check_colpermute_crsmatrix=========="
453  <<endl;
454  }
455 
456  int NumMyElements = 5;
457  long long NumGlobalElements = ((long long)NumMyElements)*NumProc;
458 
459  Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
460 
461  long long* p = new long long[NumMyElements];
462  long long firstGlobalRow = ((long long)MyPID)*NumMyElements;
463 
464  if (verbose) {
465  cout << "Permutation P:"<<endl;
466  }
467 
468  int i;
469 
470  for(i=0; i<NumMyElements; ++i) {
471  long long row = firstGlobalRow+i;
472  p[i] = NumGlobalElements - row - 1;
473  if (verbose1) {
474  cout << "p["<<firstGlobalRow+i<<"]: "<<p[i]<<endl;
475  }
476  }
477 
478  Epetra_CrsMatrix A(Copy, Map, 1);
479 
480  long long col;
481  double val;
482 
483  //set up a tri-diagonal graph.
484 
485  for(i=0; i<NumMyElements; ++i) {
486  long long row = firstGlobalRow+i;
487  col = NumGlobalElements - row - 1;
488  val = 1.0*col;
489 
490  A.InsertGlobalValues(row, 1, &val, &col);
491 
492  if (col > 0) {
493  long long colm1 = col-1;
494  val = 1.0*colm1;
495  A.InsertGlobalValues(row, 1, &val, &colm1);
496  }
497 
498  if (col < NumGlobalElements-1) {
499  long long colp1 = col+1;
500  val = 1.0*colp1;
501  A.InsertGlobalValues(row, 1, &val, &colp1);
502  }
503  }
504 
505  A.FillComplete();
506 
507  if (verbose1) {
508  cout << "*************** matrix A: ********************"<<endl;
509  cout << A << endl;
510  }
511 
513 
514  bool column_permutation = true;
515  Epetra_CrsMatrix& B = P(A, column_permutation);
516 
517  if (verbose1) {
518  cout <<"************* column-permuted matrix B: ****************"<<endl;
519  cout << B << endl;
520  }
521 
522  delete [] p;
523 
524  return(0);
525 }
526 
527 //------------------------------------------------------------------------------
529  bool verbose)
530 {
531  int MyPID = Comm.MyPID();
532  int NumProc = Comm.NumProc();
533 
534  Comm.Barrier();
535  bool verbose1 = verbose;
536 
537  if (verbose) verbose = (MyPID==0);
538 
539  if (verbose) {
540  cerr << "================check_rowpermute_multivector_local=========="
541  <<endl;
542  }
543 
544  int NumMyElements = 5;
545  long long NumGlobalElements = ((long long)NumMyElements)*NumProc;
546 
547  Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
548 
549  long long* p = new long long[NumMyElements];
550  long long firstGlobalRow = ((long long)MyPID)*NumMyElements;
551 
552  //Set up a permutation that will reverse the order of all LOCAL rows. (i.e.,
553  //this test won't cause any inter-processor data movement.)
554 
555  if (verbose) {
556  cout << "Permutation P:"<<endl;
557  }
558 
559  int i;
560 
561  for(i=0; i<NumMyElements; ++i) {
562  p[i] = firstGlobalRow+NumMyElements-1-i;
563  if (verbose1) {
564  cout << "p["<<firstGlobalRow+i<<"]: "<<p[i]<<endl;
565  }
566  }
567 
568  Epetra_MultiVector v(Map, 3);
569 
570  double* v0 = v[0];
571  double* v1 = v[1];
572  double* v2 = v[2];
573 
574  for(i=0; i<NumMyElements; ++i) {
575  v0[i] = 1.0*(firstGlobalRow+i) + 0.1;
576  v1[i] = 1.0*(firstGlobalRow+i) + 0.2;
577  v2[i] = 1.0*(firstGlobalRow+i) + 0.3;
578  }
579 
580  if (verbose1) {
581  cout << "*************** MultiVector v: ********************"<<endl;
582  cout << v << endl;
583  }
584 
586 
587  Epetra_MultiVector& Pv = P(v);
588 
589  if (verbose1) {
590  cout <<"************* permuted MultiVector Pv: ****************"<<endl;
591  cout << Pv << endl;
592  }
593 
594  delete [] p;
595 
596  return(0);
597 }
598 
#define EPETRA_CHK_ERR(a)
std::string EpetraExt_Version()
int MyPID() const
virtual void Barrier() const =0
int check_colpermute_crsgraph(Epetra_Comm &Comm, bool verbose)
virtual int MyPID() const =0
int main(int argc, char **argv)
Definition: HDF5_IO.cpp:67
int check_rowpermute_crsmatrix_global_diagonal(Epetra_Comm &Comm, bool verbose)
int check_rowpermute_crsmatrix_local_diagonal(Epetra_Comm &Comm, bool verbose)
virtual int NumProc() const =0
int check_rowpermute_crsgraph_local_diagonal(Epetra_Comm &Comm, bool verbose)
int Broadcast(double *MyVals, int Count, int Root) const
int check_rowpermute_multivector_local(Epetra_Comm &Comm, bool verbose)
int check_colpermute_crsmatrix(Epetra_Comm &Comm, bool verbose)