test/Map/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 // Epetra_Map (and Epetra_LocalMap) Test routine
00030 #include "Epetra_Time.h"
00031 #include "Epetra_Map.h"
00032 #include "Epetra_LocalMap.h"
00033 #ifdef EPETRA_MPI
00034 #include "Epetra_MpiComm.h"
00035 #include <mpi.h>
00036 #else
00037 #include "Epetra_SerialComm.h"
00038 #endif
00039 #include "checkmap.h"
00040 #include "../epetra_test_err.h"
00041 #include "Epetra_Version.h"
00042 
00043 int checkMapDataClass(Epetra_Comm& Comm, int verbose);
00044 int checkLocalMapDataClass(Epetra_Comm& Comm, int verbose);
00045 
00046 int main(int argc, char *argv[]) {
00047 
00048   int ierr=0, returnierr=0;
00049 
00050 #ifdef EPETRA_MPI
00051   MPI_Init(&argc,&argv);
00052   Epetra_MpiComm Comm(MPI_COMM_WORLD);
00053 #else
00054   Epetra_SerialComm Comm;
00055 #endif
00056 
00057   bool verbose = false;
00058 
00059   // Check if we should print results to standard out
00060   if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
00061 
00062 
00063   if (!verbose) {
00064     Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
00065   }
00066   int MyPID = Comm.MyPID();
00067   int NumProc = Comm.NumProc();
00068 
00069   if (verbose && MyPID==0)
00070     cout << Epetra_Version() << endl << endl;
00071 
00072   if (verbose) cout << Comm << endl;
00073 
00074   bool verbose1 = verbose;
00075   if (verbose) verbose = (MyPID==0);
00076 
00077   int NumMyElements = 10000;
00078   int NumMyElements1 = NumMyElements; // Used for local map
00079   int NumGlobalElements = NumMyElements*NumProc+EPETRA_MIN(NumProc,3);
00080   if (MyPID < 3) NumMyElements++;
00081   int IndexBase = 0;
00082   bool DistributedGlobal = (NumGlobalElements>NumMyElements);
00083   
00084   Epetra_Map* Map;
00085 
00086   // Test exceptions
00087 
00088   if (verbose)
00089     cout << "*******************************************************************************************" << endl
00090    << "        Testing Exceptions (Expect error messages if EPETRA_NO_ERROR_REPORTS is not defined" << endl
00091    << "*******************************************************************************************" << endl
00092    << endl << endl;
00093 
00094   try {
00095     if (verbose) cout << "Checking Epetra_Map(-2, IndexBase, Comm)" << endl;
00096     Map = new Epetra_Map(-2, IndexBase, Comm);
00097   }
00098   catch (int Error) {
00099     if (Error!=-1) {
00100       if (Error!=0) {
00101   EPETRA_TEST_ERR(Error,returnierr);
00102   if (verbose) cout << "Error code should be -1" << endl;
00103       }
00104       else {
00105   cout << "Error code = " << Error << "Should be -1" << endl;
00106   returnierr+=1;
00107       }
00108     }
00109     else if (verbose) cout << "Checked OK\n\n" << endl;
00110   }
00111 
00112   try {
00113     if (verbose) cout << "Checking Epetra_Map(2, 3, IndexBase, Comm)" << endl;
00114     Map = new Epetra_Map(2, 3, IndexBase, Comm);
00115   }
00116   catch (int Error) {
00117     if (Error!=-4) {
00118       if (Error!=0) {
00119   EPETRA_TEST_ERR(Error,returnierr);
00120   if (verbose) cout << "Error code should be -4" << endl;
00121       }
00122       else {
00123   cout << "Error code = " << Error << "Should be -4" << endl;
00124   returnierr+=1;
00125       }
00126     }
00127     else if (verbose) cout << "Checked OK\n\n" << endl;
00128   }
00129 
00130   if (verbose) cerr << flush;
00131   if (verbose) cout << flush;
00132   Comm.Barrier();
00133   if (verbose)
00134     cout << endl << endl
00135       << "*******************************************************************************************" << endl
00136       << "        Testing valid constructor now......................................................" << endl
00137       << "*******************************************************************************************" << endl
00138       << endl << endl;
00139 
00140   // Test Epetra-defined uniform linear distribution constructor
00141   Map = new Epetra_Map(NumGlobalElements, IndexBase, Comm);
00142   if (verbose) cout << "Checking Epetra_Map(NumGlobalElements, IndexBase, Comm)" << endl;
00143   ierr = checkmap(*Map, NumGlobalElements, NumMyElements, 0, 
00144       IndexBase, Comm, DistributedGlobal);
00145 
00146   EPETRA_TEST_ERR(ierr,returnierr);
00147   if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
00148 
00149   delete Map;
00150 
00151   // Test User-defined linear distribution constructor
00152   Map = new Epetra_Map(NumGlobalElements, NumMyElements, IndexBase, Comm);
00153 
00154   if (verbose) cout << "Checking Epetra_Map(NumGlobalElements, NumMyElements, IndexBase, Comm)" << endl;
00155   ierr = checkmap(*Map, NumGlobalElements, NumMyElements, 0, 
00156       IndexBase, Comm, DistributedGlobal);
00157 
00158   EPETRA_TEST_ERR(ierr,returnierr);
00159   if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
00160   delete Map;
00161 
00162   // Test User-defined arbitrary distribution constructor
00163   // Generate Global Element List.  Do in reverse for fun!
00164 
00165   int * MyGlobalElements = new int[NumMyElements];
00166   int MaxMyGID = (Comm.MyPID()+1)*NumMyElements-1+IndexBase;
00167   if (Comm.MyPID()>2) MaxMyGID+=3;
00168   for (int i = 0; i<NumMyElements; i++) MyGlobalElements[i] = MaxMyGID-i;
00169 
00170   Map = new Epetra_Map(NumGlobalElements, NumMyElements, MyGlobalElements, 
00171                        IndexBase, Comm);
00172   if (verbose) cout << "Checking Epetra_Map(NumGlobalElements, NumMyElements, MyGlobalElements,  IndexBase, Comm)" << endl;
00173   ierr = checkmap(*Map, NumGlobalElements, NumMyElements, MyGlobalElements, 
00174                   IndexBase, Comm, DistributedGlobal);
00175 
00176   EPETRA_TEST_ERR(ierr,returnierr);
00177   if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
00178   // Test Copy constructor
00179   Epetra_Map* Map1 = new Epetra_Map(*Map);
00180 
00181   // Test SameAs() method
00182   bool same = Map1->SameAs(*Map);
00183   EPETRA_TEST_ERR(!(same==true),ierr);// should return true since Map1 is a copy of Map
00184 
00185   Epetra_BlockMap* Map2 = new Epetra_Map(NumGlobalElements, NumMyElements, MyGlobalElements,  IndexBase, Comm);
00186   same = Map2->SameAs(*Map);
00187   EPETRA_TEST_ERR(!(same==true),ierr); // Map and Map2 were created with the same sets of parameters
00188   delete Map2;
00189 
00190   // now test SameAs() on a map that is different
00191 
00192   Map2 =  new Epetra_Map(NumGlobalElements, NumMyElements, MyGlobalElements, IndexBase-1, Comm);
00193   same = Map2->SameAs(*Map);
00194   EPETRA_TEST_ERR(!(same==false),ierr); // IndexBases are different
00195   delete Map2;
00196 
00197   // Back to testing copy constructor
00198   if (verbose) cout << "Checking Epetra_Map(*Map)" << endl;
00199   ierr = checkmap(*Map1, NumGlobalElements, NumMyElements, MyGlobalElements, 
00200       IndexBase, Comm, DistributedGlobal);
00201 
00202   EPETRA_TEST_ERR(ierr,returnierr);
00203   if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
00204   Epetra_Map* SmallMap = 0;
00205   if (verbose1) {
00206     // Build a small map for test cout.  Use 10 elements from current map
00207     int* MyEls = Map->MyGlobalElements();
00208     int IndBase = Map->IndexBase();
00209     int MyLen = EPETRA_MIN(10+Comm.MyPID(),Map->NumMyElements());
00210     SmallMap = new Epetra_Map(-1, MyLen, MyEls, IndBase, Comm);
00211   }
00212 
00213   delete [] MyGlobalElements;
00214   delete Map;
00215   delete Map1;
00216 
00217   // Test reference-counting in Epetra_Map
00218   if (verbose) cout << "Checking Epetra_Map reference counting" << endl;
00219   ierr = checkMapDataClass(Comm, verbose);
00220   EPETRA_TEST_ERR(ierr,returnierr);
00221   if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
00222 
00223   // Test LocalMap constructor
00224   Epetra_LocalMap* LocalMap = new Epetra_LocalMap(NumMyElements1, IndexBase, Comm);
00225   if (verbose) cout << "Checking Epetra_LocalMap(NumMyElements1, IndexBase, Comm)" << endl;
00226   ierr = checkmap(*LocalMap, NumMyElements1, NumMyElements1, 0, IndexBase, Comm, false);
00227 
00228   EPETRA_TEST_ERR(ierr,returnierr);
00229   if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
00230   // Test Copy constructor
00231   Epetra_LocalMap* LocalMap1 = new Epetra_LocalMap(*LocalMap);
00232   if (verbose) cout << "Checking Epetra_LocalMap(*LocalMap)" << endl;
00233   ierr = checkmap(*LocalMap1, NumMyElements1, NumMyElements1, 0, IndexBase, Comm, false);
00234 
00235   EPETRA_TEST_ERR(ierr,returnierr);
00236   if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
00237   delete LocalMap1;
00238   delete LocalMap;
00239 
00240   // Test reference-counting in Epetra_LocalMap
00241   if (verbose) cout << "Checking Epetra_LocalMap reference counting" << endl;
00242   ierr = checkLocalMapDataClass(Comm, verbose);
00243   EPETRA_TEST_ERR(ierr,returnierr);
00244   if (verbose && ierr==0) cout << "Checked OK\n\n" <<endl;
00245 
00246   // Test output
00247   if (verbose1) {
00248     if (verbose) cout << "Test ostream << operator" << endl << flush;
00249     cout << *SmallMap;
00250     delete SmallMap;
00251   }
00252 
00253 #ifdef EPETRA_MPI
00254   MPI_Finalize();
00255 #endif
00256 
00257   return returnierr;
00258 }
00259 
00260 int checkMapDataClass(Epetra_Comm& Comm, int verbose) {
00261   int returnierr = 0;
00262   int NumGlobalElements = 1000;
00263   int IndexBase = 0;
00264 
00265   Epetra_Map m1(NumGlobalElements, IndexBase, Comm);
00266   int m1count = m1.ReferenceCount();
00267   const Epetra_BlockMapData* m1addr = m1.DataPtr();
00268   EPETRA_TEST_ERR(!(m1count==1),returnierr); // count should be 1
00269   if(verbose) cout << "Default constructor. \nm1= " << m1count << "  " << m1addr << endl;
00270   
00271   Epetra_Map* m2 = new Epetra_Map(m1);
00272   int m2count = m2->ReferenceCount();
00273   const Epetra_BlockMapData* m2addr = m2->DataPtr();
00274   int m1countold = m1count;
00275   m1count = m1.ReferenceCount();
00276   EPETRA_TEST_ERR(!(m2count==m1count && m1count==(m1countold+1)),returnierr); // both counts should be 2
00277   EPETRA_TEST_ERR(!(m1addr==m2addr),returnierr); // addresses should be same
00278   if(verbose) cout << "Copy constructor. \nm1= " << m1count << "  " << m1addr 
00279                    << "\nm2= " << m2count << "  " << m2addr << endl;
00280 
00281   delete m2;
00282   m1countold = m1count;
00283   m1count = m1.ReferenceCount();
00284   EPETRA_TEST_ERR(!(m1count == m1countold-1), returnierr); // count should have decremented (to 1)
00285   EPETRA_TEST_ERR(!(m1addr == m1.DataPtr()),returnierr); // m1addr should be unchanged
00286   if(verbose) cout << "m2 destroyed. \nm1= " << m1count << "  " << m1addr << endl;
00287 
00288   { // inside of braces to test stack deallocation.
00289     if(verbose) cout << "Assignment operator, post construction" << endl;
00290     Epetra_Map m3(NumGlobalElements, IndexBase+1, Comm);
00291     int m3count = m3.ReferenceCount();
00292     const Epetra_BlockMapData* m3addr = m3.DataPtr();
00293     EPETRA_TEST_ERR(!(m3count==1),returnierr); // m3count should be 1 initially
00294     EPETRA_TEST_ERR(!(m1addr!=m3addr),returnierr); // m1 and m3 should have different ptr addresses
00295     if(verbose) cout << "Prior to assignment: \nm1= " << m1count << "  " << m1addr 
00296                      << "\nm3= " << m3count << "  " << m3addr << endl;
00297     m3 = m1;
00298     m3count = m3.ReferenceCount();
00299     m3addr = m3.DataPtr();
00300     m1countold = m1count;
00301     m1count = m1.ReferenceCount();
00302     EPETRA_TEST_ERR(!(m3count==m1count && m1count==m1countold+1),returnierr); // both counts should be 2
00303     EPETRA_TEST_ERR(!(m1addr==m3addr),returnierr); // addresses should be same
00304     if(verbose) cout << "After assignment: \nm1= " << m1count << "  " << m1addr 
00305                      << "\nm3= " << m3count << "  " << m3addr << endl;
00306   }
00307   m1countold = m1count;
00308   m1count = m1.ReferenceCount();
00309   EPETRA_TEST_ERR(!(m1count==m1countold-1), returnierr); // count should have decremented (to 1)
00310   EPETRA_TEST_ERR(!(m1addr== m1.DataPtr()),returnierr); // m1addr should be unchanged
00311   if(verbose) cout << "m3 destroyed. \nm1= " << m1count << "  " << m1addr << endl;
00312 
00313   return(returnierr);
00314 }
00315 
00316 int checkLocalMapDataClass(Epetra_Comm& Comm, int verbose) {
00317   int returnierr = 0;
00318   int NumMyElements = 100;
00319   int IndexBase = 0;
00320 
00321   Epetra_LocalMap m1(NumMyElements, IndexBase, Comm);
00322   int m1count = m1.ReferenceCount();
00323   const Epetra_BlockMapData* m1addr = m1.DataPtr();
00324   EPETRA_TEST_ERR(!(m1count==1),returnierr); // count should be 1
00325   if(verbose) cout << "Default constructor. \nm1= " << m1count << "  " << m1addr << endl;
00326   
00327   Epetra_LocalMap* m2 = new Epetra_LocalMap(m1);
00328   int m2count = m2->ReferenceCount();
00329   const Epetra_BlockMapData* m2addr = m2->DataPtr();
00330   int m1countold = m1count;
00331   m1count = m1.ReferenceCount();
00332   EPETRA_TEST_ERR(!(m2count==m1count && m1count==(m1countold+1)),returnierr); // both counts should be 2
00333   EPETRA_TEST_ERR(!(m1addr==m2addr),returnierr); // addresses should be same
00334   if(verbose) cout << "Copy constructor. \nm1= " << m1count << "  " << m1addr 
00335                    << "\nm2= " << m2count << "  " << m2addr << endl;
00336 
00337   delete m2;
00338   m1countold = m1count;
00339   m1count = m1.ReferenceCount();
00340   EPETRA_TEST_ERR(!(m1count == m1countold-1), returnierr); // count should have decremented (to 1)
00341   EPETRA_TEST_ERR(!(m1addr == m1.DataPtr()),returnierr); // m1addr should be unchanged
00342   if(verbose) cout << "m2 destroyed. \nm1= " << m1count << "  " << m1addr << endl;
00343 
00344   { // inside of braces to test stack deallocation.
00345     if(verbose) cout << "Assignment operator, post construction" << endl;
00346     Epetra_LocalMap m3(NumMyElements, IndexBase+1, Comm);
00347     int m3count = m3.ReferenceCount();
00348     const Epetra_BlockMapData* m3addr = m3.DataPtr();
00349     EPETRA_TEST_ERR(!(m3count==1),returnierr); // m3count should be 1 initially
00350     EPETRA_TEST_ERR(!(m1addr!=m3addr),returnierr); // m1 and m3 should have different ptr addresses
00351     if(verbose) cout << "Prior to assignment: \nm1= " << m1count << "  " << m1addr 
00352                      << "\nm3= " << m3count << "  " << m3addr << endl;
00353     m3 = m1;
00354     m3count = m3.ReferenceCount();
00355     m3addr = m3.DataPtr(); // cast int* to int
00356     m1countold = m1count;
00357     m1count = m1.ReferenceCount();
00358     EPETRA_TEST_ERR(!(m3count==m1count && m1count==m1countold+1),returnierr); // both counts should be 2
00359     EPETRA_TEST_ERR(!(m1addr==m3addr),returnierr); // addresses should be same
00360     if(verbose) cout << "After assignment: \nm1= " << m1count << "  " << m1addr 
00361                      << "\nm3= " << m3count << "  " << m3addr << endl;
00362   }
00363   m1countold = m1count;
00364   m1count = m1.ReferenceCount();
00365   EPETRA_TEST_ERR(!(m1count==m1countold-1), returnierr); // count should have decremented (to 1)
00366   EPETRA_TEST_ERR(!(m1addr== m1.DataPtr()),returnierr); // m1addr should be unchanged
00367   if(verbose) cout << "m3 destroyed. \nm1= " << m1count << "  " << m1addr << endl;
00368 
00369   return(returnierr);
00370 }

Generated on Thu Sep 18 12:37:56 2008 for Epetra Package Browser (Single Doxygen Collection) by doxygen 1.3.9.1