Epetra Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
test/ImportExport_LL/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 #include "Epetra_Map.h"
44 #include "Epetra_Time.h"
45 #include "Epetra_CrsMatrix.h"
46 #include "Epetra_Vector.h"
47 #include "Epetra_LongLongVector.h"
48 #include "Epetra_Import.h"
49 #include "Epetra_Export.h"
50 #include "Epetra_OffsetIndex.h"
51 #ifdef EPETRA_MPI
52 #include "Epetra_MpiComm.h"
53 #include "mpi.h"
54 #else
55 #include "Epetra_SerialComm.h"
56 #endif
57 #ifndef __cplusplus
58 #define __cplusplus
59 #endif
60 #include "../epetra_test_err.h"
61 #include "Epetra_Version.h"
62 
66 
67 int main(int argc, char *argv[])
68 {
69  int ierr = 0, i, j, forierr = 0;
70 
71 #ifdef EPETRA_MPI
72  // Initialize MPI
73  MPI_Init(&argc,&argv);
74  Epetra_MpiComm Comm( MPI_COMM_WORLD );
75 #else
76  Epetra_SerialComm Comm;
77 #endif
78 
79  bool verbose = false;
80 
81  // Check if we should print results to standard out
82  if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
83 
84 
85 
86 
87  //char tmp;
88  //if (Comm.MyPID()==0) cout << "Press any key to continue..."<< endl;
89  //if (Comm.MyPID()==0) cin >> tmp;
90  //Comm.Barrier();
91 
92  Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
93  int MyPID = Comm.MyPID();
94  int NumProc = Comm.NumProc();
95 
96  if (verbose && MyPID==0)
97  cout << Epetra_Version() << endl << endl;
98 
99  if (verbose) cout << "Processor "<<MyPID<<" of "<< NumProc
100  << " is alive."<<endl;
101 
102  // Redefine verbose to only print on PE 0
103  if (verbose && Comm.MyPID()!=0) verbose = false;
104 
105  int NumMyEquations = 20;
106  long long NumGlobalEquations = NumMyEquations*NumProc+EPETRA_MIN(NumProc,3);
107  if (MyPID < 3) NumMyEquations++;
108  // Construct a Source Map that puts approximately the same Number of equations on each processor in
109  // uniform global ordering
110 
111  Epetra_Map SourceMap(NumGlobalEquations, NumMyEquations, 0LL, Comm);
112 
113  // Get update list and number of local equations from newly created Map
114  int NumMyElements = SourceMap.NumMyElements();
115  long long * SourceMyGlobalElements = new long long[NumMyElements];
116  SourceMap.MyGlobalElements(SourceMyGlobalElements);
117 
118  // Construct a Target Map that will contain:
119  // some unchanged elements (relative to the soure map),
120  // some permuted elements
121  // some off-processor elements
122  Epetra_Vector RandVec(SourceMap);
123  RandVec.Random(); // This creates a vector of random numbers between negative one and one.
124 
125  long long *TargetMyGlobalElements = new long long[NumMyElements];
126 
127  long long MinGID = SourceMap.MinMyGID64();
128  for (i=0; i< NumMyEquations/2; i++) TargetMyGlobalElements[i] = i + MinGID; // Half will be the same...
129  for (i=NumMyEquations/2; i<NumMyEquations; i++) {
130  int index = abs((int)(((double) (NumGlobalEquations-1) ) * RandVec[i]));
131  TargetMyGlobalElements[i] = EPETRA_MIN(NumGlobalEquations-1,(long long) EPETRA_MAX(0,index));
132  }
133 
134  int NumSameIDs = 0;
135  int NumPermutedIDs = 0;
136  int NumRemoteIDs = 0;
137  bool StillContiguous = true;
138  for (i=0; i < NumMyEquations; i++) {
139  if (SourceMyGlobalElements[i]==TargetMyGlobalElements[i] && StillContiguous)
140  NumSameIDs++;
141  else if (SourceMap.MyGID(TargetMyGlobalElements[i])) {
142  StillContiguous = false;
143  NumPermutedIDs++;
144  }
145  else {
146  StillContiguous = false;
147  NumRemoteIDs++;
148  }
149  }
150  EPETRA_TEST_ERR(!(NumMyEquations==NumSameIDs+NumPermutedIDs+NumRemoteIDs),ierr);
151 
152  Epetra_Map TargetMap((long long) -1, NumMyElements, TargetMyGlobalElements, 0LL, Comm);
153 
154  // Create a multivector whose elements are GlobalID * (column number +1)
155 
156  int NumVectors = 3;
157  Epetra_MultiVector SourceMultiVector(SourceMap, NumVectors);
158  for (j=0; j < NumVectors; j++)
159  for (i=0; i < NumMyElements; i++)
160  SourceMultiVector[j][i] = (double) SourceMyGlobalElements[i]*(j+1);
161 
162  // Create a target multivector that we will fill using an Import
163 
164  Epetra_MultiVector TargetMultiVector(TargetMap, NumVectors);
165 
166  Epetra_Import Importer(TargetMap, SourceMap);
167 
168  EPETRA_TEST_ERR(!(TargetMultiVector.Import(SourceMultiVector, Importer, Insert)==0),ierr);
169 
170  // Test Target against expected values
171  forierr = 0;
172  for (j=0; j < NumVectors; j++)
173  for (i=0; i < NumMyElements; i++) {
174  if (TargetMultiVector[j][i]!= (double) TargetMyGlobalElements[i]*(j+1))
175  cout << "TargetMultiVector["<<i<<"]["<<j<<"] = " << TargetMultiVector[j][i]
176  << " TargetMyGlobalElements[i]*(j+1) = " << TargetMyGlobalElements[i]*(j+1) << endl;
177  forierr += !(TargetMultiVector[j][i]== (double) TargetMyGlobalElements[i]*(j+1));
178  }
179  EPETRA_TEST_ERR(forierr,ierr);
180 
181  if (verbose) cout << "MultiVector Import using Importer Check OK" << endl << endl;
182 
183 
185 
186  // Now use Importer to do an export
187 
188  Epetra_Vector TargetVector(SourceMap);
189  Epetra_Vector ExpectedTarget(SourceMap);
190  Epetra_Vector SourceVector(TargetMap);
191 
192  NumSameIDs = Importer.NumSameIDs();
193  int NumPermuteIDs = Importer.NumPermuteIDs();
194  int NumExportIDs = Importer.NumExportIDs();
195  int *PermuteFromLIDs = Importer.PermuteFromLIDs();
196  int *ExportLIDs = Importer.ExportLIDs();
197  int *ExportPIDs = Importer.ExportPIDs();
198 
199  for (i=0; i < NumSameIDs; i++) ExpectedTarget[i] = (double) (MyPID+1);
200  for (i=0; i < NumPermuteIDs; i++) ExpectedTarget[PermuteFromLIDs[i]] =
201  (double) (MyPID+1);
202  for (i=0; i < NumExportIDs; i++) ExpectedTarget[ExportLIDs[i]] +=
203  (double) (ExportPIDs[i]+1);
204 
205  for (i=0; i < NumMyElements; i++) SourceVector[i] = (double) (MyPID+1);
206 
207  EPETRA_TEST_ERR(!(TargetVector.Export(SourceVector, Importer, Add)==0),ierr);
208 
209  forierr = 0;
210  for (i=0; i < NumMyElements; i++) {
211  if (TargetVector[i]!= ExpectedTarget[i])
212  cout << " TargetVector["<<i<<"] = " << TargetVector[i]
213  << " ExpectedTarget["<<i<<"] = " << ExpectedTarget[i] << " on PE " << MyPID << endl;
214  forierr += !(TargetVector[i]== ExpectedTarget[i]);
215  }
216  EPETRA_TEST_ERR(forierr,ierr);
217 
218  if (verbose) cout << "Vector Export using Importer Check OK" << endl << endl;
219 
221  // Now use Importer to create a reverse exporter
222  TargetVector.PutScalar(0.0);
223  Epetra_Export ReversedImport(Importer);
224 
225  EPETRA_TEST_ERR(!(TargetVector.Export(SourceVector, ReversedImport, Add)==0),ierr);
226 
227  forierr = 0;
228  for (i=0; i < NumMyElements; i++) {
229  if (TargetVector[i]!= ExpectedTarget[i])
230  cout << " TargetVector["<<i<<"] = " << TargetVector[i]
231  << " ExpectedTarget["<<i<<"] = " << ExpectedTarget[i] << " on PE " << MyPID << endl;
232  forierr += !(TargetVector[i]== ExpectedTarget[i]);
233  }
234  EPETRA_TEST_ERR(forierr,ierr);
235 
236  if (verbose) cout << "Vector Export using Reversed Importer Check OK" << endl << endl;
237 
239  // Now use Exporter to create a reverse importer
240  TargetVector.PutScalar(0.0);
241  Epetra_Import ReversedExport(ReversedImport);
242 
243  EPETRA_TEST_ERR(!(TargetVector.Export(SourceVector, ReversedExport, Add)==0),ierr);
244 
245  forierr = 0;
246  for (i=0; i < NumMyElements; i++) {
247  if (TargetVector[i]!= ExpectedTarget[i])
248  cout << " TargetVector["<<i<<"] = " << TargetVector[i]
249  << " ExpectedTarget["<<i<<"] = " << ExpectedTarget[i] << " on PE " << MyPID << endl;
250  forierr += !(TargetVector[i]== ExpectedTarget[i]);
251  }
252  EPETRA_TEST_ERR(forierr,ierr);
253 
254  if (verbose) cout << "Vector Export using Reversed Exporter Check OK" << endl << endl;
255 
256 
258  // Build a tridiagonal system two ways:
259  // 1) From "standard" matrix view where equations are uniquely owned.
260  // 2) From 1D PDE view where nodes (equations) between processors are shared and partial contributions are done
261  // in parallel, then merged together at the end of the construction process.
262  //
264 
265 
266 
267  // Construct a Standard Map that puts approximately the same number of equations on each processor in
268  // uniform global ordering
269 
270  Epetra_Map StandardMap(NumGlobalEquations, NumMyEquations, 0LL, Comm);
271 
272  // Get update list and number of local equations from newly created Map
273  NumMyElements = StandardMap.NumMyElements();
274  long long * StandardMyGlobalElements = new long long[NumMyElements];
275  StandardMap.MyGlobalElements(StandardMyGlobalElements);
276 
277 
278  // Create a standard Epetra_CrsGraph
279 
280  Epetra_CrsGraph StandardGraph(Copy, StandardMap, 3);
281  EPETRA_TEST_ERR(StandardGraph.IndicesAreGlobal(),ierr);
282  EPETRA_TEST_ERR(StandardGraph.IndicesAreLocal(),ierr);
283 
284  // Add rows one-at-a-time
285  // Need some vectors to help
286  // Off diagonal Values will always be -1
287 
288 
289  long long *Indices = new long long[2];
290  int NumEntries;
291 
292  forierr = 0;
293  for (i=0; i<NumMyEquations; i++)
294  {
295  if (StandardMyGlobalElements[i]==0)
296  {
297  Indices[0] = 1;
298  NumEntries = 1;
299  }
300  else if (StandardMyGlobalElements[i] == NumGlobalEquations-1)
301  {
302  Indices[0] = NumGlobalEquations-2;
303  NumEntries = 1;
304  }
305  else
306  {
307  Indices[0] = StandardMyGlobalElements[i]-1;
308  Indices[1] = StandardMyGlobalElements[i]+1;
309  NumEntries = 2;
310  }
311  forierr += !(StandardGraph.InsertGlobalIndices(StandardMyGlobalElements[i], NumEntries, Indices)==0);
312  forierr += !(StandardGraph.InsertGlobalIndices(StandardMyGlobalElements[i], 1, StandardMyGlobalElements+i)==0); // Put in the diagonal entry
313  }
314  EPETRA_TEST_ERR(forierr,ierr);
315 
316  // Finish up
317  EPETRA_TEST_ERR(!(StandardGraph.IndicesAreGlobal()),ierr);
318  EPETRA_TEST_ERR(!(StandardGraph.FillComplete()==0),ierr);
319  EPETRA_TEST_ERR(!(StandardGraph.IndicesAreLocal()),ierr);
320  EPETRA_TEST_ERR(StandardGraph.StorageOptimized(),ierr);
321  StandardGraph.OptimizeStorage();
322  EPETRA_TEST_ERR(!(StandardGraph.StorageOptimized()),ierr);
323  EPETRA_TEST_ERR(StandardGraph.UpperTriangular(),ierr);
324  EPETRA_TEST_ERR(StandardGraph.LowerTriangular(),ierr);
325 
326  // Create Epetra_CrsMatrix using the just-built graph
327 
328  Epetra_CrsMatrix StandardMatrix(Copy, StandardGraph);
329  EPETRA_TEST_ERR(StandardMatrix.IndicesAreGlobal(),ierr);
330  EPETRA_TEST_ERR(!(StandardMatrix.IndicesAreLocal()),ierr);
331 
332  // Add rows one-at-a-time
333  // Need some vectors to help
334  // Off diagonal Values will always be -1
335 
336 
337  double *Values = new double[2];
338  Values[0] = -1.0; Values[1] = -1.0;
339  double two = 2.0;
340 
341  forierr = 0;
342  for (i=0; i<NumMyEquations; i++)
343  {
344  if (StandardMyGlobalElements[i]==0)
345  {
346  Indices[0] = 1;
347  NumEntries = 1;
348  }
349  else if (StandardMyGlobalElements[i] == NumGlobalEquations-1)
350  {
351  Indices[0] = NumGlobalEquations-2;
352  NumEntries = 1;
353  }
354  else
355  {
356  Indices[0] = StandardMyGlobalElements[i]-1;
357  Indices[1] = StandardMyGlobalElements[i]+1;
358  NumEntries = 2;
359  }
360  forierr += !(StandardMatrix.ReplaceGlobalValues(StandardMyGlobalElements[i], NumEntries, Values, Indices)==0);
361  // Put in the diagonal entry
362  forierr += !(StandardMatrix.ReplaceGlobalValues(StandardMyGlobalElements[i], 1, &two, StandardMyGlobalElements+i)==0);
363  }
364  EPETRA_TEST_ERR(forierr,ierr);
365 
366  // Finish up
367  EPETRA_TEST_ERR(!(StandardMatrix.IndicesAreLocal()),ierr);
368  EPETRA_TEST_ERR(!(StandardMatrix.FillComplete()==0),ierr);
369  EPETRA_TEST_ERR(!(StandardMatrix.IndicesAreLocal()),ierr);
370  // EPETRA_TEST_ERR((StandardMatrix.StorageOptimized()),ierr);
371  EPETRA_TEST_ERR((StandardMatrix.OptimizeStorage()),ierr);
372  EPETRA_TEST_ERR(!(StandardMatrix.StorageOptimized()),ierr);
373  EPETRA_TEST_ERR(StandardMatrix.UpperTriangular(),ierr);
374  EPETRA_TEST_ERR(StandardMatrix.LowerTriangular(),ierr);
375 
376  // Construct an Overlapped Map of StandardMap that include the endpoints from two neighboring processors.
377 
378  int OverlapNumMyElements;
379  long long OverlapMinMyGID;
380 
381  OverlapNumMyElements = NumMyElements + 1;
382  if (MyPID==0) OverlapNumMyElements--;
383 
384  if (MyPID==0) OverlapMinMyGID = StandardMap.MinMyGID64();
385  else OverlapMinMyGID = StandardMap.MinMyGID64()-1;
386 
387  long long * OverlapMyGlobalElements = new long long[OverlapNumMyElements];
388 
389  for (i=0; i< OverlapNumMyElements; i++) OverlapMyGlobalElements[i] = OverlapMinMyGID + i;
390 
391  Epetra_Map OverlapMap((long long) -1, OverlapNumMyElements, OverlapMyGlobalElements, 0LL, Comm);
392 
393  // Create the Overlap Epetra_Matrix
394 
395  Epetra_CrsMatrix OverlapMatrix(Copy, OverlapMap, 4);
396  EPETRA_TEST_ERR(OverlapMatrix.IndicesAreGlobal(),ierr);
397  EPETRA_TEST_ERR(OverlapMatrix.IndicesAreLocal(),ierr);
398 
399  // Add matrix element one cell at a time.
400  // Each cell does an incoming and outgoing flux calculation
401 
402 
403  double pos_one = 1.0;
404  double neg_one = -1.0;
405 
406  forierr = 0;
407  for (i=0; i<OverlapNumMyElements; i++)
408  {
409  long long node_left = OverlapMyGlobalElements[i]-1;
410  long long node_center = node_left + 1;
411  long long node_right = node_left + 2;
412  if (i>0) {
413  if (node_left>-1)
414  forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &neg_one, &node_left)==0);
415  forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0);
416  }
417  if (i<OverlapNumMyElements-1) {
418  forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0);
419  if (node_right<NumGlobalEquations)
420  forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &neg_one, &node_right)==0);
421  }
422  }
423  EPETRA_TEST_ERR(forierr,ierr);
424 
425  // Handle endpoints
426  if (MyPID==0) {
427  long long node_center = 0;
428  EPETRA_TEST_ERR(!(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0),ierr);
429  }
430  if (MyPID==NumProc-1) {
431  long long node_center = OverlapMyGlobalElements[OverlapNumMyElements-1];
432  EPETRA_TEST_ERR(!(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0),ierr);
433  }
434 
435  EPETRA_TEST_ERR(!(OverlapMatrix.FillComplete()==0),ierr);
436 
437  // Make a gathered matrix from OverlapMatrix. It should be identical to StandardMatrix
438 
439  Epetra_CrsMatrix GatheredMatrix(Copy, StandardGraph);
440  Epetra_Export Exporter(OverlapMap, StandardMap);
441  EPETRA_TEST_ERR(!(GatheredMatrix.Export(OverlapMatrix, Exporter, Add)==0),ierr);
442  EPETRA_TEST_ERR(!(GatheredMatrix.FillComplete()==0),ierr);
443 
444  // Check if entries of StandardMatrix and GatheredMatrix are identical
445 
446  int StandardNumEntries, GatheredNumEntries;
447  int * StandardIndices, * GatheredIndices;
448  double * StandardValues, * GatheredValues;
449 
450  int StandardNumMyNonzeros = StandardMatrix.NumMyNonzeros();
451  int GatheredNumMyNonzeros = GatheredMatrix.NumMyNonzeros();
452  EPETRA_TEST_ERR(!(StandardNumMyNonzeros==GatheredNumMyNonzeros),ierr);
453 
454  int StandardNumMyRows = StandardMatrix.NumMyRows();
455  int GatheredNumMyRows = GatheredMatrix.NumMyRows();
456  EPETRA_TEST_ERR(!(StandardNumMyRows==GatheredNumMyRows),ierr);
457 
458  forierr = 0;
459  for (i=0; i< StandardNumMyRows; i++)
460  {
461  forierr += !(StandardMatrix.ExtractMyRowView(i, StandardNumEntries, StandardValues, StandardIndices)==0);
462  forierr += !(GatheredMatrix.ExtractMyRowView(i, GatheredNumEntries, GatheredValues, GatheredIndices)==0);
463  forierr += !(StandardNumEntries==GatheredNumEntries);
464  for (j=0; j < StandardNumEntries; j++) {
465  //if (StandardIndices[j]!=GatheredIndices[j])
466  // cout << "MyPID = " << MyPID << " i = " << i << " StandardIndices[" << j << "] = " << StandardIndices[j]
467  // << " GatheredIndices[" << j << "] = " << GatheredIndices[j] << endl;
468  //if (StandardValues[j]!=GatheredValues[j])
469  //cout << "MyPID = " << MyPID << " i = " << i << " StandardValues[" << j << "] = " << StandardValues[j]
470  // << " GatheredValues[" << j << "] = " << GatheredValues[j] << endl;
471  forierr += !(StandardIndices[j]==GatheredIndices[j]);
472  forierr += !(StandardValues[j]==GatheredValues[j]);
473  }
474  }
475  EPETRA_TEST_ERR(forierr,ierr);
476 
477  if (verbose) cout << "Matrix Export Check OK" << endl << endl;
478 
479  //Do Again with use of Epetra_OffsetIndex object for speed
480  Epetra_OffsetIndex OffsetIndex( OverlapMatrix.Graph(), GatheredMatrix.Graph(), Exporter );
481  EPETRA_TEST_ERR(!(GatheredMatrix.Export(OverlapMatrix, Exporter, Add)==0),ierr);
482 
483  if (verbose) cout << "Optimized Matrix Export Check OK" << endl << endl;
484 
485  bool passed;
486  Epetra_LongLongVector v1(StandardMap); v1.PutValue(2);
487  Epetra_LongLongVector v2(StandardMap); v2.PutValue(3);
488 
489  Epetra_Export identExporter(StandardMap,StandardMap); // Identity exporter
490  EPETRA_TEST_ERR(!(v2.Export(v1, identExporter, Insert)==0),ierr);
491  passed = (v2.MinValue()==2);
492  EPETRA_TEST_ERR(!passed,ierr);
493 
494  v1.PutValue(1);
495  Epetra_Import identImporter(StandardMap,StandardMap); // Identity importer
496  EPETRA_TEST_ERR(!(v2.Import(v1, identExporter, Insert)==0),ierr);
497  passed = passed && (v2.MaxValue()==1);
498  EPETRA_TEST_ERR(!passed,ierr);
499 
500  if (verbose) {
501  if (passed) cout << "Identity Import/Export Check OK" << endl << endl;
502  else cout << "Identity Import/Export Check Failed" << endl << endl;
503  }
504 
505  int NumSubMapElements = StandardMap.NumMyElements()/2;
506  int SubStart = Comm.MyPID();
507  NumSubMapElements = EPETRA_MIN(NumSubMapElements,StandardMap.NumMyElements()-SubStart);
508  Epetra_Map SubMap((long long) -1, NumSubMapElements, StandardMyGlobalElements+SubStart, 0LL, Comm);
509 
510  Epetra_LongLongVector v3(View, SubMap, SubMap.MyGlobalElements64()); // Fill v3 with GID values for variety
511  Epetra_Export subExporter(SubMap, StandardMap); // Export to a subset of indices of standard map
512  EPETRA_TEST_ERR(!(v2.Export(v3,subExporter,Insert)==0),ierr);
513 
514  forierr = 0;
515  for (i=0; i<SubMap.NumMyElements(); i++) {
516  int i1 = StandardMap.LID(SubMap.GID64(i));
517  forierr += !(v3[i]==v2[i1]);
518  }
519  EPETRA_TEST_ERR(forierr,ierr);
520 
521  Epetra_Import subImporter(StandardMap, SubMap); // Import to a subset of indices of standard map
522  EPETRA_TEST_ERR(!(v1.Import(v3,subImporter,Insert)==0),ierr);
523 
524  for (i=0; i<SubMap.NumMyElements(); i++) {
525  int i1 = StandardMap.LID(SubMap.GID64(i));
526  forierr += !(v3[i]==v1[i1]);
527  }
528  EPETRA_TEST_ERR(forierr,ierr);
529 
530  if (verbose) {
531  if (forierr==0) cout << "SubMap Import/Export Check OK" << endl << endl;
532  else cout << "SubMap Import/Export Check Failed" << endl << endl;
533  }
534 
535 #ifdef DOESNT_WORK_IN_PARALLEL
536  forierr = special_submap_import_test(Comm);
537  EPETRA_TEST_ERR(forierr, ierr);
538 
539  if (verbose) {
540  if (forierr==0) cout << "Special SubMap Import Check OK" << endl << endl;
541  else cout << "Special SubMap Import Check Failed" << endl << endl;
542  }
543 #endif
544 
545  forierr = alternate_import_constructor_test(Comm);
546  EPETRA_TEST_ERR(forierr, ierr);
547 
548  if (verbose) {
549  if (forierr==0) cout << "Alternative Import Constructor Check OK" << endl << endl;
550  else cout << "Alternative Import Constructor Check Failed" << endl << endl;
551  }
552 
553  // Release all objects
554 
555  delete [] SourceMyGlobalElements;
556  delete [] TargetMyGlobalElements;
557  delete [] OverlapMyGlobalElements;
558  delete [] StandardMyGlobalElements;
559 
560  delete [] Values;
561  delete [] Indices;
562 
563 #ifdef EPETRA_MPI
564  MPI_Finalize() ;
565 #endif
566 
567 /* end main
568 */
569 return ierr ;
570 }
571 
573 {
574  int localProc = Comm.MyPID();
575 
576  //set up ids_source and ids_target such that ids_source are only
577  //a subset of ids_target, and furthermore that ids_target are ordered
578  //such that the LIDs don't match up. In other words, even if gid 2 does
579  //exist in both ids_source and ids_target, it will correspond to different
580  //LIDs on at least 1 proc.
581  //
582  //This is to test a certain bug-fix in Epetra_Import where the 'RemoteLIDs'
583  //array wasn't being calculated correctly on all procs.
584 
585  long long ids_source[1];
586  ids_source[0] = localProc*2+2;
587 
588  long long ids_target[3];
589  ids_target[0] = localProc*2+2;
590  ids_target[1] = localProc*2+1;
591  ids_target[2] = localProc*2+0;
592 
593  Epetra_Map map_source((long long) -1, 1, &ids_source[0], 0LL, Comm);
594  Epetra_Map map_target((long long) -1, 3, &ids_target[0], 0LL, Comm);
595 
596  Epetra_Import importer(map_target, map_source);
597 
598  Epetra_LongLongVector vec_source(map_source);
599  Epetra_LongLongVector vec_target(map_target);
600 
601  vec_target.PutValue(0);
602 
603  //set vec_source's contents so that entry[i] == GID[i].
604  long long* GIDs = map_source.MyGlobalElements64();
605  for(int i=0; i<map_source.NumMyElements(); ++i) {
606  vec_source[i] = GIDs[i];
607  }
608 
609  //Import vec_source into vec_target. This should result in the contents
610  //of vec_target remaining 0 for the entries that don't exist in vec_source,
611  //and other entries should be equal to the corresponding GID in the map.
612 
613  vec_target.Import(vec_source, importer, Insert);
614 
615  GIDs = map_target.MyGlobalElements64();
616  int test_failed = 0;
617 
618  //the test passes if the i-th entry in vec_target equals either 0 or
619  //GIDs[i].
620  for(int i=0; i<vec_target.MyLength(); ++i) {
621  if (vec_target[i] != GIDs[i] && vec_target[i] != 0) test_failed = 1;
622  }
623 
624  int global_result;
625  Comm.MaxAll(&test_failed, &global_result, 1);
626 
627  //If test didn't fail on any procs, global_result should be 0.
628  //If test failed on any proc, global_result should be 1.
629  return global_result;
630 }
631 
633 {
634  int localProc = Comm.MyPID();
635 
636 
637  long long ids_source[1];
638  ids_source[0] = localProc*2+2;
639 
640  long long ids_target[3];
641  ids_target[0] = localProc*2+2;
642  ids_target[1] = localProc*2+1;
643  ids_target[2] = localProc*2+0;
644 
645  Epetra_Map map_source((long long) -1, 1, &ids_source[0], 0LL, Comm);
646  Epetra_Map map_target((long long) -1, 3, &ids_target[0], 0LL, Comm);
647 
648  Epetra_Import importer(map_target, map_source);
649 
650  Epetra_LongLongVector vec_source(map_source);
651  Epetra_LongLongVector vec_target(map_target);
652 
653  vec_target.PutValue(0);
654 
655  //set vec_source's contents so that entry[i] == GID[i].
656  long long* GIDs = map_source.MyGlobalElements64();
657  for(int i=0; i<map_source.NumMyElements(); ++i) {
658  vec_source[i] = GIDs[i];
659  }
660 
661  //Import vec_source into vec_target. This should result in the contents
662  //of vec_target remaining 0 for the entries that don't exist in vec_source,
663  //and other entries should be equal to the corresponding GID in the map.
664 
665  vec_target.Import(vec_source, importer, Insert);
666 
667  GIDs = map_target.MyGlobalElements64();
668  int test_failed = 0;
669 
670  //the test passes if the i-th entry in vec_target equals either 0 or
671  //GIDs[i].
672  for(int i=0; i<vec_target.MyLength(); ++i) {
673  if (vec_target[i] != GIDs[i] && vec_target[i] != 0) test_failed = 1;
674  }
675 
676  int global_result;
677  Comm.MaxAll(&test_failed, &global_result, 1);
678 
679  //If test didn't fail on any procs, global_result should be 0.
680  //If test failed on any proc, global_result should be 1.
681  return global_result;
682 }
683 
684 
685 int test_import_gid(const char * name,Epetra_LongLongVector & Source, Epetra_LongLongVector & Target, const Epetra_Import & Import){
686  int i;
687  bool test_passed=true;
688 
689  // Setup
690  for(i=0; i<Source.MyLength(); i++)
691  Source[i] = Source.Map().GID64(i);
692  Target.PutValue(0);
693 
694  // Import
695  Target.Import(Source,Import,Add);
696 
697  // Test
698  for(i=0; i<Target.MyLength(); i++){
699  if(Target[i] != Target.Map().GID64(i)) test_passed=false;
700  }
701 
702  if(!test_passed){
703  printf("[%d] test_import_gid %s failed: ",Source.Map().Comm().MyPID(),name);
704  for(i=0; i<Target.MyLength(); i++)
705  printf("%2lld(%2lld) ",Target[i],Target.Map().GID64(i));
706  printf("\n");
707  fflush(stdout);
708  }
709 
710  return !test_passed;
711 }
712 
713 
714 
716  int rv=0;
717  int nodes_per_proc=10;
718  int numprocs = Comm.NumProc();
719  int mypid = Comm.MyPID();
720 
721  // Only run if we have multiple procs & MPI
722  if(numprocs==0) return 0;
723 #ifndef HAVE_MPI
724  return 0;
725 #endif
726 
727  // Build Map 1 - linear
728  Epetra_Map Map1((long long)-1,nodes_per_proc,(long long)0,Comm);
729 
730  // Build Map 2 - mod striped
731  std::vector<long long> MyGIDs(nodes_per_proc);
732  for(int i=0; i<nodes_per_proc; i++)
733  MyGIDs[i] = (mypid*nodes_per_proc + i) % numprocs;
734  Epetra_Map Map2((long long)-1,nodes_per_proc,&MyGIDs[0],(long long)0,Comm);
735 
736  // For testing
737  Epetra_LongLongVector Source(Map1), Target(Map2);
738 
739 
740  // Build Import 1 - normal
741  Epetra_Import Import1(Map2,Map1);
742  rv = rv|| test_import_gid("Alt test: 2 map constructor",Source,Target, Import1);
743 
744  // Build Import 2 - no-comm constructor
745  int Nremote=Import1.NumRemoteIDs();
746  const int * RemoteLIDs = Import1.RemoteLIDs();
747  std::vector<int> RemotePIDs(Nremote+1); // I hate you, stl vector....
748  std::vector<int> AllPIDs;
749  Epetra_Util::GetPids(Import1,AllPIDs,true);
750 
751  for(int i=0; i<Nremote; i++) {
752  RemotePIDs[i]=AllPIDs[RemoteLIDs[i]];
753  }
754  Epetra_Import Import2(Import1.TargetMap(),Import1.SourceMap(),Nremote,&RemotePIDs[0],Import1.NumExportIDs(),Import1.ExportLIDs(),Import1.ExportPIDs());
755 
756  rv = rv || test_import_gid("Alt test: no comm constructor",Source,Target,Import2);
757 
758 
759  // Build Import 3 - Remotes only
760  Epetra_Import Import3(Import1.TargetMap(),Import1.SourceMap(),Nremote,&RemotePIDs[0]);
761  rv = rv || test_import_gid("Alt test: remote only constructor",Source,Target, Import3);
762 
763 
764  return rv;
765 }
long long MinMyGID64() const
Epetra_MultiVector: A class for constructing and using dense multi-vectors, vectors and matrices in p...
int combine_mode_test(Epetra_Comm &Comm)
Epetra_Map: A class for partitioning vectors and matrices.
Definition: Epetra_Map.h:127
int Random()
Set multi-vector values to random numbers.
int FillComplete()
Tranform to local index space. Perform other operations to allow optimal matrix operations.
Epetra_LongLongVector: A class for constructing and using dense integer vectors on a parallel compute...
int NumSameIDs() const
Returns the number of elements that are identical between the source and target maps, up to the first different ID.
Epetra_OffsetIndex: This class builds index for efficient mapping of data from one Epetra_CrsGraph ba...
int special_submap_import_test(Epetra_Comm &Comm)
#define EPETRA_TEST_ERR(a, b)
int MyGlobalElements(int *MyGlobalElementList) const
Puts list of global elements on this processor into the user-provided array.
int PutScalar(double ScalarConstant)
Initialize all values in a multi-vector with constant value.
virtual int InsertGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)
Insert a list of elements in a given global row of the matrix.
static void SetTracebackMode(int TracebackModeValue)
Set the value of the Epetra_Object error traceback report mode.
int * ExportPIDs() const
List of processors to which elements will be sent, ExportLIDs() [i] will be sent to processor ExportP...
long long MinValue()
Find minimum value.
Epetra_Export: This class builds an export object for efficient exporting of off-processor elements...
Definition: Epetra_Export.h:70
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_Import: This class builds an import object for efficient importing of off-processor elements...
Definition: Epetra_Import.h:71
bool IndicesAreLocal() const
If matrix indices has been transformed to local, this query returns true, otherwise it returns false...
static int GetPids(const Epetra_Import &Importer, std::vector< int > &pids, bool use_minus_one_for_local)
Epetra_Util GetPids function.
int * ExportLIDs() const
List of elements that will be sent to other processors.
virtual int MyPID() const =0
Return my process ID.
int FillComplete(bool OptimizeDataStorage=true)
Signal that data entry is complete. Perform transformations to local index space. ...
int MyLength() const
Returns the local vector length on the calling processor of vectors in the multi-vector.
virtual int MaxAll(double *PartialMaxs, double *GlobalMaxs, int Count) const =0
Epetra_Comm Global Max function.
int NumExportIDs() const
Returns the number of elements that must be sent by the calling processor to other processors...
int ExtractMyRowView(int MyRow, int &NumEntries, double *&Values, int *&Indices) const
Returns a view of the specified local row values via pointers to internal data.
int NumMyElements() const
Number of elements on the calling processor.
long long GID64(int LID) const
int Export(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
Exports an Epetra_DistObject using the Epetra_Import object.
int alternate_import_constructor_test(Epetra_Comm &Comm)
int * PermuteFromLIDs() const
List of elements in the source map that are permuted.
Epetra_Comm: The Epetra Communication Abstract Base Class.
Definition: Epetra_Comm.h:81
long long MaxValue()
Find maximum value.
int OptimizeStorage()
Make consecutive row index sections contiguous, minimize internal storage used for constructing graph...
int NumMyRows() const
Returns the number of matrix rows owned by the calling processor.
bool IndicesAreGlobal() const
If column indices are in global range, this query returns true, otherwise it returns false...
int LID(int GID) const
Returns local ID of global ID, return -1 if not found on this processor.
bool MyGID(int GID_in) const
Returns true if the GID passed in belongs to the calling processor in this map, otherwise returns fal...
int NumProc() const
Returns total number of processes (always returns 1 for SerialComm).
int InsertGlobalIndices(int GlobalRow, int NumIndices, int *Indices)
Enter a list of elements in a specified global row of the graph.
int PutValue(long long Value)
Set all elements of the vector to Value.
const Epetra_Comm & Comm() const
Access function for Epetra_Comm communicator.
Epetra_SerialComm: The Epetra Serial Communication Class.
bool IndicesAreGlobal() const
If matrix indices has not been transformed to local, this query returns true, otherwise it returns fa...
int NumPermuteIDs() const
Returns the number of elements that are local to the calling processor, but not part of the first Num...
bool UpperTriangular() const
If graph is upper triangular in local index space, this query returns true, otherwise it returns fals...
bool StorageOptimized() const
If OptimizeStorage() has been called, this query returns true, otherwise it returns false...
virtual int NumProc() const =0
Returns total number of processes.
Epetra_CrsMatrix: A class for constructing and using real-valued double-precision sparse compressed r...
int test_import_gid(const char *name, Epetra_IntVector &Source, Epetra_IntVector &Target, const Epetra_Import &Import)
int main(int argc, char *argv[])
const Epetra_CrsGraph & Graph() const
Returns a reference to the Epetra_CrsGraph object associated with this matrix.
bool LowerTriangular() const
If graph is lower triangular in local index space, this query returns true, otherwise it returns fals...
int Import(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
Imports an Epetra_DistObject using the Epetra_Import object.
#define EPETRA_MAX(x, y)
Epetra_CrsGraph: A class for constructing and using sparse compressed row graphs. ...
const Epetra_BlockMap & Map() const
Returns the address of the Epetra_BlockMap for this multi-vector.
bool IndicesAreLocal() const
If column indices are in local range, this query returns true, otherwise it returns false...
int NumMyNonzeros() const
Returns the number of nonzero entries in the calling processor&#39;s portion of the matrix.