Epetra Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
test/Vector/cxx_main.cpp
Go to the documentation of this file.
1 //@HEADER
2 // ************************************************************************
3 //
4 // Epetra: 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 
43 // Epetra_Vector Test routine
44 
45 #include "Epetra_Time.h"
46 #include "Epetra_BlockMap.h"
47 #include "Epetra_Vector.h"
48 #include "BuildTestProblems.h"
49 #include "ExecuteTestProblems.h"
50 #include "../epetra_test_err.h"
51 #include "Epetra_Version.h"
52 
53 #ifdef EPETRA_MPI
54 # include "Epetra_MpiComm.h"
55 # include <mpi.h>
56 #else
57 # include "Epetra_SerialComm.h"
58 #endif
59 
60 #ifdef HAVE_EPETRA_TEUCHOS
61 # include "Teuchos_VerboseObject.hpp"
62 #endif
63 
64 
65 int main(int argc, char *argv[]) {
66 
67  int ierr = 0, i;
68 
69 #ifdef EPETRA_MPI
70 
71  // Initialize MPI
72 
73  MPI_Init(&argc,&argv);
74  int rank; // My process ID
75 
76  MPI_Comm_rank(MPI_COMM_WORLD, &rank);
77  Epetra_MpiComm Comm(MPI_COMM_WORLD);
78 
79 #else
80 
81  int rank = 0;
82  Epetra_SerialComm Comm;
83 
84 #endif
85 
86 #ifdef HAVE_EPETRA_TEUCHOS
88  fancyOut = Teuchos::VerboseObjectBase::getDefaultOStream();
89  if (Comm.NumProc() > 1 ) {
90  fancyOut->setShowProcRank(true);
91  fancyOut->setOutputToRootOnly(-1);
92  }
93  std::ostream &out = *fancyOut;
94 #else
95  std::ostream &out = std::cout;
96 #endif
97 
98  Comm.SetTracebackMode(0); // This should shut down any error tracing
99  bool verbose = false;
100 
101  // Check if we should print results to standard out
102  if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
103 
104  // char tmp;
105  // if (rank==0) out << "Press any key to continue..."<< endl;
106  // if (rank==0) cin >> tmp;
107  // Comm.Barrier();
108 
109  int MyPID = Comm.MyPID();
110  int NumProc = Comm.NumProc();
111 
112  if (verbose && MyPID==0)
113  out << Epetra_Version() << endl << endl;
114 
115  if (verbose) out << Comm <<endl;
116 
117  bool verbose1 = verbose;
118 
119  // Redefine verbose to only print on PE 0
120  if (verbose && rank!=0) verbose = false;
121 
122  int NumMyElements = 10000;
123  int NumMyElements1 = NumMyElements; // Needed for localmap
124  int NumGlobalElements = NumMyElements*NumProc+EPETRA_MIN(NumProc,3);
125  if (MyPID < 3) NumMyElements++;
126  int IndexBase = 0;
127  int ElementSize = 7;
128 
129  // Test LocalMap constructor
130  // and Petra-defined uniform linear distribution constructor
131 
132  if (verbose) out << "\n*********************************************************" << endl;
133  if (verbose) out << "Checking Epetra_LocalMap(NumMyElements1, IndexBase, Comm)" << endl;
134  if (verbose) out << " and Epetra_BlockMap(NumGlobalElements, ElementSize, IndexBase, Comm)" << endl;
135  if (verbose) out << "*********************************************************" << endl;
136 
137  Epetra_LocalMap *LocalMap = new Epetra_LocalMap(NumMyElements1, IndexBase,
138  Comm);
139  Epetra_BlockMap * BlockMap = new Epetra_BlockMap(NumGlobalElements, ElementSize, IndexBase, Comm);
140  EPETRA_TEST_ERR(VectorTests(*BlockMap, verbose),ierr);
141 
142  EPETRA_TEST_ERR(MatrixTests(*BlockMap, *LocalMap, verbose),ierr);
143 
144  delete BlockMap;
145 
146  // Test User-defined linear distribution constructor
147 
148  if (verbose) out << "\n*********************************************************" << endl;
149  if (verbose) out << "Checking Epetra_BlockMap(NumGlobalElements, NumMyElements, ElementSize, IndexBase, Comm)" << endl;
150  if (verbose) out << "*********************************************************" << endl;
151 
152  BlockMap = new Epetra_BlockMap(NumGlobalElements, NumMyElements, ElementSize, IndexBase, Comm);
153 
154  EPETRA_TEST_ERR(VectorTests(*BlockMap, verbose),ierr);
155 
156  EPETRA_TEST_ERR(MatrixTests(*BlockMap, *LocalMap, verbose),ierr);
157 
158  delete BlockMap;
159 
160  // Test User-defined arbitrary distribution constructor
161  // Generate Global Element List. Do in reverse for fun!
162 
163  int * MyGlobalElements = new int[NumMyElements];
164  int MaxMyGID = (Comm.MyPID()+1)*NumMyElements-1+IndexBase;
165  if (Comm.MyPID()>2) MaxMyGID+=3;
166  for (i = 0; i<NumMyElements; i++) MyGlobalElements[i] = MaxMyGID-i;
167 
168  if (verbose) out << "\n*********************************************************" << endl;
169  if (verbose) out << "Checking Epetra_BlockMap(NumGlobalElements, NumMyElements, MyGlobalElements, ElementSize, IndexBase, Comm)" << endl;
170  if (verbose) out << "*********************************************************" << endl;
171 
172  BlockMap = new Epetra_BlockMap(NumGlobalElements, NumMyElements, MyGlobalElements, ElementSize,
173  IndexBase, Comm);
174  EPETRA_TEST_ERR(VectorTests(*BlockMap, verbose),ierr);
175 
176  EPETRA_TEST_ERR(MatrixTests(*BlockMap, *LocalMap, verbose),ierr);
177 
178  delete BlockMap;
179 
180  int * ElementSizeList = new int[NumMyElements];
181  int NumMyEquations = 0;
182  int NumGlobalEquations = 0;
183  for (i = 0; i<NumMyElements; i++)
184  {
185  ElementSizeList[i] = i%6+2; // blocksizes go from 2 to 7
186  NumMyEquations += ElementSizeList[i];
187  }
188  ElementSize = 7; // Set to maximum for use in checkmap
189  NumGlobalEquations = Comm.NumProc()*NumMyEquations;
190 
191  // Adjust NumGlobalEquations based on processor ID
192  if (Comm.NumProc() > 3)
193  {
194  if (Comm.MyPID()>2)
195  NumGlobalEquations += 3*((NumMyElements)%6+2);
196  else
197  NumGlobalEquations -= (Comm.NumProc()-3)*((NumMyElements-1)%6+2);
198  }
199 
200  if (verbose) out << "\n*********************************************************" << endl;
201  if (verbose) out << "Checking Epetra_BlockMap(NumGlobalElements, NumMyElements, MyGlobalElements, ElementSizeList, IndexBase, Comm)" << endl;
202  if (verbose) out << "*********************************************************" << endl;
203 
204  BlockMap = new Epetra_BlockMap(NumGlobalElements, NumMyElements, MyGlobalElements, ElementSizeList,
205  IndexBase, Comm);
206  EPETRA_TEST_ERR(VectorTests(*BlockMap, verbose),ierr);
207 
208  EPETRA_TEST_ERR(MatrixTests(*BlockMap, *LocalMap, verbose),ierr);
209 
210  // Test Copy constructor
211 
212  if (verbose) out << "\n*********************************************************" << endl;
213  if (verbose) out << "Checking Epetra_BlockMap(*BlockMap)" << endl;
214  if (verbose) out << "*********************************************************" << endl;
215 
216  Epetra_BlockMap * BlockMap1 = new Epetra_BlockMap(*BlockMap);
217 
218  EPETRA_TEST_ERR(VectorTests(*BlockMap, verbose),ierr);
219 
220  EPETRA_TEST_ERR(MatrixTests(*BlockMap, *LocalMap, verbose),ierr);
221 
222  delete [] ElementSizeList;
223  delete [] MyGlobalElements;
224  delete BlockMap;
225  delete BlockMap1;
226 
227 
228  // Test Petra-defined uniform linear distribution constructor
229 
230  if (verbose) out << "\n*********************************************************" << endl;
231  if (verbose) out << "Checking Epetra_Map(NumGlobalElements, IndexBase, Comm)" << endl;
232  if (verbose) out << "*********************************************************" << endl;
233 
234  Epetra_Map * Map = new Epetra_Map(NumGlobalElements, IndexBase, Comm);
235  EPETRA_TEST_ERR(VectorTests(*Map, verbose),ierr);
236 
237  EPETRA_TEST_ERR(MatrixTests(*Map, *LocalMap, verbose),ierr);
238 
239  delete Map;
240 
241  // Test User-defined linear distribution constructor
242 
243  if (verbose) out << "\n*********************************************************" << endl;
244  if (verbose) out << "Checking Epetra_Map(NumGlobalElements, NumMyElements, IndexBase, Comm)" << endl;
245  if (verbose) out << "*********************************************************" << endl;
246 
247  Map = new Epetra_Map(NumGlobalElements, NumMyElements, IndexBase, Comm);
248 
249  EPETRA_TEST_ERR(VectorTests(*Map, verbose),ierr);
250 
251  EPETRA_TEST_ERR(MatrixTests(*Map, *LocalMap, verbose),ierr);
252 
253  delete Map;
254 
255  // Test User-defined arbitrary distribution constructor
256  // Generate Global Element List. Do in reverse for fun!
257 
258  MyGlobalElements = new int[NumMyElements];
259  MaxMyGID = (Comm.MyPID()+1)*NumMyElements-1+IndexBase;
260  if (Comm.MyPID()>2) MaxMyGID+=3;
261  for (i = 0; i<NumMyElements; i++) MyGlobalElements[i] = MaxMyGID-i;
262 
263  if (verbose) out << "\n*********************************************************" << endl;
264  if (verbose) out << "Checking Epetra_Map(NumGlobalElements, NumMyElements, MyGlobalElements, IndexBase, Comm)" << endl;
265  if (verbose) out << "*********************************************************" << endl;
266 
267  Map = new Epetra_Map(NumGlobalElements, NumMyElements, MyGlobalElements,
268  IndexBase, Comm);
269  EPETRA_TEST_ERR(VectorTests(*Map, verbose),ierr);
270 
271  EPETRA_TEST_ERR(MatrixTests(*Map, *LocalMap, verbose),ierr);
272 
273  // Test Copy constructor
274 
275  if (verbose) out << "\n*********************************************************" << endl;
276  if (verbose) out << "Checking Epetra_Map(*Map)" << endl;
277  if (verbose) out << "*********************************************************" << endl;
278 
279  Epetra_Map Map1(*Map);
280 
281  EPETRA_TEST_ERR(VectorTests(*Map, verbose),ierr);
282 
283  EPETRA_TEST_ERR(MatrixTests(*Map, *LocalMap, verbose),ierr);
284 
285  delete [] MyGlobalElements;
286  delete Map;
287 
288  if (verbose1)
289  {
290  // Test Vector MFLOPS for 2D Dot Product
291  int M = 1;
292  int K = 1000000;
293  Epetra_Map Map2(-1, K, IndexBase, Comm);
294  Epetra_LocalMap Map3(M, IndexBase, Comm);
295 
296  Epetra_Vector A(Map2);A.Random();
297  Epetra_Vector B(Map2);B.Random();
298  Epetra_Vector C(Map3);C.Random();
299 
300  // Test Epetra_Vector label
301  const char* VecLabel = A.Label();
302  const char* VecLabel1 = "Epetra::Vector";
303  if (verbose) out << endl << endl <<"This should say " << VecLabel1 << ": " << VecLabel << endl << endl << endl;
304  EPETRA_TEST_ERR(strcmp(VecLabel1,VecLabel),ierr);
305  if (verbose) out << "Testing Assignment operator" << endl;
306 
307  double tmp1 = 1.00001* (double) (MyPID+1);
308  double tmp2 = tmp1;
309  A[1] = tmp1;
310  tmp2 = A[1];
311  out << "On PE "<< MyPID << " A[1] should equal = " << tmp1;
312  if (tmp1==tmp2) out << " and it does!" << endl;
313  else out << " but it equals " << tmp2;
314 
315  Comm.Barrier();
316 
317  if (verbose) out << endl << endl << "Testing MFLOPs" << endl;
318  Epetra_Flops counter;
319  C.SetFlopCounter(counter);
320  Epetra_Time mytimer(Comm);
321  C.Multiply('T', 'N', 0.5, A, B, 0.0);
322  double Multiply_time = mytimer.ElapsedTime();
323  double Multiply_flops = C.Flops();
324  if (verbose) out << "\n\nTotal FLOPs = " << Multiply_flops << endl;
325  if (verbose) out << "Total Time = " << Multiply_time << endl;
326  if (verbose) out << "MFLOPs = " << Multiply_flops/Multiply_time/1000000.0 << endl;
327 
328  Comm.Barrier();
329 
330  // Test Vector ostream operator with Petra-defined uniform linear distribution constructor
331  // and a small vector
332 
333  Epetra_Map Map4(100, IndexBase, Comm);
334  double * Dp = new double[100];
335  for (i=0; i<100; i++)
336  Dp[i] = i;
337  Epetra_Vector D(View, Map4,Dp);
338 
339  if (verbose) out << "\n\nTesting ostream operator: Multivector should be 100-by-2 and print i,j indices"
340  << endl << endl;
341  out << D << endl;
342 
343  if (verbose) out << "Traceback Mode value = " << D.GetTracebackMode() << endl;
344  delete [] Dp;
345  }
346 
347  delete LocalMap;
348 
349 #ifdef EPETRA_MPI
350  MPI_Finalize();
351 #endif
352 
353  return ierr;
354 
355 }
356 
void SetFlopCounter(const Epetra_Flops &FlopCounter_in)
Set the internal Epetra_Flops() pointer.
Epetra_Map: A class for partitioning vectors and matrices.
Definition: Epetra_Map.h:119
int Random()
Set multi-vector values to random numbers.
#define EPETRA_TEST_ERR(a, b)
basic_FancyOStream & setShowProcRank(const bool showProcRank)
double ElapsedTime(void) const
Epetra_Time elapsed time function.
static void SetTracebackMode(int TracebackModeValue)
Set the value of the Epetra_Object error traceback report mode.
int MatrixTests(const Epetra_BlockMap &Map, const Epetra_LocalMap &LocalMap, int NumVectors, bool verbose)
Epetra_Vector: A class for constructing and using dense vectors on a parallel computer.
#define EPETRA_MIN(x, y)
int MyPID() const
Return my process ID.
Epetra_MpiComm: The Epetra MPI Communication Class.
std::string Epetra_Version()
Epetra_Time: The Epetra Timing Class.
Definition: Epetra_Time.h:75
virtual const char * Label() const
Epetra_Object Label access funtion.
Epetra_BlockMap: A class for partitioning block element vectors and matrices.
int NumProc() const
Returns total number of processes (always returns 1 for SerialComm).
basic_FancyOStream & setOutputToRootOnly(const int rootRank)
Epetra_SerialComm: The Epetra Serial Communication Class.
Epetra_Flops: The Epetra Floating Point Operations Class.
Definition: Epetra_Flops.h:58
int Multiply(char TransA, char TransB, double ScalarAB, const Epetra_MultiVector &A, const Epetra_MultiVector &B, double ScalarThis)
Matrix-Matrix multiplication, this = ScalarThis*this + ScalarAB*A*B.
void Barrier() const
Epetra_SerialComm Barrier function.
static int GetTracebackMode()
Get the value of the Epetra_Object error report mode.
double Flops() const
Returns the number of floating point operations with this multi-vector.
int main(int argc, char *argv[])
int VectorTests(const Epetra_BlockMap &Map, bool verbose)
Epetra_LocalMap: A class for replicating vectors and matrices across multiple processors.