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

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