Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
Bug7758.cpp
Go to the documentation of this file.
1/*
2// @HEADER
3// ***********************************************************************
4//
5// Epetra: Linear Algebra Services Package
6// Copyright (2008) Sandia Corporation
7//
8// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9// the U.S. Government retains certain rights in this software.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// ************************************************************************
39// @HEADER
40*/
41
42// Creates vectors with different maps; tests results of export into them.
43// Documents the behavior of Epetra Add in Epetra_MultiVector for some
44// common (and a few uncommon) use cases.
45// Same tests as tpetra/core/test/MultiVector/Bug7758.cpp.
46
47#include "Epetra_Map.h"
48#include "Epetra_Vector.h"
49#include "Epetra_Export.h"
50
51#ifdef EPETRA_MPI
52#include <mpi.h>
53#include "Epetra_MpiComm.h"
54#else
55#include "Epetra_SerialComm.h"
56#endif
57
60{
61 int me = comm.MyPID();
62 int np = comm.NumProc();
63 int ierr = 0;
64
65 using vector_t = Epetra_Vector;
66 using map_t = Epetra_Map;
67
68 const int nGlobalEntries = 8 * np;
69 const double tgtScalar = 100. * (me+1);
70 const double srcScalar = 2.;
71
72 // Default one-to-one linear block map in Trilinos
73
74 map_t defaultMap(nGlobalEntries, 0, comm);
75
76 std::cout << me << " DEFAULT MAP" << std::endl;
77 defaultMap.Print(std::cout);
78
79 // Create vectors; see what the result is with CombineMode=ADD
80
81 vector_t defaultVecTgt(defaultMap);
82 defaultVecTgt.PutScalar(tgtScalar);
83
84 vector_t defaultVecSrc(defaultMap);
85 defaultVecSrc.PutScalar(srcScalar);
86
87 // Export Default-to-default: should be a copy of src to tgt
88
89 Epetra_Export defaultToDefault(defaultMap, defaultMap);
90 defaultVecTgt.Export(defaultVecSrc, defaultToDefault, Add);
91
92 std::cout << me << " DEFAULT TO DEFAULT " << std::endl;
93 defaultVecTgt.Print(std::cout);
94
95 // Check result; all vector entries should be srcScalar
96 for (int i = 0; i < defaultVecTgt.MyLength(); i++)
97 if (defaultVecTgt[i] != srcScalar) ierr++;
98
99 if (ierr > 0)
100 std::cout << "TEST FAILED: DEFAULT-TO-DEFAULT-EPETRA TEST HAD " << ierr
101 << " FAILURES ON RANK " << me << std::endl;
102
103 int gerr;
104 comm.SumAll(&ierr, &gerr, 1);
105
106 return gerr;
107}
108
111{
112 int me = comm.MyPID();
113 int np = comm.NumProc();
114 int ierr = 0;
115
116 using vector_t = Epetra_Vector;
117 using map_t = Epetra_Map;
118
119 const int nGlobalEntries = 8 * np;
120 const double tgtScalar = 100. * (me+1);
121 const double srcScalar = 2.;
122 std::vector<int> myEntries(nGlobalEntries);
123
124 // Default one-to-one linear block map in Trilinos
125
126 map_t defaultMap(nGlobalEntries, 0, comm);
127
128 std::cout << me << " DEFAULT MAP" << std::endl;
129 defaultMap.Print(std::cout);
130
131 // One-to-one cyclic map: deal out entries like cards
132
133 int nMyEntries = 0;
134 for (int i = 0; i < nGlobalEntries; i++) {
135 if (i % np == me) {
136 myEntries[nMyEntries++] = i;
137 }
138 }
139
140 int dummy = -1;
141 map_t cyclicMap(dummy, nMyEntries, &myEntries[0], 0, comm);
142
143 std::cout << me << " CYCLIC MAP" << std::endl;
144 cyclicMap.Print(std::cout);
145
146 // Create vectors; see what the result is with CombineMode=ADD
147
148 vector_t defaultVecTgt(defaultMap);
149 defaultVecTgt.PutScalar(tgtScalar);
150
151 vector_t cyclicVecSrc(cyclicMap);
152 cyclicVecSrc.PutScalar(srcScalar);
153
154 // Export Cyclic-to-default
155
156 Epetra_Export cyclicToDefault(cyclicMap, defaultMap);
157 defaultVecTgt.Export(cyclicVecSrc, cyclicToDefault, Add);
158
159 std::cout << me << " CYCLIC TO DEFAULT " << std::endl;
160 defaultVecTgt.Print(std::cout);
161
162 // Check result
163 for (int i = 0; i < defaultVecTgt.MyLength(); i++) {
164 if (cyclicMap.LID(defaultMap.GID(i)) != -1) {
165 // element is in both cyclic (source) and default (target) map;
166 // initial target value was overwritten
167 if (defaultVecTgt[i] != srcScalar) ierr++;
168 }
169 else {
170 // element is in only default (target) map;
171 // initial target value was not overwritten
172 if (defaultVecTgt[i] != tgtScalar + srcScalar) ierr++;
173 }
174 }
175 if (ierr > 0)
176 std::cout << "TEST FAILED: CYCLIC-TO-DEFAULT-EPETRA TEST HAD " << ierr
177 << " FAILURES ON RANK " << me << std::endl;
178
179 int gerr;
180 comm.SumAll(&ierr, &gerr, 1);
181
182 return gerr;
183}
184
187{
188 int me = comm.MyPID();
189 int np = comm.NumProc();
190 int ierr = 0;
191
192 using vector_t = Epetra_Vector;
193 using map_t = Epetra_Map;
194
195 if (np > 1) { // Need more than one proc to avoid duplicate entries in maps
196
197 const int nGlobalEntries = 8 * np;
198 const double tgtScalar = 100. * (me+1);
199 const double srcScalar = 2.;
200 std::vector<int> myEntries(nGlobalEntries);
201
202 // Default one-to-one linear block map in Trilinos
203
204 map_t defaultMap(nGlobalEntries, 0, comm);
205
206 std::cout << me << " DEFAULT MAP" << std::endl;
207 defaultMap.Print(std::cout);
208
209 // Overlap map; some entries are stored on two procs
210 int nMyEntries = 0;
211 for (int i = 0; i < defaultMap.NumMyElements()/2; i++) {
212 myEntries[nMyEntries++] = defaultMap.GID(i);
213 }
214 for (int i = 0; i < defaultMap.NumMyElements(); i++) {
215 myEntries[nMyEntries++] =
216 (defaultMap.MaxMyGID() + 1 + i) % nGlobalEntries;
217 }
218
219 int dummy = -1;
220 map_t overlapMap(dummy, nMyEntries, &myEntries[0], 0, comm);
221
222 std::cout << me << " OVERLAP MAP" << std::endl;
223 overlapMap.Print(std::cout);
224
225 // Create vectors; see what the result is with CombineMode=ADD
226
227 vector_t defaultVecTgt(defaultMap);
228 defaultVecTgt.PutScalar(tgtScalar);
229
230 vector_t overlapVecSrc(overlapMap);
231 overlapVecSrc.PutScalar(srcScalar);
232
233 // Export Overlap-to-default
234
235 Epetra_Export overlapToDefault(overlapMap, defaultMap);
236 defaultVecTgt.Export(overlapVecSrc, overlapToDefault, Add);
237
238 std::cout << me << " OVERLAP TO DEFAULT " << std::endl;
239 defaultVecTgt.Print(std::cout);
240
241 for (int i = 0; i < defaultVecTgt.MyLength()/2; i++)
242 if (defaultVecTgt[i] != srcScalar + srcScalar) ierr++; // overlapped
243 for (int i = defaultVecTgt.MyLength()/2;
244 i < defaultVecTgt.MyLength(); i++)
245 if (defaultVecTgt[i] != tgtScalar + srcScalar) ierr++; // not overlapped
246 if (ierr > 0)
247 std::cout << "TEST FAILED: OVERLAP-TO-DEFAULT-EPETRA TEST HAD " << ierr
248 << " FAILURES ON RANK " << me << std::endl;
249 }
250
251 int gerr;
252 comm.SumAll(&ierr, &gerr, 1);
253
254 return gerr;
255}
256
258
260{
261 // Test case showing behavior when target map is all on processor zero.
262 // In the source map, even numbered entries are on even numbered processors;
263 // odd numbered entreis are on odd numbered processors.
264 // In copyAndPermute, even numbered entries are copied from processor zero's
265 // source vector to the target vector, and odd numbered entries are unchanged.
266 // Then received values are added to the target vector. The result is that
267 // odd entries include the initial target values in their sum, while the
268 // even entries are not included.
269
270 int me = comm.MyPID();
271 int np = comm.NumProc();
272 int ierr = 0;
273
274 using vector_t = Epetra_Vector;
275 using map_t = Epetra_Map;
276
277 const int nGlobalEntries = 8 * np;
278 const double tgtScalar = 100. * (me+1);
279 const double srcScalar = 2.;
280 std::vector<int> myEntries(nGlobalEntries);
281
282 // Odd entries given to odd procs; even entries given to even procs
283 int nMyEntries = 0;
284 for (int i = 0; i < nGlobalEntries; i++) {
285 if (int(i % 2) == (me % 2)) {
286 myEntries[nMyEntries++] = i;
287 }
288 }
289
290 int dummy = -1;
291 map_t oddEvenMap(dummy, nMyEntries, &myEntries[0], 0, comm);
292
293 std::cout << me << " ODDEVEN MAP" << std::endl;
294 oddEvenMap.Print(std::cout);
295
296 // Map with all entries on one processor
297
298 dummy = -1;
299 int nSerialEntries = (me == 0 ? nGlobalEntries : 0);
300 map_t serialMap(dummy, nSerialEntries, 0, comm);
301
302 std::cout << me << " SERIAL MAP" << std::endl;
303 serialMap.Print(std::cout);
304
305 // Create vectors; see what the result is with CombineMode=ADD
306
307 vector_t oddEvenVecSrc(oddEvenMap);
308 oddEvenVecSrc.PutScalar(srcScalar);
309
310 vector_t serialVecTgt(serialMap);
311 serialVecTgt.PutScalar(tgtScalar);
312
313 // Export oddEven-to-serial
314
315 Epetra_Export oddEvenToSerial(oddEvenMap, serialMap);
316 serialVecTgt.Export(oddEvenVecSrc, oddEvenToSerial, Add);
317
318 std::cout << me << " ODDEVEN TO SERIAL " << std::endl;
319 serialVecTgt.Print(std::cout);
320
321 // Check result; all vector entries should be srcScalar
322
323 for (int i = 0; i < serialVecTgt.MyLength(); i++) {
324 double nCopies = double(((np+1) / 2) - ((i % 2 == 1) && (np % 2 == 1)));
325 if (serialVecTgt[i] != (i%2 ? tgtScalar : 0.) + srcScalar * nCopies)
326 ierr++;
327 }
328 if (ierr > 0)
329 std::cout << "TEST FAILED: ODDEVEN-TO-SERIAL TEST HAD " << ierr
330 << " FAILURES ON RANK " << me << std::endl;
331
332 int gerr;
333 comm.SumAll(&ierr, &gerr, 1);
334
335 return gerr;
336}
337
340{
341 int me = comm.MyPID();
342 int np = comm.NumProc();
343 int ierr = 0;
344
345 using vector_t = Epetra_Vector;
346 using map_t = Epetra_Map;
347
348 if (np > 1) { // Need more than one proc to avoid duplicate entries in maps
349
350 const int nGlobalEntries = 8 * np;
351 const double tgtScalar = 100. * (me+1);
352 const double srcScalar = 2.;
353 std::vector<int> myEntries(nGlobalEntries);
354
355 // Default one-to-one linear block map in Trilinos
356
357 map_t defaultMap(nGlobalEntries, 0, comm);
358
359 std::cout << me << " DEFAULT MAP" << std::endl;
360 defaultMap.Print(std::cout);
361
362 // Superset map; some entries are stored on two procs
363 int nMyEntries = 0;
364 for (int i = 0; i < defaultMap.NumMyElements(); i++) {
365 myEntries[nMyEntries++] = defaultMap.GID(i);
366 }
367 for (int i = 0; i < defaultMap.NumMyElements()/2; i++) {
368 myEntries[nMyEntries++] =
369 (defaultMap.MaxMyGID() + 1 + i) % nGlobalEntries;
370 }
371
372 int dummy = -1;
373 map_t supersetMap(dummy, nMyEntries, &myEntries[0], 0, comm);
374
375 std::cout << me << " SUPERSET MAP" << std::endl;
376 supersetMap.Print(std::cout);
377
378 // Create vectors; see what the result is with CombineMode=ADD
379
380 vector_t defaultVecTgt(defaultMap);
381 defaultVecTgt.PutScalar(tgtScalar);
382
383 vector_t supersetVecSrc(supersetMap);
384 supersetVecSrc.PutScalar(srcScalar);
385
386 // Export Superset-to-default
387
388 Epetra_Export supersetToDefault(supersetMap, defaultMap);
389 defaultVecTgt.Export(supersetVecSrc, supersetToDefault, Add);
390
391 std::cout << me << " SUPERSET TO DEFAULT " << std::endl;
392 defaultVecTgt.Print(std::cout);
393
394 for (int i = 0; i < defaultVecTgt.MyLength()/2; i++)
395 if (defaultVecTgt[i] != 2 * srcScalar) ierr++; // overlapped
396 for (int i = defaultVecTgt.MyLength()/2;
397 i < defaultVecTgt.MyLength(); i++)
398 if (defaultVecTgt[i] != srcScalar) ierr++; // not overlapped
399 if (ierr > 0)
400 std::cout << "TEST FAILED: SUPERSET-TO-DEFAULT-EPETRA TEST HAD " << ierr
401 << " FAILURES ON RANK " << me << std::endl;
402 }
403
404 int gerr;
405 comm.SumAll(&ierr, &gerr, 1);
406
407 return gerr;
408}
409
412{
413 // This use case is similar to matrix assembly case in which user
414 // has a map of owned entries and a map of shared entries, with no
415 // overlap between the maps. In this case, received shared entries are
416 // added to the owned entries' values, as copyAndPermute is never called.
417
418 int me = comm.MyPID();
419 int np = comm.NumProc();
420 int ierr = 0;
421
422 if (np > 1) { // Need more than one proc to avoid duplicate entries in maps
423 using vector_t = Epetra_Vector;
424 using map_t = Epetra_Map;
425
426 const int nGlobalEntries = 8 * np;
427 const double tgtScalar = 100. * (me+1);
428 const double srcScalar = 2.;
429 std::vector<int> myEntries(nGlobalEntries);
430
431 // Default one-to-one linear block map in Trilinos
432
433 map_t defaultMap(nGlobalEntries, 0, comm);
434
435 std::cout << me << " DEFAULT MAP" << std::endl;
436 defaultMap.Print(std::cout);
437
438 // Map with no sames or permutes
439 int nMyEntries = 0;
440 for (int i = 0; i < defaultMap.NumMyElements(); i++) {
441 myEntries[nMyEntries++] =
442 (defaultMap.MaxMyGID() + 1 + i) % nGlobalEntries;
443 }
444
445 int dummy = -1;
446 map_t noSamesMap(dummy, nMyEntries, &myEntries[0], 0, comm);
447
448 std::cout << me << " NOSAMES MAP" << std::endl;
449 noSamesMap.Print(std::cout);
450
451 // Create vectors; see what the result is with CombineMode=ADD
452
453 vector_t defaultVecTgt(defaultMap);
454 defaultVecTgt.PutScalar(tgtScalar);
455
456 vector_t noSamesVecSrc(noSamesMap);
457 noSamesVecSrc.PutScalar(srcScalar);
458
459 // Export Overlap-to-default
460
461 Epetra_Export noSamesToDefault(noSamesMap, defaultMap);
462 defaultVecTgt.Export(noSamesVecSrc, noSamesToDefault, Add);
463
464 std::cout << me << " NOSAMES TO DEFAULT " << std::endl;
465 defaultVecTgt.Print(std::cout);
466
467 for (int i = 0; i < defaultVecTgt.MyLength(); i++)
468 if (defaultVecTgt[i] != tgtScalar + srcScalar) ierr++;
469 if (ierr > 0)
470 std::cout << "TEST FAILED: NOSAMES-TO-DEFAULT TEST HAD " << ierr
471 << " FAILURES ON RANK " << me << std::endl;
472 }
473
474 int gerr;
475 comm.SumAll(&ierr, &gerr, 1);
476
477 return gerr;
478}
479
481int main(int argc, char **argv)
482{
483#ifdef EPETRA_MPI
484 // Initialize MPI
485 MPI_Init(&argc,&argv);
486 Epetra_MpiComm Comm(MPI_COMM_WORLD);
487#else
489#endif
490
491 int gerr = 0;
492
493 gerr += DefaultToDefaultEpetra(Comm);
494 gerr += CyclicToDefaultEpetra(Comm);
495 gerr += OverlapToDefaultEpetra(Comm);
496 gerr += OddEvenToSerialEpetra(Comm);
497 gerr += SupersetToDefaultEpetra(Comm);
498 gerr += NoSamesToDefaultEpetra(Comm);
499
500 if (Comm.MyPID() == 0) {
501 if (gerr > 0) std::cout << "TEST FAILED" << std::endl;
502 else std::cout << "TEST PASSED" << std::endl;
503 }
504
505#ifdef EPETRA_MPI
506 MPI_Finalize();
507#endif
508
509 return gerr;
510
511}
int SupersetToDefaultEpetra(const Epetra_Comm &comm)
Definition: Bug7758.cpp:339
int OverlapToDefaultEpetra(const Epetra_Comm &comm)
Definition: Bug7758.cpp:186
int main(int argc, char **argv)
Definition: Bug7758.cpp:481
int DefaultToDefaultEpetra(const Epetra_Comm &comm)
Definition: Bug7758.cpp:59
int CyclicToDefaultEpetra(const Epetra_Comm &comm)
Definition: Bug7758.cpp:110
int NoSamesToDefaultEpetra(const Epetra_Comm &comm)
Definition: Bug7758.cpp:411
int OddEvenToSerialEpetra(const Epetra_Comm &comm)
Definition: Bug7758.cpp:259
Epetra_Comm: The Epetra Communication Abstract Base Class.
Definition: Epetra_Comm.h:73
virtual int NumProc() const =0
Returns total number of processes.
virtual int MyPID() const =0
Return my process ID.
virtual int SumAll(double *PartialSums, double *GlobalSums, int Count) const =0
Epetra_Comm Global Sum function.
Epetra_Export: This class builds an export object for efficient exporting of off-processor elements.
Definition: Epetra_Export.h:62
Epetra_Map: A class for partitioning vectors and matrices.
Definition: Epetra_Map.h:119
Epetra_MpiComm: The Epetra MPI Communication Class.
int MyPID() const
Return my process ID.
Epetra_SerialComm: The Epetra Serial Communication Class.
Epetra_Vector: A class for constructing and using dense vectors on a parallel computer.