example/CrsSingletonFilterExample/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_Comm.h"
00030 #include "Epetra_Map.h"
00031 #include "Epetra_Time.h"
00032 #include "Epetra_BlockMap.h"
00033 #include "Epetra_MultiVector.h"
00034 #include "Epetra_Vector.h"
00035 #include "Epetra_Export.h"
00036 #include "Epetra_LinearProblem.h"
00037 #include "Epetra_CrsSingletonFilter.h"
00038 
00039 #include "Epetra_VbrMatrix.h"
00040 #include "Epetra_CrsMatrix.h"
00041 #ifdef EPETRA_MPI
00042 #include "mpi.h"
00043 #include "Epetra_MpiComm.h"
00044 #else
00045 #include "Epetra_SerialComm.h"
00046 #endif
00047 #include "Trilinos_Util.h"
00048 #include "Ifpack_IlukGraph.h"
00049 #include "Ifpack_CrsRiluk.h"
00050 
00051 void BiCGSTAB(Epetra_CrsMatrix &A, Epetra_Vector &x, Epetra_Vector &b, 
00052         Ifpack_CrsRiluk *M, 
00053         int Maxiter, double Tolerance, 
00054         double *residual, bool verbose);
00055 int Statistics(const Epetra_CrsSingletonFilter & filter);
00056 int main(int argc, char *argv[]) {
00057 
00058 #ifdef EPETRA_MPI
00059   MPI_Init(&argc,&argv);
00060   Epetra_MpiComm Comm (MPI_COMM_WORLD);
00061 #else
00062   Epetra_SerialComm Comm;
00063 #endif
00064 
00065   cout << Comm << endl;
00066 
00067   int MyPID = Comm.MyPID();
00068 
00069   bool verbose = false;
00070   bool verbose1 = true;
00071   if (MyPID==0) verbose = true;
00072 
00073   if(argc < 2 && verbose) {
00074     cerr << "Usage: " << argv[0] 
00075    << " HB_filename [level_fill [level_overlap [absolute_threshold [ relative_threshold]]]]" << endl
00076    << "where:" << endl
00077    << "HB_filename        - filename and path of a Harwell-Boeing data set" << endl
00078    << "level_fill         - The amount of fill to use for ILU(k) preconditioner (default 0)" << endl
00079    << "level_overlap      - The amount of overlap used for overlapping Schwarz subdomains (default 0)" << endl
00080    << "absolute_threshold - The minimum value to place on the diagonal prior to factorization (default 0.0)" << endl
00081    << "relative_threshold - The relative amount to perturb the diagonal prior to factorization (default 1.0)" << endl << endl
00082    << "To specify a non-default value for one of these parameters, you must specify all" << endl
00083    << " preceding values but not any subsequent parameters. Example:" << endl
00084    << "ifpackHpcSerialMsr.exe mymatrix.hpc 1  - loads mymatrix.hpc, uses level fill of one, all other values are defaults" << endl
00085    << endl;
00086     return(1);
00087 
00088   }
00089 
00090   // Uncomment the next three lines to debug in mpi mode
00091   //int tmp;
00092   //if (MyPID==0) cin >> tmp;
00093   //Comm.Barrier();
00094 
00095   Epetra_Map * readMap;
00096   Epetra_CrsMatrix * readA; 
00097   Epetra_Vector * readx; 
00098   Epetra_Vector * readb;
00099   Epetra_Vector * readxexact;
00100    
00101   // Call routine to read in HB problem
00102   Trilinos_Util_ReadHb2Epetra(argv[1], Comm, readMap, readA, readx, readb, readxexact);
00103 
00104   // Create uniform distributed map
00105   Epetra_Map map(readMap->NumGlobalElements(), 0, Comm);
00106 
00107   // Create Exporter to distribute read-in matrix and vectors
00108 
00109   Epetra_Export exporter(*readMap, map);
00110   Epetra_CrsMatrix A(Copy, map, 0);
00111   Epetra_Vector x(map);
00112   Epetra_Vector b(map);
00113   Epetra_Vector xexact(map);
00114 
00115   Epetra_Time FillTimer(Comm);
00116   x.Export(*readx, exporter, Add);
00117   b.Export(*readb, exporter, Add);
00118   xexact.Export(*readxexact, exporter, Add);
00119   Comm.Barrier();
00120   double vectorRedistributeTime = FillTimer.ElapsedTime();
00121   A.Export(*readA, exporter, Add);
00122   Comm.Barrier();
00123   double matrixRedistributeTime = FillTimer.ElapsedTime() - vectorRedistributeTime;
00124   assert(A.FillComplete()==0);    
00125   Comm.Barrier();
00126   double fillCompleteTime = FillTimer.ElapsedTime() - matrixRedistributeTime;
00127   if (Comm.MyPID()==0)  {
00128     cout << "\n\n****************************************************" << endl;
00129     cout << "\n Vector redistribute  time (sec) = " << vectorRedistributeTime<< endl;
00130     cout << "    Matrix redistribute time (sec) = " << matrixRedistributeTime << endl;
00131     cout << "    Transform to Local  time (sec) = " << fillCompleteTime << endl<< endl;
00132   }
00133   Epetra_Vector tmp1(*readMap);
00134   Epetra_Vector tmp2(map);
00135   readA->Multiply(false, *readxexact, tmp1);
00136 
00137   A.Multiply(false, xexact, tmp2);
00138   double residual;
00139   tmp1.Norm2(&residual);
00140   if (verbose) cout << "Norm of Ax from file            = " << residual << endl;
00141   tmp2.Norm2(&residual);
00142   if (verbose) cout << "Norm of Ax after redistribution = " << residual << endl << endl << endl;
00143 
00144   //cout << "A from file = " << *readA << endl << endl << endl;
00145 
00146   //cout << "A after dist = " << A << endl << endl << endl;
00147 
00148   delete readA;
00149   delete readx;
00150   delete readb;
00151   delete readxexact;
00152   delete readMap;
00153 
00154   Comm.Barrier();
00155 
00156   bool smallProblem = false;
00157   if (A.RowMap().NumGlobalElements()<100) smallProblem = true;
00158 
00159   if (smallProblem)
00160     cout << "Original Matrix = " << endl << A   << endl;
00161 
00162   x.PutScalar(0.0);
00163 
00164   Epetra_LinearProblem FullProblem(&A, &x, &b);
00165   double normb, norma;
00166   b.NormInf(&normb);
00167   norma = A.NormInf();
00168   if (verbose)
00169     cout << "Inf norm of Original Matrix = " << norma << endl
00170    << "Inf norm of Original RHS    = " << normb << endl;
00171   
00172   Epetra_Time ReductionTimer(Comm);
00173   Epetra_CrsSingletonFilter SingletonFilter;
00174   Comm.Barrier();
00175   double reduceInitTime = ReductionTimer.ElapsedTime();
00176   SingletonFilter.Analyze(&A);
00177   Comm.Barrier();
00178   double reduceAnalyzeTime = ReductionTimer.ElapsedTime() - reduceInitTime;
00179 
00180   if (SingletonFilter.SingletonsDetected())
00181     cout << "Singletons found" << endl;
00182   else {
00183     cout << "Singletons not found" << endl;
00184     exit(1);
00185   }
00186   SingletonFilter.ConstructReducedProblem(&FullProblem);
00187   Comm.Barrier();
00188   double reduceConstructTime = ReductionTimer.ElapsedTime() - reduceInitTime;
00189 
00190   double totalReduceTime = ReductionTimer.ElapsedTime();
00191 
00192   if (verbose)
00193     cout << "\n\n****************************************************" << endl
00194    << "    Reduction init  time (sec)           = " << reduceInitTime<< endl
00195    << "    Reduction Analyze time (sec)         = " << reduceAnalyzeTime << endl
00196    << "    Construct Reduced Problem time (sec) = " << reduceConstructTime << endl
00197    << "    Reduction Total time (sec)           = " << totalReduceTime << endl<< endl;
00198 
00199   Statistics(SingletonFilter);
00200 
00201   Epetra_LinearProblem * ReducedProblem = SingletonFilter.ReducedProblem();
00202 
00203   Epetra_CrsMatrix * Ap = dynamic_cast<Epetra_CrsMatrix *>(ReducedProblem->GetMatrix());
00204   Epetra_Vector * bp = (*ReducedProblem->GetRHS())(0);
00205   Epetra_Vector * xp = (*ReducedProblem->GetLHS())(0);
00206   
00207 
00208   if (smallProblem)
00209     cout << " Reduced Matrix = " << endl << *Ap << endl
00210    << " LHS before sol = " << endl << *xp << endl
00211    << " RHS            = " << endl << *bp << endl;
00212 
00213   // Construct ILU preconditioner
00214 
00215   double elapsed_time, total_flops, MFLOPs;
00216   Epetra_Time timer(Comm);
00217 
00218   int LevelFill = 0;
00219   if (argc > 2)  LevelFill = atoi(argv[2]);
00220   if (verbose) cout << "Using Level Fill = " << LevelFill << endl;
00221   int Overlap = 0;
00222   if (argc > 3) Overlap = atoi(argv[3]);
00223   if (verbose) cout << "Using Level Overlap = " << Overlap << endl;
00224   double Athresh = 0.0;
00225   if (argc > 4) Athresh = atof(argv[4]);
00226   if (verbose) cout << "Using Absolute Threshold Value of = " << Athresh << endl;
00227 
00228   double Rthresh = 1.0;
00229   if (argc > 5) Rthresh = atof(argv[5]);
00230   if (verbose) cout << "Using Relative Threshold Value of = " << Rthresh << endl;
00231 
00232   Ifpack_IlukGraph * IlukGraph = 0;
00233   Ifpack_CrsRiluk * ILUK = 0;
00234 
00235   if (LevelFill>-1) {
00236     elapsed_time = timer.ElapsedTime();
00237     IlukGraph = new Ifpack_IlukGraph(Ap->Graph(), LevelFill, Overlap);
00238     assert(IlukGraph->ConstructFilledGraph()==0);
00239     elapsed_time = timer.ElapsedTime() - elapsed_time;
00240     if (verbose) cout << "Time to construct ILUK graph = " << elapsed_time << endl;
00241 
00242 
00243     Epetra_Flops fact_counter;
00244   
00245     elapsed_time = timer.ElapsedTime();
00246     ILUK = new Ifpack_CrsRiluk(*IlukGraph);
00247     ILUK->SetFlopCounter(fact_counter);
00248     ILUK->SetAbsoluteThreshold(Athresh);
00249     ILUK->SetRelativeThreshold(Rthresh);
00250     //assert(ILUK->InitValues()==0);
00251     int initerr = ILUK->InitValues(*Ap);
00252     if (initerr!=0) {
00253       cout << endl << Comm << endl << "  InitValues error = " << initerr;
00254       if (initerr==1) cout << "  Zero diagonal found, warning error only";
00255       cout << endl << endl;
00256     }
00257     assert(ILUK->Factor()==0);
00258     elapsed_time = timer.ElapsedTime() - elapsed_time;
00259     total_flops = ILUK->Flops();
00260     MFLOPs = total_flops/elapsed_time/1000000.0;
00261     if (verbose) cout << "Time to compute preconditioner values = " 
00262         << elapsed_time << endl
00263         << "MFLOPS for Factorization = " << MFLOPs << endl;
00264     //cout << *ILUK << endl;
00265   double Condest;
00266   ILUK->Condest(false, Condest);
00267 
00268   if (verbose) cout << "Condition number estimate for this preconditioner = " << Condest << endl;
00269   }
00270   int Maxiter = 100;
00271   double Tolerance = 1.0E-8;
00272 
00273   Epetra_Flops counter;
00274   Ap->SetFlopCounter(counter);
00275   xp->SetFlopCounter(*Ap);
00276   bp->SetFlopCounter(*Ap);
00277   if (ILUK!=0) ILUK->SetFlopCounter(*Ap);
00278 
00279   elapsed_time = timer.ElapsedTime();
00280 
00281   double normreducedb, normreduceda;
00282   bp->NormInf(&normreducedb);
00283   normreduceda = Ap->NormInf();
00284   if (verbose) 
00285     cout << "Inf norm of Reduced Matrix = " << normreduceda << endl
00286    << "Inf norm of Reduced RHS    = " << normreducedb << endl;
00287 
00288   BiCGSTAB(*Ap, *xp, *bp, ILUK, Maxiter, Tolerance, &residual, verbose);
00289 
00290   elapsed_time = timer.ElapsedTime() - elapsed_time;
00291   total_flops = counter.Flops();
00292   MFLOPs = total_flops/elapsed_time/1000000.0;
00293   if (verbose) cout << "Time to compute solution = " 
00294         << elapsed_time << endl
00295         << "Number of operations in solve = " << total_flops << endl
00296         << "MFLOPS for Solve = " << MFLOPs<< endl << endl;
00297 
00298   SingletonFilter.ComputeFullSolution();
00299 
00300   if (smallProblem)
00301   cout << " Reduced LHS after sol = " << endl << *xp << endl
00302        << " Full    LHS after sol = " << endl << x << endl
00303        << " Full  Exact LHS         = " << endl << xexact << endl;
00304 
00305   Epetra_Vector resid(x);
00306 
00307   resid.Update(1.0, x, -1.0, xexact, 0.0); // resid = xcomp - xexact
00308 
00309   resid.Norm2(&residual);
00310   double normx, normxexact;
00311   x.Norm2(&normx);
00312   xexact.Norm2(&normxexact);
00313 
00314   if (verbose) 
00315     cout << "2-norm of computed solution                               = " << normx << endl
00316    << "2-norm of exact solution                                  = " << normxexact << endl
00317    << "2-norm of difference between computed and exact solution  = " << residual << endl;
00318     
00319   if (verbose1 && residual>1.0e-5) {
00320     if (verbose)
00321       cout << "Difference between computed and exact solution appears large..." << endl
00322      << "Computing norm of A times this difference.  If this norm is small, then matrix is singular"
00323      << endl;
00324     Epetra_Vector bdiff(b);
00325     assert(A.Multiply(false, resid, bdiff)==0);
00326     assert(bdiff.Norm2(&residual)==0);
00327     if (verbose) 
00328       cout << "2-norm of A times difference between computed and exact solution  = " << residual << endl;
00329     
00330   }
00331   if (verbose) 
00332     cout << "********************************************************" << endl
00333    << "              Solving again with 2*Ax=2*b" << endl
00334    << "********************************************************" << endl;
00335 
00336   A.Scale(1.0); // A = 2*A
00337   b.Scale(1.0); // b = 2*b
00338   x.PutScalar(0.0);
00339   b.NormInf(&normb);
00340   norma = A.NormInf();
00341   if (verbose)
00342     cout << "Inf norm of Original Matrix = " << norma << endl
00343    << "Inf norm of Original RHS    = " << normb << endl;
00344   double updateReducedProblemTime = ReductionTimer.ElapsedTime();
00345   SingletonFilter.UpdateReducedProblem(&FullProblem);
00346   Comm.Barrier();
00347   updateReducedProblemTime = ReductionTimer.ElapsedTime() - updateReducedProblemTime;
00348   if (verbose)
00349     cout << "\n\n****************************************************" << endl
00350    << "    Update Reduced Problem time (sec)           = " << updateReducedProblemTime<< endl
00351    << "****************************************************" << endl;
00352   Statistics(SingletonFilter);
00353 
00354   if (LevelFill>-1) {
00355 
00356     Epetra_Flops fact_counter;
00357   
00358     elapsed_time = timer.ElapsedTime();
00359 
00360     int initerr = ILUK->InitValues(*Ap);
00361     if (initerr!=0) {
00362       cout << endl << Comm << endl << "  InitValues error = " << initerr;
00363       if (initerr==1) cout << "  Zero diagonal found, warning error only";
00364       cout << endl << endl;
00365     }
00366     assert(ILUK->Factor()==0);
00367     elapsed_time = timer.ElapsedTime() - elapsed_time;
00368     total_flops = ILUK->Flops();
00369     MFLOPs = total_flops/elapsed_time/1000000.0;
00370     if (verbose) cout << "Time to compute preconditioner values = " 
00371         << elapsed_time << endl
00372         << "MFLOPS for Factorization = " << MFLOPs << endl;
00373     double Condest;
00374     ILUK->Condest(false, Condest);
00375     
00376     if (verbose) cout << "Condition number estimate for this preconditioner = " << Condest << endl;
00377   }
00378   bp->NormInf(&normreducedb);
00379   normreduceda = Ap->NormInf();
00380   if (verbose) 
00381     cout << "Inf norm of Reduced Matrix = " << normreduceda << endl
00382    << "Inf norm of Reduced RHS    = " << normreducedb << endl;
00383 
00384   BiCGSTAB(*Ap, *xp, *bp, ILUK, Maxiter, Tolerance, &residual, verbose);
00385 
00386   elapsed_time = timer.ElapsedTime() - elapsed_time;
00387   total_flops = counter.Flops();
00388   MFLOPs = total_flops/elapsed_time/1000000.0;
00389   if (verbose) cout << "Time to compute solution = " 
00390         << elapsed_time << endl
00391         << "Number of operations in solve = " << total_flops << endl
00392         << "MFLOPS for Solve = " << MFLOPs<< endl << endl;
00393 
00394   SingletonFilter.ComputeFullSolution();
00395 
00396   if (smallProblem)
00397   cout << " Reduced LHS after sol = " << endl << *xp << endl
00398        << " Full    LHS after sol = " << endl << x << endl
00399        << " Full  Exact LHS         = " << endl << xexact << endl;
00400 
00401   resid.Update(1.0, x, -1.0, xexact, 0.0); // resid = xcomp - xexact
00402 
00403   resid.Norm2(&residual);
00404   x.Norm2(&normx);
00405   xexact.Norm2(&normxexact);
00406 
00407   if (verbose) 
00408     cout << "2-norm of computed solution                               = " << normx << endl
00409    << "2-norm of exact solution                                  = " << normxexact << endl
00410    << "2-norm of difference between computed and exact solution  = " << residual << endl;
00411     
00412   if (verbose1 && residual>1.0e-5) {
00413     if (verbose)
00414       cout << "Difference between computed and exact solution appears large..." << endl
00415      << "Computing norm of A times this difference.  If this norm is small, then matrix is singular"
00416      << endl;
00417     Epetra_Vector bdiff(b);
00418     assert(A.Multiply(false, resid, bdiff)==0);
00419     assert(bdiff.Norm2(&residual)==0);
00420     if (verbose) 
00421       cout << "2-norm of A times difference between computed and exact solution  = " << residual << endl;
00422     
00423   }
00424  
00425 
00426 
00427   if (ILUK!=0) delete ILUK;
00428   if (IlukGraph!=0) delete IlukGraph;
00429                
00430 #ifdef EPETRA_MPI
00431   MPI_Finalize() ;
00432 #endif
00433 
00434 return 0 ;
00435 }
00436 void BiCGSTAB(Epetra_CrsMatrix &A, 
00437         Epetra_Vector &x, 
00438         Epetra_Vector &b, 
00439         Ifpack_CrsRiluk *M, 
00440         int Maxiter, 
00441         double Tolerance, 
00442         double *residual, bool verbose) {
00443 
00444   // Allocate vectors needed for iterations
00445   Epetra_Vector phat(x.Map()); phat.SetFlopCounter(x);
00446   Epetra_Vector p(x.Map()); p.SetFlopCounter(x);
00447   Epetra_Vector shat(x.Map()); shat.SetFlopCounter(x);
00448   Epetra_Vector s(x.Map()); s.SetFlopCounter(x);
00449   Epetra_Vector r(b.Map()); r.SetFlopCounter(x);
00450   Epetra_Vector rtilde(x.Map()); rtilde.PutScalar(1.0); rtilde.SetFlopCounter(x);
00451   Epetra_Vector v(x.Map()); v.SetFlopCounter(x);
00452   
00453 
00454   A.Multiply(false, x, r); // r = A*x
00455 
00456   r.Update(1.0, b, -1.0); // r = b - A*x
00457 
00458   double r_norm, b_norm, scaled_r_norm, rhon, rhonm1 = 1.0;
00459   double alpha = 1.0, omega = 1.0, sigma;
00460   double omega_num, omega_den;
00461   r.Norm2(&r_norm);
00462   b.Norm2(&b_norm);
00463   scaled_r_norm = r_norm/b_norm;
00464   r.Dot(rtilde,&rhon);
00465 
00466   if (verbose) cout << "Initial residual = " << r_norm
00467         << " Scaled residual = " << scaled_r_norm << endl;
00468 
00469 
00470   for (int i=0; i<Maxiter; i++) { // Main iteration loop   
00471 
00472     double beta = (rhon/rhonm1) * (alpha/omega);
00473     rhonm1 = rhon;
00474 
00475     /* p    = r + beta*(p - omega*v)       */
00476     /* phat = M^-1 p                       */
00477     /* v    = A phat                       */
00478 
00479     double dtemp = - beta*omega;
00480 
00481     p.Update(1.0, r, dtemp, v, beta);
00482     if (M==0) 
00483       phat.Scale(1.0, p);
00484     else
00485       M->Solve(false, p, phat);
00486     A.Multiply(false, phat, v);
00487 
00488     
00489     rtilde.Dot(v,&sigma);
00490     alpha = rhon/sigma;    
00491 
00492     /* s = r - alpha*v                     */
00493     /* shat = M^-1 s                       */
00494     /* r = A shat (r is a tmp here for t ) */
00495 
00496     s.Update(-alpha, v, 1.0, r, 0.0);
00497     if (M==0) 
00498       shat.Scale(1.0, s);
00499     else
00500       M->Solve(false, s, shat);
00501     A.Multiply(false, shat, r);
00502 
00503     r.Dot(s, &omega_num);
00504     r.Dot(r, &omega_den);
00505     omega = omega_num / omega_den;
00506 
00507     /* x = x + alpha*phat + omega*shat */
00508     /* r = s - omega*r */
00509 
00510     x.Update(alpha, phat, omega, shat, 1.0);
00511     r.Update(1.0, s, -omega); 
00512     
00513     r.Norm2(&r_norm);
00514     scaled_r_norm = r_norm/b_norm;
00515     r.Dot(rtilde,&rhon);
00516 
00517     if (verbose) cout << "Iter "<< i << " residual = " << r_norm
00518           << " Scaled residual = " << scaled_r_norm << endl;
00519 
00520     if (r_norm < Tolerance) break;
00521   }
00522   return;
00523 }
00524 //==============================================================================
00525 int Statistics(const Epetra_CrsSingletonFilter & filter) {
00526 
00527   // Create local variables of some of the stats we need because we don't know if the
00528   // method calls are collective or not for the particular implementation of the Row Matrix
00529   int fmaxentries;
00530   int maxentries = filter.FullMatrix()->MaxNumEntries();
00531   filter.FullMatrix()->RowMatrixRowMap().Comm().SumAll(&maxentries, &fmaxentries, 1);
00532   int fnrows = filter.FullMatrix()->NumGlobalRows();
00533   int fnnz = filter.FullMatrix()->NumGlobalNonzeros();
00534   if (filter.FullMatrix()->RowMatrixRowMap().Comm().MyPID()!=0) return(0);
00535 
00536   cout << "Full System characteristics:" << endl << endl
00537        << "  Dimension                             = " << fnrows << endl
00538        << "  Number of nonzeros                    = " <<fnnz << endl
00539        << "  Maximum Number of Row Entries         = " << fmaxentries << endl << endl
00540        << "Reduced System characteristics:" << endl << endl
00541        << "  Dimension                             = " << filter.ReducedMatrix()->NumGlobalRows() << endl
00542        << "  Number of nonzeros                    = " << filter.ReducedMatrix()->NumGlobalNonzeros() << endl
00543        << "  Maximum Number of Row Entries         = " << filter.ReducedMatrix()->GlobalMaxNumEntries() << endl << endl
00544        << "Singleton information: " << endl
00545        << "  Total number of singletons            = " << filter.NumSingletons() << endl
00546        << "  Number of rows with single entries    = " << filter.NumRowSingletons() << endl
00547        << "  Number of columns with single entries " << endl
00548        << "    (that were not already counted as " << endl
00549        << "     row singletons)                    = " << filter.NumColSingletons() << endl << endl
00550        << "Ratios: " << endl
00551        << "  Percent reduction in dimension        = " << filter.RatioOfDimensions()*100.0 << endl
00552        << "  Percent reduction in nonzero count    = " << filter.RatioOfNonzeros()*100.0 << endl << endl;
00553 
00554   return(0);
00555 
00556 }

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