test/Comm/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_Comm Test routine
00030 #include "../epetra_test_err.h"
00031 #include "Epetra_Time.h"
00032 #include "Epetra_Util.h"
00033 #include "Epetra_Distributor.h"
00034 #include "Epetra_SerialComm.h"
00035 #include "Epetra_IntSerialDenseVector.h"
00036 #include "Epetra_Version.h"
00037 
00038 #ifdef EPETRA_MPI
00039 #include <mpi.h>
00040 #include "Epetra_MpiComm.h"
00041 
00042 int checkMpiDataClass(bool verbose);
00043 #endif
00044 
00045 int checkSerialDataClass(bool verbose);
00046 int checkCommMethods(Epetra_Comm& petracomm,
00047                      bool verbose, bool verbose1,
00048                      int& NumProc, int& rank);
00049 int checkRankAndSize(Epetra_Comm& petracomm, bool verbose, int rank, int size);
00050 void checkBarrier(Epetra_Comm& petracomm, bool verbose, int rank);
00051 
00052 int checkDistributor(Epetra_Distributor* distr,
00053                      Epetra_Comm& Comm);
00054 
00055 int main(int argc, char* argv[]) {
00056   bool verbose = false;  // used to set verbose false on non-root processors
00057   bool verbose1 = false; // user's command-line argument
00058   // Check if we should print results to standard out
00059   if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose1 = true;
00060 
00061   int ierr = 0;
00062   int returnierr = 0;
00063   int size = 1;
00064   int rank = 0;
00065 
00066   if (verbose1)
00067     cout << Epetra_Version() << endl << endl;
00068 
00069   // Test Epetra_SerialComm
00070   if(verbose1) cout << "Testing Epetra_SerialComm..." << endl;
00071   Epetra_SerialComm serialcomm;
00072   if (verbose1) cout << serialcomm << endl;
00073   ierr = checkRankAndSize(serialcomm, verbose1, rank, size);
00074   EPETRA_TEST_ERR(ierr,returnierr);
00075   if (verbose1 && ierr==0) cout << "Checked OK\n\n" <<endl;
00076   // method testing
00077   int numProc = serialcomm.NumProc();
00078   ierr = checkCommMethods(serialcomm, verbose, verbose1, numProc, rank);
00079   EPETRA_TEST_ERR(ierr,returnierr);
00080   if (verbose1 && ierr==0) cout << "Checked OK\n\n" <<endl;
00081   // clone
00082   if(verbose1) cout << "SerialComm Clone.." << endl;
00083   Epetra_Comm* cloned_serialcomm = serialcomm.Clone();
00084   ierr = checkCommMethods(*cloned_serialcomm, verbose, verbose1, numProc, rank);
00085   delete cloned_serialcomm;
00086   EPETRA_TEST_ERR(ierr,returnierr);
00087   if (verbose1 && ierr==0) cout << "Checked OK\n\n" <<endl;
00088   // check inner data class
00089   ierr = checkSerialDataClass(verbose1);
00090   EPETRA_TEST_ERR(ierr,returnierr);
00091   if (verbose1 && ierr==0) cout << "Checked OK\n\n" <<endl;
00092 
00093   Epetra_Distributor* serialdistr = serialcomm.CreateDistributor();
00094   ierr = checkDistributor(serialdistr, serialcomm);
00095   delete serialdistr;
00096   EPETRA_TEST_ERR(ierr, returnierr);
00097 
00098   // Test Epetra_MpiComm
00099 #ifdef EPETRA_MPI
00100   // Initialize MPI
00101   if(verbose1) cout << "Testing Epetra_MpiComm..." << endl;
00102   MPI_Init(&argc,&argv);
00103   MPI_Comm_size(MPI_COMM_WORLD, &size);
00104   MPI_Comm_rank(MPI_COMM_WORLD, &rank);
00105   Epetra_MpiComm petracomm( MPI_COMM_WORLD );
00106   ierr = checkRankAndSize(petracomm, verbose1, rank, size);
00107   EPETRA_TEST_ERR(ierr,returnierr);
00108   if (verbose1 && ierr==0) cout << "Checked OK\n\n" <<endl;
00109   MPI_Comm MPIComm1 = petracomm.Comm();
00110   int size1, rank1;
00111   MPI_Comm_size(MPIComm1, &size1);
00112   MPI_Comm_rank(MPIComm1, &rank1);
00113   if (verbose1) cout << petracomm <<  ".  Using MPI_Comm from Petra_Comm:"
00114                            << " Processor " << rank1 << " of " << size1
00115                            << " (should be the same)." << endl;
00116   EPETRA_TEST_ERR(!(rank1==rank),ierr);
00117   EPETRA_TEST_ERR(!(size1==size),ierr);
00118   checkBarrier(petracomm, verbose1, rank);
00119 
00120   // method testing
00121   numProc = petracomm.NumProc();
00122   ierr = checkCommMethods(petracomm, verbose, verbose1, numProc, rank);
00123   EPETRA_TEST_ERR(ierr,returnierr);
00124   if (verbose1 && ierr==0) cout << "Checked OK\n\n" <<endl;
00125 
00126   // clone
00127   if(verbose1) cout << "MpiComm Clone.." << endl;
00128   Epetra_Comm* cloned_mpicomm = petracomm.Clone();
00129   ierr = checkCommMethods(*cloned_mpicomm, verbose, verbose1, numProc, rank);
00130   delete cloned_mpicomm;
00131   EPETRA_TEST_ERR(ierr,returnierr);
00132   if (verbose1 && ierr==0) cout << "Checked OK\n\n" <<endl;
00133 
00134   // check inner data class
00135   petracomm.Barrier();
00136   ierr = checkMpiDataClass(verbose1);
00137   EPETRA_TEST_ERR(ierr,returnierr);
00138   if (verbose1 && ierr==0) cout << "Checked OK\n\n" <<endl;
00139 
00140   Epetra_Distributor* plldistr = petracomm.CreateDistributor();
00141   ierr = checkDistributor(plldistr, petracomm);
00142   delete plldistr;
00143   EPETRA_TEST_ERR(ierr, returnierr);
00144 
00145   petracomm.Barrier();
00146   MPI_Finalize();
00147 #endif
00148 
00149   return(returnierr);
00150 }
00151 
00152 //=============================================================================
00153 void checkBarrier(Epetra_Comm& petracomm, bool verbose, int rank) {
00154   // Do some timing to test barrier function
00155   int MyPID = petracomm.MyPID();
00156   Epetra_Time before_barrier(petracomm);
00157   Epetra_Time after_barrier(petracomm);
00158   // Give each processor rank+1 amount of work
00159   // Time before barrier should increase roughly linearly
00160   // Time after barrier should be same for all processors
00161   double sum = 0.0;
00162   for (int j=0; j<rank+1; j++)
00163     for (int i=0; i<1000000; i++) 
00164       sum += ((double )rand())/((double) RAND_MAX);
00165   sum /= rank+1;
00166   if (verbose) cout << "Processor " << MyPID
00167                     << " Time to reach barrier: "
00168                     << before_barrier.ElapsedTime() << endl;
00169   petracomm.Barrier();
00170   if (verbose) cout << "Processor " << MyPID << " Sum result  "
00171                     << sum << " Time to get beyond barrier: "
00172                     << after_barrier.ElapsedTime() << endl;
00173  
00174   petracomm.Barrier();
00175 }
00176 
00177 //=============================================================================
00178 int checkRankAndSize(Epetra_Comm& petracomm, bool verbose, int rank, int size) {
00179   int ierr = 0;
00180   //if(verbose) cout << "CRS Breakpoint 1" << endl;
00181   int MyPID = petracomm.MyPID();
00182   //if(verbose) cout << "CRS Breakpoint 2" << endl;
00183   int NumProc = petracomm.NumProc();
00184   EPETRA_TEST_ERR(!(MyPID==rank),ierr);
00185   EPETRA_TEST_ERR(!(NumProc==size),ierr);
00186   petracomm.Barrier();
00187   return(ierr);
00188 }
00189 
00190 //=============================================================================
00191 int checkCommMethods(Epetra_Comm& petracomm, bool verbose, bool verbose1, int& NumProc, int& rank) {
00192   int i,j;
00193   int forierr = 0;
00194   int ierr = 0;
00195 
00196   verbose = (petracomm.MyPID() == 0);  // Turn verbose on; 
00197                                        // it is always false in main.
00198 
00199   // Some vars needed for the following tests
00200   int count = 4;
00201   int* iInputs = new int[count]; // General array for int type tests
00202   for (i=0; i<count; i++)
00203     iInputs[i] = 10*(i + rank - 2) + rank; 
00204     // if these values are changed, the expected maxs, mins, sums, etc must also change.  
00205     //NOTE: Broadcst() does not use these values.  The lines that need to be changed are located
00206     //in the "Values for ****** tests" sections directly below.
00207 
00208   double* dInputs = new double[count]; // General array for double type tests
00209   for (i=0; i<count; i++)
00210     dInputs[i] = pow(2.0,i-rank); 
00211     // if these values are changed, the expected maxs, mins, sums, etc must also change.  
00212     //NOTE: Broadcst() does not use these values.  The lines that need to be changed are located
00213     //in the "Values for ****** tests" sections directly below.
00214 
00215 
00216   // Values for Broadcast tests
00217   int* iVals = new int[count];
00218   if (rank == 0) {
00219     for (i=0; i<count; i++)
00220        iVals[i] = i; // if these values are changed, the values in iBVals must also be changed
00221   }
00222   
00223   int* iBVals = new int[count]; // Values to be checked against the values broadcast to the non root processors
00224   for (i=0; i<count; i++)
00225     iBVals[i] = i; // if these values are changed, the values in iVals must also be changed
00226   double* dVals = new double[count];
00227   if (rank == 0) {
00228      for (i=0; i<count; i++)
00229        dVals[i] = double(i); // if these values are changed, the values in dBVals must also be changed
00230   }    
00231   double* dBVals = new double[count]; // Values to be checked against the values broadcast to the non root processors
00232   for (i=0; i<count; i++)
00233     dBVals[i] = i; // if these values are changed, the values in dVals must also be changed
00234 
00235   const char *cConst = "Heidi, do you want a cookie?";
00236   int cCount = strlen(cConst)+1;
00237   char* cVals = new char[cCount];
00238   if (rank == 0) {
00239      strcpy(cVals, cConst);        // if these values are changed, 
00240      cVals[cCount-1] = '\0';       // the values in cBVals must also be changed
00241   }
00242   char* cBVals = new char[cCount]; // Values to be checked against the values 
00243                                    // broadcast to the non root processors
00244   strcpy(cBVals, cConst);          // if these values are changed, 
00245   cBVals[cCount-1] = '\0';         // the values in cVals must also be changed
00246 
00247   // Values for MaxAll tests
00248   int* iMyGlobalMaxs = new int[count];
00249   for (i=0; i<count; i++)
00250     iMyGlobalMaxs[i]=10 * (i + NumProc-1 -2) +  NumProc-1; // if these values are changed, iInput must be changed 
00251                                                            //as well as all other values dependent on iInput
00252   double* dMyGlobalMaxs = new double[count];
00253   for (i=0; i<count; i++)
00254     dMyGlobalMaxs[i]= pow(2.0,i); //if these values are changed, dInput must be changed 
00255                                   //as well as all other values dependent on dInput
00256 
00257 
00258   // Values for MinAll tests
00259   int* iMyGlobalMins = new int[count];
00260   for (i=0; i<count; i++)
00261     iMyGlobalMins[i]= 10 * (i - 2); //if these values are changed, iInput must be changed 
00262                                     //as well as all other values dependent on iInput
00263   double* dMyGlobalMins = new double[count];
00264   for (i=0; i<count; i++)
00265     dMyGlobalMins[i]= pow(2.0,i-(NumProc-1)); //if these values are changed, dInput must be changed 
00266                                               //as well as all other values dependent on dInput
00267 
00268 
00269   // Values for SumAll tests
00270   int* iMyGlobalSums = new int[count];
00271   for (i=0; i<count; i++){
00272     iMyGlobalSums[i]=0;
00273     for (j=0; j<NumProc; j++)
00274       iMyGlobalSums[i] += 10*(i+j-2) + j;// if these values are changed, iInput must be changed                                          
00275   }                                      //as well as all other values dependent on iInput
00276 
00277   double* dMyGlobalSums = new double[count];
00278   for (i=0; i<count; i++){
00279     dMyGlobalSums[i]=0;
00280     for (j=0; j<NumProc; j++)
00281       dMyGlobalSums[i] += pow(2.0,i-j);// if these values are changed, dInput must be changed 
00282   }                                    //as well as all other values dependent on dInput
00283 
00284 
00285   // Values for ScanSum tests
00286   int* iMyScanSums = new int[count];
00287   for (i=0; i<count; i++)
00288     iMyScanSums[i] = int((rank+1)*(10*(2*i+rank-4)+rank)*.5);// if these values are changed, 
00289                                                              //iInput must be changed as well as 
00290                                                              //all other values dependent on iInput
00291   double* dMyScanSums = new double[count];
00292   for (i=0; i<count; i++) {
00293     dMyScanSums[i] = 0;
00294     for (j=0; j<=rank; j++)
00295       dMyScanSums[i] += pow(2.0,i-j); //if these values are changed, dInput must be changed 
00296   }                                   //as well as all other values dependent on dInput
00297 
00298 
00299   // Values for Gather tests
00300   int totalVals = count*NumProc;
00301   int* iMyOrderedVals = new int[totalVals];
00302   double* dMyOrderedVals = new double[totalVals];
00303   int k=0;
00304   for (j=0; j<NumProc; j++) {
00305     for (i=0; i<count; i++) {
00306       iMyOrderedVals[k] = 10*(i + j - 2) + j;; // if these values are changed, iInput must be changed 
00307                                                //as well as all other values dependent on iInput
00308       dMyOrderedVals[k] = pow(2.0,i-j); // if these values are changed, dInput must be changed 
00309                                         //as well as all other values dependent on dInput
00310       k++;
00311     }
00312   }
00313   petracomm.Barrier();
00314 
00315 
00316   // Method testing section
00317   // Test the Broadcast functions
00318   EPETRA_TEST_ERR(petracomm.Broadcast(iVals,count,0),ierr);
00319   if (verbose1) {
00320     if (rank == 0) 
00321       cout << "The values on the root processor are: ";
00322     else
00323       cout << "The values on processor " << rank << " are: ";
00324     for (i=0; i<count; i++) 
00325       cout << iVals[i] << " ";
00326     cout << endl;
00327   }
00328   // ierr = 0; need to track errors the whole way through the file - this line of code seems like a bad idea 
00329   forierr = 0;
00330   for (i=0; i<count; i++)
00331     forierr += !(iVals[i] == iBVals[i]); // otherwise Broadcast didn't occur properly
00332   EPETRA_TEST_ERR(forierr,ierr);
00333   delete [] iVals;
00334   delete [] iBVals;
00335   petracomm.Barrier();
00336   if (verbose) cout << endl << "Broadcast (type int) test passed!" << endl << endl;// If test gets to here the test passed, 
00337                                                                                    //only output on one node
00338   petracomm.Barrier();
00339 
00340   EPETRA_TEST_ERR(petracomm.Broadcast(dVals,count,0),ierr);
00341   if (verbose1) {
00342     if (rank == 0)
00343       cout << "The values on the root processor are: ";
00344     else
00345       cout << "The values on processor " << rank << " are: ";
00346     for (i=0; i<count; i++) 
00347       cout << dVals[i] << " ";
00348     cout << endl;
00349   }
00350   forierr = 0;
00351   for (i=0; i<count; i++)
00352     forierr += !(dVals[i] == dBVals[i]); // otherwise Broadcast didn't occur properly
00353   EPETRA_TEST_ERR(forierr,ierr);
00354   delete [] dVals;
00355   delete [] dBVals;
00356   petracomm.Barrier();
00357   if (verbose) cout << endl << "Broadcast (type double) test passed!" << endl << endl;// If test gets to here the test passed, 
00358                                                                                       //only output on one node
00359   petracomm.Barrier();
00360 
00361   EPETRA_TEST_ERR(petracomm.Broadcast(cVals,cCount,0),ierr);
00362   if (verbose1) {
00363     if (rank == 0)
00364       cout << "The values on the root processor are: " << cVals << endl;
00365     else
00366       cout << "The values on processor " << rank << " are: " << cVals << endl;
00367   }
00368   forierr = 0;
00369   for (i=0; i<cCount; i++)
00370     forierr += !(cVals[i] == cBVals[i]); // otherwise Broadcast didn't work.
00371   EPETRA_TEST_ERR(forierr,ierr);
00372   delete [] cVals;
00373   delete [] cBVals;
00374   petracomm.Barrier();
00375   // If test gets to here the test passed, 
00376   if (verbose) 
00377     cout << endl << "Broadcast (type char) test passed!" << endl << endl;
00378                                                                                       //only output on one node
00379   petracomm.Barrier();
00380 
00381  // Test the MaxAll functions
00382   int* iGlobalMaxs = new int[count];
00383   if (verbose1) {
00384     cout << "The values on processor " << rank << " are: ";
00385     for (i=0; i<count; i++) 
00386       cout << iInputs[i] << " ";
00387     cout << endl;
00388   }
00389   EPETRA_TEST_ERR(petracomm.MaxAll(iInputs,iGlobalMaxs,count),ierr);
00390   petracomm.Barrier();
00391   
00392   if (verbose1) {
00393     cout << "The max values according to processor " << rank << " are: ";
00394     for (i=0; i<count; i++) 
00395       cout << iGlobalMaxs[i] << " ";
00396     cout << endl;
00397   }
00398   forierr = 0;
00399   for (i=0; i<count; i++) 
00400     forierr += !(iMyGlobalMaxs[i] == iGlobalMaxs[i]);
00401   EPETRA_TEST_ERR(forierr,ierr);
00402 
00403   delete [] iGlobalMaxs;
00404   delete [] iMyGlobalMaxs;
00405   petracomm.Barrier();
00406   if (verbose) cout << endl << "MaxAll (type int) test passed!" << endl << endl;// If test gets to here the test passed, 
00407                                                                                 //only output on one node
00408   petracomm.Barrier();
00409 
00410   double* dGlobalMaxs = new double[count];
00411   if (verbose1) {
00412     cout << "The values on processor " << rank << " are: ";
00413     for (i=0; i<count; i++) 
00414       cout << dInputs[i] << " ";
00415     cout << endl;
00416   }
00417   EPETRA_TEST_ERR(petracomm.MaxAll(dInputs,dGlobalMaxs,count),ierr);
00418   petracomm.Barrier();
00419   
00420   if (verbose1) {
00421     cout << "The max values according to processor " << rank << " are: ";
00422     for (i=0; i<count; i++) 
00423       cout << dGlobalMaxs[i] << " ";
00424     cout << endl;
00425   }
00426   forierr = 0;
00427   for (i=0; i<count; i++)
00428     forierr += !(Epetra_Util::Chop(dMyGlobalMaxs[i] - dGlobalMaxs[i]) == 0);
00429   EPETRA_TEST_ERR(forierr,ierr);
00430   delete [] dGlobalMaxs;
00431   delete [] dMyGlobalMaxs;
00432   petracomm.Barrier();
00433   if (verbose) cout << endl << "MaxAll (type double) test passed!" << endl << endl;// If test gets to here the test passed, 
00434                                                                                    //only output on one node
00435   petracomm.Barrier();
00436 
00437 
00438  // Test the MinAll functions
00439   int* iGlobalMins = new int[count];
00440   if (verbose1) {
00441     cout << "The values on processor " << rank << " are: ";
00442     for (i=0; i<count; i++) 
00443       cout << iInputs[i] << " ";
00444     cout << endl;
00445   }
00446   EPETRA_TEST_ERR(petracomm.MinAll(iInputs,iGlobalMins,count),ierr);
00447   petracomm.Barrier();
00448 
00449   if (verbose1) {
00450     cout << "The min values according to processor " << rank << " are: ";
00451     for (i=0; i<count; i++) 
00452       cout << iGlobalMins[i] << " ";
00453     cout << endl;
00454   }
00455   forierr = 0;
00456   for (i=0; i<count; i++) 
00457     forierr += !(iMyGlobalMins[i] == iGlobalMins[i]); // otherwise calculated min is wrong
00458   EPETRA_TEST_ERR(forierr,ierr);
00459   delete [] iGlobalMins;
00460   delete [] iMyGlobalMins;
00461   petracomm.Barrier();
00462   if (verbose) cout << endl << "MinAll (type int) test passed!" << endl << endl;// If test gets to here the test passed, 
00463                                                                                 //only output on one node
00464   petracomm.Barrier();
00465 
00466   double* dGlobalMins = new double[count];
00467   if (verbose1) {
00468     cout << "The values on processor " << rank << " are: ";
00469     for (i=0; i<count; i++) 
00470       cout << dInputs[i] << " ";
00471     cout << endl;
00472   }
00473   EPETRA_TEST_ERR(petracomm.MinAll(dInputs,dGlobalMins,count),ierr);
00474   petracomm.Barrier();
00475 
00476   if (verbose1) {
00477     cout << "The min values according to processor " << rank << " are: ";
00478     for (i=0; i<count; i++) 
00479       cout << dGlobalMins[i] << " ";
00480     cout << endl;
00481   }
00482   forierr = 0;
00483   for (i=0; i<count; i++)
00484     forierr += !(Epetra_Util::Chop(dMyGlobalMins[i] - dGlobalMins[i]) == 0); // otherwise calculated min is wrong
00485   EPETRA_TEST_ERR(forierr,ierr);
00486   delete [] dGlobalMins;
00487   delete [] dMyGlobalMins;
00488   petracomm.Barrier();
00489   if (verbose) cout << endl << "MinAll (type double) test passed!" << endl << endl;// If test gets to here the test passed, 
00490                                                                                    //only output on one node
00491   petracomm.Barrier();
00492 
00493 
00494  // Test the SumAll functions
00495   int* iGlobalSums = new int[count];
00496   if (verbose1) {
00497     cout << "The values on processor " << rank << " are: ";
00498     for (i=0; i<count; i++) 
00499       cout << iInputs[i] << " ";
00500     cout << endl;
00501   }
00502   EPETRA_TEST_ERR(petracomm.SumAll(iInputs,iGlobalSums,count),ierr);
00503   petracomm.Barrier();
00504 
00505   if (verbose1) {
00506     cout << "The sums of all values according to processor " << rank << " are: ";
00507     for (i=0; i<count; i++) 
00508       cout << iGlobalSums[i] << " ";
00509     cout << endl;
00510   }
00511   forierr = 0;
00512   for (i=0; i<count; i++)
00513     forierr += !(iMyGlobalSums[i] == iGlobalSums[i]); // otherwise calculated sum is wrong
00514   EPETRA_TEST_ERR(forierr,ierr);
00515   delete [] iGlobalSums;
00516   delete [] iMyGlobalSums;
00517   petracomm.Barrier();
00518   if (verbose) cout << endl << "SumAll (type int) test passed!" << endl << endl;// If test gets to here the test passed, 
00519                                                                                 //only output on one node
00520   petracomm.Barrier();
00521 
00522   double* dGlobalSums = new double[count];
00523   if (verbose1) {
00524     cout << "The values on processor " << rank << " are: ";
00525     for (i=0; i<count; i++) 
00526       cout << dInputs[i] << " ";
00527     cout << endl;
00528   }
00529   EPETRA_TEST_ERR(petracomm.SumAll(dInputs,dGlobalSums,count),ierr);
00530   petracomm.Barrier();
00531 
00532   if (verbose1) {
00533     cout << "The sums of all values according to processor " << rank << " are: ";
00534     for (i=0; i<count; i++) 
00535       cout << dGlobalSums[i] << " ";
00536     cout << endl;
00537   }
00538   forierr = 0;
00539   for (i=0; i<count; i++)
00540     forierr += !(Epetra_Util::Chop(dMyGlobalSums[i] - dGlobalSums[i]) == 0); // otherwise calculated sum is wrong
00541   EPETRA_TEST_ERR(forierr,ierr);
00542 
00543   delete [] dGlobalSums;
00544   delete [] dMyGlobalSums;
00545   petracomm.Barrier();
00546   if (verbose) cout << endl << "SumAll (type double) test passed!" << endl << endl;// If test gets to here the test passed, 
00547                                                                                    //only output on one node
00548   petracomm.Barrier();
00549 
00550 
00551  // Test the ScanSum functions
00552   int* iScanSums = new int[count];
00553   if (verbose1) {
00554     cout << "The values on processor " << rank << " are: ";
00555     for (i=0; i<count; i++) 
00556       cout << iInputs[i] << " ";
00557     cout << endl;
00558   }
00559   
00560   EPETRA_TEST_ERR(petracomm.ScanSum(iInputs,iScanSums,count),ierr);
00561   petracomm.Barrier();
00562 
00563   if (verbose1) {
00564     cout << "The sums of all values on processors 0 - " << rank << " are: ";
00565     for (i=0; i<count; i++) {
00566       cout << iScanSums[i] << " ";
00567     }
00568     cout << endl;
00569   }
00570   forierr = 0;
00571   for (i=0; i<count; i++)
00572     forierr += !(iMyScanSums[i] == iScanSums[i]);
00573   EPETRA_TEST_ERR(forierr,ierr);
00574   delete [] iScanSums;
00575   delete [] iMyScanSums;
00576   petracomm.Barrier();
00577   if (verbose) cout << endl << "ScanSum (type int) test passed!" << endl << endl;// If test gets to here the test passed, 
00578                                                                                  //only output on one node
00579   petracomm.Barrier();
00580 
00581   double* dScanSums = new double[count];
00582   if (verbose1) {
00583     cout << "The values on processor " << rank << " are: ";
00584     for (i=0; i<count; i++)
00585       cout << dInputs[i] << " ";
00586     cout << endl;
00587   }
00588 
00589   EPETRA_TEST_ERR(petracomm.ScanSum(dInputs,dScanSums,count),ierr);
00590   petracomm.Barrier();
00591 
00592   if (verbose1) {
00593     cout << "The sums of all values on processors 0 - " << rank << " are: ";
00594     for (i=0; i<count; i++) {
00595       cout << dScanSums[i] << " ";
00596     }
00597     cout << endl;
00598   }
00599   forierr = 0;
00600   for (i=0; i<count; i++)
00601     forierr += !(Epetra_Util::Chop(dMyScanSums[i] - dScanSums[i])== 0);
00602   EPETRA_TEST_ERR(forierr,ierr);
00603   delete [] dScanSums;
00604   delete [] dMyScanSums;
00605   petracomm.Barrier();
00606   if (verbose) cout << endl << "ScanSum (type double) test passed!" << endl << endl;// If test gets to here the test passed, 
00607                                                                                     //only output on one node
00608   petracomm.Barrier();
00609 
00610 
00611  // Test the Gather functions
00612   int* iOrderedVals = new int[totalVals];
00613   if (verbose1) {
00614     cout << "The values on processor " << rank << " are: ";
00615     for (i=0; i<count; i++) 
00616       cout << iInputs[i] << " ";
00617     cout << endl;
00618   }
00619   
00620   EPETRA_TEST_ERR(petracomm.GatherAll(iInputs,iOrderedVals,count),ierr);
00621   petracomm.Barrier();
00622 
00623   if (verbose1) {
00624     cout << "The combined list of all values from all processors according to processor " << rank << " is: ";
00625     for (i=0; i<totalVals; i++) {
00626       cout << iOrderedVals[i] << " ";
00627     }
00628     cout << endl;
00629   }
00630   forierr = 0;
00631   for (i=0; i<totalVals; i++)
00632     forierr += !(iMyOrderedVals[i] == iOrderedVals[i]);
00633   EPETRA_TEST_ERR(forierr,ierr);
00634   delete [] iOrderedVals;
00635   delete [] iMyOrderedVals;
00636   petracomm.Barrier();
00637   if (verbose) cout << endl << "GatherAll (type int) test passed!" << endl << endl;// If test gets to here the test passed, 
00638                                                                                    //only output on one node
00639   petracomm.Barrier();
00640 
00641   double* dOrderedVals = new double[totalVals];
00642   if (verbose1) {
00643     cout << "The values on processor " << rank << " are: ";
00644     for (i=0; i<count; i++) 
00645       cout << dInputs[i] << " ";
00646     cout << endl;
00647   }
00648   
00649   EPETRA_TEST_ERR(petracomm.GatherAll(dInputs,dOrderedVals,count),ierr);
00650   petracomm.Barrier();
00651 
00652   if (verbose1) {
00653     cout << "The combined list of all values from all processors according to processor " << rank << " is: ";
00654     for (i=0; i<totalVals; i++) {
00655       cout << dOrderedVals[i] << " ";
00656     }
00657     cout << endl;
00658   }
00659   forierr = 0;
00660   for (i=0; i<totalVals; i++)
00661     forierr += !(Epetra_Util::Chop(dMyOrderedVals[i] - dOrderedVals[i]) == 0);
00662   EPETRA_TEST_ERR(forierr,ierr);
00663   delete [] dOrderedVals;
00664   delete [] dMyOrderedVals;
00665   petracomm.Barrier();
00666   if (verbose) cout << endl << "GatherAll (type double) test passed!" << endl << endl;// If test gets to here the test passed, 
00667                                                                                       //only output on one node
00668   petracomm.Barrier();
00669 
00670   delete[] dInputs;
00671   delete[] iInputs;
00672 
00673   return(ierr);
00674 }
00675 
00676 //=============================================================================
00677 int checkSerialDataClass(bool verbose) {
00678   int ierr = 0;
00679   if(verbose) cout << "Testing Reference Counting... ";               
00680   Epetra_SerialComm c1;
00681   int c1count = c1.ReferenceCount();
00682   const Epetra_SerialCommData* c1addr = c1.DataPtr();
00683   EPETRA_TEST_ERR(!(c1count==1),ierr); // count should be 1
00684   if(verbose) cout << "Default constructor. \nc1= " << c1count << "  " << c1addr << endl;
00685 
00686   Epetra_SerialComm* c2 = new Epetra_SerialComm(c1);
00687   int c2count = c2->ReferenceCount();
00688   const Epetra_SerialCommData* c2addr = c2->DataPtr();
00689   int c1countold = c1count;
00690   c1count = c1.ReferenceCount();
00691   EPETRA_TEST_ERR(!(c2count==c1count && c1count==(c1countold+1)),ierr); // both counts should be 2
00692   EPETRA_TEST_ERR(!(c1addr==c2addr),ierr); // addresses should be same
00693   if(verbose) cout << "Copy constructor(heap). \nc1= " << c1count << "  " << c1addr 
00694                     << "\nc2= " << c2count << "  " << c2addr << endl;
00695   delete c2;
00696   c1countold = c1count;
00697   c1count = c1.ReferenceCount();
00698   EPETRA_TEST_ERR(!(c1count==c1countold-1),ierr); // count should have decremented (to 1)
00699   EPETRA_TEST_ERR(!(c1addr==c1.DataPtr()),ierr); // c1addr should be unchanged
00700   if(verbose) cout << "c2 Deleted. \nc1= " << c1count << "  " << c1addr << endl;
00701   { // inside own set of brackets so that c2a will be automatically at end of brackets
00702     // so that we can test to make sure objects on the stack deallocate correctly
00703     Epetra_SerialComm c2a(c1);
00704     c2count = c2a.ReferenceCount();
00705     c2addr = c2a.DataPtr();
00706     c1countold = c1count;
00707     c1count = c1.ReferenceCount();
00708     EPETRA_TEST_ERR(!(c2count==c1count && c1count==c1countold+1),ierr); // both counts should be 2
00709     EPETRA_TEST_ERR(!(c1addr==c2addr),ierr); // addresses should be same
00710     if(verbose) cout << "Copy constructor(stack). \nc1= " << c1count << "  " << c1addr 
00711                       << "\nc2a= " << c2count << "  " << c2addr << endl;
00712   }
00713   c1countold = c1count;
00714   c1count = c1.ReferenceCount();
00715   EPETRA_TEST_ERR(!(c1count==c1countold-1),ierr); // count should have decremented (to 1)
00716   EPETRA_TEST_ERR(!(c1addr==c1.DataPtr()),ierr); // c1addr should be unchanged
00717   if(verbose) cout << "c2a Destroyed. \nc1= " << c1count << "  " << c1addr << endl;
00718   if(verbose) cout << "Assignment operator, post construction" << endl;
00719   Epetra_SerialComm c3;
00720   int c3count = c3.ReferenceCount();
00721   const Epetra_SerialCommData* c3addr = c3.DataPtr();
00722   EPETRA_TEST_ERR(!(c3count==1),ierr); // c3count should be 1 initially
00723   EPETRA_TEST_ERR(!(c1addr!=c3addr),ierr); // c1 and c3 should have different ptr addresses
00724   if(verbose)cout << "Prior to assignment: \nc1=" << c1count << "  " << c1addr 
00725                    << "\nc3=" << c3count << "  " << c3addr << endl;
00726   c3 = c1;
00727   c3count = c3.ReferenceCount();
00728   c3addr = c3.DataPtr();
00729   c1countold = c1count;
00730   c1count = c1.ReferenceCount();
00731   EPETRA_TEST_ERR(!(c3count==c1count && c1count==c1countold+1),ierr); // both counts should be 2
00732   EPETRA_TEST_ERR(!(c1addr==c3addr),ierr); // addresses should be same
00733   if(verbose)cout << "After assignment: \nc1=" << c1count << "  " << c1addr 
00734                    << "\nc3=" << c3count << "  " << c3addr << endl;
00735   return(ierr);
00736 }
00737 
00738 //=============================================================================
00739 #ifdef EPETRA_MPI
00740 int checkMpiDataClass(bool verbose) {
00741   int ierr = 0;
00742   if(verbose) cout << "Testing Reference Counting... ";               
00743   Epetra_MpiComm c1( MPI_COMM_WORLD );
00744   int c1count = c1.ReferenceCount();
00745   const Epetra_MpiCommData* c1addr = c1.DataPtr();
00746   EPETRA_TEST_ERR(!(c1count==1),ierr); // count should be 1
00747   if(verbose) cout << "Default constructor. \nc1= " << c1count << "  " << c1addr << endl;
00748 
00749   Epetra_MpiComm* c2 = new Epetra_MpiComm(c1);
00750   int c2count = c2->ReferenceCount();
00751   const Epetra_MpiCommData* c2addr = c2->DataPtr();
00752   int c1countold = c1count;
00753   c1count = c1.ReferenceCount();
00754   EPETRA_TEST_ERR(!(c2count==c1count && c1count==(c1countold+1)),ierr); // both counts should be 2
00755   EPETRA_TEST_ERR(!(c1addr==c2addr),ierr); // addresses should be same
00756   if(verbose) cout << "Copy constructor(heap). \nc1= " << c1count << "  " << c1addr 
00757                     << "\nc2= " << c2count << "  " << c2addr << endl;
00758   delete c2;
00759   c1countold = c1count;
00760   c1count = c1.ReferenceCount();
00761   EPETRA_TEST_ERR(!(c1count==c1countold-1),ierr); // count should have decremented (to 1)
00762   EPETRA_TEST_ERR(!(c1addr==c1.DataPtr()),ierr); // c1addr should be unchanged
00763   if(verbose) cout << "c2 Deleted. \nc1= " << c1count << "  " << c1addr << endl;
00764   { // inside own set of brackets so that c2a will be automatically at end of brackets
00765     // so that we can test to make sure objects on the stack deallocate correctly
00766     Epetra_MpiComm c2a(c1);
00767     c2count = c2a.ReferenceCount();
00768     c2addr = c2a.DataPtr();
00769     c1countold = c1count;
00770     c1count = c1.ReferenceCount();
00771     EPETRA_TEST_ERR(!(c2count==c1count && c1count==c1countold+1),ierr); // both counts should be 2
00772     EPETRA_TEST_ERR(!(c1addr==c2addr),ierr); // addresses should be same
00773     if(verbose) cout << "Copy constructor(stack). \nc1= " << c1count << "  " << c1addr 
00774                       << "\nc2a= " << c2count << "  " << c2addr << endl;
00775   }
00776   c1countold = c1count;
00777   c1count = c1.ReferenceCount();
00778   EPETRA_TEST_ERR(!(c1count==c1countold-1),ierr); // count should have decremented (to 1)
00779   EPETRA_TEST_ERR(!(c1addr==c1.DataPtr()),ierr); // c1addr should be unchanged
00780   if(verbose) cout << "c2a Destroyed. \nc1= " << c1count << "  " << c1addr << endl;
00781   if(verbose) cout << "Assignment operator, post construction" << endl;
00782   Epetra_MpiComm c3( MPI_COMM_WORLD );
00783   int c3count = c3.ReferenceCount();
00784   const Epetra_MpiCommData* c3addr = c3.DataPtr();
00785   EPETRA_TEST_ERR(!(c3count==1),ierr); // c3count should be 1 initially
00786   EPETRA_TEST_ERR(!(c1addr!=c3addr),ierr); // c1 and c3 should have different ptr addresses
00787   if(verbose)cout << "Prior to assignment: \nc1=" << c1count << "  " << c1addr 
00788                    << "\nc3=" << c3count << "  " << c3addr << endl;
00789   c3 = c1;
00790   c3count = c3.ReferenceCount();
00791   c3addr = c3.DataPtr();
00792   c1countold = c1count;
00793   c1count = c1.ReferenceCount();
00794   EPETRA_TEST_ERR(!(c3count==c1count && c1count==c1countold+1),ierr); // both counts should be 2
00795   EPETRA_TEST_ERR(!(c1addr==c3addr),ierr); // addresses should be same
00796   if(verbose)cout << "After assignment: \nc1=" << c1count << "  " << c1addr 
00797                    << "\nc3=" << c3count << "  " << c3addr << endl;
00798   return(ierr);
00799 }
00800 #endif
00801 
00802 int checkDistributor(Epetra_Distributor* distr,
00803                      Epetra_Comm& Comm)
00804 {
00805   int numprocs = Comm.NumProc();
00806 
00807   int numExportIDs = numprocs;
00808   int* exportPIDs = new int[numExportIDs];
00809   for(int p=0; p<numExportIDs; ++p) {
00810     exportPIDs[p] = p;
00811   }
00812 
00813   bool deterministic = true;
00814   int numRemoteIDs = 0;
00815 
00816   int err = distr->CreateFromSends(numExportIDs, exportPIDs,
00817                                    deterministic, numRemoteIDs);
00818 
00819   //numRemoteIDs should equal numExportIDs.
00820 
00821   int returnValue = numRemoteIDs == numExportIDs ? 0 : -99;
00822 
00823   delete [] exportPIDs;
00824 
00825   if (returnValue + err != 0) {
00826     return(returnValue + err);
00827   }
00828 
00829   int* exportIDs = new int[numExportIDs];
00830   for(int i=0; i<numExportIDs; ++i) {
00831     exportIDs[i] = i+1;
00832   }
00833 
00834   int len_imports = 0;
00835   char* imports = NULL;
00836 
00837   err = distr->Do((char*)exportIDs, sizeof(int),
00838                   len_imports, imports);
00839 
00840   delete [] exportIDs;
00841 
00842   if (len_imports > 0) {
00843     delete [] imports;
00844   }
00845 
00846   return(err);
00847 }
00848 
00849 /*
00850   end of file cxx_main.cpp
00851 */

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