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   // Some vars needed for the following tests
00197   int count = 4;
00198   int* iInputs = new int[count]; // General array for int type tests
00199   for (i=0; i<count; i++)
00200     iInputs[i] = 10*(i + rank - 2) + rank; 
00201     // if these values are changed, the expected maxs, mins, sums, etc must also change.  
00202     //NOTE: Broadcst() does not use these values.  The lines that need to be changed are located
00203     //in the "Values for ****** tests" sections directly below.
00204 
00205   double* dInputs = new double[count]; // General array for double type tests
00206   for (i=0; i<count; i++)
00207     dInputs[i] = pow(2.0,i-rank); 
00208     // if these values are changed, the expected maxs, mins, sums, etc must also change.  
00209     //NOTE: Broadcst() does not use these values.  The lines that need to be changed are located
00210     //in the "Values for ****** tests" sections directly below.
00211 
00212 
00213   // Values for Broadcast tests
00214   int* iVals = new int[count];
00215   if (rank == 0) {
00216     for (i=0; i<count; i++)
00217        iVals[i] = i; // if these values are changed, the values in iBVals must also be changed
00218   }
00219   
00220   int* iBVals = new int[count]; // Values to be checked against the values broadcast to the non root processors
00221   for (i=0; i<count; i++)
00222     iBVals[i] = i; // if these values are changed, the values in iVals must also be changed
00223   double* dVals = new double[count];
00224   if (rank == 0) {
00225      for (i=0; i<count; i++)
00226        dVals[i] = double(i); // if these values are changed, the values in dBVals must also be changed
00227   }    
00228   double* dBVals = new double[count]; // Values to be checked against the values broadcast to the non root processors
00229   for (i=0; i<count; i++)
00230     dBVals[i] = i; // if these values are changed, the values in dVals must also be changed
00231 
00232 
00233   // Values for MaxAll tests
00234   int* iMyGlobalMaxs = new int[count];
00235   for (i=0; i<count; i++)
00236     iMyGlobalMaxs[i]=10 * (i + NumProc-1 -2) +  NumProc-1; // if these values are changed, iInput must be changed 
00237                                                            //as well as all other values dependent on iInput
00238   double* dMyGlobalMaxs = new double[count];
00239   for (i=0; i<count; i++)
00240     dMyGlobalMaxs[i]= pow(2.0,i); //if these values are changed, dInput must be changed 
00241                                   //as well as all other values dependent on dInput
00242 
00243 
00244   // Values for MinAll tests
00245   int* iMyGlobalMins = new int[count];
00246   for (i=0; i<count; i++)
00247     iMyGlobalMins[i]= 10 * (i - 2); //if these values are changed, iInput must be changed 
00248                                     //as well as all other values dependent on iInput
00249   double* dMyGlobalMins = new double[count];
00250   for (i=0; i<count; i++)
00251     dMyGlobalMins[i]= pow(2.0,i-(NumProc-1)); //if these values are changed, dInput must be changed 
00252                                               //as well as all other values dependent on dInput
00253 
00254 
00255   // Values for SumAll tests
00256   int* iMyGlobalSums = new int[count];
00257   for (i=0; i<count; i++){
00258     iMyGlobalSums[i]=0;
00259     for (j=0; j<NumProc; j++)
00260       iMyGlobalSums[i] += 10*(i+j-2) + j;// if these values are changed, iInput must be changed                                          
00261   }                                      //as well as all other values dependent on iInput
00262 
00263   double* dMyGlobalSums = new double[count];
00264   for (i=0; i<count; i++){
00265     dMyGlobalSums[i]=0;
00266     for (j=0; j<NumProc; j++)
00267       dMyGlobalSums[i] += pow(2.0,i-j);// if these values are changed, dInput must be changed 
00268   }                                    //as well as all other values dependent on dInput
00269 
00270 
00271   // Values for ScanSum tests
00272   int* iMyScanSums = new int[count];
00273   for (i=0; i<count; i++)
00274     iMyScanSums[i] = int((rank+1)*(10*(2*i+rank-4)+rank)*.5);// if these values are changed, 
00275                                                              //iInput must be changed as well as 
00276                                                              //all other values dependent on iInput
00277   double* dMyScanSums = new double[count];
00278   for (i=0; i<count; i++) {
00279     dMyScanSums[i] = 0;
00280     for (j=0; j<=rank; j++)
00281       dMyScanSums[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 Gather tests
00286   int totalVals = count*NumProc;
00287   int* iMyOrderedVals = new int[totalVals];
00288   double* dMyOrderedVals = new double[totalVals];
00289   int k=0;
00290   for (j=0; j<NumProc; j++) {
00291     for (i=0; i<count; i++) {
00292       iMyOrderedVals[k] = 10*(i + j - 2) + j;; // if these values are changed, iInput must be changed 
00293                                                //as well as all other values dependent on iInput
00294       dMyOrderedVals[k] = pow(2.0,i-j); // if these values are changed, dInput must be changed 
00295                                         //as well as all other values dependent on dInput
00296       k++;
00297     }
00298   }
00299   petracomm.Barrier();
00300 
00301 
00302   // Method testing section
00303   // Test the Broadcast functions
00304   EPETRA_TEST_ERR(petracomm.Broadcast(iVals,count,0),ierr);
00305   if (verbose1) {
00306     if (rank == 0) 
00307       cout << "The values on the root processor are: ";
00308     else
00309       cout << "The values on processor " << rank << " are: ";
00310     for (i=0; i<count; i++) 
00311       cout << iVals[i] << " ";
00312     cout << endl;
00313   }
00314   // ierr = 0; need to track errors the whole way through the file - this line of code seems like a bad idea 
00315   forierr = 0;
00316   for (i=0; i<count; i++)
00317     forierr += !(iVals[i] == iBVals[i]); // otherwise Broadcast didn't occur properly
00318   EPETRA_TEST_ERR(forierr,ierr);
00319   delete [] iVals;
00320   delete [] iBVals;
00321   petracomm.Barrier();
00322   if (verbose) cout << endl << "Broadcast (type int) test passed!" << endl << endl;// If test gets to here the test passed, 
00323                                                                                    //only output on one node
00324   petracomm.Barrier();
00325 
00326   EPETRA_TEST_ERR(petracomm.Broadcast(dVals,count,0),ierr);
00327   if (verbose1) {
00328     if (rank == 0)
00329       cout << "The values on the root processor are: ";
00330     else
00331       cout << "The values on processor " << rank << " are: ";
00332     for (i=0; i<count; i++) 
00333       cout << dVals[i] << " ";
00334     cout << endl;
00335   }
00336   forierr = 0;
00337   for (i=0; i<count; i++)
00338     forierr += !(dVals[i] == dBVals[i]); // otherwise Broadcast didn't occur properly
00339   EPETRA_TEST_ERR(forierr,ierr);
00340   delete [] dVals;
00341   delete [] dBVals;
00342   petracomm.Barrier();
00343   if (verbose) cout << endl << "Broadcast (type double) test passed!" << endl << endl;// If test gets to here the test passed, 
00344                                                                                       //only output on one node
00345   petracomm.Barrier();
00346 
00347 
00348  // Test the MaxAll functions
00349   int* iGlobalMaxs = new int[count];
00350   if (verbose1) {
00351     cout << "The values on processor " << rank << " are: ";
00352     for (i=0; i<count; i++) 
00353       cout << iInputs[i] << " ";
00354     cout << endl;
00355   }
00356   EPETRA_TEST_ERR(petracomm.MaxAll(iInputs,iGlobalMaxs,count),ierr);
00357   petracomm.Barrier();
00358   
00359   if (verbose1) {
00360     cout << "The max values according to processor " << rank << " are: ";
00361     for (i=0; i<count; i++) 
00362       cout << iGlobalMaxs[i] << " ";
00363     cout << endl;
00364   }
00365   forierr = 0;
00366   for (i=0; i<count; i++) 
00367     forierr += !(iMyGlobalMaxs[i] == iGlobalMaxs[i]);
00368   EPETRA_TEST_ERR(forierr,ierr);
00369 
00370   delete [] iGlobalMaxs;
00371   delete [] iMyGlobalMaxs;
00372   petracomm.Barrier();
00373   if (verbose) cout << endl << "MaxAll (type int) test passed!" << endl << endl;// If test gets to here the test passed, 
00374                                                                                 //only output on one node
00375   petracomm.Barrier();
00376 
00377   double* dGlobalMaxs = new double[count];
00378   if (verbose1) {
00379     cout << "The values on processor " << rank << " are: ";
00380     for (i=0; i<count; i++) 
00381       cout << dInputs[i] << " ";
00382     cout << endl;
00383   }
00384   EPETRA_TEST_ERR(petracomm.MaxAll(dInputs,dGlobalMaxs,count),ierr);
00385   petracomm.Barrier();
00386   
00387   if (verbose1) {
00388     cout << "The max values according to processor " << rank << " are: ";
00389     for (i=0; i<count; i++) 
00390       cout << dGlobalMaxs[i] << " ";
00391     cout << endl;
00392   }
00393   forierr = 0;
00394   for (i=0; i<count; i++)
00395     forierr += !(Epetra_Util::Chop(dMyGlobalMaxs[i] - dGlobalMaxs[i]) == 0);
00396   EPETRA_TEST_ERR(forierr,ierr);
00397   delete [] dGlobalMaxs;
00398   delete [] dMyGlobalMaxs;
00399   petracomm.Barrier();
00400   if (verbose) cout << endl << "MaxAll (type double) test passed!" << endl << endl;// If test gets to here the test passed, 
00401                                                                                    //only output on one node
00402   petracomm.Barrier();
00403 
00404 
00405  // Test the MinAll functions
00406   int* iGlobalMins = new int[count];
00407   if (verbose1) {
00408     cout << "The values on processor " << rank << " are: ";
00409     for (i=0; i<count; i++) 
00410       cout << iInputs[i] << " ";
00411     cout << endl;
00412   }
00413   EPETRA_TEST_ERR(petracomm.MinAll(iInputs,iGlobalMins,count),ierr);
00414   petracomm.Barrier();
00415 
00416   if (verbose1) {
00417     cout << "The min values according to processor " << rank << " are: ";
00418     for (i=0; i<count; i++) 
00419       cout << iGlobalMins[i] << " ";
00420     cout << endl;
00421   }
00422   forierr = 0;
00423   for (i=0; i<count; i++) 
00424     forierr += !(iMyGlobalMins[i] == iGlobalMins[i]); // otherwise calculated min is wrong
00425   EPETRA_TEST_ERR(forierr,ierr);
00426   delete [] iGlobalMins;
00427   delete [] iMyGlobalMins;
00428   petracomm.Barrier();
00429   if (verbose) cout << endl << "MinAll (type int) test passed!" << endl << endl;// If test gets to here the test passed, 
00430                                                                                 //only output on one node
00431   petracomm.Barrier();
00432 
00433   double* dGlobalMins = new double[count];
00434   if (verbose1) {
00435     cout << "The values on processor " << rank << " are: ";
00436     for (i=0; i<count; i++) 
00437       cout << dInputs[i] << " ";
00438     cout << endl;
00439   }
00440   EPETRA_TEST_ERR(petracomm.MinAll(dInputs,dGlobalMins,count),ierr);
00441   petracomm.Barrier();
00442 
00443   if (verbose1) {
00444     cout << "The min values according to processor " << rank << " are: ";
00445     for (i=0; i<count; i++) 
00446       cout << dGlobalMins[i] << " ";
00447     cout << endl;
00448   }
00449   forierr = 0;
00450   for (i=0; i<count; i++)
00451     forierr += !(Epetra_Util::Chop(dMyGlobalMins[i] - dGlobalMins[i]) == 0); // otherwise calculated min is wrong
00452   EPETRA_TEST_ERR(forierr,ierr);
00453   delete [] dGlobalMins;
00454   delete [] dMyGlobalMins;
00455   petracomm.Barrier();
00456   if (verbose) cout << endl << "MinAll (type double) test passed!" << endl << endl;// If test gets to here the test passed, 
00457                                                                                    //only output on one node
00458   petracomm.Barrier();
00459 
00460 
00461  // Test the SumAll functions
00462   int* iGlobalSums = new int[count];
00463   if (verbose1) {
00464     cout << "The values on processor " << rank << " are: ";
00465     for (i=0; i<count; i++) 
00466       cout << iInputs[i] << " ";
00467     cout << endl;
00468   }
00469   EPETRA_TEST_ERR(petracomm.SumAll(iInputs,iGlobalSums,count),ierr);
00470   petracomm.Barrier();
00471 
00472   if (verbose1) {
00473     cout << "The sums of all values according to processor " << rank << " are: ";
00474     for (i=0; i<count; i++) 
00475       cout << iGlobalSums[i] << " ";
00476     cout << endl;
00477   }
00478   forierr = 0;
00479   for (i=0; i<count; i++)
00480     forierr += !(iMyGlobalSums[i] == iGlobalSums[i]); // otherwise calculated sum is wrong
00481   EPETRA_TEST_ERR(forierr,ierr);
00482   delete [] iGlobalSums;
00483   delete [] iMyGlobalSums;
00484   petracomm.Barrier();
00485   if (verbose) cout << endl << "SumAll (type int) test passed!" << endl << endl;// If test gets to here the test passed, 
00486                                                                                 //only output on one node
00487   petracomm.Barrier();
00488 
00489   double* dGlobalSums = new double[count];
00490   if (verbose1) {
00491     cout << "The values on processor " << rank << " are: ";
00492     for (i=0; i<count; i++) 
00493       cout << dInputs[i] << " ";
00494     cout << endl;
00495   }
00496   EPETRA_TEST_ERR(petracomm.SumAll(dInputs,dGlobalSums,count),ierr);
00497   petracomm.Barrier();
00498 
00499   if (verbose1) {
00500     cout << "The sums of all values according to processor " << rank << " are: ";
00501     for (i=0; i<count; i++) 
00502       cout << dGlobalSums[i] << " ";
00503     cout << endl;
00504   }
00505   forierr = 0;
00506   for (i=0; i<count; i++)
00507     forierr += !(Epetra_Util::Chop(dMyGlobalSums[i] - dGlobalSums[i]) == 0); // otherwise calculated sum is wrong
00508   EPETRA_TEST_ERR(forierr,ierr);
00509 
00510   delete [] dGlobalSums;
00511   delete [] dMyGlobalSums;
00512   petracomm.Barrier();
00513   if (verbose) cout << endl << "SumAll (type double) test passed!" << endl << endl;// If test gets to here the test passed, 
00514                                                                                    //only output on one node
00515   petracomm.Barrier();
00516 
00517 
00518  // Test the ScanSum functions
00519   int* iScanSums = new int[count];
00520   if (verbose1) {
00521     cout << "The values on processor " << rank << " are: ";
00522     for (i=0; i<count; i++) 
00523       cout << iInputs[i] << " ";
00524     cout << endl;
00525   }
00526   
00527   EPETRA_TEST_ERR(petracomm.ScanSum(iInputs,iScanSums,count),ierr);
00528   petracomm.Barrier();
00529 
00530   if (verbose1) {
00531     cout << "The sums of all values on processors 0 - " << rank << " are: ";
00532     for (i=0; i<count; i++) {
00533       cout << iScanSums[i] << " ";
00534     }
00535     cout << endl;
00536   }
00537   forierr = 0;
00538   for (i=0; i<count; i++)
00539     forierr += !(iMyScanSums[i] == iScanSums[i]);
00540   EPETRA_TEST_ERR(forierr,ierr);
00541   delete [] iScanSums;
00542   delete [] iMyScanSums;
00543   petracomm.Barrier();
00544   if (verbose) cout << endl << "ScanSum (type int) test passed!" << endl << endl;// If test gets to here the test passed, 
00545                                                                                  //only output on one node
00546   petracomm.Barrier();
00547 
00548   double* dScanSums = new double[count];
00549   if (verbose1) {
00550     cout << "The values on processor " << rank << " are: ";
00551     for (i=0; i<count; i++)
00552       cout << dInputs[i] << " ";
00553     cout << endl;
00554   }
00555 
00556   EPETRA_TEST_ERR(petracomm.ScanSum(dInputs,dScanSums,count),ierr);
00557   petracomm.Barrier();
00558 
00559   if (verbose1) {
00560     cout << "The sums of all values on processors 0 - " << rank << " are: ";
00561     for (i=0; i<count; i++) {
00562       cout << dScanSums[i] << " ";
00563     }
00564     cout << endl;
00565   }
00566   forierr = 0;
00567   for (i=0; i<count; i++)
00568     forierr += !(Epetra_Util::Chop(dMyScanSums[i] - dScanSums[i])== 0);
00569   EPETRA_TEST_ERR(forierr,ierr);
00570   delete [] dScanSums;
00571   delete [] dMyScanSums;
00572   petracomm.Barrier();
00573   if (verbose) cout << endl << "ScanSum (type double) test passed!" << endl << endl;// If test gets to here the test passed, 
00574                                                                                     //only output on one node
00575   petracomm.Barrier();
00576 
00577 
00578  // Test the Gather functions
00579   int* iOrderedVals = new int[totalVals];
00580   if (verbose1) {
00581     cout << "The values on processor " << rank << " are: ";
00582     for (i=0; i<count; i++) 
00583       cout << iInputs[i] << " ";
00584     cout << endl;
00585   }
00586   
00587   EPETRA_TEST_ERR(petracomm.GatherAll(iInputs,iOrderedVals,count),ierr);
00588   petracomm.Barrier();
00589 
00590   if (verbose1) {
00591     cout << "The combined list of all values from all processors according to processor " << rank << " is: ";
00592     for (i=0; i<totalVals; i++) {
00593       cout << iOrderedVals[i] << " ";
00594     }
00595     cout << endl;
00596   }
00597   forierr = 0;
00598   for (i=0; i<totalVals; i++)
00599     forierr += !(iMyOrderedVals[i] == iOrderedVals[i]);
00600   EPETRA_TEST_ERR(forierr,ierr);
00601   delete [] iOrderedVals;
00602   delete [] iMyOrderedVals;
00603   petracomm.Barrier();
00604   if (verbose) cout << endl << "GatherAll (type int) test passed!" << endl << endl;// If test gets to here the test passed, 
00605                                                                                    //only output on one node
00606   petracomm.Barrier();
00607 
00608   double* dOrderedVals = new double[totalVals];
00609   if (verbose1) {
00610     cout << "The values on processor " << rank << " are: ";
00611     for (i=0; i<count; i++) 
00612       cout << dInputs[i] << " ";
00613     cout << endl;
00614   }
00615   
00616   EPETRA_TEST_ERR(petracomm.GatherAll(dInputs,dOrderedVals,count),ierr);
00617   petracomm.Barrier();
00618 
00619   if (verbose1) {
00620     cout << "The combined list of all values from all processors according to processor " << rank << " is: ";
00621     for (i=0; i<totalVals; i++) {
00622       cout << dOrderedVals[i] << " ";
00623     }
00624     cout << endl;
00625   }
00626   forierr = 0;
00627   for (i=0; i<totalVals; i++)
00628     forierr += !(Epetra_Util::Chop(dMyOrderedVals[i] - dOrderedVals[i]) == 0);
00629   EPETRA_TEST_ERR(forierr,ierr);
00630   delete [] dOrderedVals;
00631   delete [] dMyOrderedVals;
00632   petracomm.Barrier();
00633   if (verbose) cout << endl << "GatherAll (type double) test passed!" << endl << endl;// If test gets to here the test passed, 
00634                                                                                       //only output on one node
00635   petracomm.Barrier();
00636 
00637   delete[] dInputs;
00638   delete[] iInputs;
00639 
00640   return(ierr);
00641 }
00642 
00643 //=============================================================================
00644 int checkSerialDataClass(bool verbose) {
00645   int ierr = 0;
00646   if(verbose) cout << "Testing Reference Counting... ";               
00647   Epetra_SerialComm c1;
00648   int c1count = c1.ReferenceCount();
00649   const Epetra_SerialCommData* c1addr = c1.DataPtr();
00650   EPETRA_TEST_ERR(!(c1count==1),ierr); // count should be 1
00651   if(verbose) cout << "Default constructor. \nc1= " << c1count << "  " << c1addr << endl;
00652 
00653   Epetra_SerialComm* c2 = new Epetra_SerialComm(c1);
00654   int c2count = c2->ReferenceCount();
00655   const Epetra_SerialCommData* c2addr = c2->DataPtr();
00656   int c1countold = c1count;
00657   c1count = c1.ReferenceCount();
00658   EPETRA_TEST_ERR(!(c2count==c1count && c1count==(c1countold+1)),ierr); // both counts should be 2
00659   EPETRA_TEST_ERR(!(c1addr==c2addr),ierr); // addresses should be same
00660   if(verbose) cout << "Copy constructor(heap). \nc1= " << c1count << "  " << c1addr 
00661                     << "\nc2= " << c2count << "  " << c2addr << endl;
00662   delete c2;
00663   c1countold = c1count;
00664   c1count = c1.ReferenceCount();
00665   EPETRA_TEST_ERR(!(c1count==c1countold-1),ierr); // count should have decremented (to 1)
00666   EPETRA_TEST_ERR(!(c1addr==c1.DataPtr()),ierr); // c1addr should be unchanged
00667   if(verbose) cout << "c2 Deleted. \nc1= " << c1count << "  " << c1addr << endl;
00668   { // inside own set of brackets so that c2a will be automatically at end of brackets
00669     // so that we can test to make sure objects on the stack deallocate correctly
00670     Epetra_SerialComm c2a(c1);
00671     c2count = c2a.ReferenceCount();
00672     c2addr = c2a.DataPtr();
00673     c1countold = c1count;
00674     c1count = c1.ReferenceCount();
00675     EPETRA_TEST_ERR(!(c2count==c1count && c1count==c1countold+1),ierr); // both counts should be 2
00676     EPETRA_TEST_ERR(!(c1addr==c2addr),ierr); // addresses should be same
00677     if(verbose) cout << "Copy constructor(stack). \nc1= " << c1count << "  " << c1addr 
00678                       << "\nc2a= " << c2count << "  " << c2addr << endl;
00679   }
00680   c1countold = c1count;
00681   c1count = c1.ReferenceCount();
00682   EPETRA_TEST_ERR(!(c1count==c1countold-1),ierr); // count should have decremented (to 1)
00683   EPETRA_TEST_ERR(!(c1addr==c1.DataPtr()),ierr); // c1addr should be unchanged
00684   if(verbose) cout << "c2a Destroyed. \nc1= " << c1count << "  " << c1addr << endl;
00685   if(verbose) cout << "Assignment operator, post construction" << endl;
00686   Epetra_SerialComm c3;
00687   int c3count = c3.ReferenceCount();
00688   const Epetra_SerialCommData* c3addr = c3.DataPtr();
00689   EPETRA_TEST_ERR(!(c3count==1),ierr); // c3count should be 1 initially
00690   EPETRA_TEST_ERR(!(c1addr!=c3addr),ierr); // c1 and c3 should have different ptr addresses
00691   if(verbose)cout << "Prior to assignment: \nc1=" << c1count << "  " << c1addr 
00692                    << "\nc3=" << c3count << "  " << c3addr << endl;
00693   c3 = c1;
00694   c3count = c3.ReferenceCount();
00695   c3addr = c3.DataPtr();
00696   c1countold = c1count;
00697   c1count = c1.ReferenceCount();
00698   EPETRA_TEST_ERR(!(c3count==c1count && c1count==c1countold+1),ierr); // both counts should be 2
00699   EPETRA_TEST_ERR(!(c1addr==c3addr),ierr); // addresses should be same
00700   if(verbose)cout << "After assignment: \nc1=" << c1count << "  " << c1addr 
00701                    << "\nc3=" << c3count << "  " << c3addr << endl;
00702   return(ierr);
00703 }
00704 
00705 //=============================================================================
00706 #ifdef EPETRA_MPI
00707 int checkMpiDataClass(bool verbose) {
00708   int ierr = 0;
00709   if(verbose) cout << "Testing Reference Counting... ";               
00710   Epetra_MpiComm c1( MPI_COMM_WORLD );
00711   int c1count = c1.ReferenceCount();
00712   const Epetra_MpiCommData* c1addr = c1.DataPtr();
00713   EPETRA_TEST_ERR(!(c1count==1),ierr); // count should be 1
00714   if(verbose) cout << "Default constructor. \nc1= " << c1count << "  " << c1addr << endl;
00715 
00716   Epetra_MpiComm* c2 = new Epetra_MpiComm(c1);
00717   int c2count = c2->ReferenceCount();
00718   const Epetra_MpiCommData* c2addr = c2->DataPtr();
00719   int c1countold = c1count;
00720   c1count = c1.ReferenceCount();
00721   EPETRA_TEST_ERR(!(c2count==c1count && c1count==(c1countold+1)),ierr); // both counts should be 2
00722   EPETRA_TEST_ERR(!(c1addr==c2addr),ierr); // addresses should be same
00723   if(verbose) cout << "Copy constructor(heap). \nc1= " << c1count << "  " << c1addr 
00724                     << "\nc2= " << c2count << "  " << c2addr << endl;
00725   delete c2;
00726   c1countold = c1count;
00727   c1count = c1.ReferenceCount();
00728   EPETRA_TEST_ERR(!(c1count==c1countold-1),ierr); // count should have decremented (to 1)
00729   EPETRA_TEST_ERR(!(c1addr==c1.DataPtr()),ierr); // c1addr should be unchanged
00730   if(verbose) cout << "c2 Deleted. \nc1= " << c1count << "  " << c1addr << endl;
00731   { // inside own set of brackets so that c2a will be automatically at end of brackets
00732     // so that we can test to make sure objects on the stack deallocate correctly
00733     Epetra_MpiComm c2a(c1);
00734     c2count = c2a.ReferenceCount();
00735     c2addr = c2a.DataPtr();
00736     c1countold = c1count;
00737     c1count = c1.ReferenceCount();
00738     EPETRA_TEST_ERR(!(c2count==c1count && c1count==c1countold+1),ierr); // both counts should be 2
00739     EPETRA_TEST_ERR(!(c1addr==c2addr),ierr); // addresses should be same
00740     if(verbose) cout << "Copy constructor(stack). \nc1= " << c1count << "  " << c1addr 
00741                       << "\nc2a= " << c2count << "  " << c2addr << endl;
00742   }
00743   c1countold = c1count;
00744   c1count = c1.ReferenceCount();
00745   EPETRA_TEST_ERR(!(c1count==c1countold-1),ierr); // count should have decremented (to 1)
00746   EPETRA_TEST_ERR(!(c1addr==c1.DataPtr()),ierr); // c1addr should be unchanged
00747   if(verbose) cout << "c2a Destroyed. \nc1= " << c1count << "  " << c1addr << endl;
00748   if(verbose) cout << "Assignment operator, post construction" << endl;
00749   Epetra_MpiComm c3( MPI_COMM_WORLD );
00750   int c3count = c3.ReferenceCount();
00751   const Epetra_MpiCommData* c3addr = c3.DataPtr();
00752   EPETRA_TEST_ERR(!(c3count==1),ierr); // c3count should be 1 initially
00753   EPETRA_TEST_ERR(!(c1addr!=c3addr),ierr); // c1 and c3 should have different ptr addresses
00754   if(verbose)cout << "Prior to assignment: \nc1=" << c1count << "  " << c1addr 
00755                    << "\nc3=" << c3count << "  " << c3addr << endl;
00756   c3 = c1;
00757   c3count = c3.ReferenceCount();
00758   c3addr = c3.DataPtr();
00759   c1countold = c1count;
00760   c1count = c1.ReferenceCount();
00761   EPETRA_TEST_ERR(!(c3count==c1count && c1count==c1countold+1),ierr); // both counts should be 2
00762   EPETRA_TEST_ERR(!(c1addr==c3addr),ierr); // addresses should be same
00763   if(verbose)cout << "After assignment: \nc1=" << c1count << "  " << c1addr 
00764                    << "\nc3=" << c3count << "  " << c3addr << endl;
00765   return(ierr);
00766 }
00767 #endif
00768 
00769 int checkDistributor(Epetra_Distributor* distr,
00770                      Epetra_Comm& Comm)
00771 {
00772   int numprocs = Comm.NumProc();
00773 
00774   int numExportIDs = numprocs;
00775   int* exportPIDs = new int[numExportIDs];
00776   for(int p=0; p<numExportIDs; ++p) {
00777     exportPIDs[p] = p;
00778   }
00779 
00780   bool deterministic = true;
00781   int numRemoteIDs = 0;
00782 
00783   int err = distr->CreateFromSends(numExportIDs, exportPIDs,
00784                                    deterministic, numRemoteIDs);
00785 
00786   //numRemoteIDs should equal numExportIDs.
00787 
00788   int returnValue = numRemoteIDs == numExportIDs ? 0 : -99;
00789 
00790   delete [] exportPIDs;
00791 
00792   if (returnValue + err != 0) {
00793     return(returnValue + err);
00794   }
00795 
00796   int* exportIDs = new int[numExportIDs];
00797   for(int i=0; i<numExportIDs; ++i) {
00798     exportIDs[i] = i+1;
00799   }
00800 
00801   int len_imports = 0;
00802   char* imports = NULL;
00803 
00804   err = distr->Do((char*)exportIDs, sizeof(int),
00805                   len_imports, imports);
00806 
00807   delete [] exportIDs;
00808 
00809   if (len_imports > 0) {
00810     delete [] imports;
00811   }
00812 
00813   return(err);
00814 }
00815 
00816 /*
00817   end of file cxx_main.cpp
00818 */

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