EpetraExt Package Browser (Single Doxygen Collection) Development
test/hypre/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 "EpetraExt_HypreIJMatrix.h"
00033 #include "HYPRE_IJ_mv.h"
00034 #include "hypre_Helpers.hpp"
00035 #include "Epetra_Vector.h"
00036 #include "Epetra_Flops.h"
00037 #ifdef EPETRA_MPI
00038 #include "Epetra_MpiComm.h"
00039 #include "mpi.h"
00040 #else
00041 #include "Epetra_SerialComm.h"
00042 #endif
00043 #include "../epetra_test_err.h"
00044 #include "Epetra_Version.h"
00045 #include <vector>
00046 #include <algorithm>
00047 #include <string>
00048 #include <iostream>
00049 #include <fstream>
00050 
00051 // prototypes
00052 
00053 int checkValues( double x, double y, string message = "", bool verbose = false) { 
00054   if (fabs((x-y)/x) > 0.01) {
00055     return(1); 
00056     if (verbose) cout << "********** " << message << " check failed.********** " << endl;
00057   }
00058   else {
00059     if (verbose) cout << message << " check OK." << endl;    
00060     return(0);
00061   }
00062 }
00063 
00064 int checkMultiVectors( Epetra_MultiVector & X, Epetra_MultiVector & Y, string message = "", bool verbose = false) {
00065   int numVectors = X.NumVectors();
00066   int length = Y.MyLength();
00067   int badvalue = 0;
00068   int globalbadvalue = 0;
00069   for (int j=0; j<numVectors; j++)
00070     for (int i=0; i< length; i++) 
00071       if (checkValues(X[j][i], Y[j][i])==1) badvalue = 1;
00072   X.Map().Comm().MaxAll(&badvalue, &globalbadvalue, 1);
00073 
00074   if (verbose) {
00075     if (globalbadvalue==0) cout << message << " check OK." << endl;
00076     else cout << "********* " << message << " check failed.********** " << endl;
00077   }
00078   return(globalbadvalue);
00079 }
00080 
00081 int check(Epetra_RowMatrix & A, Epetra_RowMatrix & B, bool verbose);
00082 
00083 int powerMethodTests(Epetra_RowMatrix & A, Epetra_RowMatrix & JadA, Epetra_Map & Map, 
00084          Epetra_Vector & q, Epetra_Vector & z, Epetra_Vector & resid, bool verbose);
00085 
00086 int power_method(bool TransA, Epetra_RowMatrix& A, 
00087      Epetra_Vector& q,
00088      Epetra_Vector& z0, 
00089      Epetra_Vector& resid, 
00090      double * lambda, int niters, double tolerance,
00091      bool verbose);
00092 
00093 int main(int argc, char *argv[])
00094 {
00095   int ierr = 0, i, forierr = 0;
00096 #ifdef EPETRA_MPI
00097 
00098   // Initialize MPI
00099 
00100   MPI_Init(&argc,&argv);
00101   int rank; // My process ID
00102 
00103   MPI_Comm_rank(MPI_COMM_WORLD, &rank);
00104   Epetra_MpiComm Comm( MPI_COMM_WORLD );
00105 
00106 #else
00107 
00108   int rank = 0;
00109   Epetra_SerialComm Comm;
00110 
00111 #endif
00112 
00113   bool verbose = false;
00114 
00115   // Check if we should print results to standard out
00116   if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
00117 
00118   int verbose_int = verbose ? 1 : 0;
00119   Comm.Broadcast(&verbose_int, 1, 0);
00120   verbose = verbose_int==1 ? true : false;
00121 
00122 
00123   //  char tmp;
00124   //  if (rank==0) cout << "Press any key to continue..."<< endl;
00125   //  if (rank==0) cin >> tmp;
00126   //  Comm.Barrier();
00127 
00128   Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
00129   int MyPID = Comm.MyPID();
00130   int NumProc = Comm.NumProc();
00131 
00132   if(verbose && MyPID==0)
00133     cout << Epetra_Version() << endl << endl;
00134 
00135   if (verbose) cout << "Processor "<<MyPID<<" of "<< NumProc
00136         << " is alive."<<endl;
00137 
00138   // Redefine verbose to only print on PE 0
00139   if(verbose && rank!=0) 
00140     verbose = false;
00141 
00142   int NumMyEquations = 10000;
00143   int NumGlobalEquations = (NumMyEquations * NumProc) + EPETRA_MIN(NumProc,3);
00144   if(MyPID < 3) 
00145     NumMyEquations++;
00146 
00147   // Construct a Map that puts approximately the same Number of equations on each processor
00148 
00149   Epetra_Map Map(NumGlobalEquations, NumMyEquations, 0, Comm);
00150   
00151   // Get update list and number of local equations from newly created Map
00152   vector<int> MyGlobalElements(Map.NumMyElements());
00153   Map.MyGlobalElements(&MyGlobalElements[0]);
00154 
00155   // Create an integer vector NumNz that is used to build the Petra Matrix.
00156   // NumNz[i] is the Number of OFF-DIAGONAL term for the ith global equation on this processor
00157 
00158   vector<int> NumNz(NumMyEquations);
00159 
00160   // We are building a tridiagonal matrix where each row has (-1 2 -1)
00161   // So we need 2 off-diagonal terms (except for the first and last equation)
00162 
00163   for(i = 0; i < NumMyEquations; i++)
00164     if((MyGlobalElements[i] == 0) || (MyGlobalElements[i] == NumGlobalEquations - 1))
00165       NumNz[i] = 1;
00166     else
00167       NumNz[i] = 2;
00168 
00169   // Create a Epetra_Matrix
00170 
00171   Epetra_CrsMatrix A(Copy, Map, &NumNz[0]);
00172   EPETRA_TEST_ERR(A.IndicesAreGlobal(),ierr);
00173   EPETRA_TEST_ERR(A.IndicesAreLocal(),ierr);
00174   
00175   // Add  rows one-at-a-time
00176   // Need some vectors to help
00177   // Off diagonal Values will always be -1
00178 
00179 
00180   vector<double> Values(2);
00181   Values[0] = -1.0; 
00182   Values[1] = -1.0;
00183   vector<int> Indices(2);
00184   double two = 2.0;
00185   int NumEntries;
00186 
00187   forierr = 0;
00188   for(i = 0; i < NumMyEquations; i++) {
00189     if(MyGlobalElements[i] == 0) {
00190       Indices[0] = 1;
00191       NumEntries = 1;
00192     }
00193     else if (MyGlobalElements[i] == NumGlobalEquations-1) {
00194       Indices[0] = NumGlobalEquations-2;
00195       NumEntries = 1;
00196     }
00197     else {
00198       Indices[0] = MyGlobalElements[i]-1;
00199       Indices[1] = MyGlobalElements[i]+1;
00200       NumEntries = 2;
00201     }
00202     forierr += !(A.InsertGlobalValues(MyGlobalElements[i], NumEntries, &Values[0], &Indices[0])==0);
00203     forierr += !(A.InsertGlobalValues(MyGlobalElements[i], 1, &two, &MyGlobalElements[i])>0); // Put in the diagonal entry
00204   }
00205   EPETRA_TEST_ERR(forierr,ierr);
00206 
00207   // Finish up
00208   A.FillComplete();
00209   A.OptimizeStorage();
00210 
00211   HYPRE_IJMatrix Matrix;
00212   int ilower = Map.MinMyGID();
00213   int iupper = Map.MaxMyGID();
00214 
00215   //printf("Proc[%d], ilower = %d, iupper = %d.\n", MyPID, ilower, iupper);
00216   HYPRE_IJMatrixCreate(MPI_COMM_WORLD, ilower, iupper, ilower, iupper, &Matrix);
00217   HYPRE_IJMatrixSetObjectType(Matrix, HYPRE_PARCSR);
00218   HYPRE_IJMatrixInitialize(Matrix);
00219   
00220   for(i = 0; i < A.NumMyRows(); i++){
00221     int numElements;
00222     A.NumMyRowEntries(i, numElements);
00223     vector<int> my_indices; my_indices.resize(numElements);
00224     vector<double> my_values; my_values.resize(numElements);
00225     int numEntries;
00226     A.ExtractMyRowCopy(i, numElements, numEntries, &my_values[0], &my_indices[0]);
00227     for(int j = 0; j < numEntries; j++) {
00228       my_indices[j] = A.GCID(my_indices[j]);
00229     }
00230     int GlobalRow[1];
00231     GlobalRow[0] = A.GRID(i);
00232     HYPRE_IJMatrixSetValues(Matrix, 1, &numEntries, GlobalRow, &my_indices[0], &my_values[0]);
00233   } 
00234   HYPRE_IJMatrixAssemble(Matrix);
00235 
00236   EpetraExt_HypreIJMatrix JadA(Matrix);
00237  
00238   JadA.SetMaps(JadA.RowMatrixRowMap(), A.RowMatrixColMap());
00239   // Create vectors for Power method
00240 
00241   Epetra_Vector q(Map);
00242   Epetra_Vector z(Map); z.Random();
00243   Epetra_Vector resid(Map);
00244 
00245   Epetra_Flops flopcounter;
00246   A.SetFlopCounter(flopcounter);
00247   q.SetFlopCounter(A);
00248   z.SetFlopCounter(A);
00249   resid.SetFlopCounter(A);
00250   JadA.SetFlopCounter(A);
00251 
00252   if (verbose) cout << "=======================================" << endl
00253         << "Testing Jad using CrsMatrix as input..." << endl
00254         << "=======================================" << endl;
00255 
00256   A.ResetFlops();
00257   powerMethodTests(A, JadA, Map, q, z, resid, verbose);
00258 
00259 #ifdef EPETRA_MPI
00260   MPI_Finalize() ;
00261 #endif
00262 
00263 return ierr ;
00264 }
00265 
00266 int powerMethodTests(Epetra_RowMatrix & A, Epetra_RowMatrix & JadA, Epetra_Map & Map, 
00267          Epetra_Vector & q, Epetra_Vector & z, Epetra_Vector & resid, bool verbose) {
00268 
00269   // variable needed for iteration
00270   double lambda = 0.0;
00271   // int niters = 10000;
00272   int niters = 300;
00273   double tolerance = 1.0e-2;
00274   int ierr = 0;
00275 
00277   
00278   // Iterate
00279 
00280   Epetra_Time timer(Map.Comm());
00281   
00282   double startTime = timer.ElapsedTime();
00283   EPETRA_TEST_ERR(power_method(false, A, q, z, resid, &lambda, niters, tolerance, verbose),ierr);
00284   double elapsed_time = timer.ElapsedTime() - startTime;
00285   double total_flops = q.Flops();
00286   double MFLOPs = total_flops/elapsed_time/1000000.0;
00287   double lambdaref = lambda;
00288   double flopsref = total_flops;
00289 
00290   if (verbose) 
00291     cout << "\n\nTotal MFLOPs for reference first solve = " << MFLOPs << endl
00292       <<     "Total FLOPS                            = " <<total_flops <<endl<<endl;
00293 
00294   lambda = 0.0;
00295   startTime = timer.ElapsedTime();
00296   
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   
00404   EPETRA_TEST_ERR(!A.RowMatrixColMap().SameAs(B.RowMatrixColMap()),ierr);
00405   EPETRA_TEST_ERR(!A.RowMatrixRowMap().SameAs(B.RowMatrixRowMap()),ierr);
00406   EPETRA_TEST_ERR(!A.UpperTriangular()==B.UpperTriangular(),ierr);
00407   EPETRA_TEST_ERR(!A.UseTranspose()==B.UseTranspose(),ierr);
00408 
00409   int NumVectors = 5;
00410   { // No transpose case
00411     Epetra_MultiVector X(A.OperatorDomainMap(), NumVectors);
00412     Epetra_MultiVector YA1(A.OperatorRangeMap(), NumVectors);
00413     Epetra_MultiVector YA2(YA1);
00414     Epetra_MultiVector YB1(YA1);
00415     Epetra_MultiVector YB2(YA1);
00416     X.Random();
00417     
00418     bool transA = false;
00419     A.SetUseTranspose(transA);
00420     B.SetUseTranspose(transA);
00421     A.Apply(X,YA1);
00422     A.Multiply(transA, X, YA2);
00423     EPETRA_TEST_ERR(checkMultiVectors(YA1,YA2,"A Multiply and A Apply", verbose),ierr);
00424     B.Apply(X,YB1);
00425     EPETRA_TEST_ERR(checkMultiVectors(YA1,YB1,"A Multiply and B Multiply", verbose),ierr);
00426     B.Multiply(transA, X, YB2);
00427     EPETRA_TEST_ERR(checkMultiVectors(YA1,YB2,"A Multiply and B Apply", verbose), ierr);
00428     
00429   }
00430   {// transpose case
00431     Epetra_MultiVector X(A.OperatorRangeMap(), NumVectors);
00432     Epetra_MultiVector YA1(A.OperatorDomainMap(), NumVectors);
00433     Epetra_MultiVector YA2(YA1);
00434     Epetra_MultiVector YB1(YA1);
00435     Epetra_MultiVector YB2(YA1);
00436     X.Random();
00437     
00438     bool transA = true;
00439     A.SetUseTranspose(transA);
00440     B.SetUseTranspose(transA);
00441     A.Apply(X,YA1);
00442     A.Multiply(transA, X, YA2);
00443     EPETRA_TEST_ERR(checkMultiVectors(YA1,YA2, "A Multiply and A Apply (transpose)", verbose),ierr);
00444     B.Apply(X,YB1);
00445     EPETRA_TEST_ERR(checkMultiVectors(YA1,YB1, "A Multiply and B Multiply (transpose)", verbose),ierr);
00446     B.Multiply(transA, X,YB2);
00447     EPETRA_TEST_ERR(checkMultiVectors(YA1,YB2, "A Multiply and B Apply (transpose)", verbose),ierr);
00448     
00449   }
00450 
00451   Epetra_Vector diagA(A.RowMatrixRowMap());
00452   EPETRA_TEST_ERR(A.ExtractDiagonalCopy(diagA),ierr);
00453   Epetra_Vector diagB(B.RowMatrixRowMap());
00454   EPETRA_TEST_ERR(B.ExtractDiagonalCopy(diagB),ierr);
00455   EPETRA_TEST_ERR(checkMultiVectors(diagA,diagB, "ExtractDiagonalCopy", verbose),ierr);
00456 
00457   Epetra_Vector rowA(A.RowMatrixRowMap());
00458   EPETRA_TEST_ERR(A.InvRowSums(rowA),ierr);
00459   Epetra_Vector rowB(B.RowMatrixRowMap());
00460   EPETRA_TEST_ERR(B.InvRowSums(rowB),ierr)
00461   EPETRA_TEST_ERR(checkMultiVectors(rowA,rowB, "InvRowSums", verbose),ierr);
00462 
00463   Epetra_Vector colA(A.RowMatrixColMap());
00464   EPETRA_TEST_ERR(A.InvColSums(colA),ierr);
00465   Epetra_Vector colB(B.RowMatrixColMap());
00466   EPETRA_TEST_ERR(B.InvColSums(colB),ierr);
00467   EPETRA_TEST_ERR(checkMultiVectors(colA,colB, "InvColSums", verbose),ierr);
00468 
00469   EPETRA_TEST_ERR(checkValues(A.NormInf(), B.NormInf(), "NormInf before scaling", verbose), ierr);
00470   EPETRA_TEST_ERR(checkValues(A.NormOne(), B.NormOne(), "NormOne before scaling", verbose),ierr);
00471 
00472   //EPETRA_TEST_ERR(A.RightScale(colA),ierr);  
00473   //EPETRA_TEST_ERR(B.RightScale(colB),ierr);
00474 
00475 
00476   EPETRA_TEST_ERR(A.LeftScale(rowA),ierr);
00477   EPETRA_TEST_ERR(B.LeftScale(rowB),ierr);
00478 
00479 
00480   EPETRA_TEST_ERR(checkValues(A.NormInf(), B.NormInf(), "NormInf after scaling", verbose), ierr);
00481   EPETRA_TEST_ERR(checkValues(A.NormOne(), B.NormOne(), "NormOne after scaling", verbose),ierr);
00482 
00483   vector<double> valuesA(A.MaxNumEntries());
00484   vector<int> indicesA(A.MaxNumEntries());  
00485   vector<double> valuesB(B.MaxNumEntries());
00486   vector<int> indicesB(B.MaxNumEntries());
00487   return(0);
00488   for (int i=0; i<A.NumMyRows(); i++) {
00489     int nA, nB;
00490     EPETRA_TEST_ERR(A.ExtractMyRowCopy(i, A.MaxNumEntries(), nA, &valuesA[0], &indicesA[0]),ierr); 
00491     EPETRA_TEST_ERR(B.ExtractMyRowCopy(i, B.MaxNumEntries(), nB, &valuesB[0], &indicesB[0]),ierr);
00492     EPETRA_TEST_ERR(!nA==nB,ierr);
00493     for (int j=0; j<nA; j++) {
00494       double curVal = valuesA[j];
00495       int curIndex = indicesA[j];
00496       bool notfound = true;
00497       int jj = 0;
00498       while (notfound && jj< nB) {
00499   if (!checkValues(curVal, valuesB[jj])) notfound = false;
00500   jj++;
00501       }
00502       EPETRA_TEST_ERR(notfound, ierr);
00503       vector<int>::iterator p = find(indicesB.begin(),indicesB.end(),curIndex);  // find curIndex in indicesB
00504       EPETRA_TEST_ERR(p==indicesB.end(), ierr);
00505     }
00506 
00507   }
00508   if (verbose) cout << "RowMatrix Methods check OK" << endl;
00509 
00510   return (ierr);
00511 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines