Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
test/Map/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_Map (and Epetra_LocalMap) Test routine
44#include "Epetra_Time.h"
45#include "Epetra_Map.h"
46#include "Epetra_LocalMap.h"
47#ifdef EPETRA_MPI
48#include "Epetra_MpiComm.h"
49#include <mpi.h>
50#else
51#include "Epetra_SerialComm.h"
52#endif
53#include "checkmap.h"
54#include "../epetra_test_err.h"
55#include "Epetra_Version.h"
56
57int checkMapDataClass(Epetra_Comm& Comm, int verbose);
58int checkLocalMapDataClass(Epetra_Comm& Comm, int verbose);
59
60int main(int argc, char *argv[]) {
61
62 int ierr=0, returnierr=0;
63
64#ifdef EPETRA_MPI
65 MPI_Init(&argc,&argv);
66 Epetra_MpiComm Comm(MPI_COMM_WORLD);
67#else
69#endif
70
71 bool verbose = false;
72
73 // Check if we should print results to standard out
74 if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
75
76
77 if (!verbose) {
78 Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
79 }
80 int MyPID = Comm.MyPID();
81 int NumProc = Comm.NumProc();
82
83 if (verbose && MyPID==0)
84 cout << Epetra_Version() << endl << endl;
85
86 if (verbose) cout << Comm << endl;
87
88 bool verbose1 = verbose;
89 if (verbose) verbose = (MyPID==0);
90
91 int NumMyElements = 10000;
92 int NumMyElements1 = NumMyElements; // Used for local map
93 int NumGlobalElements = NumMyElements*NumProc+EPETRA_MIN(NumProc,3);
94 if (MyPID < 3) NumMyElements++;
95 int IndexBase = 0;
96 bool DistributedGlobal = (NumGlobalElements>NumMyElements);
97
98 Epetra_Map* Map;
99
100 // Test exceptions
101
102 if (verbose)
103 cout << "*******************************************************************************************" << endl
104 << " Testing Exceptions (Expect error messages if EPETRA_NO_ERROR_REPORTS is not defined" << endl
105 << "*******************************************************************************************" << endl
106 << endl << endl;
107
108 try {
109 if (verbose) cout << "Checking Epetra_Map(-2, IndexBase, Comm)" << endl;
110 Map = new Epetra_Map(-2, IndexBase, Comm);
111 }
112 catch (int Error) {
113 if (Error!=-1) {
114 if (Error!=0) {
115 EPETRA_TEST_ERR(Error,returnierr);
116 if (verbose) cout << "Error code should be -1" << endl;
117 }
118 else {
119 cout << "Error code = " << Error << "Should be -1" << endl;
120 returnierr+=1;
121 }
122 }
123 else if (verbose) cout << "Checked OK\n\n" << endl;
124 }
125
126 try {
127 if (verbose) cout << "Checking Epetra_Map(2, 3, IndexBase, Comm)" << endl;
128 Map = new Epetra_Map(2, 3, IndexBase, Comm);
129 }
130 catch (int Error) {
131 if (Error!=-4) {
132 if (Error!=0) {
133 EPETRA_TEST_ERR(Error,returnierr);
134 if (verbose) cout << "Error code should be -4" << endl;
135 }
136 else {
137 cout << "Error code = " << Error << "Should be -4" << endl;
138 returnierr+=1;
139 }
140 }
141 else if (verbose) cout << "Checked OK\n\n" << endl;
142 }
143
144 if (verbose) cerr << flush;
145 if (verbose) cout << flush;
146 Comm.Barrier();
147 if (verbose)
148 cout << endl << endl
149 << "*******************************************************************************************" << endl
150 << " Testing valid constructor now......................................................" << endl
151 << "*******************************************************************************************" << endl
152 << endl << endl;
153
154 // Test Epetra-defined uniform linear distribution constructor
155 Map = new Epetra_Map(NumGlobalElements, IndexBase, Comm);
156 if (verbose) cout << "Checking Epetra_Map(NumGlobalElements, IndexBase, Comm)" << endl;
157 ierr = checkmap(*Map, NumGlobalElements, NumMyElements, 0,
158 IndexBase, Comm, DistributedGlobal);
159
160 EPETRA_TEST_ERR(ierr,returnierr);
161 if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
162
163 delete Map;
164
165 // Test User-defined linear distribution constructor
166 Map = new Epetra_Map(NumGlobalElements, NumMyElements, IndexBase, Comm);
167
168 if (verbose) cout << "Checking Epetra_Map(NumGlobalElements, NumMyElements, IndexBase, Comm)" << endl;
169 ierr = checkmap(*Map, NumGlobalElements, NumMyElements, 0,
170 IndexBase, Comm, DistributedGlobal);
171
172 EPETRA_TEST_ERR(ierr,returnierr);
173 if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
174 delete Map;
175
176 // Test User-defined arbitrary distribution constructor
177 // Generate Global Element List. Do in reverse for fun!
178
179 int * MyGlobalElements = new int[NumMyElements];
180 int MaxMyGID = (Comm.MyPID()+1)*NumMyElements-1+IndexBase;
181 if (Comm.MyPID()>2) MaxMyGID+=3;
182 for (int i = 0; i<NumMyElements; i++) MyGlobalElements[i] = MaxMyGID-i;
183
184 Map = new Epetra_Map(NumGlobalElements, NumMyElements, MyGlobalElements,
185 IndexBase, Comm);
186 if (verbose) cout << "Checking Epetra_Map(NumGlobalElements, NumMyElements, MyGlobalElements, IndexBase, Comm)" << endl;
187 ierr = checkmap(*Map, NumGlobalElements, NumMyElements, MyGlobalElements,
188 IndexBase, Comm, DistributedGlobal);
189
190 EPETRA_TEST_ERR(ierr,returnierr);
191 if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
192 // Test Copy constructor
193 Epetra_Map* Map1 = new Epetra_Map(*Map);
194
195 // Test SameAs() method
196 bool same = Map1->SameAs(*Map);
197 EPETRA_TEST_ERR(!(same==true),ierr);// should return true since Map1 is a copy of Map
198
199 Epetra_BlockMap* Map2 = new Epetra_Map(NumGlobalElements, NumMyElements, MyGlobalElements, IndexBase, Comm);
200 same = Map2->SameAs(*Map);
201 EPETRA_TEST_ERR(!(same==true),ierr); // Map and Map2 were created with the same sets of parameters
202 delete Map2;
203
204 // now test SameAs() on a map that is different
205
206 Map2 = new Epetra_Map(NumGlobalElements, NumMyElements, MyGlobalElements, IndexBase-1, Comm);
207 same = Map2->SameAs(*Map);
208 EPETRA_TEST_ERR(!(same==false),ierr); // IndexBases are different
209 delete Map2;
210
211 // Back to testing copy constructor
212 if (verbose) cout << "Checking Epetra_Map(*Map)" << endl;
213 ierr = checkmap(*Map1, NumGlobalElements, NumMyElements, MyGlobalElements,
214 IndexBase, Comm, DistributedGlobal);
215
216 EPETRA_TEST_ERR(ierr,returnierr);
217 if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
218 Epetra_Map* SmallMap = 0;
219 if (verbose1) {
220 // Build a small map for test cout. Use 10 elements from current map
221 int* MyEls = Map->MyGlobalElements();
222 int IndBase = Map->IndexBase();
223 int MyLen = EPETRA_MIN(10+Comm.MyPID(),Map->NumMyElements());
224 SmallMap = new Epetra_Map(-1, MyLen, MyEls, IndBase, Comm);
225 }
226
227 delete [] MyGlobalElements;
228 delete Map;
229 delete Map1;
230
231 // Test reference-counting in Epetra_Map
232 if (verbose) cout << "Checking Epetra_Map reference counting" << endl;
233 ierr = checkMapDataClass(Comm, verbose);
234 EPETRA_TEST_ERR(ierr,returnierr);
235 if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
236
237 // Test LocalMap constructor
238 Epetra_LocalMap* LocalMap = new Epetra_LocalMap(NumMyElements1, IndexBase, Comm);
239 if (verbose) cout << "Checking Epetra_LocalMap(NumMyElements1, IndexBase, Comm)" << endl;
240 ierr = checkmap(*LocalMap, NumMyElements1, NumMyElements1, 0, IndexBase, Comm, false);
241
242 EPETRA_TEST_ERR(ierr,returnierr);
243 if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
244 // Test Copy constructor
245 Epetra_LocalMap* LocalMap1 = new Epetra_LocalMap(*LocalMap);
246 if (verbose) cout << "Checking Epetra_LocalMap(*LocalMap)" << endl;
247 ierr = checkmap(*LocalMap1, NumMyElements1, NumMyElements1, 0, IndexBase, Comm, false);
248
249 EPETRA_TEST_ERR(ierr,returnierr);
250 if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
251 delete LocalMap1;
252 delete LocalMap;
253
254 // Test reference-counting in Epetra_LocalMap
255 if (verbose) cout << "Checking Epetra_LocalMap reference counting" << endl;
256 ierr = checkLocalMapDataClass(Comm, verbose);
257 EPETRA_TEST_ERR(ierr,returnierr);
258 if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
259
260 // Test output
261 if (verbose1) {
262 if (verbose) cout << "Test ostream << operator" << endl << flush;
263 cout << *SmallMap;
264 delete SmallMap;
265 }
266
267#ifdef EPETRA_MPI
268 MPI_Finalize();
269#endif
270
271 return returnierr;
272}
273
274int checkMapDataClass(Epetra_Comm& Comm, int verbose) {
275 int returnierr = 0;
276 int NumGlobalElements = 1000;
277 int IndexBase = 0;
278
279 Epetra_Map m1(NumGlobalElements, IndexBase, Comm);
280 int m1count = m1.ReferenceCount();
281 const Epetra_BlockMapData* m1addr = m1.DataPtr();
282 EPETRA_TEST_ERR(!(m1count==1),returnierr); // count should be 1
283 if(verbose) cout << "Default constructor. \nm1= " << m1count << " " << m1addr << endl;
284
285 Epetra_Map* m2 = new Epetra_Map(m1);
286 int m2count = m2->ReferenceCount();
287 const Epetra_BlockMapData* m2addr = m2->DataPtr();
288 int m1countold = m1count;
289 m1count = m1.ReferenceCount();
290 EPETRA_TEST_ERR(!(m2count==m1count && m1count==(m1countold+1)),returnierr); // both counts should be 2
291 EPETRA_TEST_ERR(!(m1addr==m2addr),returnierr); // addresses should be same
292 if(verbose) cout << "Copy constructor. \nm1= " << m1count << " " << m1addr
293 << "\nm2= " << m2count << " " << m2addr << endl;
294
295 delete m2;
296 m1countold = m1count;
297 m1count = m1.ReferenceCount();
298 EPETRA_TEST_ERR(!(m1count == m1countold-1), returnierr); // count should have decremented (to 1)
299 EPETRA_TEST_ERR(!(m1addr == m1.DataPtr()),returnierr); // m1addr should be unchanged
300 if(verbose) cout << "m2 destroyed. \nm1= " << m1count << " " << m1addr << endl;
301
302 { // inside of braces to test stack deallocation.
303 if(verbose) cout << "Assignment operator, post construction" << endl;
304 Epetra_Map m3(NumGlobalElements, IndexBase+1, Comm);
305 int m3count = m3.ReferenceCount();
306 const Epetra_BlockMapData* m3addr = m3.DataPtr();
307 EPETRA_TEST_ERR(!(m3count==1),returnierr); // m3count should be 1 initially
308 EPETRA_TEST_ERR(!(m1addr!=m3addr),returnierr); // m1 and m3 should have different ptr addresses
309 if(verbose) cout << "Prior to assignment: \nm1= " << m1count << " " << m1addr
310 << "\nm3= " << m3count << " " << m3addr << endl;
311 m3 = m1;
312 m3count = m3.ReferenceCount();
313 m3addr = m3.DataPtr();
314 m1countold = m1count;
315 m1count = m1.ReferenceCount();
316 EPETRA_TEST_ERR(!(m3count==m1count && m1count==m1countold+1),returnierr); // both counts should be 2
317 EPETRA_TEST_ERR(!(m1addr==m3addr),returnierr); // addresses should be same
318 if(verbose) cout << "After assignment: \nm1= " << m1count << " " << m1addr
319 << "\nm3= " << m3count << " " << m3addr << endl;
320 }
321 m1countold = m1count;
322 m1count = m1.ReferenceCount();
323 EPETRA_TEST_ERR(!(m1count==m1countold-1), returnierr); // count should have decremented (to 1)
324 EPETRA_TEST_ERR(!(m1addr== m1.DataPtr()),returnierr); // m1addr should be unchanged
325 if(verbose) cout << "m3 destroyed. \nm1= " << m1count << " " << m1addr << endl;
326
327 return(returnierr);
328}
329
330int checkLocalMapDataClass(Epetra_Comm& Comm, int verbose) {
331 int returnierr = 0;
332 int NumMyElements = 100;
333 int IndexBase = 0;
334
335 Epetra_LocalMap m1(NumMyElements, IndexBase, Comm);
336 int m1count = m1.ReferenceCount();
337 const Epetra_BlockMapData* m1addr = m1.DataPtr();
338 EPETRA_TEST_ERR(!(m1count==1),returnierr); // count should be 1
339 if(verbose) cout << "Default constructor. \nm1= " << m1count << " " << m1addr << endl;
340
341 Epetra_LocalMap* m2 = new Epetra_LocalMap(m1);
342 int m2count = m2->ReferenceCount();
343 const Epetra_BlockMapData* m2addr = m2->DataPtr();
344 int m1countold = m1count;
345 m1count = m1.ReferenceCount();
346 EPETRA_TEST_ERR(!(m2count==m1count && m1count==(m1countold+1)),returnierr); // both counts should be 2
347 EPETRA_TEST_ERR(!(m1addr==m2addr),returnierr); // addresses should be same
348 if(verbose) cout << "Copy constructor. \nm1= " << m1count << " " << m1addr
349 << "\nm2= " << m2count << " " << m2addr << endl;
350
351 delete m2;
352 m1countold = m1count;
353 m1count = m1.ReferenceCount();
354 EPETRA_TEST_ERR(!(m1count == m1countold-1), returnierr); // count should have decremented (to 1)
355 EPETRA_TEST_ERR(!(m1addr == m1.DataPtr()),returnierr); // m1addr should be unchanged
356 if(verbose) cout << "m2 destroyed. \nm1= " << m1count << " " << m1addr << endl;
357
358 { // inside of braces to test stack deallocation.
359 if(verbose) cout << "Assignment operator, post construction" << endl;
360 Epetra_LocalMap m3(NumMyElements, IndexBase+1, Comm);
361 int m3count = m3.ReferenceCount();
362 const Epetra_BlockMapData* m3addr = m3.DataPtr();
363 EPETRA_TEST_ERR(!(m3count==1),returnierr); // m3count should be 1 initially
364 EPETRA_TEST_ERR(!(m1addr!=m3addr),returnierr); // m1 and m3 should have different ptr addresses
365 if(verbose) cout << "Prior to assignment: \nm1= " << m1count << " " << m1addr
366 << "\nm3= " << m3count << " " << m3addr << endl;
367 m3 = m1;
368 m3count = m3.ReferenceCount();
369 m3addr = m3.DataPtr(); // cast int* to int
370 m1countold = m1count;
371 m1count = m1.ReferenceCount();
372 EPETRA_TEST_ERR(!(m3count==m1count && m1count==m1countold+1),returnierr); // both counts should be 2
373 EPETRA_TEST_ERR(!(m1addr==m3addr),returnierr); // addresses should be same
374 if(verbose) cout << "After assignment: \nm1= " << m1count << " " << m1addr
375 << "\nm3= " << m3count << " " << m3addr << endl;
376 }
377 m1countold = m1count;
378 m1count = m1.ReferenceCount();
379 EPETRA_TEST_ERR(!(m1count==m1countold-1), returnierr); // count should have decremented (to 1)
380 EPETRA_TEST_ERR(!(m1addr== m1.DataPtr()),returnierr); // m1addr should be unchanged
381 if(verbose) cout << "m3 destroyed. \nm1= " << m1count << " " << m1addr << endl;
382
383 return(returnierr);
384}
int checkmap(Epetra_BlockMap &Map, int NumGlobalElements, int NumMyElements, int *MyGlobalElements, int ElementSize, int *ElementSizeList, int NumGlobalPoints, int NumMyPoints, int IndexBase, Epetra_Comm &Comm, bool DistributedGlobal, bool IsOneToOne)
#define EPETRA_MIN(x, y)
std::string Epetra_Version()
Epetra_BlockMapData: The Epetra BlockMap Data Class.
Epetra_BlockMap: A class for partitioning block element vectors and matrices.
int MyGlobalElements(int *MyGlobalElementList) const
Puts list of global elements on this processor into the user-provided array.
int ReferenceCount() const
Returns the reference count of BlockMapData.
int IndexBase() const
Index base for this map.
bool SameAs(const Epetra_BlockMap &Map) const
Returns true if this and Map are identical maps.
const Epetra_BlockMapData * DataPtr() const
Returns a pointer to the BlockMapData instance this BlockMap uses.
int NumMyElements() const
Number of elements on the calling processor.
Epetra_Comm: The Epetra Communication Abstract Base Class.
Definition: Epetra_Comm.h:73
Epetra_LocalMap: A class for replicating vectors and matrices across multiple processors.
Epetra_Map: A class for partitioning vectors and matrices.
Definition: Epetra_Map.h:119
Epetra_MpiComm: The Epetra MPI Communication Class.
void Barrier() const
Epetra_MpiComm Barrier function.
int NumProc() const
Returns total number of processes.
int MyPID() const
Return my process ID.
static void SetTracebackMode(int TracebackModeValue)
Set the value of the Epetra_Object error traceback report mode.
Epetra_SerialComm: The Epetra Serial Communication Class.
#define EPETRA_TEST_ERR(a, b)
int main(int argc, char *argv[])
int checkLocalMapDataClass(Epetra_Comm &Comm, int verbose)
int checkMapDataClass(Epetra_Comm &Comm, int verbose)