Epetra Package Browser (Single Doxygen Collection) Development
test/ImportExport_LL/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_LongLongVector.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 int alternate_import_constructor_test(Epetra_Comm& Comm);
00066 
00067 int main(int argc, char *argv[])
00068 {
00069   int ierr = 0, i, j, forierr = 0;
00070 
00071 #ifdef EPETRA_MPI
00072   // Initialize MPI
00073   MPI_Init(&argc,&argv);
00074   Epetra_MpiComm Comm( MPI_COMM_WORLD );
00075 #else
00076   Epetra_SerialComm Comm;
00077 #endif
00078 
00079   bool verbose = false;
00080 
00081   // Check if we should print results to standard out
00082   if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
00083 
00084 
00085 
00086 
00087   //char tmp;
00088   //if (Comm.MyPID()==0) cout << "Press any key to continue..."<< endl;
00089   //if (Comm.MyPID()==0) cin >> tmp;
00090   //Comm.Barrier();
00091 
00092   Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
00093   int MyPID = Comm.MyPID();
00094   int NumProc = Comm.NumProc();
00095 
00096   if (verbose && MyPID==0)
00097     cout << Epetra_Version() << endl << endl;
00098 
00099   if (verbose) cout << "Processor "<<MyPID<<" of "<< NumProc
00100               << " is alive."<<endl;
00101 
00102   // Redefine verbose to only print on PE 0
00103   if (verbose && Comm.MyPID()!=0) verbose = false;
00104 
00105   int NumMyEquations = 20;
00106   long long NumGlobalEquations = NumMyEquations*NumProc+EPETRA_MIN(NumProc,3);
00107   if (MyPID < 3) NumMyEquations++;
00108   // Construct a Source Map that puts approximately the same Number of equations on each processor in 
00109   // uniform global ordering
00110 
00111   Epetra_Map SourceMap(NumGlobalEquations, NumMyEquations, 0LL, Comm);
00112   
00113   // Get update list and number of local equations from newly created Map
00114   int NumMyElements = SourceMap.NumMyElements();
00115   long long * SourceMyGlobalElements = new long long[NumMyElements];
00116   SourceMap.MyGlobalElements(SourceMyGlobalElements);
00117 
00118   // Construct a Target Map that will contain:
00119   //  some unchanged elements (relative to the soure map),
00120   //  some permuted elements
00121   //  some off-processor elements
00122   Epetra_Vector RandVec(SourceMap);
00123   RandVec.Random(); // This creates a vector of random numbers between negative one and one.
00124 
00125   long long *TargetMyGlobalElements = new long long[NumMyElements];
00126 
00127   long long MinGID = SourceMap.MinMyGID64();
00128   for (i=0; i< NumMyEquations/2; i++) TargetMyGlobalElements[i] = i + MinGID; // Half will be the same...
00129   for (i=NumMyEquations/2; i<NumMyEquations; i++) {
00130     int index = abs((int)(((double) (NumGlobalEquations-1) ) * RandVec[i]));
00131     TargetMyGlobalElements[i] = EPETRA_MIN(NumGlobalEquations-1,(long long) EPETRA_MAX(0,index));
00132   }
00133 
00134   int NumSameIDs = 0;
00135   int NumPermutedIDs = 0;
00136   int NumRemoteIDs = 0;
00137   bool StillContiguous = true;
00138   for (i=0; i < NumMyEquations; i++) {
00139     if (SourceMyGlobalElements[i]==TargetMyGlobalElements[i] && StillContiguous)
00140       NumSameIDs++;
00141     else if (SourceMap.MyGID(TargetMyGlobalElements[i])) {
00142       StillContiguous = false;
00143       NumPermutedIDs++;
00144     }
00145     else {
00146       StillContiguous = false;
00147       NumRemoteIDs++;
00148     }
00149   }
00150   EPETRA_TEST_ERR(!(NumMyEquations==NumSameIDs+NumPermutedIDs+NumRemoteIDs),ierr);
00151 
00152   Epetra_Map TargetMap((long long) -1, NumMyElements, TargetMyGlobalElements, 0LL, Comm);
00153 
00154   // Create a multivector whose elements are GlobalID * (column number +1)
00155 
00156   int NumVectors = 3;
00157   Epetra_MultiVector SourceMultiVector(SourceMap, NumVectors);
00158   for (j=0; j < NumVectors; j++)
00159     for (i=0; i < NumMyElements; i++)
00160       SourceMultiVector[j][i] = (double) SourceMyGlobalElements[i]*(j+1);
00161 
00162   // Create a target multivector that we will fill using an Import
00163 
00164   Epetra_MultiVector TargetMultiVector(TargetMap, NumVectors);
00165 
00166   Epetra_Import Importer(TargetMap, SourceMap);
00167 
00168   EPETRA_TEST_ERR(!(TargetMultiVector.Import(SourceMultiVector, Importer, Insert)==0),ierr);
00169 
00170   // Test Target against expected values
00171   forierr = 0;
00172   for (j=0; j < NumVectors; j++)
00173     for (i=0; i < NumMyElements; i++) {
00174       if (TargetMultiVector[j][i]!= (double) TargetMyGlobalElements[i]*(j+1))
00175   cout << "TargetMultiVector["<<i<<"]["<<j<<"] = " << TargetMultiVector[j][i] 
00176        <<  "  TargetMyGlobalElements[i]*(j+1) = " <<  TargetMyGlobalElements[i]*(j+1) << endl;
00177       forierr += !(TargetMultiVector[j][i]== (double) TargetMyGlobalElements[i]*(j+1));
00178     }
00179   EPETRA_TEST_ERR(forierr,ierr);
00180 
00181   if (verbose) cout << "MultiVector Import using Importer Check OK" << endl << endl;
00182 
00183 
00185 
00186   // Now use Importer to do an export
00187 
00188   Epetra_Vector TargetVector(SourceMap);
00189   Epetra_Vector ExpectedTarget(SourceMap);
00190   Epetra_Vector SourceVector(TargetMap);
00191 
00192   NumSameIDs = Importer.NumSameIDs();
00193   int NumPermuteIDs = Importer.NumPermuteIDs();
00194   int NumExportIDs = Importer.NumExportIDs();
00195   int *PermuteFromLIDs = Importer.PermuteFromLIDs();
00196   int *ExportLIDs = Importer.ExportLIDs();
00197   int *ExportPIDs = Importer.ExportPIDs();
00198 
00199   for (i=0; i < NumSameIDs; i++) ExpectedTarget[i] = (double) (MyPID+1);
00200   for (i=0; i < NumPermuteIDs; i++) ExpectedTarget[PermuteFromLIDs[i]] = 
00201               (double) (MyPID+1);
00202   for (i=0; i < NumExportIDs; i++) ExpectedTarget[ExportLIDs[i]] += 
00203              (double) (ExportPIDs[i]+1);
00204 
00205   for (i=0; i < NumMyElements; i++) SourceVector[i] =  (double) (MyPID+1);
00206 
00207   EPETRA_TEST_ERR(!(TargetVector.Export(SourceVector, Importer, Add)==0),ierr);
00208 
00209   forierr = 0;
00210   for (i=0; i < NumMyElements; i++) {
00211     if (TargetVector[i]!= ExpectedTarget[i])
00212       cout <<  "     TargetVector["<<i<<"] = " << TargetVector[i] 
00213      <<  "   ExpectedTarget["<<i<<"] = " <<  ExpectedTarget[i] << " on PE " << MyPID << endl;
00214     forierr += !(TargetVector[i]== ExpectedTarget[i]);
00215   }
00216   EPETRA_TEST_ERR(forierr,ierr);
00217 
00218   if (verbose) cout << "Vector Export using Importer Check OK" << endl << endl;
00219 
00221   // Now use Importer to create a reverse exporter
00222   TargetVector.PutScalar(0.0);
00223   Epetra_Export ReversedImport(Importer);
00224 
00225   EPETRA_TEST_ERR(!(TargetVector.Export(SourceVector, ReversedImport, Add)==0),ierr);
00226 
00227   forierr = 0;
00228   for (i=0; i < NumMyElements; i++) {
00229     if (TargetVector[i]!= ExpectedTarget[i])
00230       cout <<  "     TargetVector["<<i<<"] = " << TargetVector[i] 
00231      <<  "   ExpectedTarget["<<i<<"] = " <<  ExpectedTarget[i] << " on PE " << MyPID << endl;
00232     forierr += !(TargetVector[i]== ExpectedTarget[i]);
00233   }
00234   EPETRA_TEST_ERR(forierr,ierr);
00235 
00236   if (verbose) cout << "Vector Export using Reversed Importer Check OK" << endl << endl;
00237 
00239   // Now use Exporter to create a reverse importer
00240   TargetVector.PutScalar(0.0);
00241   Epetra_Import ReversedExport(ReversedImport);
00242 
00243   EPETRA_TEST_ERR(!(TargetVector.Export(SourceVector, ReversedExport, Add)==0),ierr);
00244 
00245   forierr = 0;
00246   for (i=0; i < NumMyElements; i++) {
00247     if (TargetVector[i]!= ExpectedTarget[i])
00248       cout <<  "     TargetVector["<<i<<"] = " << TargetVector[i] 
00249      <<  "   ExpectedTarget["<<i<<"] = " <<  ExpectedTarget[i] << " on PE " << MyPID << endl;
00250     forierr += !(TargetVector[i]== ExpectedTarget[i]);
00251   }
00252   EPETRA_TEST_ERR(forierr,ierr);
00253 
00254   if (verbose) cout << "Vector Export using Reversed Exporter Check OK" << endl << endl;
00255 
00256 
00258   //  Build a tridiagonal system two ways:
00259   //  1) From "standard" matrix view where equations are uniquely owned.
00260   //  2) From 1D PDE view where nodes (equations) between processors are shared and partial contributions are done
00261   //     in parallel, then merged together at the end of the construction process.
00262   //
00264 
00265 
00266 
00267   // Construct a Standard Map that puts approximately the same number of equations on each processor in 
00268   // uniform global ordering
00269 
00270   Epetra_Map StandardMap(NumGlobalEquations, NumMyEquations, 0LL, Comm);
00271   
00272   // Get update list and number of local equations from newly created Map
00273   NumMyElements = StandardMap.NumMyElements();
00274   long long * StandardMyGlobalElements = new long long[NumMyElements];
00275   StandardMap.MyGlobalElements(StandardMyGlobalElements);
00276 
00277 
00278   // Create a standard Epetra_CrsGraph
00279 
00280   Epetra_CrsGraph StandardGraph(Copy, StandardMap, 3);
00281   EPETRA_TEST_ERR(StandardGraph.IndicesAreGlobal(),ierr);
00282   EPETRA_TEST_ERR(StandardGraph.IndicesAreLocal(),ierr);
00283   
00284   // Add  rows one-at-a-time
00285   // Need some vectors to help
00286   // Off diagonal Values will always be -1
00287 
00288 
00289   long long *Indices = new long long[2];
00290   int NumEntries;
00291   
00292   forierr = 0;
00293   for (i=0; i<NumMyEquations; i++)
00294     {
00295     if (StandardMyGlobalElements[i]==0)
00296       {
00297   Indices[0] = 1;
00298   NumEntries = 1;
00299       }
00300     else if (StandardMyGlobalElements[i] == NumGlobalEquations-1)
00301       {
00302   Indices[0] = NumGlobalEquations-2;
00303   NumEntries = 1;
00304       }
00305     else
00306       {
00307   Indices[0] = StandardMyGlobalElements[i]-1;
00308   Indices[1] = StandardMyGlobalElements[i]+1;
00309   NumEntries = 2;
00310       }
00311     forierr += !(StandardGraph.InsertGlobalIndices(StandardMyGlobalElements[i], NumEntries, Indices)==0);
00312     forierr += !(StandardGraph.InsertGlobalIndices(StandardMyGlobalElements[i], 1, StandardMyGlobalElements+i)==0); // Put in the diagonal entry
00313     }
00314   EPETRA_TEST_ERR(forierr,ierr);
00315 
00316   // Finish up
00317   EPETRA_TEST_ERR(!(StandardGraph.IndicesAreGlobal()),ierr);
00318   EPETRA_TEST_ERR(!(StandardGraph.FillComplete()==0),ierr);
00319   EPETRA_TEST_ERR(!(StandardGraph.IndicesAreLocal()),ierr);
00320   EPETRA_TEST_ERR(StandardGraph.StorageOptimized(),ierr);
00321   StandardGraph.OptimizeStorage();
00322   EPETRA_TEST_ERR(!(StandardGraph.StorageOptimized()),ierr);
00323   EPETRA_TEST_ERR(StandardGraph.UpperTriangular(),ierr);
00324   EPETRA_TEST_ERR(StandardGraph.LowerTriangular(),ierr);
00325 
00326   // Create Epetra_CrsMatrix using the just-built graph
00327 
00328   Epetra_CrsMatrix StandardMatrix(Copy, StandardGraph);
00329   EPETRA_TEST_ERR(StandardMatrix.IndicesAreGlobal(),ierr);
00330   EPETRA_TEST_ERR(!(StandardMatrix.IndicesAreLocal()),ierr);
00331   
00332   // Add  rows one-at-a-time
00333   // Need some vectors to help
00334   // Off diagonal Values will always be -1
00335 
00336 
00337   double *Values = new double[2];
00338   Values[0] = -1.0; Values[1] = -1.0;
00339   double two = 2.0;
00340   
00341   forierr = 0;
00342   for (i=0; i<NumMyEquations; i++)
00343     {
00344     if (StandardMyGlobalElements[i]==0)
00345       {
00346   Indices[0] = 1;
00347   NumEntries = 1;
00348       }
00349     else if (StandardMyGlobalElements[i] == NumGlobalEquations-1)
00350       {
00351   Indices[0] = NumGlobalEquations-2;
00352   NumEntries = 1;
00353       }
00354     else
00355       {
00356   Indices[0] = StandardMyGlobalElements[i]-1;
00357   Indices[1] = StandardMyGlobalElements[i]+1;
00358   NumEntries = 2;
00359       }
00360     forierr += !(StandardMatrix.ReplaceGlobalValues(StandardMyGlobalElements[i], NumEntries, Values, Indices)==0);
00361     // Put in the diagonal entry
00362     forierr += !(StandardMatrix.ReplaceGlobalValues(StandardMyGlobalElements[i], 1, &two, StandardMyGlobalElements+i)==0); 
00363     }
00364   EPETRA_TEST_ERR(forierr,ierr);
00365 
00366   // Finish up
00367   EPETRA_TEST_ERR(!(StandardMatrix.IndicesAreLocal()),ierr);
00368   EPETRA_TEST_ERR(!(StandardMatrix.FillComplete()==0),ierr);
00369   EPETRA_TEST_ERR(!(StandardMatrix.IndicesAreLocal()),ierr);
00370   //  EPETRA_TEST_ERR((StandardMatrix.StorageOptimized()),ierr);
00371   EPETRA_TEST_ERR((StandardMatrix.OptimizeStorage()),ierr);
00372   EPETRA_TEST_ERR(!(StandardMatrix.StorageOptimized()),ierr);
00373   EPETRA_TEST_ERR(StandardMatrix.UpperTriangular(),ierr);
00374   EPETRA_TEST_ERR(StandardMatrix.LowerTriangular(),ierr);
00375 
00376   // Construct an Overlapped Map of StandardMap that include the endpoints from two neighboring processors.
00377 
00378   int OverlapNumMyElements;
00379   long long OverlapMinMyGID;
00380 
00381   OverlapNumMyElements = NumMyElements + 1;
00382   if (MyPID==0) OverlapNumMyElements--;
00383 
00384   if (MyPID==0) OverlapMinMyGID = StandardMap.MinMyGID64();
00385   else OverlapMinMyGID = StandardMap.MinMyGID64()-1;
00386 
00387   long long * OverlapMyGlobalElements = new long long[OverlapNumMyElements];
00388 
00389   for (i=0; i< OverlapNumMyElements; i++) OverlapMyGlobalElements[i] = OverlapMinMyGID + i;
00390 
00391   Epetra_Map OverlapMap((long long) -1, OverlapNumMyElements, OverlapMyGlobalElements, 0LL, Comm);
00392 
00393   // Create the Overlap Epetra_Matrix
00394 
00395   Epetra_CrsMatrix OverlapMatrix(Copy, OverlapMap, 4);
00396   EPETRA_TEST_ERR(OverlapMatrix.IndicesAreGlobal(),ierr);
00397   EPETRA_TEST_ERR(OverlapMatrix.IndicesAreLocal(),ierr);
00398   
00399   // Add  matrix element one cell at a time.
00400   // Each cell does an incoming and outgoing flux calculation
00401 
00402 
00403   double pos_one = 1.0;
00404   double neg_one = -1.0;
00405 
00406   forierr = 0;
00407   for (i=0; i<OverlapNumMyElements; i++)
00408     {
00409       long long node_left = OverlapMyGlobalElements[i]-1;
00410       long long node_center = node_left + 1;
00411       long long node_right = node_left + 2;
00412       if (i>0) {
00413   if (node_left>-1)
00414     forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &neg_one, &node_left)==0);
00415   forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0);
00416       }
00417       if (i<OverlapNumMyElements-1) {
00418   forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0);
00419   if (node_right<NumGlobalEquations) 
00420     forierr += !(OverlapMatrix.InsertGlobalValues(node_center, 1, &neg_one, &node_right)==0);
00421       }
00422     }
00423   EPETRA_TEST_ERR(forierr,ierr);
00424 
00425   // Handle endpoints
00426   if (MyPID==0) {
00427     long long node_center = 0;
00428     EPETRA_TEST_ERR(!(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0),ierr);
00429   }
00430   if (MyPID==NumProc-1) {
00431     long long node_center = OverlapMyGlobalElements[OverlapNumMyElements-1];
00432     EPETRA_TEST_ERR(!(OverlapMatrix.InsertGlobalValues(node_center, 1, &pos_one, &node_center)==0),ierr);
00433   }
00434     
00435   EPETRA_TEST_ERR(!(OverlapMatrix.FillComplete()==0),ierr);
00436 
00437   // Make a gathered matrix from OverlapMatrix.  It should be identical to StandardMatrix
00438 
00439   Epetra_CrsMatrix GatheredMatrix(Copy, StandardGraph);
00440   Epetra_Export Exporter(OverlapMap, StandardMap);
00441   EPETRA_TEST_ERR(!(GatheredMatrix.Export(OverlapMatrix, Exporter, Add)==0),ierr);
00442   EPETRA_TEST_ERR(!(GatheredMatrix.FillComplete()==0),ierr);
00443 
00444   // Check if entries of StandardMatrix and GatheredMatrix are identical
00445 
00446   int StandardNumEntries, GatheredNumEntries;
00447   int * StandardIndices, * GatheredIndices;
00448   double * StandardValues, * GatheredValues;
00449 
00450   int StandardNumMyNonzeros = StandardMatrix.NumMyNonzeros();
00451   int GatheredNumMyNonzeros = GatheredMatrix.NumMyNonzeros();
00452   EPETRA_TEST_ERR(!(StandardNumMyNonzeros==GatheredNumMyNonzeros),ierr);
00453 
00454   int StandardNumMyRows = StandardMatrix.NumMyRows();
00455   int GatheredNumMyRows = GatheredMatrix.NumMyRows();
00456   EPETRA_TEST_ERR(!(StandardNumMyRows==GatheredNumMyRows),ierr);
00457 
00458   forierr = 0;
00459   for (i=0; i< StandardNumMyRows; i++)
00460     {
00461       forierr += !(StandardMatrix.ExtractMyRowView(i, StandardNumEntries, StandardValues, StandardIndices)==0);
00462       forierr += !(GatheredMatrix.ExtractMyRowView(i, GatheredNumEntries, GatheredValues, GatheredIndices)==0);
00463       forierr += !(StandardNumEntries==GatheredNumEntries);
00464       for (j=0; j < StandardNumEntries; j++) {
00465   //if (StandardIndices[j]!=GatheredIndices[j])
00466   // cout << "MyPID = " << MyPID << " i = " << i << "   StandardIndices[" << j << "] = " << StandardIndices[j] 
00467   //      << "   GatheredIndices[" << j << "] = " << GatheredIndices[j] << endl;
00468   //if (StandardValues[j]!=GatheredValues[j])
00469   //cout << "MyPID = " << MyPID << " i = " << i << "    StandardValues[" << j << "] = " <<  StandardValues[j] 
00470   //     << "    GatheredValues[" << j << "] = " <<  GatheredValues[j] << endl;
00471   forierr += !(StandardIndices[j]==GatheredIndices[j]);
00472   forierr += !(StandardValues[j]==GatheredValues[j]);
00473       }
00474     }
00475   EPETRA_TEST_ERR(forierr,ierr);
00476 
00477   if (verbose) cout << "Matrix Export Check OK" << endl << endl;
00478 
00479   //Do Again with use of Epetra_OffsetIndex object for speed
00480   Epetra_OffsetIndex OffsetIndex( OverlapMatrix.Graph(), GatheredMatrix.Graph(), Exporter );
00481   EPETRA_TEST_ERR(!(GatheredMatrix.Export(OverlapMatrix, Exporter, Add)==0),ierr);
00482 
00483   if (verbose) cout << "Optimized Matrix Export Check OK" << endl << endl;
00484 
00485   bool passed;
00486   Epetra_LongLongVector v1(StandardMap); v1.PutValue(2);
00487   Epetra_LongLongVector v2(StandardMap); v2.PutValue(3);
00488 
00489   Epetra_Export identExporter(StandardMap,StandardMap); // Identity exporter
00490   EPETRA_TEST_ERR(!(v2.Export(v1, identExporter, Insert)==0),ierr);
00491   passed = (v2.MinValue()==2);
00492   EPETRA_TEST_ERR(!passed,ierr);
00493 
00494   v1.PutValue(1);
00495   Epetra_Import identImporter(StandardMap,StandardMap); // Identity importer
00496   EPETRA_TEST_ERR(!(v2.Import(v1, identExporter, Insert)==0),ierr);
00497   passed = passed && (v2.MaxValue()==1);
00498   EPETRA_TEST_ERR(!passed,ierr);
00499 
00500   if (verbose) {
00501     if (passed) cout << "Identity Import/Export Check OK" << endl << endl;
00502     else cout << "Identity Import/Export Check Failed" << endl << endl;
00503   }
00504 
00505   int NumSubMapElements = StandardMap.NumMyElements()/2;
00506   int SubStart = Comm.MyPID();
00507   NumSubMapElements = EPETRA_MIN(NumSubMapElements,StandardMap.NumMyElements()-SubStart);
00508   Epetra_Map SubMap((long long) -1, NumSubMapElements, StandardMyGlobalElements+SubStart, 0LL, Comm);
00509 
00510   Epetra_LongLongVector v3(View, SubMap, SubMap.MyGlobalElements64()); // Fill v3 with GID values for variety
00511   Epetra_Export subExporter(SubMap, StandardMap); // Export to a subset of indices of standard map
00512   EPETRA_TEST_ERR(!(v2.Export(v3,subExporter,Insert)==0),ierr);
00513 
00514   forierr = 0;
00515   for (i=0; i<SubMap.NumMyElements(); i++) {
00516     int i1 = StandardMap.LID(SubMap.GID64(i));
00517     forierr += !(v3[i]==v2[i1]);
00518   }
00519   EPETRA_TEST_ERR(forierr,ierr);
00520 
00521   Epetra_Import subImporter(StandardMap, SubMap); // Import to a subset of indices of standard map
00522   EPETRA_TEST_ERR(!(v1.Import(v3,subImporter,Insert)==0),ierr);
00523 
00524   for (i=0; i<SubMap.NumMyElements(); i++) {
00525     int i1 = StandardMap.LID(SubMap.GID64(i));
00526     forierr += !(v3[i]==v1[i1]);
00527   }
00528   EPETRA_TEST_ERR(forierr,ierr);
00529 
00530   if (verbose) {
00531     if (forierr==0) cout << "SubMap Import/Export Check OK" << endl << endl;
00532     else cout << "SubMap Import/Export Check Failed" << endl << endl;
00533   }
00534 
00535 #ifdef DOESNT_WORK_IN_PARALLEL
00536   forierr = special_submap_import_test(Comm);
00537   EPETRA_TEST_ERR(forierr, ierr);
00538 
00539   if (verbose) {
00540     if (forierr==0) cout << "Special SubMap Import Check OK" << endl << endl;
00541     else cout << "Special SubMap Import Check Failed" << endl << endl;
00542   }
00543 #endif
00544 
00545   forierr =  alternate_import_constructor_test(Comm);
00546   EPETRA_TEST_ERR(forierr, ierr);
00547 
00548   if (verbose) {
00549     if (forierr==0) cout << "Alternative Import Constructor Check OK" << endl << endl;
00550     else cout << "Alternative Import Constructor Check Failed" << endl << endl;
00551   }
00552 
00553   // Release all objects
00554 
00555   delete [] SourceMyGlobalElements;
00556   delete [] TargetMyGlobalElements;
00557   delete [] OverlapMyGlobalElements;
00558   delete [] StandardMyGlobalElements;
00559 
00560   delete [] Values;
00561   delete [] Indices;
00562 
00563 #ifdef EPETRA_MPI
00564   MPI_Finalize() ;
00565 #endif
00566 
00567 /* end main
00568 */
00569 return ierr ;
00570 }
00571 
00572 int special_submap_import_test(Epetra_Comm& Comm)
00573 {
00574   int localProc = Comm.MyPID();
00575 
00576   //set up ids_source and ids_target such that ids_source are only
00577   //a subset of ids_target, and furthermore that ids_target are ordered
00578   //such that the LIDs don't match up. In other words, even if gid 2 does
00579   //exist in both ids_source and ids_target, it will correspond to different
00580   //LIDs on at least 1 proc.
00581   //
00582   //This is to test a certain bug-fix in Epetra_Import where the 'RemoteLIDs'
00583   //array wasn't being calculated correctly on all procs.
00584 
00585   long long ids_source[1];
00586   ids_source[0] = localProc*2+2;
00587 
00588   long long ids_target[3];
00589   ids_target[0] = localProc*2+2;
00590   ids_target[1] = localProc*2+1;
00591   ids_target[2] = localProc*2+0;
00592 
00593   Epetra_Map map_source((long long) -1, 1, &ids_source[0], 0LL, Comm);
00594   Epetra_Map map_target((long long) -1, 3, &ids_target[0], 0LL, Comm);
00595 
00596   Epetra_Import importer(map_target, map_source);
00597 
00598   Epetra_LongLongVector vec_source(map_source);
00599   Epetra_LongLongVector vec_target(map_target);
00600 
00601   vec_target.PutValue(0);
00602 
00603   //set vec_source's contents so that entry[i] == GID[i].
00604   long long* GIDs = map_source.MyGlobalElements64();
00605   for(int i=0; i<map_source.NumMyElements(); ++i) {
00606     vec_source[i] = GIDs[i];
00607   }
00608 
00609   //Import vec_source into vec_target. This should result in the contents
00610   //of vec_target remaining 0 for the entries that don't exist in vec_source,
00611   //and other entries should be equal to the corresponding GID in the map.
00612 
00613   vec_target.Import(vec_source, importer, Insert);
00614 
00615   GIDs = map_target.MyGlobalElements64();
00616   int test_failed = 0;
00617 
00618   //the test passes if the i-th entry in vec_target equals either 0 or
00619   //GIDs[i].
00620   for(int i=0; i<vec_target.MyLength(); ++i) {
00621     if (vec_target[i] != GIDs[i] && vec_target[i] != 0) test_failed = 1;
00622   }
00623 
00624   int global_result;
00625   Comm.MaxAll(&test_failed, &global_result, 1);
00626 
00627   //If test didn't fail on any procs, global_result should be 0.
00628   //If test failed on any proc, global_result should be 1.
00629   return global_result;
00630 }
00631 
00632 int combine_mode_test(Epetra_Comm& Comm)
00633 {
00634   int localProc = Comm.MyPID();
00635 
00636 
00637   long long ids_source[1];
00638   ids_source[0] = localProc*2+2;
00639 
00640   long long ids_target[3];
00641   ids_target[0] = localProc*2+2;
00642   ids_target[1] = localProc*2+1;
00643   ids_target[2] = localProc*2+0;
00644 
00645   Epetra_Map map_source((long long) -1, 1, &ids_source[0], 0LL, Comm);
00646   Epetra_Map map_target((long long) -1, 3, &ids_target[0], 0LL, Comm);
00647 
00648   Epetra_Import importer(map_target, map_source);
00649 
00650   Epetra_LongLongVector vec_source(map_source);
00651   Epetra_LongLongVector vec_target(map_target);
00652 
00653   vec_target.PutValue(0);
00654 
00655   //set vec_source's contents so that entry[i] == GID[i].
00656   long long* GIDs = map_source.MyGlobalElements64();
00657   for(int i=0; i<map_source.NumMyElements(); ++i) {
00658     vec_source[i] = GIDs[i];
00659   }
00660 
00661   //Import vec_source into vec_target. This should result in the contents
00662   //of vec_target remaining 0 for the entries that don't exist in vec_source,
00663   //and other entries should be equal to the corresponding GID in the map.
00664 
00665   vec_target.Import(vec_source, importer, Insert);
00666 
00667   GIDs = map_target.MyGlobalElements64();
00668   int test_failed = 0;
00669 
00670   //the test passes if the i-th entry in vec_target equals either 0 or
00671   //GIDs[i].
00672   for(int i=0; i<vec_target.MyLength(); ++i) {
00673     if (vec_target[i] != GIDs[i] && vec_target[i] != 0) test_failed = 1;
00674   }
00675 
00676   int global_result;
00677   Comm.MaxAll(&test_failed, &global_result, 1);
00678 
00679   //If test didn't fail on any procs, global_result should be 0.
00680   //If test failed on any proc, global_result should be 1.
00681   return global_result;
00682 }
00683 
00684 
00685 int test_import_gid(const char * name,Epetra_LongLongVector & Source, Epetra_LongLongVector & Target, const Epetra_Import & Import){
00686   int i;
00687   bool test_passed=true;
00688 
00689   // Setup
00690   for(i=0; i<Source.MyLength(); i++)
00691     Source[i] = Source.Map().GID64(i);
00692   Target.PutValue(0);
00693 
00694   // Import
00695   Target.Import(Source,Import,Add);
00696 
00697   // Test
00698   for(i=0; i<Target.MyLength(); i++){
00699     if(Target[i] != Target.Map().GID64(i)) test_passed=false;
00700   }
00701 
00702   if(!test_passed){
00703     printf("[%d] test_import_gid %s failed: ",Source.Map().Comm().MyPID(),name);
00704     for(i=0; i<Target.MyLength(); i++)
00705       printf("%2lld(%2lld) ",Target[i],Target.Map().GID64(i));
00706     printf("\n");
00707     fflush(stdout);
00708   }
00709 
00710   return !test_passed;
00711 }
00712 
00713 
00714 
00715 int alternate_import_constructor_test(Epetra_Comm& Comm) {
00716   int rv=0;
00717   int nodes_per_proc=10;
00718   int numprocs = Comm.NumProc();
00719   int mypid    = Comm.MyPID();
00720 
00721   // Only run if we have multiple procs & MPI
00722   if(numprocs==0) return 0;
00723 #ifndef HAVE_MPI
00724   return 0;
00725 #endif
00726 
00727   // Build Map 1 - linear
00728   Epetra_Map Map1((long long)-1,nodes_per_proc,(long long)0,Comm);
00729   
00730   // Build Map 2 - mod striped
00731   std::vector<long long> MyGIDs(nodes_per_proc);
00732   for(int i=0; i<nodes_per_proc; i++)
00733     MyGIDs[i] = (mypid*nodes_per_proc + i) % numprocs;
00734   Epetra_Map Map2((long long)-1,nodes_per_proc,&MyGIDs[0],(long long)0,Comm);
00735 
00736   // For testing
00737   Epetra_LongLongVector Source(Map1), Target(Map2);
00738 
00739 
00740   // Build Import 1 - normal
00741   Epetra_Import Import1(Map2,Map1);
00742   rv = rv|| test_import_gid("Alt test: 2 map constructor",Source,Target, Import1);
00743 
00744   // Build Import 2 - no-comm constructor
00745   int Nremote=Import1.NumRemoteIDs();
00746   const int * RemoteLIDs = Import1.RemoteLIDs();
00747   std::vector<int> RemotePIDs(Nremote+1); // I hate you, stl vector....
00748   std::vector<int> AllPIDs;  
00749   Epetra_Util::GetPids(Import1,AllPIDs,true);
00750 
00751   for(int i=0; i<Nremote; i++) {
00752     RemotePIDs[i]=AllPIDs[RemoteLIDs[i]];
00753   }
00754   Epetra_Import Import2(Import1.TargetMap(),Import1.SourceMap(),Nremote,&RemotePIDs[0],Import1.NumExportIDs(),Import1.ExportLIDs(),Import1.ExportPIDs());
00755 
00756   rv = rv || test_import_gid("Alt test: no comm constructor",Source,Target,Import2);
00757 
00758 
00759   // Build Import 3 - Remotes only
00760   Epetra_Import Import3(Import1.TargetMap(),Import1.SourceMap(),Nremote,&RemotePIDs[0]);
00761   rv = rv || test_import_gid("Alt test: remote only constructor",Source,Target, Import3);
00762 
00763 
00764   return rv;
00765 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines