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  return(0);
195 }
196 
197 //------------------------------------------------------------------------------
199  bool verbose)
200 {
201  int MyPID = Comm.MyPID();
202  int NumProc = Comm.NumProc();
203 
204  Comm.Barrier();
205  bool verbose1 = verbose;
206 
207  if (verbose) verbose = (MyPID==0);
208 
209  if (verbose) {
210  cerr << "================check_rowpermute_crsgraph_local_diagonal=========="
211  <<endl;
212  }
213 
214  int NumMyElements = 5;
215  long long NumGlobalElements = ((long long)NumMyElements)*NumProc;
216 
217  Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
218 
219  long long* p = new long long[NumMyElements];
220  long long firstGlobalRow = ((long long)MyPID)*NumMyElements;
221 
222  //Set up a permutation that will reverse the order of all LOCAL rows. (i.e.,
223  //this test won't cause any inter-processor data movement.)
224 
225  if (verbose) {
226  cout << "Permutation P:"<<endl;
227  }
228 
229  int i;
230 
231  for(i=0; i<NumMyElements; ++i) {
232  p[i] = firstGlobalRow+NumMyElements-1-i;
233  if (verbose1) {
234  cout << "p["<<firstGlobalRow+i<<"]: "<<p[i]<<endl;
235  }
236  }
237 
238  Epetra_CrsGraph Agrph(Copy, Map, 1);
239 
240  long long col;
241 
242  //set up a diagonal graph. It's diagonal because that's the easiest
243  //to fill and to examine output before and after permutation...
244 
245  for(i=0; i<NumMyElements; ++i) {
246  long long row = firstGlobalRow+i;
247  col = row;
248 
249  Agrph.InsertGlobalIndices(row, 1, &col);
250  }
251 
252  Agrph.FillComplete();
253 
254  if (verbose1) {
255  cout << "*************** graph Agrph: ********************"<<endl;
256  cout << Agrph << endl;
257  }
258 
260 
261  Epetra_CrsGraph& Bgrph = P(Agrph);
262 
263  if (verbose1) {
264  cout <<"************* permuted graph Bgrph: ****************"<<endl;
265  cout << Bgrph << endl;
266  }
267 
268  return(0);
269 }
270 
271 //------------------------------------------------------------------------------
273  bool verbose)
274 {
275  int MyPID = Comm.MyPID();
276  int NumProc = Comm.NumProc();
277 
278  Comm.Barrier();
279  bool verbose1 = verbose;
280 
281  if (verbose) verbose = (MyPID==0);
282 
283  if (verbose) {
284  cerr << "================check_colpermute_crsgraph=========="
285  <<endl;
286  }
287 
288  int NumMyElements = 5;
289  long long NumGlobalElements = ((long long)NumMyElements)*NumProc;
290 
291  Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
292 
293  long long* p = new long long[NumMyElements];
294  long long firstGlobalRow = ((long long)MyPID)*NumMyElements;
295 
296  if (verbose) {
297  cout << "Permutation P:"<<endl;
298  }
299 
300  int i;
301 
302  for(i=0; i<NumMyElements; ++i) {
303  long long row = firstGlobalRow+i;
304  p[i] = NumGlobalElements - row - 1;
305  if (verbose1) {
306  cout << "p["<<firstGlobalRow+i<<"]: "<<p[i]<<endl;
307  }
308  }
309 
310  Epetra_CrsGraph Agrph(Copy, Map, 1);
311 
312  long long col;
313 
314  //set up a tri-diagonal graph.
315 
316  for(i=0; i<NumMyElements; ++i) {
317  long long row = firstGlobalRow+i;
318  col = NumGlobalElements - row - 1;
319 
320  Agrph.InsertGlobalIndices(row, 1, &col);
321 
322  if (col > 0) {
323  long long colm1 = col-1;
324  Agrph.InsertGlobalIndices(row, 1, &colm1);
325  }
326 
327  if (col < NumGlobalElements-1) {
328  long long colp1 = col+1;
329  Agrph.InsertGlobalIndices(row, 1, &colp1);
330  }
331  }
332 
333  Agrph.FillComplete();
334 
335  if (verbose1) {
336  cout << "*************** graph Agrph: ********************"<<endl;
337  cout << Agrph << endl;
338  }
339 
341 
342  bool column_permutation = true;
343  Epetra_CrsGraph& Bgrph = P(Agrph, column_permutation);
344 
345  if (verbose1) {
346  cout <<"************* column-permuted graph Bgrph: ****************"<<endl;
347  cout << Bgrph << endl;
348  }
349 
350  delete [] p;
351 
352  return(0);
353 }
354 
355 //-------------------------------------------------------------------------------
357  bool verbose)
358 {
359  int MyPID = Comm.MyPID();
360  int NumProc = Comm.NumProc();
361 
362  Comm.Barrier();
363  bool verbose1 = verbose;
364 
365  if (verbose) verbose = (MyPID==0);
366 
367  if (verbose) {
368  cerr << "================check_rowpermute_crsmatrix_global_diagonal=========="
369  <<endl;
370  }
371 
372  int NumMyElements = 5;
373  long long NumGlobalElements = ((long long)NumMyElements)*NumProc;
374 
375  Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
376 
377  long long* p = new long long[NumMyElements];
378  long long firstGlobalRow = ((long long)MyPID)*NumMyElements;
379 
380  //Now set up a permutation that will GLOBALLY reverse the order of all rows.
381  //(i.e., if there are multiple processors, there will be inter-processor
382  //data movement as rows are migrated.)
383 
384  int i;
385 
386  Epetra_CrsMatrix A(Copy, Map, 1);
387 
388  long long col;
389  double val;
390 
391  //set up a diagonal matrix A. It's diagonal because that's the easiest
392  //to fill and to examine output before and after permutation...
393 
394  for(i=0; i<NumMyElements; ++i) {
395  long long row = firstGlobalRow+i;
396  val = 1.0*row;
397  col = row;
398 
399  A.InsertGlobalValues(row, 1, &val, &col);
400  }
401 
402  A.FillComplete();
403 
404  if (verbose1) {
405  cout << "******************* matrix A: ****************************"<<endl;
406  cout << A << endl;
407  }
408 
409  if (verbose) {
410  cout << "Permutation P:"<<endl;
411  }
412 
413  for(i=0; i<NumMyElements; ++i) {
414  long long globalrow = NumGlobalElements-(firstGlobalRow+i)-1;
415  p[i] = globalrow;
416  if (verbose1) {
417  cout << "p["<<firstGlobalRow+i<<"]: "<<p[i]<<endl;
418  }
419  }
420 
421  EpetraExt::Permutation64<Epetra_CrsMatrix> Pglobal(Copy, Map, p);
422 
423  Epetra_CrsMatrix& Bglobal = Pglobal(A);
424 
425  if (verbose1) {
426  cout << "******************* permuted matrix Bglobal: *******************" <<endl;
427  cout << Bglobal << endl;
428  }
429 
430  return(0);
431 }
432 
433 //-------------------------------------------------------------------------------
435  bool verbose)
436 {
437  int MyPID = Comm.MyPID();
438  int NumProc = Comm.NumProc();
439 
440  Comm.Barrier();
441  bool verbose1 = verbose;
442 
443  if (verbose) verbose = (MyPID==0);
444 
445  if (verbose) {
446  cerr << "================check_colpermute_crsmatrix=========="
447  <<endl;
448  }
449 
450  int NumMyElements = 5;
451  long long NumGlobalElements = ((long long)NumMyElements)*NumProc;
452 
453  Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
454 
455  long long* p = new long long[NumMyElements];
456  long long firstGlobalRow = ((long long)MyPID)*NumMyElements;
457 
458  if (verbose) {
459  cout << "Permutation P:"<<endl;
460  }
461 
462  int i;
463 
464  for(i=0; i<NumMyElements; ++i) {
465  long long row = firstGlobalRow+i;
466  p[i] = NumGlobalElements - row - 1;
467  if (verbose1) {
468  cout << "p["<<firstGlobalRow+i<<"]: "<<p[i]<<endl;
469  }
470  }
471 
472  Epetra_CrsMatrix A(Copy, Map, 1);
473 
474  long long col;
475  double val;
476 
477  //set up a tri-diagonal graph.
478 
479  for(i=0; i<NumMyElements; ++i) {
480  long long row = firstGlobalRow+i;
481  col = NumGlobalElements - row - 1;
482  val = 1.0*col;
483 
484  A.InsertGlobalValues(row, 1, &val, &col);
485 
486  if (col > 0) {
487  long long colm1 = col-1;
488  val = 1.0*colm1;
489  A.InsertGlobalValues(row, 1, &val, &colm1);
490  }
491 
492  if (col < NumGlobalElements-1) {
493  long long colp1 = col+1;
494  val = 1.0*colp1;
495  A.InsertGlobalValues(row, 1, &val, &colp1);
496  }
497  }
498 
499  A.FillComplete();
500 
501  if (verbose1) {
502  cout << "*************** matrix A: ********************"<<endl;
503  cout << A << endl;
504  }
505 
507 
508  bool column_permutation = true;
509  Epetra_CrsMatrix& B = P(A, column_permutation);
510 
511  if (verbose1) {
512  cout <<"************* column-permuted matrix B: ****************"<<endl;
513  cout << B << endl;
514  }
515 
516  delete [] p;
517 
518  return(0);
519 }
520 
521 //------------------------------------------------------------------------------
523  bool verbose)
524 {
525  int MyPID = Comm.MyPID();
526  int NumProc = Comm.NumProc();
527 
528  Comm.Barrier();
529  bool verbose1 = verbose;
530 
531  if (verbose) verbose = (MyPID==0);
532 
533  if (verbose) {
534  cerr << "================check_rowpermute_multivector_local=========="
535  <<endl;
536  }
537 
538  int NumMyElements = 5;
539  long long NumGlobalElements = ((long long)NumMyElements)*NumProc;
540 
541  Epetra_Map Map(NumGlobalElements, NumMyElements, 0, Comm);
542 
543  long long* p = new long long[NumMyElements];
544  long long firstGlobalRow = ((long long)MyPID)*NumMyElements;
545 
546  //Set up a permutation that will reverse the order of all LOCAL rows. (i.e.,
547  //this test won't cause any inter-processor data movement.)
548 
549  if (verbose) {
550  cout << "Permutation P:"<<endl;
551  }
552 
553  int i;
554 
555  for(i=0; i<NumMyElements; ++i) {
556  p[i] = firstGlobalRow+NumMyElements-1-i;
557  if (verbose1) {
558  cout << "p["<<firstGlobalRow+i<<"]: "<<p[i]<<endl;
559  }
560  }
561 
562  Epetra_MultiVector v(Map, 3);
563 
564  double* v0 = v[0];
565  double* v1 = v[1];
566  double* v2 = v[2];
567 
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;
572  }
573 
574  if (verbose1) {
575  cout << "*************** MultiVector v: ********************"<<endl;
576  cout << v << endl;
577  }
578 
580 
581  Epetra_MultiVector& Pv = P(v);
582 
583  if (verbose1) {
584  cout <<"************* permuted MultiVector Pv: ****************"<<endl;
585  cout << Pv << endl;
586  }
587 
588  return(0);
589 }
590 
#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)