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