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