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 2001 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 
00065 int main(int argc, char *argv[])
00066 {
00067   int ierr = 0, i, j, forierr = 0;
00068 
00069 #ifdef EPETRA_MPI
00070   // Initialize MPI
00071   MPI_Init(&argc,&argv);
00072   Epetra_MpiComm Comm( MPI_COMM_WORLD );
00073 #else
00074   Epetra_SerialComm Comm;
00075 #endif
00076 
00077   bool verbose = false;
00078 
00079   // Check if we should print results to standard out
00080   if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
00081 
00082 
00083 
00084 
00085   //char tmp;
00086   //if (Comm.MyPID()==0) cout << "Press any key to continue..."<< endl;
00087   //if (Comm.MyPID()==0) cin >> tmp;
00088   //Comm.Barrier();
00089 
00090   Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
00091   int MyPID = Comm.MyPID();
00092   int NumProc = Comm.NumProc();
00093 
00094   if (verbose && MyPID==0)
00095     cout << Epetra_Version() << endl << endl;
00096 
00097   if (verbose) cout << "Processor "<<MyPID<<" of "<< NumProc
00098               << " is alive."<<endl;
00099 
00100   // Redefine verbose to only print on PE 0
00101   if (verbose && Comm.MyPID()!=0) verbose = false;
00102 
00103   int NumMyEquations = 20;
00104   int NumGlobalEquations = NumMyEquations*NumProc+EPETRA_MIN(NumProc,3);
00105   if (MyPID < 3) NumMyEquations++;
00106   // Construct a Source Map that puts approximately the same Number of equations on each processor in 
00107   // uniform global ordering
00108 
00109   Epetra_Map SourceMap(NumGlobalEquations, NumMyEquations, 0, Comm);
00110   
00111   // Get update list and number of local equations from newly created Map
00112   int NumMyElements = SourceMap.NumMyElements();
00113   int * SourceMyGlobalElements = new int[NumMyElements];
00114   SourceMap.MyGlobalElements(SourceMyGlobalElements);
00115 
00116   // Construct a Target Map that will contain:
00117   //  some unchanged elements (relative to the soure map),
00118   //  some permuted elements
00119   //  some off-processor elements
00120   Epetra_Vector RandVec(SourceMap);
00121   RandVec.Random(); // This creates a vector of random numbers between negative one and one.
00122 
00123   int *TargetMyGlobalElements = new int[NumMyElements];
00124 
00125   int MinGID = SourceMap.MinMyGID();
00126   for (i=0; i< NumMyEquations/2; i++) TargetMyGlobalElements[i] = i + MinGID; // Half will be the same...
00127   for (i=NumMyEquations/2; i<NumMyEquations; i++) {
00128     int index = abs((int)(((double) (NumGlobalEquations-1) ) * RandVec[i]));
00129     TargetMyGlobalElements[i] = EPETRA_MIN(NumGlobalEquations-1,EPETRA_MAX(0,index));
00130   }
00131 
00132   int NumSameIDs = 0;
00133   int NumPermutedIDs = 0;
00134   int NumRemoteIDs = 0;
00135   bool StillContiguous = true;
00136   for (i=0; i < NumMyEquations; i++) {
00137     if (SourceMyGlobalElements[i]==TargetMyGlobalElements[i] && StillContiguous)
00138       NumSameIDs++;
00139     else if (SourceMap.MyGID(TargetMyGlobalElements[i])) {
00140       StillContiguous = false;
00141       NumPermutedIDs++;
00142     }
00143     else {
00144       StillContiguous = false;
00145       NumRemoteIDs++;
00146     }
00147   }
00148   EPETRA_TEST_ERR(!(NumMyEquations==NumSameIDs+NumPermutedIDs+NumRemoteIDs),ierr);
00149 
00150   Epetra_Map TargetMap(-1, NumMyElements, TargetMyGlobalElements, 0, Comm);
00151 
00152   // Create a multivector whose elements are GlobalID * (column number +1)
00153 
00154   int NumVectors = 3;
00155   Epetra_MultiVector SourceMultiVector(SourceMap, NumVectors);
00156   for (j=0; j < NumVectors; j++)
00157     for (i=0; i < NumMyElements; i++)
00158       SourceMultiVector[j][i] = (double) SourceMyGlobalElements[i]*(j+1);
00159 
00160   // Create a target multivector that we will fill using an Import
00161 
00162   Epetra_MultiVector TargetMultiVector(TargetMap, NumVectors);
00163 
00164   Epetra_Import Importer(TargetMap, SourceMap);
00165 
00166   EPETRA_TEST_ERR(!(TargetMultiVector.Import(SourceMultiVector, Importer, Insert)==0),ierr);
00167 
00168   // Test Target against expected values
00169   forierr = 0;
00170   for (j=0; j < NumVectors; j++)
00171     for (i=0; i < NumMyElements; i++) {
00172       if (TargetMultiVector[j][i]!= (double) TargetMyGlobalElements[i]*(j+1))
00173   cout << "TargetMultiVector["<<i<<"]["<<j<<"] = " << TargetMultiVector[j][i] 
00174        <<  "  TargetMyGlobalElements[i]*(j+1) = " <<  TargetMyGlobalElements[i]*(j+1) << endl;
00175       forierr += !(TargetMultiVector[j][i]== (double) TargetMyGlobalElements[i]*(j+1));
00176     }
00177   EPETRA_TEST_ERR(forierr,ierr);
00178 
00179   if (verbose) cout << "MultiVector Import using Importer Check OK" << endl << endl;
00180 
00181 
00183 
00184   // Now use Importer to do an export
00185 
00186   Epetra_Vector TargetVector(SourceMap);
00187   Epetra_Vector ExpectedTarget(SourceMap);
00188   Epetra_Vector SourceVector(TargetMap);
00189 
00190   NumSameIDs = Importer.NumSameIDs();
00191   int NumPermuteIDs = Importer.NumPermuteIDs();
00192   int NumExportIDs = Importer.NumExportIDs();
00193   int *PermuteFromLIDs = Importer.PermuteFromLIDs();
00194   int *ExportLIDs = Importer.ExportLIDs();
00195   int *ExportPIDs = Importer.ExportPIDs();
00196 
00197   for (i=0; i < NumSameIDs; i++) ExpectedTarget[i] = (double) (MyPID+1);
00198   for (i=0; i < NumPermuteIDs; i++) ExpectedTarget[PermuteFromLIDs[i]] = 
00199               (double) (MyPID+1);
00200   for (i=0; i < NumExportIDs; i++) ExpectedTarget[ExportLIDs[i]] += 
00201              (double) (ExportPIDs[i]+1);
00202 
00203   for (i=0; i < NumMyElements; i++) SourceVector[i] =  (double) (MyPID+1);
00204 
00205   EPETRA_TEST_ERR(!(TargetVector.Export(SourceVector, Importer, Add)==0),ierr);
00206 
00207   forierr = 0;
00208   for (i=0; i < NumMyElements; i++) {
00209     if (TargetVector[i]!= ExpectedTarget[i])
00210       cout <<  "     TargetVector["<<i<<"] = " << TargetVector[i] 
00211      <<  "   ExpectedTarget["<<i<<"] = " <<  ExpectedTarget[i] << " on PE " << MyPID << endl;
00212     forierr += !(TargetVector[i]== ExpectedTarget[i]);
00213   }
00214   EPETRA_TEST_ERR(forierr,ierr);
00215 
00216   if (verbose) cout << "Vector Export using Importer Check OK" << endl << endl;
00217 
00218 
00219 
00221   //  Build a tridiagonal system two ways:
00222   //  1) From "standard" matrix view where equations are uniquely owned.
00223   //  2) From 1D PDE view where nodes (equations) between processors are shared and partial contributions are done
00224   //     in parallel, then merged together at the end of the construction process.
00225   //
00227 
00228 
00229 
00230   // Construct a Standard Map that puts approximately the same number of equations on each processor in 
00231   // uniform global ordering
00232 
00233   Epetra_Map StandardMap(NumGlobalEquations, NumMyEquations, 0, Comm);
00234   
00235   // Get update list and number of local equations from newly created Map
00236   NumMyElements = StandardMap.NumMyElements();
00237   int * StandardMyGlobalElements = new int[NumMyElements];
00238   StandardMap.MyGlobalElements(StandardMyGlobalElements);
00239 
00240 
00241   // Create a standard Epetra_CrsGraph
00242 
00243   Epetra_CrsGraph StandardGraph(Copy, StandardMap, 3);
00244   EPETRA_TEST_ERR(StandardGraph.IndicesAreGlobal(),ierr);
00245   EPETRA_TEST_ERR(StandardGraph.IndicesAreLocal(),ierr);
00246   
00247   // Add  rows one-at-a-time
00248   // Need some vectors to help
00249   // Off diagonal Values will always be -1
00250 
00251 
00252   int *Indices = new int[2];
00253   int NumEntries;
00254   
00255   forierr = 0;
00256   for (i=0; i<NumMyEquations; i++)
00257     {
00258     if (StandardMyGlobalElements[i]==0)
00259       {
00260   Indices[0] = 1;
00261   NumEntries = 1;
00262       }
00263     else if (StandardMyGlobalElements[i] == NumGlobalEquations-1)
00264       {
00265   Indices[0] = NumGlobalEquations-2;
00266   NumEntries = 1;
00267       }
00268     else
00269       {
00270   Indices[0] = StandardMyGlobalElements[i]-1;
00271   Indices[1] = StandardMyGlobalElements[i]+1;
00272   NumEntries = 2;
00273       }
00274     forierr += !(StandardGraph.InsertGlobalIndices(StandardMyGlobalElements[i], NumEntries, Indices)==0);
00275     forierr += !(StandardGraph.InsertGlobalIndices(StandardMyGlobalElements[i], 1, StandardMyGlobalElements+i)==0); // Put in the diagonal entry
00276     }
00277   EPETRA_TEST_ERR(forierr,ierr);
00278 
00279   // Finish up
00280   EPETRA_TEST_ERR(!(StandardGraph.IndicesAreGlobal()),ierr);
00281   EPETRA_TEST_ERR(!(StandardGraph.FillComplete()==0),ierr);
00282   EPETRA_TEST_ERR(!(StandardGraph.IndicesAreLocal()),ierr);
00283   EPETRA_TEST_ERR(StandardGraph.StorageOptimized(),ierr);
00284   StandardGraph.OptimizeStorage();
00285   EPETRA_TEST_ERR(!(StandardGraph.StorageOptimized()),ierr);
00286   EPETRA_TEST_ERR(StandardGraph.UpperTriangular(),ierr);
00287   EPETRA_TEST_ERR(StandardGraph.LowerTriangular(),ierr);
00288 
00289   // Create Epetra_CrsMatrix using the just-built graph
00290 
00291   Epetra_CrsMatrix StandardMatrix(Copy, StandardGraph);
00292   EPETRA_TEST_ERR(StandardMatrix.IndicesAreGlobal(),ierr);
00293   EPETRA_TEST_ERR(!(StandardMatrix.IndicesAreLocal()),ierr);
00294   
00295   // Add  rows one-at-a-time
00296   // Need some vectors to help
00297   // Off diagonal Values will always be -1
00298 
00299 
00300   double *Values = new double[2];
00301   Values[0] = -1.0; Values[1] = -1.0;
00302   double two = 2.0;
00303   
00304   forierr = 0;
00305   for (i=0; i<NumMyEquations; i++)
00306     {
00307     if (StandardMyGlobalElements[i]==0)
00308       {
00309   Indices[0] = 1;
00310   NumEntries = 1;
00311       }
00312     else if (StandardMyGlobalElements[i] == NumGlobalEquations-1)
00313       {
00314   Indices[0] = NumGlobalEquations-2;
00315   NumEntries = 1;
00316       }
00317     else
00318       {
00319   Indices[0] = StandardMyGlobalElements[i]-1;
00320   Indices[1] = StandardMyGlobalElements[i]+1;
00321   NumEntries = 2;
00322       }
00323     forierr += !(StandardMatrix.ReplaceGlobalValues(StandardMyGlobalElements[i], NumEntries, Values, Indices)==0);
00324     // Put in the diagonal entry
00325     forierr += !(StandardMatrix.ReplaceGlobalValues(StandardMyGlobalElements[i], 1, &two, StandardMyGlobalElements+i)==0); 
00326     }
00327   EPETRA_TEST_ERR(forierr,ierr);
00328 
00329   // Finish up
00330   EPETRA_TEST_ERR(!(StandardMatrix.IndicesAreLocal()),ierr);
00331   EPETRA_TEST_ERR(!(StandardMatrix.FillComplete()==0),ierr);
00332   EPETRA_TEST_ERR(!(StandardMatrix.IndicesAreLocal()),ierr);
00333   EPETRA_TEST_ERR((StandardMatrix.StorageOptimized()),ierr);
00334   EPETRA_TEST_ERR((StandardMatrix.OptimizeStorage()),ierr);
00335   EPETRA_TEST_ERR(!(StandardMatrix.StorageOptimized()),ierr);
00336   EPETRA_TEST_ERR(StandardMatrix.UpperTriangular(),ierr);
00337   EPETRA_TEST_ERR(StandardMatrix.LowerTriangular(),ierr);
00338 
00339   // Construct an Overlapped Map of StandardMap that include the endpoints from two neighboring processors.
00340 
00341   int OverlapNumMyElements;
00342   int OverlapMinMyGID;
00343 
00344   OverlapNumMyElements = NumMyElements + 1;
00345   if (MyPID==0) OverlapNumMyElements--;
00346 
00347   if (MyPID==0) OverlapMinMyGID = StandardMap.MinMyGID();
00348   else OverlapMinMyGID = StandardMap.MinMyGID()-1;
00349 
00350   int * OverlapMyGlobalElements = new int[OverlapNumMyElements];
00351 
00352   for (i=0; i< OverlapNumMyElements; i++) OverlapMyGlobalElements[i] = OverlapMinMyGID + i;
00353 
00354   Epetra_Map OverlapMap(-1, OverlapNumMyElements, OverlapMyGlobalElements, 0, Comm);
00355 
00356   // Create the Overlap Epetra_Matrix
00357 
00358   Epetra_CrsMatrix OverlapMatrix(Copy, OverlapMap, 4);
00359   EPETRA_TEST_ERR(OverlapMatrix.IndicesAreGlobal(),ierr);
00360   EPETRA_TEST_ERR(OverlapMatrix.IndicesAreLocal(),ierr);
00361   
00362   // Add  matrix element one cell at a time.
00363   // Each cell does an incoming and outgoing flux calculation
00364 
00365 
00366   double pos_one = 1.0;
00367   double neg_one = -1.0;
00368 
00369   forierr = 0;
00370   for (i=0; i<OverlapNumMyElements; i++)
00371     {
00372       int node_left = OverlapMyGlobalElements[i]-1;
00373       int node_center = node_left + 1;
00374       int node_right = node_left + 2;
00375       if (i>0) {
00376   if (node_left>-1)
00377     forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &neg_one, &node_left)==0);
00378   forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0);
00379       }
00380       if (i<OverlapNumMyElements-1) {
00381   forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0);
00382   if (node_right<NumGlobalEquations) 
00383     forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &neg_one, &node_right)==0);
00384       }
00385     }
00386   EPETRA_TEST_ERR(forierr,ierr);
00387 
00388   // Handle endpoints
00389   if (MyPID==0) {
00390     int node_center = 0;
00391     EPETRA_TEST_ERR(!(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0),ierr);
00392   }
00393   if (MyPID==NumProc-1) {
00394     int node_center = OverlapMyGlobalElements[OverlapNumMyElements-1];
00395     EPETRA_TEST_ERR(!(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0),ierr);
00396   }
00397     
00398   EPETRA_TEST_ERR(!(OverlapMatrix.FillComplete()==0),ierr);
00399 
00400   // Make a gathered matrix from OverlapMatrix.  It should be identical to StandardMatrix
00401 
00402   Epetra_CrsMatrix GatheredMatrix(Copy, StandardGraph);
00403   Epetra_Export Exporter(OverlapMap, StandardMap);
00404   EPETRA_TEST_ERR(!(GatheredMatrix.Export(OverlapMatrix, Exporter, Add)==0),ierr);
00405   EPETRA_TEST_ERR(!(GatheredMatrix.FillComplete()==0),ierr);
00406 
00407   // Check if entries of StandardMatrix and GatheredMatrix are identical
00408 
00409   int StandardNumEntries, GatheredNumEntries;
00410   int * StandardIndices, * GatheredIndices;
00411   double * StandardValues, * GatheredValues;
00412 
00413   int StandardNumMyNonzeros = StandardMatrix.NumMyNonzeros();
00414   int GatheredNumMyNonzeros = GatheredMatrix.NumMyNonzeros();
00415   EPETRA_TEST_ERR(!(StandardNumMyNonzeros==GatheredNumMyNonzeros),ierr);
00416 
00417   int StandardNumMyRows = StandardMatrix.NumMyRows();
00418   int GatheredNumMyRows = GatheredMatrix.NumMyRows();
00419   EPETRA_TEST_ERR(!(StandardNumMyRows==GatheredNumMyRows),ierr);
00420 
00421   forierr = 0;
00422   for (i=0; i< StandardNumMyRows; i++)
00423     {
00424       forierr += !(StandardMatrix.ExtractMyRowView(i, StandardNumEntries, StandardValues, StandardIndices)==0);
00425       forierr += !(GatheredMatrix.ExtractMyRowView(i, GatheredNumEntries, GatheredValues, GatheredIndices)==0);
00426       forierr += !(StandardNumEntries==GatheredNumEntries);
00427       for (j=0; j < StandardNumEntries; j++) {
00428   //if (StandardIndices[j]!=GatheredIndices[j])
00429   // cout << "MyPID = " << MyPID << " i = " << i << "   StandardIndices[" << j << "] = " << StandardIndices[j] 
00430   //      << "   GatheredIndices[" << j << "] = " << GatheredIndices[j] << endl;
00431   //if (StandardValues[j]!=GatheredValues[j])
00432   //cout << "MyPID = " << MyPID << " i = " << i << "    StandardValues[" << j << "] = " <<  StandardValues[j] 
00433   //     << "    GatheredValues[" << j << "] = " <<  GatheredValues[j] << endl;
00434   forierr += !(StandardIndices[j]==GatheredIndices[j]);
00435   forierr += !(StandardValues[j]==GatheredValues[j]);
00436       }
00437     }
00438   EPETRA_TEST_ERR(forierr,ierr);
00439 
00440   if (verbose) cout << "Matrix Export Check OK" << endl << endl;
00441 
00442   //Do Again with use of Epetra_OffsetIndex object for speed
00443   Epetra_OffsetIndex OffsetIndex( OverlapMatrix.Graph(), GatheredMatrix.Graph(), Exporter );
00444   EPETRA_TEST_ERR(!(GatheredMatrix.Export(OverlapMatrix, Exporter, Add)==0),ierr);
00445 
00446   if (verbose) cout << "Optimized Matrix Export Check OK" << endl << endl;
00447 
00448   bool passed;
00449   Epetra_IntVector v1(StandardMap); v1.PutValue(2);
00450   Epetra_IntVector v2(StandardMap); v2.PutValue(3);
00451 
00452   Epetra_Export identExporter(StandardMap,StandardMap); // Identity exporter
00453   EPETRA_TEST_ERR(!(v2.Export(v1, identExporter, Insert)==0),ierr);
00454   passed = (v2.MinValue()==2);
00455   EPETRA_TEST_ERR(!passed,ierr);
00456 
00457   v1.PutValue(1);
00458   Epetra_Import identImporter(StandardMap,StandardMap); // Identity importer
00459   EPETRA_TEST_ERR(!(v2.Import(v1, identExporter, Insert)==0),ierr);
00460   passed = passed && (v2.MaxValue()==1);
00461   EPETRA_TEST_ERR(!passed,ierr);
00462 
00463   if (verbose) {
00464     if (passed) cout << "Identity Import/Export Check OK" << endl << endl;
00465     else cout << "Identity Import/Export Check Failed" << endl << endl;
00466   }
00467 
00468   int NumSubMapElements = StandardMap.NumMyElements()/2;
00469   int SubStart = Comm.MyPID();
00470   NumSubMapElements = EPETRA_MIN(NumSubMapElements,StandardMap.NumMyElements()-SubStart);
00471   Epetra_Map SubMap(-1, NumSubMapElements, StandardMyGlobalElements+SubStart, 0, Comm);
00472 
00473   Epetra_IntVector v3(View, SubMap, SubMap.MyGlobalElements()); // Fill v3 with GID values for variety
00474   Epetra_Export subExporter(SubMap, StandardMap); // Export to a subset of indices of standard map
00475   EPETRA_TEST_ERR(!(v2.Export(v3,subExporter,Insert)==0),ierr);
00476 
00477   forierr = 0;
00478   for (i=0; i<SubMap.NumMyElements(); i++) {
00479     int i1 = StandardMap.LID(SubMap.GID(i));
00480     forierr += !(v3[i]==v2[i1]);
00481   }
00482   EPETRA_TEST_ERR(forierr,ierr);
00483 
00484   Epetra_Import subImporter(StandardMap, SubMap); // Import to a subset of indices of standard map
00485   EPETRA_TEST_ERR(!(v1.Import(v3,subImporter,Insert)==0),ierr);
00486 
00487   for (i=0; i<SubMap.NumMyElements(); i++) {
00488     int i1 = StandardMap.LID(SubMap.GID(i));
00489     forierr += !(v3[i]==v1[i1]);
00490   }
00491   EPETRA_TEST_ERR(forierr,ierr);
00492 
00493   if (verbose) {
00494     if (forierr==0) cout << "SubMap Import/Export Check OK" << endl << endl;
00495     else cout << "SubMap Import/Export Check Failed" << endl << endl;
00496   }
00497 
00498   forierr = special_submap_import_test(Comm);
00499   EPETRA_TEST_ERR(forierr, ierr);
00500 
00501   if (verbose) {
00502     if (forierr==0) cout << "Special SubMap Import Check OK" << endl << endl;
00503     else cout << "Special SubMap Import Check Failed" << endl << endl;
00504   }
00505 
00506   // Release all objects
00507 
00508   delete [] SourceMyGlobalElements;
00509   delete [] TargetMyGlobalElements;
00510   delete [] OverlapMyGlobalElements;
00511   delete [] StandardMyGlobalElements;
00512 
00513   delete [] Values;
00514   delete [] Indices;
00515 
00516 #ifdef EPETRA_MPI
00517   MPI_Finalize() ;
00518 #endif
00519 
00520 /* end main
00521 */
00522 return ierr ;
00523 }
00524 
00525 int special_submap_import_test(Epetra_Comm& Comm)
00526 {
00527   int localProc = Comm.MyPID();
00528 
00529   //set up ids_source and ids_target such that ids_source are only
00530   //a subset of ids_target, and furthermore that ids_target are ordered
00531   //such that the LIDs don't match up. In other words, even if gid 2 does
00532   //exist in both ids_source and ids_target, it will correspond to different
00533   //LIDs on at least 1 proc.
00534   //
00535   //This is to test a certain bug-fix in Epetra_Import where the 'RemoteLIDs'
00536   //array wasn't being calculated correctly on all procs.
00537 
00538   int ids_source[1];
00539   ids_source[0] = localProc*2+2;
00540 
00541   int ids_target[3];
00542   ids_target[0] = localProc*2+2;
00543   ids_target[1] = localProc*2+1;
00544   ids_target[2] = localProc*2+0;
00545 
00546   Epetra_Map map_source(-1, 1, &ids_source[0], 0, Comm);
00547   Epetra_Map map_target(-1, 3, &ids_target[0], 0, Comm);
00548 
00549   Epetra_Import importer(map_target, map_source);
00550 
00551   Epetra_IntVector vec_source(map_source);
00552   Epetra_IntVector vec_target(map_target);
00553 
00554   vec_target.PutValue(0);
00555 
00556   //set vec_source's contents so that entry[i] == GID[i].
00557   int* GIDs = map_source.MyGlobalElements();
00558   for(int i=0; i<map_source.NumMyElements(); ++i) {
00559     vec_source[i] = GIDs[i];
00560   }
00561 
00562   //Import vec_source into vec_target. This should result in the contents
00563   //of vec_target remaining 0 for the entries that don't exist in vec_source,
00564   //and other entries should be equal to the corresponding GID in the map.
00565 
00566   vec_target.Import(vec_source, importer, Insert);
00567 
00568   GIDs = map_target.MyGlobalElements();
00569   int test_failed = 0;
00570 
00571   //the test passes if the i-th entry in vec_target equals either 0 or
00572   //GIDs[i].
00573   for(int i=0; i<vec_target.MyLength(); ++i) {
00574     if (vec_target[i] != GIDs[i] && vec_target[i] != 0) test_failed = 1;
00575   }
00576 
00577   int global_result;
00578   Comm.MaxAll(&test_failed, &global_result, 1);
00579 
00580   //If test didn't fail on any procs, global_result should be 0.
00581   //If test failed on any proc, global_result should be 1.
00582   return global_result;
00583 }
00584 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines