test/CrsGraph/cxx_main.cpp

Go to the documentation of this file.
00001 //@HEADER
00002 // ************************************************************************
00003 // 
00004 //               Epetra: Linear Algebra Services Package 
00005 //                 Copyright (2001) Sandia Corporation
00006 // 
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00025 // 
00026 // ************************************************************************
00027 //@HEADER
00028 
00029 #include "Epetra_CrsGraph.h"
00030 #include "Epetra_Map.h"
00031 #ifdef EPETRA_MPI
00032 #include "Epetra_MpiComm.h"
00033 #include <mpi.h>
00034 #else
00035 #include "Epetra_SerialComm.h"
00036 #endif
00037 #include "../epetra_test_err.h"
00038 #include "Epetra_Version.h"
00039 
00040 // Prototype
00041 int check(Epetra_CrsGraph& A, int NumMyRows1, int NumGlobalRows1, int NumMyNonzeros1,
00042           int NumGlobalNonzeros1, int* MyGlobalElements, bool verbose);
00043 
00044 int checkSharedOwnership(Epetra_Comm& Comm, bool verbose);
00045 int checkCopyAndAssignment(Epetra_Comm& Comm, bool verbose);
00046 
00047 int main(int argc, char *argv[]) {
00048   int ierr = 0;
00049   int i;
00050   int forierr = 0;
00051   int* Indices;
00052   bool debug = true;
00053 
00054 #ifdef EPETRA_MPI
00055 
00056   // Initialize MPI
00057 
00058   MPI_Init(&argc,&argv);
00059   int rank; // My process ID
00060 
00061   MPI_Comm_rank(MPI_COMM_WORLD, &rank);
00062   Epetra_MpiComm Comm( MPI_COMM_WORLD );
00063 
00064 #else
00065 
00066   int rank = 0;
00067   Epetra_SerialComm Comm;
00068 
00069 #endif
00070 
00071   bool verbose = false;
00072 
00073   // Check if we should print results to standard out
00074   if(argc > 1) {
00075     if(argv[1][0]=='-' && argv[1][1]=='v') {
00076       verbose = true;
00077     }
00078   }
00079 
00080   //char tmp;
00081   //if (rank==0) cout << "Press any key to continue..."<< endl;
00082   //if (rank==0) cin >> tmp;
00083   //Comm.Barrier();
00084 
00085   Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
00086   int MyPID = Comm.MyPID();
00087   int NumProc = Comm.NumProc();
00088 
00089   if(verbose && MyPID==0)
00090     cout << Epetra_Version() << endl << endl;
00091 
00092   if(verbose) cout << "Processor "<<MyPID<<" of "<< NumProc << " is alive." << endl;
00093 
00094   bool verbose1 = verbose;
00095 
00096   // Redefine verbose to only print on PE 0
00097   if(verbose && rank != 0) 
00098     verbose = false;
00099 
00100   int NumMyEquations = 5;
00101   int NumGlobalEquations = NumMyEquations*NumProc+EPETRA_MIN(NumProc,3);
00102   if(MyPID < 3) 
00103     NumMyEquations++;
00104 
00105   // Construct a Map that puts approximately the same Number of equations on each processor
00106 
00107   Epetra_Map& Map = *new Epetra_Map(NumGlobalEquations, NumMyEquations, 0, Comm);
00108   
00109   // Get update list and number of local equations from newly created Map
00110   int* MyGlobalElements = new int[Map.NumMyElements()];
00111   Map.MyGlobalElements(MyGlobalElements);
00112 
00113   // Create an integer vector NumNz that is used to build the Petra Matrix.
00114   // NumNz[i] is the Number of OFF-DIAGONAL term for the ith global equation on this processor
00115 
00116   int* NumNz = new int[NumMyEquations];
00117 
00118   // We are building a tridiagonal matrix where each row has (-1 2 -1)
00119   // So we need 2 off-diagonal terms (except for the first and last equation)
00120 
00121   for(i=0; i<NumMyEquations; i++)
00122     if(MyGlobalElements[i]==0 || MyGlobalElements[i] == NumGlobalEquations-1)
00123       NumNz[i] = 1;
00124     else
00125       NumNz[i] = 2;
00126 
00127   // Create a Epetra_CrsGraph
00128 
00129   Epetra_CrsGraph& A = *new Epetra_CrsGraph(Copy, Map, NumNz);
00130   EPETRA_TEST_ERR(A.IndicesAreGlobal(),ierr);
00131   EPETRA_TEST_ERR(A.IndicesAreLocal(),ierr);
00132   
00133   // Add  rows one-at-a-time
00134   // Need some vectors to help
00135   // Off diagonal Values will always be -1
00136 
00137   Indices = new int[2];
00138   int NumEntries;
00139 
00140   forierr = 0;
00141   for(i = 0; i < NumMyEquations; i++) {
00142     if(MyGlobalElements[i] == 0) {
00143       Indices[0] = 1;
00144       NumEntries = 1;
00145     }
00146     else if(MyGlobalElements[i] == NumGlobalEquations-1) {
00147       Indices[0] = NumGlobalEquations-2;
00148       NumEntries = 1;
00149     }
00150     else {
00151       Indices[0] = MyGlobalElements[i]-1;
00152       Indices[1] = MyGlobalElements[i]+1;
00153       NumEntries = 2;
00154     }
00155     forierr += !(A.InsertGlobalIndices(MyGlobalElements[i], NumEntries, Indices)==0);
00156     forierr += !(A.InsertGlobalIndices(MyGlobalElements[i], 1, MyGlobalElements+i)>0); // Put in the diagonal entry (should cause realloc)
00157   }
00158   EPETRA_TEST_ERR(forierr,ierr);
00159   //A.PrintGraphData(cout);
00160   delete[] Indices;
00161   
00162   // Finish up
00163   EPETRA_TEST_ERR(!(A.IndicesAreGlobal()),ierr);
00164   EPETRA_TEST_ERR(!(A.FillComplete()==0),ierr);
00165   EPETRA_TEST_ERR(!(A.IndicesAreLocal()),ierr);
00166   EPETRA_TEST_ERR(A.StorageOptimized(),ierr);
00167 
00168   A.OptimizeStorage();
00169 
00170   EPETRA_TEST_ERR(!(A.StorageOptimized()),ierr);
00171   EPETRA_TEST_ERR(A.UpperTriangular(),ierr);
00172   EPETRA_TEST_ERR(A.LowerTriangular(),ierr);
00173 
00174   if(verbose) cout << "\n*****Testing variable entry constructor\n" << endl;
00175 
00176   int NumMyNonzeros = 3 * NumMyEquations;
00177   if(A.LRID(0) >= 0) 
00178     NumMyNonzeros--; // If I own first global row, then there is one less nonzero
00179   if(A.LRID(NumGlobalEquations-1) >= 0) 
00180     NumMyNonzeros--; // If I own last global row, then there is one less nonzero
00181 
00182   EPETRA_TEST_ERR(check(A, NumMyEquations, NumGlobalEquations, NumMyNonzeros, 3*NumGlobalEquations-2, 
00183                         MyGlobalElements, verbose),ierr);
00184   forierr = 0;
00185   for(i = 0; i < NumMyEquations; i++) 
00186     forierr += !(A.NumGlobalIndices(MyGlobalElements[i])==NumNz[i]+1);
00187   EPETRA_TEST_ERR(forierr,ierr);
00188   for(i = 0; i < NumMyEquations; i++) 
00189     forierr += !(A.NumMyIndices(i)==NumNz[i]+1);
00190   EPETRA_TEST_ERR(forierr,ierr);
00191 
00192   if(verbose) cout << "NumIndices function check OK" << endl;
00193 
00194   delete &A;
00195 
00196   if(debug) Comm.Barrier();
00197 
00198   if(verbose) cout << "\n*****Testing constant entry constructor\n" << endl;
00199 
00200   Epetra_CrsGraph& AA = *new Epetra_CrsGraph(Copy, Map, 5);
00201   
00202   if(debug) Comm.Barrier();
00203 
00204   for(i = 0; i < NumMyEquations; i++) 
00205     AA.InsertGlobalIndices(MyGlobalElements[i], 1, MyGlobalElements+i);
00206 
00207   // Note:  All processors will call the following Insert routines, but only the processor
00208   //        that owns it will actually do anything
00209 
00210   int One = 1;
00211   if(AA.MyGlobalRow(0)) {
00212     EPETRA_TEST_ERR(!(AA.InsertGlobalIndices(0, 0, &One)==0),ierr);
00213   }
00214   else 
00215     EPETRA_TEST_ERR(!(AA.InsertGlobalIndices(0, 1, &One)==-2),ierr);
00216   EPETRA_TEST_ERR(!(AA.FillComplete()==0),ierr);
00217   EPETRA_TEST_ERR(AA.StorageOptimized(),ierr);
00218   EPETRA_TEST_ERR(!(AA.UpperTriangular()),ierr);
00219   EPETRA_TEST_ERR(!(AA.LowerTriangular()),ierr);
00220 
00221   if(debug) Comm.Barrier();
00222   EPETRA_TEST_ERR(check(AA, NumMyEquations, NumGlobalEquations, NumMyEquations, NumGlobalEquations, 
00223                         MyGlobalElements, verbose),ierr);
00224 
00225   if(debug) Comm.Barrier();
00226 
00227   forierr = 0;
00228   for(i = 0; i < NumMyEquations; i++) 
00229     forierr += !(AA.NumGlobalIndices(MyGlobalElements[i])==1);
00230   EPETRA_TEST_ERR(forierr,ierr);
00231 
00232   if(verbose) cout << "NumIndices function check OK" << endl;
00233 
00234   if(debug) Comm.Barrier();
00235 
00236   if(verbose) cout << "\n*****Testing copy constructor\n" << endl;
00237 
00238   Epetra_CrsGraph& B = *new Epetra_CrsGraph(AA);
00239   delete &AA;
00240 
00241   EPETRA_TEST_ERR(check(B, NumMyEquations, NumGlobalEquations, NumMyEquations, NumGlobalEquations, 
00242                         MyGlobalElements, verbose),ierr);
00243 
00244   forierr = 0;
00245   for(i = 0; i < NumMyEquations; i++) 
00246     forierr += !(B.NumGlobalIndices(MyGlobalElements[i])==1);
00247   EPETRA_TEST_ERR(forierr,ierr);
00248 
00249   if(verbose) cout << "NumIndices function check OK" << endl;
00250 
00251   if(debug) Comm.Barrier();
00252 
00253   if(verbose) cout << "\n*****Testing post construction modifications\n" << endl;
00254 
00255   EPETRA_TEST_ERR(!(B.InsertGlobalIndices(0, 1, &One)==-2),ierr);
00256   delete &B;
00257 
00258   // Release all objects
00259   delete[] MyGlobalElements;
00260   delete[] NumNz;
00261   delete &Map;
00262       
00263 
00264   if (verbose1) {
00265     // Test ostream << operator (if verbose1)
00266     // Construct a Map that puts 2 equations on each PE
00267     
00268     int NumMyElements1 = 4;
00269     int NumMyEquations1 = NumMyElements1;
00270     int NumGlobalEquations1 = NumMyEquations1*NumProc;
00271 
00272     Epetra_Map& Map1 = *new Epetra_Map(NumGlobalEquations1, NumMyElements1, 1, Comm);
00273     
00274     // Get update list and number of local equations from newly created Map
00275     int* MyGlobalElements1 = new int[Map1.NumMyElements()];
00276     Map1.MyGlobalElements(MyGlobalElements1);
00277     
00278     // Create an integer vector NumNz that is used to build the Petra Matrix.
00279     // NumNz[i] is the Number of OFF-DIAGONAL term for the ith global equation on this processor
00280     
00281     int* NumNz1 = new int[NumMyEquations1];
00282     
00283     // We are building a tridiagonal matrix where each row has (-1 2 -1)
00284     // So we need 2 off-diagonal terms (except for the first and last equation)
00285     
00286     for(i = 0; i < NumMyEquations1; i++)
00287       if(MyGlobalElements1[i]==1 || MyGlobalElements1[i] == NumGlobalEquations1)
00288         NumNz1[i] = 1;
00289       else
00290         NumNz1[i] = 2;
00291     
00292     // Create a Epetra_Graph using 1-based arithmetic
00293     
00294     Epetra_CrsGraph& A1 = *new Epetra_CrsGraph(Copy, Map1, NumNz1);
00295     
00296     // Add  rows one-at-a-time
00297     // Need some vectors to help
00298     // Off diagonal Values will always be -1
00299     
00300     
00301     int* Indices1 = new int[2];
00302     int NumEntries1;
00303 
00304     forierr = 0;
00305     for(i = 0; i < NumMyEquations1; i++) {
00306       if(MyGlobalElements1[i]==1) {
00307         Indices1[0] = 2;
00308         NumEntries1 = 1;
00309       }
00310       else if(MyGlobalElements1[i] == NumGlobalEquations1) {
00311       Indices1[0] = NumGlobalEquations1-1;
00312       NumEntries1 = 1;
00313       }
00314       else {
00315         Indices1[0] = MyGlobalElements1[i]-1;
00316         Indices1[1] = MyGlobalElements1[i]+1;
00317         NumEntries1 = 2;
00318       }
00319       forierr += !(A1.InsertGlobalIndices(MyGlobalElements1[i], NumEntries1, Indices1)==0);
00320       forierr += !(A1.InsertGlobalIndices(MyGlobalElements1[i], 1, MyGlobalElements1+i)>0); // Put in the diagonal entry
00321     }
00322     EPETRA_TEST_ERR(forierr,ierr);
00323     
00324     // Finish up
00325     EPETRA_TEST_ERR(!(A1.FillComplete()==0),ierr);
00326     
00327     if(verbose) cout << "Print out tridiagonal matrix, each part on each processor. Index base is one.\n" << endl;
00328     cout << A1 << endl;
00329     
00330   // Release all objects
00331   delete[] NumNz1;
00332   delete[] Indices1;
00333   delete[] MyGlobalElements1;
00334 
00335   delete &A1;
00336   delete &Map1;
00337   }
00338 
00339   // Test copy constructor, op=, and reference-counting
00340   int tempierr = 0;
00341   if(verbose) cout << "\n*****Checking cpy ctr, op=, and reference counting." << endl;
00342   tempierr = checkCopyAndAssignment(Comm, verbose);
00343   EPETRA_TEST_ERR(tempierr, ierr);
00344   if(verbose && (tempierr == 0)) cout << "Checked OK." << endl;
00345 
00346   // Test shared-ownership code (not implemented yet)
00347   tempierr = 0;
00348   if(verbose) cout << "\n*****Checking shared-ownership tests." << endl;
00349   tempierr = checkSharedOwnership(Comm, verbose);
00350   EPETRA_TEST_ERR(tempierr, ierr);
00351   if(verbose && (tempierr == 0)) cout << "Checked OK." << endl;
00352 
00353       
00354 #ifdef EPETRA_MPI
00355   MPI_Finalize() ;
00356 #endif
00357 /* end main */
00358   return(ierr);
00359 }
00360 
00361 //==============================================================================
00362 int checkSharedOwnership(Epetra_Comm& Comm, bool verbose) {
00363   // check to make sure each function returns 1 when it should
00364   // check to make sure each function doesn't return 1 when it shouldn't
00365   int ierr = 0;
00366 
00367   // initialize Map
00368   const int NumMyElements = 10;
00369   const int IndexBase = 0;
00370   Epetra_Map Map1(-1, NumMyElements, IndexBase, Comm);
00371   // initialize Graphs
00372   const int NumIndicesPerRow = 5;
00373   Epetra_CrsGraph * SoleOwner = new Epetra_CrsGraph(Copy, Map1, Map1, NumIndicesPerRow);
00374   Epetra_CrsGraph SharedOrig(Copy, Map1, Map1, NumIndicesPerRow);
00375   Epetra_CrsGraph SharedOwner(SharedOrig);
00376   // arrays used by Insert & Remove
00377   Epetra_IntSerialDenseVector array1(2);
00378   array1[0] = NumIndicesPerRow / 2;
00379   array1[1] = array1[0] + 1;
00380   Epetra_IntSerialDenseVector array2(NumIndicesPerRow);
00381   for(int i = 0; i < NumIndicesPerRow; i++)
00382     array2[i] = i;
00383   // output variables (declaring them here lets us comment out indiv. tests)
00384   int soleOutput, sharedOutput;
00385 
00386   // InsertMyIndices
00387   if(verbose) cout << "InsertMyIndices..." << endl;
00388   soleOutput = SoleOwner->InsertMyIndices(0, 2, array1.Values());
00389   sharedOutput = SharedOwner.InsertMyIndices(0, 2, array1.Values());
00390   EPETRA_TEST_ERR(!(soleOutput == 0), ierr);
00391   EPETRA_TEST_ERR(!(sharedOutput == 1), ierr);
00392   if(verbose && ierr > 0) cout << "soleOutput = " << soleOutput << " sharedOutput = " << sharedOutput << endl;
00393 
00394   // SortIndices
00395   //if(verbose) cout << "SortIndices..." << endl;
00396   //soleOutput = SoleOwner.SortIndices();
00397   //sharedOutput = SharedOwner.SortIndices();
00398   //EPETRA_TEST_ERR(!(soleOutput == 0), ierr);
00399   //EPETRA_TEST_ERR(!(sharedOutput == 1), ierr);
00400   //if(verbose && ierr > 0) cout << "soleOutput = " << soleOutput << " sharedOutput = " << sharedOutput << endl;
00401 
00402   // RemoveRedundantIndices
00403   //if(verbose) cout << "RemoveRedundantIndices..." << endl;
00404   //SoleOwner.InsertGlobalIndices(0, 1, array1.Values());
00405   //SharedOwner.InsertGlobalIndices(0, 1, array1.Values());
00406   //soleOutput = SoleOwner.RemoveRedundantIndices();
00407   //sharedOutput = SharedOwner.RemoveRedundantIndices();
00408   //EPETRA_TEST_ERR(!(soleOutput == 0), ierr);
00409   //EPETRA_TEST_ERR(!(sharedOutput == 1), ierr);
00410   //if(verbose && ierr > 0) cout << "soleOutput = " << soleOutput << " sharedOutput = " << sharedOutput << endl;
00411 
00412   // FillComplete (#1)
00413   if(verbose) cout << "FillComplete..." << endl;
00414   soleOutput = SoleOwner->FillComplete();
00415   sharedOutput = SharedOwner.FillComplete();
00416   EPETRA_TEST_ERR(!(soleOutput == 0), ierr);
00417   EPETRA_TEST_ERR(!(sharedOutput == 1), ierr);
00418   if(verbose && ierr > 0) cout << "soleOutput = " << soleOutput << " sharedOutput = " << sharedOutput << endl;
00419 
00420   // OptimizeStorage
00421   if(verbose) cout << "OptimizeStorage..." << endl;
00422   soleOutput = SoleOwner->OptimizeStorage();
00423   sharedOutput = SharedOwner.OptimizeStorage();
00424   EPETRA_TEST_ERR(!(soleOutput == 0), ierr);
00425   EPETRA_TEST_ERR(!(sharedOutput == 0), ierr);
00426   if(verbose && ierr > 0) cout << "soleOutput = " << soleOutput << " sharedOutput = " << sharedOutput << endl;
00427 
00428   // RemoveMyIndices (#1)
00429   if(verbose) cout << "RemoveMyIndices..." << endl;
00430   soleOutput = SoleOwner->RemoveMyIndices(0, 1, &array1[1]);
00431   sharedOutput = SharedOwner.RemoveMyIndices(0, 1, &array1[1]);
00432   EPETRA_TEST_ERR(!(soleOutput == -1), ierr);
00433   EPETRA_TEST_ERR(!(sharedOutput == -1), ierr);
00434   if(verbose && ierr > 0) cout << "soleOutput = " << soleOutput << " sharedOutput = " << sharedOutput << endl;
00435 
00436   // RemoveMyIndices (#2)
00437   if(verbose) cout << "RemoveMyIndices(#2)..." << endl;
00438   soleOutput = SoleOwner->RemoveMyIndices(0);
00439   sharedOutput = SharedOwner.RemoveMyIndices(0);
00440   EPETRA_TEST_ERR(!(soleOutput == -1), ierr);
00441   EPETRA_TEST_ERR(!(sharedOutput == -1), ierr);
00442   if(verbose && ierr > 0) cout << "soleOutput = " << soleOutput << " sharedOutput = " << sharedOutput << endl;
00443 
00444   // FillComplete (#2)
00445   if(verbose) cout << "FillComplete(#2)..." << endl;
00446   soleOutput = SoleOwner->FillComplete(SoleOwner->DomainMap(), SoleOwner->RangeMap());
00447   sharedOutput = SharedOwner.FillComplete(SharedOwner.DomainMap(), SharedOwner.RangeMap());
00448   EPETRA_TEST_ERR(!(soleOutput == 0), ierr);
00449   EPETRA_TEST_ERR(!(sharedOutput == 1), ierr);
00450   if(verbose && ierr > 0) cout << "soleOutput = " << soleOutput << " sharedOutput = " << sharedOutput << endl;
00451 
00452   {
00453     // make new Graphs so that we can insert Global instead of Local
00454     // inside of new scope so that we can use same names
00455     Epetra_CrsGraph SoleOwner(Copy, Map1, NumIndicesPerRow);
00456     Epetra_CrsGraph SharedOrig(Copy, Map1, NumIndicesPerRow);
00457     Epetra_CrsGraph SharedOwner(SharedOrig);
00458     
00459     int GlobalRow = SoleOwner.GRID(0);
00460 
00461     // InsertGlobalIndices
00462     if(verbose) cout << "InsertGlobalIndices..." << endl;
00463     soleOutput = SoleOwner.InsertGlobalIndices(GlobalRow, 2, array2.Values());
00464     sharedOutput = SharedOwner.InsertGlobalIndices(GlobalRow, 2, array2.Values());
00465     EPETRA_TEST_ERR(!(soleOutput == 0), ierr);
00466     EPETRA_TEST_ERR(!(sharedOutput == 1), ierr);
00467     if(verbose && ierr > 0) cout << "soleOutput = " << soleOutput << " sharedOutput = " << sharedOutput << endl;
00468     
00469     // RemoveGlobalIndices (#1)
00470     if(verbose) cout << "RemoveGlobalIndices..." << endl;
00471     soleOutput = SoleOwner.RemoveGlobalIndices(GlobalRow, 1, &array2[1]);
00472     sharedOutput = SharedOwner.RemoveGlobalIndices(GlobalRow, 1, &array2[1]);
00473     EPETRA_TEST_ERR(!(soleOutput == 0), ierr);
00474     EPETRA_TEST_ERR(!(sharedOutput == 1), ierr);
00475     if(verbose && ierr > 0) cout << "soleOutput = " << soleOutput << " sharedOutput = " << sharedOutput << endl;
00476     
00477     // RemoveGlobalIndices (#2)
00478     if(verbose) cout << "RemoveGlobalIndices(#2)..." << endl;
00479     soleOutput = SoleOwner.RemoveGlobalIndices(GlobalRow);
00480     sharedOutput = SharedOwner.RemoveGlobalIndices(GlobalRow);
00481     EPETRA_TEST_ERR(!(soleOutput == 0), ierr);
00482     EPETRA_TEST_ERR(!(sharedOutput == 1), ierr);
00483     if(verbose && ierr > 0) cout << "soleOutput = " << soleOutput << " sharedOutput = " << sharedOutput << endl;
00484     }
00485 
00486 
00487   // *PROT* InsertIndices
00488   // *PROT* MakeIndicesLocal
00489   delete SoleOwner;
00490   return(ierr);
00491 
00492 }
00493 
00494 //==============================================================================
00495 int checkCopyAndAssignment(Epetra_Comm& Comm, bool verbose) {
00496   int ierr = 0;
00497   // check to make sure that copy ctr and op= are working correctly
00498   // (making level 1 deep copies)
00499 
00500   // create initial Map and Graph
00501   const int NumIndicesPerRow = 10;
00502   const int NumGlobalElements = 50;
00503   const int IndexBase = 0;
00504   Epetra_Map Map1(NumGlobalElements, IndexBase, Comm);
00505   Epetra_CrsGraph Graph1(Copy, Map1, NumIndicesPerRow);
00506   int g1count = Graph1.ReferenceCount();
00507   const Epetra_CrsGraphData* g1addr = Graph1.DataPtr();
00508   EPETRA_TEST_ERR(!(g1count == 1), ierr);
00509   if(verbose) cout << "Graph1 created (def ctr). data addr = " << g1addr << " ref. count = " << g1count << endl;
00510 
00511   //test copy constructor
00512   {
00513     Epetra_CrsGraph Graph2(Graph1);
00514     int g2count = Graph2.ReferenceCount();
00515     const Epetra_CrsGraphData* g2addr = Graph2.DataPtr();
00516     EPETRA_TEST_ERR(!(g2count == g1count+1), ierr);
00517     EPETRA_TEST_ERR(!(g2addr == g1addr), ierr);
00518     if(verbose) cout << "Graph2 created (cpy ctr). data addr = " << g2addr << " ref. count = " << g2count << endl;
00519   }
00520   int g1newcount = Graph1.ReferenceCount();
00521   const Epetra_CrsGraphData* g1newaddr = Graph1.DataPtr();
00522   EPETRA_TEST_ERR(!(g1newcount == g1count), ierr);
00523   EPETRA_TEST_ERR(!(g1newaddr = g1addr), ierr);
00524   if(verbose) cout << "Graph2 destroyed. Graph1 data addr = " << g1newaddr << " ref. count = " << g1newcount << endl;
00525 
00526   //test op=
00527   Epetra_CrsGraph Graph3(Copy, Map1, NumIndicesPerRow+1);
00528   int g3count = Graph3.ReferenceCount();
00529   const Epetra_CrsGraphData* g3addr = Graph3.DataPtr();
00530   EPETRA_TEST_ERR(!(g3addr != g1addr), ierr);
00531   if(verbose) cout << "Graph3 created (op= before). data addr = " << g3addr << " ref. count = " << g3count << endl;
00532   Graph3 = Graph1;
00533   g3count = Graph3.ReferenceCount();
00534   g3addr = Graph3.DataPtr();
00535   EPETRA_TEST_ERR(!(g3count == g1count+1), ierr);
00536   EPETRA_TEST_ERR(!(g3addr == g1addr), ierr);
00537   if(verbose) cout << "Graph3 set equal to Graph1 (op= after). data addr = " << g3addr << " ref. count = " << g3count << endl;
00538 
00539   return(ierr);
00540 }
00541 
00542 //==============================================================================
00543 int check(Epetra_CrsGraph& A, int NumMyRows1, int NumGlobalRows1, int NumMyNonzeros1,
00544     int NumGlobalNonzeros1, int* MyGlobalElements, bool verbose)
00545 {  
00546 
00547   int ierr = 0;
00548   int i;
00549   int j;
00550   int forierr = 0;
00551   int NumGlobalIndices;
00552   int NumMyIndices;
00553   int* MyViewIndices;
00554   int MaxNumIndices = A.MaxNumIndices();
00555   int* MyCopyIndices = new int[MaxNumIndices];
00556   int* GlobalCopyIndices = new int[MaxNumIndices];
00557 
00558   // Test query functions
00559 
00560   int NumMyRows = A.NumMyRows();
00561   if(verbose) cout << "Number of local Rows = " << NumMyRows << endl;
00562 
00563   EPETRA_TEST_ERR(!(NumMyRows==NumMyRows1),ierr);
00564 
00565   int NumMyNonzeros = A.NumMyNonzeros();
00566   if(verbose) cout << "Number of local Nonzero entries = " << NumMyNonzeros << endl;
00567 
00568   EPETRA_TEST_ERR(!(NumMyNonzeros==NumMyNonzeros1),ierr);
00569 
00570   int NumGlobalRows = A.NumGlobalRows();
00571   if(verbose) cout << "Number of global Rows = " << NumGlobalRows << endl;
00572 
00573   EPETRA_TEST_ERR(!(NumGlobalRows==NumGlobalRows1),ierr);
00574 
00575   int NumGlobalNonzeros = A.NumGlobalNonzeros();
00576   if(verbose) cout << "Number of global Nonzero entries = " << NumGlobalNonzeros << endl;
00577 
00578   EPETRA_TEST_ERR(!(NumGlobalNonzeros==NumGlobalNonzeros1),ierr);
00579 
00580   // GlobalRowView should be illegal (since we have local indices)
00581 
00582   EPETRA_TEST_ERR(!(A.ExtractGlobalRowView(A.RowMap().MaxMyGID(), NumGlobalIndices, GlobalCopyIndices)==-2),ierr);
00583 
00584   // Other binary tests
00585 
00586   EPETRA_TEST_ERR(A.NoDiagonal(),ierr);
00587   EPETRA_TEST_ERR(!(A.Filled()),ierr);
00588   EPETRA_TEST_ERR(!(A.MyGRID(A.RowMap().MaxMyGID())),ierr);
00589   EPETRA_TEST_ERR(!(A.MyGRID(A.RowMap().MinMyGID())),ierr);
00590   EPETRA_TEST_ERR(A.MyGRID(1+A.RowMap().MaxMyGID()),ierr);
00591   EPETRA_TEST_ERR(A.MyGRID(-1+A.RowMap().MinMyGID()),ierr);
00592   EPETRA_TEST_ERR(!(A.MyLRID(0)),ierr);
00593   EPETRA_TEST_ERR(!(A.MyLRID(NumMyRows-1)),ierr);
00594   EPETRA_TEST_ERR(A.MyLRID(-1),ierr);
00595   EPETRA_TEST_ERR(A.MyLRID(NumMyRows),ierr);
00596     
00597   forierr = 0;
00598   for(i = 0; i < NumMyRows; i++) {
00599     int Row = A.GRID(i);
00600     A.ExtractGlobalRowCopy(Row, MaxNumIndices, NumGlobalIndices, GlobalCopyIndices);
00601     A.ExtractMyRowView(i, NumMyIndices, MyViewIndices);
00602     forierr += !(NumGlobalIndices==NumMyIndices);
00603     for(j = 1; j < NumMyIndices; j++) EPETRA_TEST_ERR(!(MyViewIndices[j-1]<MyViewIndices[j]),ierr);
00604     for(j = 0; j < NumGlobalIndices; j++) {
00605       forierr += !(GlobalCopyIndices[j]==A.GCID(MyViewIndices[j]));
00606       forierr += !(A.LCID(GlobalCopyIndices[j])==MyViewIndices[j]);
00607     }
00608   }
00609   EPETRA_TEST_ERR(forierr,ierr);
00610   forierr = 0;
00611   for(i = 0; i < NumMyRows; i++) {
00612     int Row = A.GRID(i);
00613     A.ExtractGlobalRowCopy(Row, MaxNumIndices, NumGlobalIndices, GlobalCopyIndices);
00614     A.ExtractMyRowCopy(i, MaxNumIndices, NumMyIndices, MyCopyIndices);
00615     forierr += !(NumGlobalIndices==NumMyIndices);
00616     for(j = 1; j < NumMyIndices; j++) 
00617       EPETRA_TEST_ERR(!(MyCopyIndices[j-1]<MyCopyIndices[j]),ierr);
00618     for(j = 0; j < NumGlobalIndices; j++) {
00619       forierr += !(GlobalCopyIndices[j]==A.GCID(MyCopyIndices[j]));
00620       forierr += !(A.LCID(GlobalCopyIndices[j])==MyCopyIndices[j]);
00621     }
00622 
00623   }
00624   EPETRA_TEST_ERR(forierr,ierr);
00625 
00626   delete[] MyCopyIndices;
00627   delete[] GlobalCopyIndices;
00628 
00629   if(verbose) cout << "Rows sorted check OK" << endl;
00630 
00631   return(ierr);
00632 }

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