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

Generated on Wed May 12 21:41:04 2010 for Epetra Package Browser (Single Doxygen Collection) by  doxygen 1.4.7