test/RowMatrix/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_Map.h"
00030 #include "Epetra_Time.h"
00031 #include "Epetra_CrsMatrix.h"
00032 #include "Epetra_JadMatrix.h"
00033 #include "Epetra_Vector.h"
00034 #include "Epetra_Flops.h"
00035 #ifdef EPETRA_MPI
00036 #include "Epetra_MpiComm.h"
00037 #include "mpi.h"
00038 #else
00039 #include "Epetra_SerialComm.h"
00040 #endif
00041 #include "../epetra_test_err.h"
00042 #include "Epetra_Version.h"
00043 #include <vector>
00044 #include <algorithm>
00045 #include <string>
00046 
00047 // prototypes
00048 
00049 int checkValues( double x, double y, string message = "", bool verbose = false) { 
00050   if (fabs((x-y)/x) > 0.01) {
00051     return(1); 
00052     if (verbose) cout << "********** " << message << " check failed.********** " << endl;
00053   }
00054   else {
00055     if (verbose) cout << message << " check OK." << endl;    
00056     return(0);
00057   }
00058 }
00059 
00060 int checkMultiVectors( Epetra_MultiVector & X, Epetra_MultiVector & Y, string message = "", bool verbose = false) {
00061   int numVectors = X.NumVectors();
00062   int length = Y.MyLength();
00063   int badvalue = 0;
00064   int globalbadvalue = 0;
00065   for (int j=0; j<numVectors; j++)
00066     for (int i=0; i< length; i++) 
00067       if (checkValues(X[j][i], Y[j][i])==1) badvalue = 1;
00068   X.Map().Comm().MaxAll(&badvalue, &globalbadvalue, 1);
00069 
00070   if (verbose) {
00071     if (globalbadvalue==0) cout << message << " check OK." << endl;
00072     else cout << "********* " << message << " check failed.********** " << endl;
00073   }
00074   return(globalbadvalue);
00075 }
00076 
00077 int check(Epetra_RowMatrix & A, Epetra_RowMatrix & B, bool verbose);
00078 
00079 int powerMethodTests(Epetra_RowMatrix & A, Epetra_RowMatrix & JadA, Epetra_Map & Map, 
00080          Epetra_Vector & q, Epetra_Vector & z, Epetra_Vector & resid, bool verbose);
00081 
00082 int power_method(bool TransA, Epetra_RowMatrix& A, 
00083      Epetra_Vector& q,
00084      Epetra_Vector& z0, 
00085      Epetra_Vector& resid, 
00086      double * lambda, int niters, double tolerance,
00087      bool verbose);
00088 
00089 int main(int argc, char *argv[])
00090 {
00091   int ierr = 0, i, forierr = 0;
00092 #ifdef EPETRA_MPI
00093 
00094   // Initialize MPI
00095 
00096   MPI_Init(&argc,&argv);
00097   int rank; // My process ID
00098 
00099   MPI_Comm_rank(MPI_COMM_WORLD, &rank);
00100   Epetra_MpiComm Comm( MPI_COMM_WORLD );
00101 
00102 #else
00103 
00104   int rank = 0;
00105   Epetra_SerialComm Comm;
00106 
00107 #endif
00108 
00109   bool verbose = false;
00110 
00111   // Check if we should print results to standard out
00112   if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
00113 
00114   int verbose_int = verbose ? 1 : 0;
00115   Comm.Broadcast(&verbose_int, 1, 0);
00116   verbose = verbose_int==1 ? true : false;
00117 
00118 
00119   //  char tmp;
00120   //  if (rank==0) cout << "Press any key to continue..."<< endl;
00121   //  if (rank==0) cin >> tmp;
00122   //  Comm.Barrier();
00123 
00124   Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
00125   int MyPID = Comm.MyPID();
00126   int NumProc = Comm.NumProc();
00127 
00128   if(verbose && MyPID==0)
00129     cout << Epetra_Version() << endl << endl;
00130 
00131   if (verbose) cout << "Processor "<<MyPID<<" of "<< NumProc
00132         << " is alive."<<endl;
00133 
00134   // Redefine verbose to only print on PE 0
00135   if(verbose && rank!=0) 
00136     verbose = false;
00137 
00138   int NumMyEquations = 10000;
00139   int NumGlobalEquations = (NumMyEquations * NumProc) + EPETRA_MIN(NumProc,3);
00140   if(MyPID < 3) 
00141     NumMyEquations++;
00142 
00143   // Construct a Map that puts approximately the same Number of equations on each processor
00144 
00145   Epetra_Map Map(NumGlobalEquations, NumMyEquations, 0, Comm);
00146   
00147   // Get update list and number of local equations from newly created Map
00148   vector<int> MyGlobalElements(Map.NumMyElements());
00149   Map.MyGlobalElements(&MyGlobalElements[0]);
00150 
00151   // Create an integer vector NumNz that is used to build the Petra Matrix.
00152   // NumNz[i] is the Number of OFF-DIAGONAL term for the ith global equation on this processor
00153 
00154   vector<int> NumNz(NumMyEquations);
00155 
00156   // We are building a tridiagonal matrix where each row has (-1 2 -1)
00157   // So we need 2 off-diagonal terms (except for the first and last equation)
00158 
00159   for(i = 0; i < NumMyEquations; i++)
00160     if((MyGlobalElements[i] == 0) || (MyGlobalElements[i] == NumGlobalEquations - 1))
00161       NumNz[i] = 1;
00162     else
00163       NumNz[i] = 2;
00164 
00165   // Create a Epetra_Matrix
00166 
00167   Epetra_CrsMatrix A(Copy, Map, &NumNz[0]);
00168   EPETRA_TEST_ERR(A.IndicesAreGlobal(),ierr);
00169   EPETRA_TEST_ERR(A.IndicesAreLocal(),ierr);
00170   
00171   // Add  rows one-at-a-time
00172   // Need some vectors to help
00173   // Off diagonal Values will always be -1
00174 
00175 
00176   vector<double> Values(2);
00177   Values[0] = -1.0; 
00178   Values[1] = -1.0;
00179   vector<int> Indices(2);
00180   double two = 2.0;
00181   int NumEntries;
00182 
00183   forierr = 0;
00184   for(i = 0; i < NumMyEquations; i++) {
00185     if(MyGlobalElements[i] == 0) {
00186       Indices[0] = 1;
00187       NumEntries = 1;
00188     }
00189     else if (MyGlobalElements[i] == NumGlobalEquations-1) {
00190       Indices[0] = NumGlobalEquations-2;
00191       NumEntries = 1;
00192     }
00193     else {
00194       Indices[0] = MyGlobalElements[i]-1;
00195       Indices[1] = MyGlobalElements[i]+1;
00196       NumEntries = 2;
00197     }
00198     forierr += !(A.InsertGlobalValues(MyGlobalElements[i], NumEntries, &Values[0], &Indices[0])==0);
00199     forierr += !(A.InsertGlobalValues(MyGlobalElements[i], 1, &two, &MyGlobalElements[i])>0); // Put in the diagonal entry
00200   }
00201   EPETRA_TEST_ERR(forierr,ierr);
00202 
00203   // Finish up
00204   A.FillComplete();
00205   A.OptimizeStorage();
00206 
00207   Epetra_JadMatrix JadA(A);
00208   Epetra_JadMatrix JadA1(A);
00209   Epetra_JadMatrix JadA2(A);
00210 
00211   // Create vectors for Power method
00212 
00213   Epetra_Vector q(Map);
00214   Epetra_Vector z(Map); z.Random();
00215   Epetra_Vector resid(Map);
00216 
00217   Epetra_Flops flopcounter;
00218   A.SetFlopCounter(flopcounter);
00219   q.SetFlopCounter(A);
00220   z.SetFlopCounter(A);
00221   resid.SetFlopCounter(A);
00222   JadA.SetFlopCounter(A);
00223   JadA1.SetFlopCounter(A);
00224   JadA2.SetFlopCounter(A);
00225   
00226 
00227   if (verbose) cout << "=======================================" << endl
00228         << "Testing Jad using CrsMatrix as input..." << endl
00229         << "=======================================" << endl;
00230 
00231   A.ResetFlops();
00232   powerMethodTests(A, JadA, Map, q, z, resid, verbose);
00233 
00234   // Increase diagonal dominance
00235 
00236   if (verbose) cout << "\n\nIncreasing the magnitude of first diagonal term and solving again\n\n"
00237         << endl;
00238 
00239   
00240   if (A.MyGlobalRow(0)) {
00241     int numvals = A.NumGlobalEntries(0);
00242     vector<double> Rowvals(numvals);
00243     vector<int> Rowinds(numvals);
00244     A.ExtractGlobalRowCopy(0, numvals, numvals, &Rowvals[0], &Rowinds[0]); // Get A[0,0]
00245 
00246     for (i=0; i<numvals; i++) if (Rowinds[i] == 0) Rowvals[i] *= 10.0;
00247     
00248     A.ReplaceGlobalValues(0, numvals, &Rowvals[0], &Rowinds[0]);
00249   }
00250   JadA.UpdateValues(A);
00251   A.ResetFlops();
00252   powerMethodTests(A, JadA, Map, q, z, resid, verbose);
00253 
00254   if (verbose) cout << "================================================================" << endl
00255               << "Testing Jad using Jad matrix as input matrix for construction..." << endl
00256               << "================================================================" << endl;
00257   JadA1.ResetFlops();
00258   powerMethodTests(JadA1, JadA2, Map, q, z, resid, verbose);
00259 
00260 #ifdef EPETRA_MPI
00261   MPI_Finalize() ;
00262 #endif
00263 
00264 return ierr ;
00265 }
00266 
00267 int powerMethodTests(Epetra_RowMatrix & A, Epetra_RowMatrix & JadA, Epetra_Map & Map, 
00268          Epetra_Vector & q, Epetra_Vector & z, Epetra_Vector & resid, bool verbose) {
00269 
00270   // variable needed for iteration
00271   double lambda = 0.0;
00272   // int niters = 10000;
00273   int niters = 300;
00274   double tolerance = 1.0e-2;
00275   int ierr = 0;
00276 
00278   
00279   // Iterate
00280 
00281   Epetra_Time timer(Map.Comm());
00282   
00283   double startTime = timer.ElapsedTime();
00284   EPETRA_TEST_ERR(power_method(false, A, q, z, resid, &lambda, niters, tolerance, verbose),ierr);
00285   double elapsed_time = timer.ElapsedTime() - startTime;
00286   double total_flops = q.Flops();
00287   double MFLOPs = total_flops/elapsed_time/1000000.0;
00288   double lambdaref = lambda;
00289   double flopsref = total_flops;
00290 
00291   if (verbose) 
00292     cout << "\n\nTotal MFLOPs for reference first solve = " << MFLOPs << endl
00293       <<     "Total FLOPS                            = " <<total_flops <<endl<<endl;
00294 
00295   lambda = 0.0;
00296   startTime = timer.ElapsedTime();
00297   EPETRA_TEST_ERR(power_method(false, JadA, q, z, resid, &lambda, niters, tolerance, verbose),ierr);
00298   elapsed_time = timer.ElapsedTime() - startTime;
00299   total_flops = q.Flops();
00300   MFLOPs = total_flops/elapsed_time/1000000.0;
00301 
00302   if (verbose) 
00303     cout << "\n\nTotal MFLOPs for candidate first solve = " << MFLOPs << endl
00304       <<     "Total FLOPS                            = " <<total_flops <<endl<<endl;
00305 
00306   EPETRA_TEST_ERR(checkValues(lambda,lambdaref," No-transpose Power Method result", verbose),ierr);
00307   EPETRA_TEST_ERR(checkValues(total_flops,flopsref," No-transpose Power Method flop count", verbose),ierr);
00308 
00310   
00311   // Solve transpose problem
00312 
00313   if (verbose) cout << "\n\nUsing transpose of matrix and solving again (should give same result).\n\n"
00314         << endl;
00315   // Iterate
00316   lambda = 0.0;
00317   startTime = timer.ElapsedTime();
00318   EPETRA_TEST_ERR(power_method(true, A, q, z, resid, &lambda, niters, tolerance, verbose),ierr);
00319   elapsed_time = timer.ElapsedTime() - startTime;
00320   total_flops = q.Flops();
00321   MFLOPs = total_flops/elapsed_time/1000000.0;
00322   lambdaref = lambda;
00323   flopsref = total_flops;
00324 
00325   if (verbose) 
00326    cout << "\n\nTotal MFLOPs for reference transpose solve = " << MFLOPs << endl
00327      <<     "Total FLOPS                                = " <<total_flops <<endl<<endl;
00328 
00329   lambda = 0.0;
00330   startTime = timer.ElapsedTime();
00331   EPETRA_TEST_ERR(power_method(true, JadA, q, z, resid, &lambda, niters, tolerance, verbose),ierr);
00332   elapsed_time = timer.ElapsedTime() - startTime;
00333   total_flops = q.Flops();
00334   MFLOPs = total_flops/elapsed_time/1000000.0;
00335 
00336   if (verbose) 
00337     cout << "\n\nTotal MFLOPs for candidate transpose solve = " << MFLOPs << endl
00338       <<     "Total FLOPS                                = " <<total_flops <<endl<<endl;
00339 
00340   EPETRA_TEST_ERR(checkValues(lambda,lambdaref,"Transpose Power Method result", verbose),ierr);
00341   EPETRA_TEST_ERR(checkValues(total_flops,flopsref,"Transpose Power Method flop count", verbose),ierr);
00342 
00343   EPETRA_TEST_ERR(check(A, JadA, verbose),ierr);
00344 
00345   return(0);
00346 }
00347 int power_method(bool TransA, Epetra_RowMatrix& A, Epetra_Vector& q, Epetra_Vector& z0, 
00348      Epetra_Vector& resid, double* lambda, int niters, double tolerance, bool verbose) 
00349 {  
00350   
00351   // Fill z with random Numbers
00352   Epetra_Vector z(z0);
00353   
00354   // variable needed for iteration
00355   double normz, residual;
00356 
00357   int ierr = 1;
00358   
00359   for(int iter = 0; iter < niters; iter++) {
00360     z.Norm2(&normz); // Compute 2-norm of z
00361     q.Scale(1.0/normz, z);
00362     A.Multiply(TransA, q, z); // Compute z = A*q // SEGFAULT HAPPENS HERE
00363     q.Dot(z, lambda); // Approximate maximum eigenvaluE
00364     if(iter%100==0 || iter+1==niters) {
00365       resid.Update(1.0, z, -(*lambda), q, 0.0); // Compute A*q - lambda*q
00366       resid.Norm2(&residual);
00367       if(verbose) cout << "Iter = " << iter << "  Lambda = " << *lambda 
00368                        << "  Residual of A*q - lambda*q = " << residual << endl;
00369     }
00370     if(residual < tolerance) {
00371       ierr = 0;
00372       break;
00373     }
00374   }
00375   return(ierr);
00376 }
00377 
00378 int check(Epetra_RowMatrix& A, Epetra_RowMatrix & B, bool verbose)  {  
00379 
00380   int ierr = 0;
00381   EPETRA_TEST_ERR(!A.Comm().NumProc()==B.Comm().NumProc(),ierr);
00382   EPETRA_TEST_ERR(!A.Comm().MyPID()==B.Comm().MyPID(),ierr);
00383   EPETRA_TEST_ERR(!A.Filled()==B.Filled(),ierr);
00384   EPETRA_TEST_ERR(!A.HasNormInf()==B.HasNormInf(),ierr);
00385   EPETRA_TEST_ERR(!A.LowerTriangular()==B.LowerTriangular(),ierr);
00386   EPETRA_TEST_ERR(!A.Map().SameAs(B.Map()),ierr);
00387   EPETRA_TEST_ERR(!A.MaxNumEntries()==B.MaxNumEntries(),ierr);
00388   EPETRA_TEST_ERR(!A.NumGlobalCols()==B.NumGlobalCols(),ierr);
00389   EPETRA_TEST_ERR(!A.NumGlobalDiagonals()==B.NumGlobalDiagonals(),ierr);
00390   EPETRA_TEST_ERR(!A.NumGlobalNonzeros()==B.NumGlobalNonzeros(),ierr);
00391   EPETRA_TEST_ERR(!A.NumGlobalRows()==B.NumGlobalRows(),ierr);
00392   EPETRA_TEST_ERR(!A.NumMyCols()==B.NumMyCols(),ierr);
00393   EPETRA_TEST_ERR(!A.NumMyDiagonals()==B.NumMyDiagonals(),ierr);
00394   EPETRA_TEST_ERR(!A.NumMyNonzeros()==B.NumMyNonzeros(),ierr);
00395   for (int i=0; i<A.NumMyRows(); i++) {
00396     int nA, nB;
00397     A.NumMyRowEntries(i,nA); B.NumMyRowEntries(i,nB);
00398     EPETRA_TEST_ERR(!nA==nB,ierr);
00399   }
00400   EPETRA_TEST_ERR(!A.NumMyRows()==B.NumMyRows(),ierr);
00401   EPETRA_TEST_ERR(!A.OperatorDomainMap().SameAs(B.OperatorDomainMap()),ierr);
00402   EPETRA_TEST_ERR(!A.OperatorRangeMap().SameAs(B.OperatorRangeMap()),ierr);
00403   EPETRA_TEST_ERR(!A.RowMatrixColMap().SameAs(B.RowMatrixColMap()),ierr);
00404   EPETRA_TEST_ERR(!A.RowMatrixRowMap().SameAs(B.RowMatrixRowMap()),ierr);
00405   EPETRA_TEST_ERR(!A.UpperTriangular()==B.UpperTriangular(),ierr);
00406   EPETRA_TEST_ERR(!A.UseTranspose()==B.UseTranspose(),ierr);
00407 
00408   int NumVectors = 5;
00409   { // No transpose case
00410     Epetra_MultiVector X(A.OperatorDomainMap(), NumVectors);
00411     Epetra_MultiVector YA1(A.OperatorRangeMap(), NumVectors);
00412     Epetra_MultiVector YA2(YA1);
00413     Epetra_MultiVector YB1(YA1);
00414     Epetra_MultiVector YB2(YA1);
00415     X.Random();
00416     
00417     bool transA = false;
00418     A.SetUseTranspose(transA);
00419     B.SetUseTranspose(transA);
00420     A.Apply(X,YA1);
00421     A.Multiply(transA, X, YA2);
00422     EPETRA_TEST_ERR(checkMultiVectors(YA1,YA2,"A Multiply and A Apply", verbose),ierr);
00423     B.Apply(X,YB1);
00424     EPETRA_TEST_ERR(checkMultiVectors(YA1,YB1,"A Multiply and B Multiply", verbose),ierr);
00425     B.Multiply(transA, X, YB2);
00426     EPETRA_TEST_ERR(checkMultiVectors(YA1,YB2,"A Multiply and B Apply", verbose), ierr);
00427     
00428   }
00429   {// transpose case
00430     Epetra_MultiVector X(A.OperatorRangeMap(), NumVectors);
00431     Epetra_MultiVector YA1(A.OperatorDomainMap(), NumVectors);
00432     Epetra_MultiVector YA2(YA1);
00433     Epetra_MultiVector YB1(YA1);
00434     Epetra_MultiVector YB2(YA1);
00435     X.Random();
00436     
00437     bool transA = true;
00438     A.SetUseTranspose(transA);
00439     B.SetUseTranspose(transA);
00440     A.Apply(X,YA1);
00441     A.Multiply(transA, X, YA2);
00442     EPETRA_TEST_ERR(checkMultiVectors(YA1,YA2, "A Multiply and A Apply (transpose)", verbose),ierr);
00443     B.Apply(X,YB1);
00444     EPETRA_TEST_ERR(checkMultiVectors(YA1,YB1, "A Multiply and B Multiply (transpose)", verbose),ierr);
00445     B.Multiply(transA, X,YB2);
00446     EPETRA_TEST_ERR(checkMultiVectors(YA1,YB2, "A Multiply and B Apply (transpose)", verbose),ierr);
00447     
00448   }
00449 
00450   Epetra_Vector diagA(A.RowMatrixRowMap());
00451   EPETRA_TEST_ERR(A.ExtractDiagonalCopy(diagA),ierr);
00452   Epetra_Vector diagB(B.RowMatrixRowMap());
00453   EPETRA_TEST_ERR(B.ExtractDiagonalCopy(diagB),ierr);
00454   EPETRA_TEST_ERR(checkMultiVectors(diagA,diagB, "ExtractDiagonalCopy", verbose),ierr);
00455 
00456   Epetra_Vector rowA(A.RowMatrixRowMap());
00457   EPETRA_TEST_ERR(A.InvRowSums(rowA),ierr);
00458   Epetra_Vector rowB(B.RowMatrixRowMap());
00459   EPETRA_TEST_ERR(B.InvRowSums(rowB),ierr)
00460   EPETRA_TEST_ERR(checkMultiVectors(rowA,rowB, "InvRowSums", verbose),ierr);
00461 
00462   Epetra_Vector colA(A.RowMatrixColMap());
00463   EPETRA_TEST_ERR(A.InvColSums(colA),ierr);
00464   Epetra_Vector colB(B.RowMatrixColMap());
00465   EPETRA_TEST_ERR(B.InvColSums(colB),ierr);
00466   EPETRA_TEST_ERR(checkMultiVectors(colA,colB, "InvColSums", verbose),ierr);
00467 
00468   EPETRA_TEST_ERR(checkValues(A.NormInf(), B.NormInf(), "NormInf before scaling", verbose), ierr);
00469   EPETRA_TEST_ERR(checkValues(A.NormOne(), B.NormOne(), "NormOne before scaling", verbose),ierr);
00470 
00471   EPETRA_TEST_ERR(A.RightScale(colA),ierr);  
00472   EPETRA_TEST_ERR(B.RightScale(colB),ierr);
00473 
00474 
00475   EPETRA_TEST_ERR(A.LeftScale(rowA),ierr);
00476   EPETRA_TEST_ERR(B.LeftScale(rowB),ierr);
00477 
00478 
00479   EPETRA_TEST_ERR(checkValues(A.NormInf(), B.NormInf(), "NormInf after scaling", verbose), ierr);
00480   EPETRA_TEST_ERR(checkValues(A.NormOne(), B.NormOne(), "NormOne after scaling", verbose),ierr);
00481 
00482   vector<double> valuesA(A.MaxNumEntries());
00483   vector<int> indicesA(A.MaxNumEntries());  
00484   vector<double> valuesB(B.MaxNumEntries());
00485   vector<int> indicesB(B.MaxNumEntries());
00486   return(0);
00487   for (int i=0; i<A.NumMyRows(); i++) {
00488     int nA, nB;
00489     EPETRA_TEST_ERR(A.ExtractMyRowCopy(i, A.MaxNumEntries(), nA, &valuesA[0], &indicesA[0]),ierr); 
00490     EPETRA_TEST_ERR(B.ExtractMyRowCopy(i, B.MaxNumEntries(), nB, &valuesB[0], &indicesB[0]),ierr);
00491     EPETRA_TEST_ERR(!nA==nB,ierr);
00492     for (int j=0; j<nA; j++) {
00493       double curVal = valuesA[j];
00494       int curIndex = indicesA[j];
00495       bool notfound = true;
00496       int jj = 0;
00497       while (notfound && jj< nB) {
00498   if (!checkValues(curVal, valuesB[jj])) notfound = false;
00499   jj++;
00500       }
00501       EPETRA_TEST_ERR(notfound, ierr);
00502       vector<int>::iterator p = find(indicesB.begin(),indicesB.end(),curIndex);  // find curIndex in indicesB
00503       EPETRA_TEST_ERR(p==indicesB.end(), ierr);
00504     }
00505 
00506   }
00507   if (verbose) cout << "RowMatrix Methods check OK" << endl;
00508 
00509   return (ierr);
00510 }

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