Epetra Package Browser (Single Doxygen Collection) Development
test/ImportExport/cxx_main.cpp
Go to the documentation of this file.
00001 //@HEADER
00002 // ************************************************************************
00003 // 
00004 //               Epetra: Linear Algebra Services Package 
00005 //                 Copyright 2011 Sandia Corporation
00006 // 
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00008 // the U.S. Government retains certain rights in this software.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00038 // 
00039 // ************************************************************************
00040 //@HEADER
00041 
00042 
00043 #include "Epetra_Map.h"
00044 #include "Epetra_Time.h"
00045 #include "Epetra_CrsMatrix.h"
00046 #include "Epetra_Vector.h"
00047 #include "Epetra_IntVector.h"
00048 #include "Epetra_Import.h"
00049 #include "Epetra_Export.h"
00050 #include "Epetra_OffsetIndex.h"
00051 #ifdef EPETRA_MPI
00052 #include "Epetra_MpiComm.h"
00053 #include "mpi.h"
00054 #else
00055 #include "Epetra_SerialComm.h"
00056 #endif
00057 #ifndef __cplusplus
00058 #define __cplusplus
00059 #endif
00060 #include "../epetra_test_err.h"
00061 #include "Epetra_Version.h"
00062 
00063 int special_submap_import_test(Epetra_Comm& Comm);
00064 int combine_mode_test(Epetra_Comm& Comm);
00065 
00066 int main(int argc, char *argv[])
00067 {
00068   int ierr = 0, i, j, forierr = 0;
00069 
00070 #ifdef EPETRA_MPI
00071   // Initialize MPI
00072   MPI_Init(&argc,&argv);
00073   Epetra_MpiComm Comm( MPI_COMM_WORLD );
00074 #else
00075   Epetra_SerialComm Comm;
00076 #endif
00077 
00078   bool verbose = false;
00079 
00080   // Check if we should print results to standard out
00081   if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
00082 
00083 
00084 
00085 
00086   //char tmp;
00087   //if (Comm.MyPID()==0) cout << "Press any key to continue..."<< endl;
00088   //if (Comm.MyPID()==0) cin >> tmp;
00089   //Comm.Barrier();
00090 
00091   Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
00092   int MyPID = Comm.MyPID();
00093   int NumProc = Comm.NumProc();
00094 
00095   if (verbose && MyPID==0)
00096     cout << Epetra_Version() << endl << endl;
00097 
00098   if (verbose) cout << "Processor "<<MyPID<<" of "<< NumProc
00099               << " is alive."<<endl;
00100 
00101   // Redefine verbose to only print on PE 0
00102   if (verbose && Comm.MyPID()!=0) verbose = false;
00103 
00104   int NumMyEquations = 20;
00105   int NumGlobalEquations = NumMyEquations*NumProc+EPETRA_MIN(NumProc,3);
00106   if (MyPID < 3) NumMyEquations++;
00107   // Construct a Source Map that puts approximately the same Number of equations on each processor in 
00108   // uniform global ordering
00109 
00110   Epetra_Map SourceMap(NumGlobalEquations, NumMyEquations, 0, Comm);
00111   
00112   // Get update list and number of local equations from newly created Map
00113   int NumMyElements = SourceMap.NumMyElements();
00114   int * SourceMyGlobalElements = new int[NumMyElements];
00115   SourceMap.MyGlobalElements(SourceMyGlobalElements);
00116 
00117   // Construct a Target Map that will contain:
00118   //  some unchanged elements (relative to the soure map),
00119   //  some permuted elements
00120   //  some off-processor elements
00121   Epetra_Vector RandVec(SourceMap);
00122   RandVec.Random(); // This creates a vector of random numbers between negative one and one.
00123 
00124   int *TargetMyGlobalElements = new int[NumMyElements];
00125 
00126   int MinGID = SourceMap.MinMyGID();
00127   for (i=0; i< NumMyEquations/2; i++) TargetMyGlobalElements[i] = i + MinGID; // Half will be the same...
00128   for (i=NumMyEquations/2; i<NumMyEquations; i++) {
00129     int index = abs((int)(((double) (NumGlobalEquations-1) ) * RandVec[i]));
00130     TargetMyGlobalElements[i] = EPETRA_MIN(NumGlobalEquations-1,EPETRA_MAX(0,index));
00131   }
00132 
00133   int NumSameIDs = 0;
00134   int NumPermutedIDs = 0;
00135   int NumRemoteIDs = 0;
00136   bool StillContiguous = true;
00137   for (i=0; i < NumMyEquations; i++) {
00138     if (SourceMyGlobalElements[i]==TargetMyGlobalElements[i] && StillContiguous)
00139       NumSameIDs++;
00140     else if (SourceMap.MyGID(TargetMyGlobalElements[i])) {
00141       StillContiguous = false;
00142       NumPermutedIDs++;
00143     }
00144     else {
00145       StillContiguous = false;
00146       NumRemoteIDs++;
00147     }
00148   }
00149   EPETRA_TEST_ERR(!(NumMyEquations==NumSameIDs+NumPermutedIDs+NumRemoteIDs),ierr);
00150 
00151   Epetra_Map TargetMap(-1, NumMyElements, TargetMyGlobalElements, 0, Comm);
00152 
00153   // Create a multivector whose elements are GlobalID * (column number +1)
00154 
00155   int NumVectors = 3;
00156   Epetra_MultiVector SourceMultiVector(SourceMap, NumVectors);
00157   for (j=0; j < NumVectors; j++)
00158     for (i=0; i < NumMyElements; i++)
00159       SourceMultiVector[j][i] = (double) SourceMyGlobalElements[i]*(j+1);
00160 
00161   // Create a target multivector that we will fill using an Import
00162 
00163   Epetra_MultiVector TargetMultiVector(TargetMap, NumVectors);
00164 
00165   Epetra_Import Importer(TargetMap, SourceMap);
00166 
00167   EPETRA_TEST_ERR(!(TargetMultiVector.Import(SourceMultiVector, Importer, Insert)==0),ierr);
00168 
00169   // Test Target against expected values
00170   forierr = 0;
00171   for (j=0; j < NumVectors; j++)
00172     for (i=0; i < NumMyElements; i++) {
00173       if (TargetMultiVector[j][i]!= (double) TargetMyGlobalElements[i]*(j+1))
00174   cout << "TargetMultiVector["<<i<<"]["<<j<<"] = " << TargetMultiVector[j][i] 
00175        <<  "  TargetMyGlobalElements[i]*(j+1) = " <<  TargetMyGlobalElements[i]*(j+1) << endl;
00176       forierr += !(TargetMultiVector[j][i]== (double) TargetMyGlobalElements[i]*(j+1));
00177     }
00178   EPETRA_TEST_ERR(forierr,ierr);
00179 
00180   if (verbose) cout << "MultiVector Import using Importer Check OK" << endl << endl;
00181 
00182 
00184 
00185   // Now use Importer to do an export
00186 
00187   Epetra_Vector TargetVector(SourceMap);
00188   Epetra_Vector ExpectedTarget(SourceMap);
00189   Epetra_Vector SourceVector(TargetMap);
00190 
00191   NumSameIDs = Importer.NumSameIDs();
00192   int NumPermuteIDs = Importer.NumPermuteIDs();
00193   int NumExportIDs = Importer.NumExportIDs();
00194   int *PermuteFromLIDs = Importer.PermuteFromLIDs();
00195   int *ExportLIDs = Importer.ExportLIDs();
00196   int *ExportPIDs = Importer.ExportPIDs();
00197 
00198   for (i=0; i < NumSameIDs; i++) ExpectedTarget[i] = (double) (MyPID+1);
00199   for (i=0; i < NumPermuteIDs; i++) ExpectedTarget[PermuteFromLIDs[i]] = 
00200               (double) (MyPID+1);
00201   for (i=0; i < NumExportIDs; i++) ExpectedTarget[ExportLIDs[i]] += 
00202              (double) (ExportPIDs[i]+1);
00203 
00204   for (i=0; i < NumMyElements; i++) SourceVector[i] =  (double) (MyPID+1);
00205 
00206   EPETRA_TEST_ERR(!(TargetVector.Export(SourceVector, Importer, Add)==0),ierr);
00207 
00208   forierr = 0;
00209   for (i=0; i < NumMyElements; i++) {
00210     if (TargetVector[i]!= ExpectedTarget[i])
00211       cout <<  "     TargetVector["<<i<<"] = " << TargetVector[i] 
00212      <<  "   ExpectedTarget["<<i<<"] = " <<  ExpectedTarget[i] << " on PE " << MyPID << endl;
00213     forierr += !(TargetVector[i]== ExpectedTarget[i]);
00214   }
00215   EPETRA_TEST_ERR(forierr,ierr);
00216 
00217   if (verbose) cout << "Vector Export using Importer Check OK" << endl << endl;
00218 
00219 
00220 
00222   //  Build a tridiagonal system two ways:
00223   //  1) From "standard" matrix view where equations are uniquely owned.
00224   //  2) From 1D PDE view where nodes (equations) between processors are shared and partial contributions are done
00225   //     in parallel, then merged together at the end of the construction process.
00226   //
00228 
00229 
00230 
00231   // Construct a Standard Map that puts approximately the same number of equations on each processor in 
00232   // uniform global ordering
00233 
00234   Epetra_Map StandardMap(NumGlobalEquations, NumMyEquations, 0, Comm);
00235   
00236   // Get update list and number of local equations from newly created Map
00237   NumMyElements = StandardMap.NumMyElements();
00238   int * StandardMyGlobalElements = new int[NumMyElements];
00239   StandardMap.MyGlobalElements(StandardMyGlobalElements);
00240 
00241 
00242   // Create a standard Epetra_CrsGraph
00243 
00244   Epetra_CrsGraph StandardGraph(Copy, StandardMap, 3);
00245   EPETRA_TEST_ERR(StandardGraph.IndicesAreGlobal(),ierr);
00246   EPETRA_TEST_ERR(StandardGraph.IndicesAreLocal(),ierr);
00247   
00248   // Add  rows one-at-a-time
00249   // Need some vectors to help
00250   // Off diagonal Values will always be -1
00251 
00252 
00253   int *Indices = new int[2];
00254   int NumEntries;
00255   
00256   forierr = 0;
00257   for (i=0; i<NumMyEquations; i++)
00258     {
00259     if (StandardMyGlobalElements[i]==0)
00260       {
00261   Indices[0] = 1;
00262   NumEntries = 1;
00263       }
00264     else if (StandardMyGlobalElements[i] == NumGlobalEquations-1)
00265       {
00266   Indices[0] = NumGlobalEquations-2;
00267   NumEntries = 1;
00268       }
00269     else
00270       {
00271   Indices[0] = StandardMyGlobalElements[i]-1;
00272   Indices[1] = StandardMyGlobalElements[i]+1;
00273   NumEntries = 2;
00274       }
00275     forierr += !(StandardGraph.InsertGlobalIndices(StandardMyGlobalElements[i], NumEntries, Indices)==0);
00276     forierr += !(StandardGraph.InsertGlobalIndices(StandardMyGlobalElements[i], 1, StandardMyGlobalElements+i)==0); // Put in the diagonal entry
00277     }
00278   EPETRA_TEST_ERR(forierr,ierr);
00279 
00280   // Finish up
00281   EPETRA_TEST_ERR(!(StandardGraph.IndicesAreGlobal()),ierr);
00282   EPETRA_TEST_ERR(!(StandardGraph.FillComplete()==0),ierr);
00283   EPETRA_TEST_ERR(!(StandardGraph.IndicesAreLocal()),ierr);
00284   EPETRA_TEST_ERR(StandardGraph.StorageOptimized(),ierr);
00285   StandardGraph.OptimizeStorage();
00286   EPETRA_TEST_ERR(!(StandardGraph.StorageOptimized()),ierr);
00287   EPETRA_TEST_ERR(StandardGraph.UpperTriangular(),ierr);
00288   EPETRA_TEST_ERR(StandardGraph.LowerTriangular(),ierr);
00289 
00290   // Create Epetra_CrsMatrix using the just-built graph
00291 
00292   Epetra_CrsMatrix StandardMatrix(Copy, StandardGraph);
00293   EPETRA_TEST_ERR(StandardMatrix.IndicesAreGlobal(),ierr);
00294   EPETRA_TEST_ERR(!(StandardMatrix.IndicesAreLocal()),ierr);
00295   
00296   // Add  rows one-at-a-time
00297   // Need some vectors to help
00298   // Off diagonal Values will always be -1
00299 
00300 
00301   double *Values = new double[2];
00302   Values[0] = -1.0; Values[1] = -1.0;
00303   double two = 2.0;
00304   
00305   forierr = 0;
00306   for (i=0; i<NumMyEquations; i++)
00307     {
00308     if (StandardMyGlobalElements[i]==0)
00309       {
00310   Indices[0] = 1;
00311   NumEntries = 1;
00312       }
00313     else if (StandardMyGlobalElements[i] == NumGlobalEquations-1)
00314       {
00315   Indices[0] = NumGlobalEquations-2;
00316   NumEntries = 1;
00317       }
00318     else
00319       {
00320   Indices[0] = StandardMyGlobalElements[i]-1;
00321   Indices[1] = StandardMyGlobalElements[i]+1;
00322   NumEntries = 2;
00323       }
00324     forierr += !(StandardMatrix.ReplaceGlobalValues(StandardMyGlobalElements[i], NumEntries, Values, Indices)==0);
00325     // Put in the diagonal entry
00326     forierr += !(StandardMatrix.ReplaceGlobalValues(StandardMyGlobalElements[i], 1, &two, StandardMyGlobalElements+i)==0); 
00327     }
00328   EPETRA_TEST_ERR(forierr,ierr);
00329 
00330   // Finish up
00331   EPETRA_TEST_ERR(!(StandardMatrix.IndicesAreLocal()),ierr);
00332   EPETRA_TEST_ERR(!(StandardMatrix.FillComplete()==0),ierr);
00333   EPETRA_TEST_ERR(!(StandardMatrix.IndicesAreLocal()),ierr);
00334   EPETRA_TEST_ERR((StandardMatrix.StorageOptimized()),ierr);
00335   EPETRA_TEST_ERR((StandardMatrix.OptimizeStorage()),ierr);
00336   EPETRA_TEST_ERR(!(StandardMatrix.StorageOptimized()),ierr);
00337   EPETRA_TEST_ERR(StandardMatrix.UpperTriangular(),ierr);
00338   EPETRA_TEST_ERR(StandardMatrix.LowerTriangular(),ierr);
00339 
00340   // Construct an Overlapped Map of StandardMap that include the endpoints from two neighboring processors.
00341 
00342   int OverlapNumMyElements;
00343   int OverlapMinMyGID;
00344 
00345   OverlapNumMyElements = NumMyElements + 1;
00346   if (MyPID==0) OverlapNumMyElements--;
00347 
00348   if (MyPID==0) OverlapMinMyGID = StandardMap.MinMyGID();
00349   else OverlapMinMyGID = StandardMap.MinMyGID()-1;
00350 
00351   int * OverlapMyGlobalElements = new int[OverlapNumMyElements];
00352 
00353   for (i=0; i< OverlapNumMyElements; i++) OverlapMyGlobalElements[i] = OverlapMinMyGID + i;
00354 
00355   Epetra_Map OverlapMap(-1, OverlapNumMyElements, OverlapMyGlobalElements, 0, Comm);
00356 
00357   // Create the Overlap Epetra_Matrix
00358 
00359   Epetra_CrsMatrix OverlapMatrix(Copy, OverlapMap, 4);
00360   EPETRA_TEST_ERR(OverlapMatrix.IndicesAreGlobal(),ierr);
00361   EPETRA_TEST_ERR(OverlapMatrix.IndicesAreLocal(),ierr);
00362   
00363   // Add  matrix element one cell at a time.
00364   // Each cell does an incoming and outgoing flux calculation
00365 
00366 
00367   double pos_one = 1.0;
00368   double neg_one = -1.0;
00369 
00370   forierr = 0;
00371   for (i=0; i<OverlapNumMyElements; i++)
00372     {
00373       int node_left = OverlapMyGlobalElements[i]-1;
00374       int node_center = node_left + 1;
00375       int node_right = node_left + 2;
00376       if (i>0) {
00377   if (node_left>-1)
00378     forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &neg_one, &node_left)==0);
00379   forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0);
00380       }
00381       if (i<OverlapNumMyElements-1) {
00382   forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0);
00383   if (node_right<NumGlobalEquations) 
00384     forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &neg_one, &node_right)==0);
00385       }
00386     }
00387   EPETRA_TEST_ERR(forierr,ierr);
00388 
00389   // Handle endpoints
00390   if (MyPID==0) {
00391     int node_center = 0;
00392     EPETRA_TEST_ERR(!(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0),ierr);
00393   }
00394   if (MyPID==NumProc-1) {
00395     int node_center = OverlapMyGlobalElements[OverlapNumMyElements-1];
00396     EPETRA_TEST_ERR(!(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0),ierr);
00397   }
00398     
00399   EPETRA_TEST_ERR(!(OverlapMatrix.FillComplete()==0),ierr);
00400 
00401   // Make a gathered matrix from OverlapMatrix.  It should be identical to StandardMatrix
00402 
00403   Epetra_CrsMatrix GatheredMatrix(Copy, StandardGraph);
00404   Epetra_Export Exporter(OverlapMap, StandardMap);
00405   EPETRA_TEST_ERR(!(GatheredMatrix.Export(OverlapMatrix, Exporter, Add)==0),ierr);
00406   EPETRA_TEST_ERR(!(GatheredMatrix.FillComplete()==0),ierr);
00407 
00408   // Check if entries of StandardMatrix and GatheredMatrix are identical
00409 
00410   int StandardNumEntries, GatheredNumEntries;
00411   int * StandardIndices, * GatheredIndices;
00412   double * StandardValues, * GatheredValues;
00413 
00414   int StandardNumMyNonzeros = StandardMatrix.NumMyNonzeros();
00415   int GatheredNumMyNonzeros = GatheredMatrix.NumMyNonzeros();
00416   EPETRA_TEST_ERR(!(StandardNumMyNonzeros==GatheredNumMyNonzeros),ierr);
00417 
00418   int StandardNumMyRows = StandardMatrix.NumMyRows();
00419   int GatheredNumMyRows = GatheredMatrix.NumMyRows();
00420   EPETRA_TEST_ERR(!(StandardNumMyRows==GatheredNumMyRows),ierr);
00421 
00422   forierr = 0;
00423   for (i=0; i< StandardNumMyRows; i++)
00424     {
00425       forierr += !(StandardMatrix.ExtractMyRowView(i, StandardNumEntries, StandardValues, StandardIndices)==0);
00426       forierr += !(GatheredMatrix.ExtractMyRowView(i, GatheredNumEntries, GatheredValues, GatheredIndices)==0);
00427       forierr += !(StandardNumEntries==GatheredNumEntries);
00428       for (j=0; j < StandardNumEntries; j++) {
00429   //if (StandardIndices[j]!=GatheredIndices[j])
00430   // cout << "MyPID = " << MyPID << " i = " << i << "   StandardIndices[" << j << "] = " << StandardIndices[j] 
00431   //      << "   GatheredIndices[" << j << "] = " << GatheredIndices[j] << endl;
00432   //if (StandardValues[j]!=GatheredValues[j])
00433   //cout << "MyPID = " << MyPID << " i = " << i << "    StandardValues[" << j << "] = " <<  StandardValues[j] 
00434   //     << "    GatheredValues[" << j << "] = " <<  GatheredValues[j] << endl;
00435   forierr += !(StandardIndices[j]==GatheredIndices[j]);
00436   forierr += !(StandardValues[j]==GatheredValues[j]);
00437       }
00438     }
00439   EPETRA_TEST_ERR(forierr,ierr);
00440 
00441   if (verbose) cout << "Matrix Export Check OK" << endl << endl;
00442 
00443   //Do Again with use of Epetra_OffsetIndex object for speed
00444   Epetra_OffsetIndex OffsetIndex( OverlapMatrix.Graph(), GatheredMatrix.Graph(), Exporter );
00445   EPETRA_TEST_ERR(!(GatheredMatrix.Export(OverlapMatrix, Exporter, Add)==0),ierr);
00446 
00447   if (verbose) cout << "Optimized Matrix Export Check OK" << endl << endl;
00448 
00449   bool passed;
00450   Epetra_IntVector v1(StandardMap); v1.PutValue(2);
00451   Epetra_IntVector v2(StandardMap); v2.PutValue(3);
00452 
00453   Epetra_Export identExporter(StandardMap,StandardMap); // Identity exporter
00454   EPETRA_TEST_ERR(!(v2.Export(v1, identExporter, Insert)==0),ierr);
00455   passed = (v2.MinValue()==2);
00456   EPETRA_TEST_ERR(!passed,ierr);
00457 
00458   v1.PutValue(1);
00459   Epetra_Import identImporter(StandardMap,StandardMap); // Identity importer
00460   EPETRA_TEST_ERR(!(v2.Import(v1, identExporter, Insert)==0),ierr);
00461   passed = passed && (v2.MaxValue()==1);
00462   EPETRA_TEST_ERR(!passed,ierr);
00463 
00464   if (verbose) {
00465     if (passed) cout << "Identity Import/Export Check OK" << endl << endl;
00466     else cout << "Identity Import/Export Check Failed" << endl << endl;
00467   }
00468 
00469   int NumSubMapElements = StandardMap.NumMyElements()/2;
00470   int SubStart = Comm.MyPID();
00471   NumSubMapElements = EPETRA_MIN(NumSubMapElements,StandardMap.NumMyElements()-SubStart);
00472   Epetra_Map SubMap(-1, NumSubMapElements, StandardMyGlobalElements+SubStart, 0, Comm);
00473 
00474   Epetra_IntVector v3(View, SubMap, SubMap.MyGlobalElements()); // Fill v3 with GID values for variety
00475   Epetra_Export subExporter(SubMap, StandardMap); // Export to a subset of indices of standard map
00476   EPETRA_TEST_ERR(!(v2.Export(v3,subExporter,Insert)==0),ierr);
00477 
00478   forierr = 0;
00479   for (i=0; i<SubMap.NumMyElements(); i++) {
00480     int i1 = StandardMap.LID(SubMap.GID(i));
00481     forierr += !(v3[i]==v2[i1]);
00482   }
00483   EPETRA_TEST_ERR(forierr,ierr);
00484 
00485   Epetra_Import subImporter(StandardMap, SubMap); // Import to a subset of indices of standard map
00486   EPETRA_TEST_ERR(!(v1.Import(v3,subImporter,Insert)==0),ierr);
00487 
00488   for (i=0; i<SubMap.NumMyElements(); i++) {
00489     int i1 = StandardMap.LID(SubMap.GID(i));
00490     forierr += !(v3[i]==v1[i1]);
00491   }
00492   EPETRA_TEST_ERR(forierr,ierr);
00493 
00494   if (verbose) {
00495     if (forierr==0) cout << "SubMap Import/Export Check OK" << endl << endl;
00496     else cout << "SubMap Import/Export Check Failed" << endl << endl;
00497   }
00498 
00499   forierr = special_submap_import_test(Comm);
00500   EPETRA_TEST_ERR(forierr, ierr);
00501 
00502   if (verbose) {
00503     if (forierr==0) cout << "Special SubMap Import Check OK" << endl << endl;
00504     else cout << "Special SubMap Import Check Failed" << endl << endl;
00505   }
00506 
00507   // Release all objects
00508 
00509   delete [] SourceMyGlobalElements;
00510   delete [] TargetMyGlobalElements;
00511   delete [] OverlapMyGlobalElements;
00512   delete [] StandardMyGlobalElements;
00513 
00514   delete [] Values;
00515   delete [] Indices;
00516 
00517 #ifdef EPETRA_MPI
00518   MPI_Finalize() ;
00519 #endif
00520 
00521 /* end main
00522 */
00523 return ierr ;
00524 }
00525 
00526 int special_submap_import_test(Epetra_Comm& Comm)
00527 {
00528   int localProc = Comm.MyPID();
00529 
00530   //set up ids_source and ids_target such that ids_source are only
00531   //a subset of ids_target, and furthermore that ids_target are ordered
00532   //such that the LIDs don't match up. In other words, even if gid 2 does
00533   //exist in both ids_source and ids_target, it will correspond to different
00534   //LIDs on at least 1 proc.
00535   //
00536   //This is to test a certain bug-fix in Epetra_Import where the 'RemoteLIDs'
00537   //array wasn't being calculated correctly on all procs.
00538 
00539   int ids_source[1];
00540   ids_source[0] = localProc*2+2;
00541 
00542   int ids_target[3];
00543   ids_target[0] = localProc*2+2;
00544   ids_target[1] = localProc*2+1;
00545   ids_target[2] = localProc*2+0;
00546 
00547   Epetra_Map map_source(-1, 1, &ids_source[0], 0, Comm);
00548   Epetra_Map map_target(-1, 3, &ids_target[0], 0, Comm);
00549 
00550   Epetra_Import importer(map_target, map_source);
00551 
00552   Epetra_IntVector vec_source(map_source);
00553   Epetra_IntVector vec_target(map_target);
00554 
00555   vec_target.PutValue(0);
00556 
00557   //set vec_source's contents so that entry[i] == GID[i].
00558   int* GIDs = map_source.MyGlobalElements();
00559   for(int i=0; i<map_source.NumMyElements(); ++i) {
00560     vec_source[i] = GIDs[i];
00561   }
00562 
00563   //Import vec_source into vec_target. This should result in the contents
00564   //of vec_target remaining 0 for the entries that don't exist in vec_source,
00565   //and other entries should be equal to the corresponding GID in the map.
00566 
00567   vec_target.Import(vec_source, importer, Insert);
00568 
00569   GIDs = map_target.MyGlobalElements();
00570   int test_failed = 0;
00571 
00572   //the test passes if the i-th entry in vec_target equals either 0 or
00573   //GIDs[i].
00574   for(int i=0; i<vec_target.MyLength(); ++i) {
00575     if (vec_target[i] != GIDs[i] && vec_target[i] != 0) test_failed = 1;
00576   }
00577 
00578   int global_result;
00579   Comm.MaxAll(&test_failed, &global_result, 1);
00580 
00581   //If test didn't fail on any procs, global_result should be 0.
00582   //If test failed on any proc, global_result should be 1.
00583   return global_result;
00584 }
00585 
00586 int combine_mode_test(Epetra_Comm& Comm)
00587 {
00588   int localProc = Comm.MyPID();
00589 
00590 
00591   int ids_source[1];
00592   ids_source[0] = localProc*2+2;
00593 
00594   int ids_target[3];
00595   ids_target[0] = localProc*2+2;
00596   ids_target[1] = localProc*2+1;
00597   ids_target[2] = localProc*2+0;
00598 
00599   Epetra_Map map_source(-1, 1, &ids_source[0], 0, Comm);
00600   Epetra_Map map_target(-1, 3, &ids_target[0], 0, Comm);
00601 
00602   Epetra_Import importer(map_target, map_source);
00603 
00604   Epetra_IntVector vec_source(map_source);
00605   Epetra_IntVector vec_target(map_target);
00606 
00607   vec_target.PutValue(0);
00608 
00609   //set vec_source's contents so that entry[i] == GID[i].
00610   int* GIDs = map_source.MyGlobalElements();
00611   for(int i=0; i<map_source.NumMyElements(); ++i) {
00612     vec_source[i] = GIDs[i];
00613   }
00614 
00615   //Import vec_source into vec_target. This should result in the contents
00616   //of vec_target remaining 0 for the entries that don't exist in vec_source,
00617   //and other entries should be equal to the corresponding GID in the map.
00618 
00619   vec_target.Import(vec_source, importer, Insert);
00620 
00621   GIDs = map_target.MyGlobalElements();
00622   int test_failed = 0;
00623 
00624   //the test passes if the i-th entry in vec_target equals either 0 or
00625   //GIDs[i].
00626   for(int i=0; i<vec_target.MyLength(); ++i) {
00627     if (vec_target[i] != GIDs[i] && vec_target[i] != 0) test_failed = 1;
00628   }
00629 
00630   int global_result;
00631   Comm.MaxAll(&test_failed, &global_result, 1);
00632 
00633   //If test didn't fail on any procs, global_result should be 0.
00634   //If test failed on any proc, global_result should be 1.
00635   return global_result;
00636 }
00637 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines