Ifpack Package Browser (Single Doxygen Collection) Development
test/BlockCheby/cxx_main.cpp
Go to the documentation of this file.
00001 // @HEADER
00002 // ***********************************************************************
00003 // 
00004 //                IFPACK
00005 //                 Copyright (2004) 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 "Ifpack_ConfigDefs.h"
00030 #if defined(HAVE_IFPACK_AZTECOO) && defined(HAVE_IFPACK_EPETRAEXT)
00031 #ifdef HAVE_MPI
00032 #include "Epetra_MpiComm.h"
00033 #else
00034 #include "Epetra_SerialComm.h"
00035 #endif
00036 #include "Epetra_CrsMatrix.h"
00037 #include "Epetra_Vector.h"
00038 #include "Epetra_LinearProblem.h"
00039 #include "Trilinos_Util_CrsMatrixGallery.h"
00040 #include "Teuchos_ParameterList.hpp"
00041 #include "Ifpack_PointRelaxation.h"
00042 #include "Ifpack_BlockRelaxation.h"
00043 #include "Ifpack_DenseContainer.h"
00044 #include "Ifpack_AdditiveSchwarz.h"
00045 #include "AztecOO.h"
00046 #include "EpetraExt_BlockDiagMatrix.h"
00047 #include "EpetraExt_PointToBlockDiagPermute.h"
00048 #include "Ifpack_Chebyshev.h"
00049 
00050 using namespace Trilinos_Util;
00051 static bool verbose = false;
00052 static bool SymmetricGallery = false;
00053 static bool Solver = AZ_gmres;
00054 
00055 //=============================================
00056 // Test the BlockDiagMatrix
00057 bool TestBlockDiagMatrix(const Epetra_Comm& Comm){
00058   const int NUM_BLOCKS=30;
00059   bool TestPassed=true;
00060   int my_blockgids[NUM_BLOCKS]; 
00061   int my_blocksizes[NUM_BLOCKS];
00062 
00063   for(int i=0;i<NUM_BLOCKS;i++){
00064     my_blockgids[i]=i+NUM_BLOCKS*Comm.MyPID();
00065     if(i<NUM_BLOCKS/3)
00066       my_blocksizes[i]=1;
00067     else if(i<2*NUM_BLOCKS/3)
00068       my_blocksizes[i]=2;    
00069     else
00070       my_blocksizes[i]=3;
00071    }
00072 
00073    
00074   // Build me a map and a DBM to go with it...
00075   Epetra_BlockMap BDMap(-1,NUM_BLOCKS,my_blockgids,my_blocksizes,0,Comm);
00076   EpetraExt_BlockDiagMatrix BMAT(BDMap,true);
00077 
00078   // Fill the matrix - This tests all three block size cases in the code, 1x1, 2x2 and larger.
00079   for(int i=0;i<BMAT.NumMyBlocks();i++){
00080     double*start=BMAT[i];
00081     if(BMAT.BlockSize(i)==1)      
00082       *start=2.0;
00083     else if(BMAT.BlockSize(i)==2){
00084       start[0]=2.0;
00085       start[1]=-1.0;
00086       start[2]=-1.0;
00087       start[3]=2.0;
00088     }
00089     else if(BMAT.BlockSize(i)==3){
00090       start[0]=4.0;
00091       start[1]=-1.0;
00092       start[2]=-1.0;
00093       start[3]=-1.0;
00094       start[4]=4.0;
00095       start[5]=-1.0;
00096       start[2]=-1.0;
00097       start[7]=-1.0;
00098       start[8]=4.0;
00099     }
00100     else
00101       exit(1);
00102   } 
00103 
00104 
00105   // Allocations for tests
00106   double norm2;
00107   Epetra_MultiVector X(BDMap,1);
00108   Epetra_MultiVector Y(BDMap,1);
00109   Epetra_MultiVector Z(BDMap,1);
00110   EpetraExt_BlockDiagMatrix BMAT_forward(BMAT);
00111   EpetraExt_BlockDiagMatrix BMAT_factor(BMAT);
00112   Teuchos::ParameterList List;
00113 
00114   //***************************
00115   // Test the Forward/Invert
00116   //***************************
00117   List.set("apply mode","invert");
00118   BMAT.SetParameters(List);
00119   X.SetSeed(24601); X.Random();  
00120   BMAT_forward.ApplyInverse(X,Y);
00121   BMAT.Compute();
00122   BMAT.ApplyInverse(Y,Z);
00123   X.Update(1.0,Z,-1.0);
00124   X.Norm2(&norm2);
00125   if(!Comm.MyPID()) cout<<"Forward/Invert Error = "<<norm2<<endl;
00126   if(norm2 > 1e-12) TestPassed=false;
00127 
00128   //***************************
00129   // Test the Forward/Factor
00130   //***************************
00131   List.set("apply mode","factor");
00132   BMAT_factor.SetParameters(List);    
00133   X.SetSeed(24601); X.Random();
00134   BMAT_forward.ApplyInverse(X,Y);
00135   BMAT_factor.Compute();
00136   BMAT_factor.ApplyInverse(Y,Z);
00137   X.Update(1.0,Z,-1.0);
00138   X.Norm2(&norm2);
00139   if(!Comm.MyPID()) cout<<"Forward/Factor Error = "<<norm2<<endl;
00140   if(norm2 > 1e-12) TestPassed=false;
00141 
00142   return TestPassed;
00143 }
00144 
00145 //=============================================
00146 //============================================= 
00147 //============================================= 
00148 void Build_Local_Contiguous_Size2_BlockMatrix(const Epetra_Comm & Comm, int NUM_ROWS,int &NumBlocks,
00149                                         int *&Blockstart_, int *&Blockids_,Epetra_CrsMatrix *&MAT){
00150   double values[3]={-1.0,2.0,-1.0};
00151   int indices[3];
00152 
00153   // Build me a CrsMatrix
00154   Epetra_Map Map(-1,NUM_ROWS,0,Comm);
00155   MAT=new Epetra_CrsMatrix(Copy,Map,2);
00156   assert(MAT->NumMyRows()%2==0);
00157   
00158   NumBlocks=MAT->NumMyRows()/2;
00159   Blockstart_=new int [NumBlocks+1];
00160   Blockids_=new int [MAT->NumMyRows()];
00161   Blockstart_[0]=0;
00162   int curr_block=0;
00163 
00164   for(int i=0;i<MAT->NumMyRows();i++){
00165     // local contiguous blocks of constant size 2
00166     int row_in_block=i%2;
00167     if(row_in_block==0){
00168       Blockstart_[curr_block]=i;      
00169       indices[0]=i; indices[1]=i+1;
00170       Blockids_[i]=i;
00171       MAT->InsertGlobalValues(Map.GID(i),2,&values[1],&indices[0]);
00172     }
00173     else if(row_in_block==1){
00174       indices[0]=i-1; indices[1]=i;      
00175       MAT->InsertGlobalValues(Map.GID(i),2,&values[0],&indices[0]);
00176       Blockids_[i]=i;
00177       curr_block++;
00178     }
00179   }
00180   Blockstart_[NumBlocks]=MAT->NumMyRows();
00181   
00182   MAT->FillComplete();
00183 }
00184 
00185 //============================================= 
00186 void Build_Local_NonContiguous_Size2_BlockMatrix(const Epetra_Comm & Comm, int NUM_ROWS,int &NumBlocks,
00187                                                  int *&Blockstart_, int *&Blockids_,Epetra_CrsMatrix *&MAT){
00188   double values[3]={-1.0,2.0,-1.0};
00189   int indices[3];
00190 
00191   // Build me a CrsMatrix
00192   Epetra_Map Map(-1,NUM_ROWS,0,Comm);
00193   MAT=new Epetra_CrsMatrix(Copy,Map,2);
00194   assert(MAT->NumMyRows()%2==0);
00195 
00196   NumBlocks=MAT->NumMyRows()/2;
00197   Blockstart_=new int [NumBlocks+1];
00198   Blockids_=new int [MAT->NumMyRows()];
00199   Blockstart_[0]=0;
00200   int curr_block=0;
00201 
00202   for(int i=0;i<MAT->NumMyRows();i++){
00203     int row_in_block=(i%2)?1:0;
00204     if(row_in_block==0){
00205       int idx=i/2;
00206       Blockstart_[curr_block]=i;      
00207       indices[0]=Map.GID(idx); indices[1]=Map.GID(idx+NumBlocks);
00208       Blockids_[i]=idx;
00209       MAT->InsertGlobalValues(Map.GID(idx),2,&values[1],&indices[0]);
00210 
00211     }
00212     else if(row_in_block==1){
00213       int idx=(i-1)/2+NumBlocks;
00214       indices[0]=Map.GID(idx-NumBlocks); indices[1]=Map.GID(idx);      
00215       MAT->InsertGlobalValues(Map.GID(idx),2,&values[0],&indices[0]);
00216       Blockids_[i]=idx;
00217       curr_block++;
00218     }
00219   }
00220   Blockstart_[NumBlocks]=MAT->NumMyRows();
00221   
00222   MAT->FillComplete();
00223 }
00224 //============================================= 
00225 void Build_NonLocal_BlockMatrix(const Epetra_Comm & Comm, int NUM_ROWS,int &NumBlocks,
00226                                 int *&Blockstart_, int *&Blockids_,Epetra_CrsMatrix *&MAT){
00227   double values[3]={-1.0,2.0,-1.0};
00228   int indices[3];
00229   int NumProcs=Comm.NumProc();
00230   
00231   // Build me a CrsMatrix
00232   Epetra_Map Map(-1,NUM_ROWS,0,Comm);
00233   MAT=new Epetra_CrsMatrix(Copy,Map,2); 
00234   int MyPID=Comm.MyPID();
00235 
00236   for(int i=0;i<MAT->NumMyRows();i++){
00237     int GID=Map.GID(i);
00238     indices[0]=GID;     
00239     if(i==0) values[1]=2.0;
00240     else values[1]=NumProcs+1;
00241     
00242     MAT->InsertGlobalValues(GID,1,&values[1],&indices[0]);    
00243 
00244     // A little off-diagonal for good measure
00245     if(NumProcs>1 && MyPID==0 && i>0){
00246       indices[1]=GID;
00247       for(int j=1;j<NumProcs;j++){
00248         indices[1]+=NUM_ROWS;//HAQ
00249         MAT->InsertGlobalValues(GID,1,&values[0],&indices[1]);
00250       }
00251     }
00252     else if(NumProcs>1 && MyPID!=0 && i>0){
00253       indices[1]=GID%NUM_ROWS;//HAQ
00254       MAT->InsertGlobalValues(GID,1,&values[0],&indices[1]);
00255     }    
00256   }
00257   MAT->FillComplete();
00258 
00259   // Build me some block structure
00260   // The first row on each proc is a solo block.  All others get blocked to
00261   // PID=0's second block  
00262   if(MyPID==0) NumBlocks=NUM_ROWS;
00263   else NumBlocks=1;  
00264   Blockstart_=new int [NumBlocks+1];
00265   Blockstart_[0]=0;
00266   int curr_idx,curr_block;
00267   
00268   if(MyPID){
00269     // PID > 0
00270     Blockids_=new int[1];
00271     Blockstart_[0]=0; Blockstart_[1]=1;
00272     Blockids_[0]=Map.GID(0);    
00273   }
00274   else{
00275     // PID 0
00276     int nnz=NumProcs*NumBlocks;
00277     Blockids_=new int[nnz+1];
00278     Blockstart_[0]=0;
00279     Blockids_[0]=Map.GID(0);
00280     curr_idx=1; curr_block=1;
00281     for(int j=1;j<NUM_ROWS;j++){
00282       Blockstart_[curr_block]=curr_idx;
00283       for(int i=0;i<NumProcs;i++){
00284         Blockids_[curr_idx]=NUM_ROWS*i+j;//FIX: THIS IS A HACK
00285         curr_idx++;
00286       }
00287       curr_block++;
00288     }
00289     Blockstart_[curr_block]=curr_idx;
00290   }
00291 }
00292 
00293 //============================================= 
00294 void Build_DiagonalStructure(const Epetra_Map &Map,int &NumBlocks,int *&Blockstart_, int *&Blockids_,bool local_ids){
00295   NumBlocks=Map.NumMyElements();
00296   Blockstart_=new int[NumBlocks+1];
00297   Blockids_=new int[NumBlocks];
00298   for(int i=0;i<NumBlocks;i++){
00299     Blockstart_[i]=i;
00300     if(local_ids) Blockids_[i]=i;
00301     else Blockids_[i]=Map.GID(i);
00302   }
00303   Blockstart_[NumBlocks]=NumBlocks;
00304 }
00305 
00306   
00307 
00308 //============================================= 
00309 //============================================= 
00310 //============================================= 
00311 double Test_PTBDP(const Epetra_CrsMatrix& MAT, int NumBlocks,int* Blockstart_,int* Blockids_,bool is_lid){
00312   // Build the block lists
00313   Teuchos::ParameterList List,Sublist;
00314   List.set("number of local blocks",NumBlocks);
00315   List.set("block start index",Blockstart_);
00316   if(is_lid) List.set("block entry lids",Blockids_);
00317   else List.set("block entry gids",Blockids_);
00318 
00319   Sublist.set("apply mode","invert");
00320   //Sublist.set("apply mode","multiply");
00321   List.set("blockdiagmatrix: list",Sublist);
00322   
00323   EpetraExt_PointToBlockDiagPermute Perm(MAT);
00324   Perm.SetParameters(List);
00325 
00326   Perm.Compute();
00327   Epetra_MultiVector X(MAT.RowMap(),1);
00328   Epetra_MultiVector Y(MAT.RowMap(),1);
00329   Epetra_MultiVector Z(MAT.RowMap(),1);
00330   X.SetSeed(24601); X.Random();  
00331 
00332   double norm2;
00333   Perm.ApplyInverse(X,Y);
00334   MAT.Apply(Y,Z);
00335   X.Update(1.0,Z,-1.0);
00336   X.Norm2(&norm2);
00337   return norm2;
00338 }
00339 
00340 
00341 //============================================= 
00342 double Test_PTBDP_C(const Epetra_CrsMatrix& MAT,int BlockSize){
00343   // Build the block lists
00344   Teuchos::ParameterList List,Sublist;
00345   List.set("contiguous block size",BlockSize);
00346 
00347   Sublist.set("apply mode","invert");
00348   //Sublist.set("apply mode","multiply");
00349   List.set("blockdiagmatrix: list",Sublist);
00350   
00351   EpetraExt_PointToBlockDiagPermute Perm(MAT);
00352   Perm.SetParameters(List);
00353 
00354   Perm.Compute();
00355   Epetra_MultiVector X(MAT.RowMap(),1);
00356   Epetra_MultiVector Y(MAT.RowMap(),1);
00357   Epetra_MultiVector Z(MAT.RowMap(),1);
00358   X.SetSeed(24601); X.Random();  
00359 
00360   double norm2;
00361   Perm.ApplyInverse(X,Y);
00362   MAT.Apply(Y,Z);
00363   X.Update(1.0,Z,-1.0);
00364   X.Norm2(&norm2);
00365   return norm2;
00366 }
00367 
00368 //============================================= 
00369 bool TestPointToBlockDiagPermute(const Epetra_Comm & Comm){
00370   const int NUM_ROWS=64;
00371   
00372   bool TestPassed=true;
00373   Epetra_CrsMatrix *MAT;
00374   int NumBlocks, *Blockstart_,*Blockids_;
00375   double norm2;
00376 
00377   // TEST #1 - Local, Contiguous
00378   Build_Local_Contiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00379   norm2=Test_PTBDP(*MAT,NumBlocks,Blockstart_,Blockids_,true);
00380   if(norm2 > 1e-12) TestPassed=false;
00381   if(!Comm.MyPID()) cout<<"P2BDP LCMAT    Error = "<<norm2<<endl;
00382   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00383 
00384   // TEST #2 - Local, Non-Contiguous
00385   Build_Local_NonContiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00386   norm2=Test_PTBDP(*MAT,NumBlocks,Blockstart_,Blockids_,true);
00387   if(norm2 > 1e-12) TestPassed=false;
00388   if(!Comm.MyPID()) cout<<"P2BDP LNCMat   Error = "<<norm2<<endl;
00389   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00390   
00391   // TEST #3 - Non-Local
00392   Build_NonLocal_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00393   norm2=Test_PTBDP(*MAT,NumBlocks,Blockstart_,Blockids_,false);
00394   if(norm2 > 1e-12) TestPassed=false;
00395   if(!Comm.MyPID()) cout<<"P2BDP NLMat    Error = "<<norm2<<endl;
00396   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00397 
00398   // TEST #4 - Local, Contiguous in ContiguousMode
00399   Build_Local_Contiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00400   norm2=Test_PTBDP_C(*MAT,2);
00401   if(norm2 > 1e-12) TestPassed=false;
00402   if(!Comm.MyPID()) cout<<"P2BDP LCMAT-C  Error = "<<norm2<<endl;
00403   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00404 
00405 
00406   return TestPassed;
00407 }
00408 
00409 
00410 //============================================= 
00411 double Test_Cheby(const Epetra_CrsMatrix& MAT, int NumBlocks,int* Blockstart_,int* Blockids_,int maxits,bool is_lid){
00412   double norm2,norm0;
00413   // Build the block lists  
00414   Teuchos::ParameterList ChebyList,List,Sublist;
00415   List.set("number of local blocks",NumBlocks);
00416   List.set("block start index",Blockstart_);
00417   if(is_lid) List.set("block entry lids",Blockids_);
00418   else List.set("block entry gids",Blockids_);
00419 
00420   Sublist.set("apply mode","invert");
00421   List.set("blockdiagmatrix: list",Sublist);
00422 
00423   ChebyList.set("chebyshev: use block mode",true);
00424   ChebyList.set("chebyshev: block list",List);
00425   ChebyList.set("chebyshev: eigenvalue autocompute ratio",30.0);//HAQ
00426   ChebyList.set("chebyshev: degree",maxits);
00427 
00428   // Build a Chebyshev
00429   Ifpack_Chebyshev Cheby(&MAT);
00430   Cheby.SetParameters(ChebyList);
00431   Cheby.Compute();
00432 
00433   Epetra_MultiVector X(MAT.RowMap(),1);
00434   Epetra_MultiVector Y(MAT.RowMap(),1);
00435   Epetra_MultiVector Z(MAT.RowMap(),1);
00436   X.SetSeed(24601); X.Random();
00437   MAT.Apply(X,Y);
00438   Y.Norm2(&norm0);
00439   
00440   Cheby.ApplyInverse(Y,Z);
00441   X.Update(1.0,Z,-1.0);
00442   X.Norm2(&norm2);
00443   return norm2 / norm0;
00444 }
00445 
00446 //============================================= 
00447 double Test_Cheby_C(const Epetra_CrsMatrix& MAT, int BlockSize,int maxits){
00448   double norm2,norm0;
00449   // Build the block lists  
00450   Teuchos::ParameterList ChebyList,List,Sublist;
00451   List.set("contiguous block size",BlockSize);
00452   Sublist.set("apply mode","invert");
00453   List.set("blockdiagmatrix: list",Sublist);
00454 
00455   ChebyList.set("chebyshev: use block mode",true);
00456   ChebyList.set("chebyshev: block list",List);
00457   ChebyList.set("chebyshev: eigenvalue autocompute ratio",30.0);//HAQ
00458   ChebyList.set("chebyshev: degree",maxits);
00459 
00460   // Build a Chebyshev
00461   Ifpack_Chebyshev Cheby(&MAT);
00462   Cheby.SetParameters(ChebyList);
00463   Cheby.Compute();
00464 
00465   Epetra_MultiVector X(MAT.RowMap(),1);
00466   Epetra_MultiVector Y(MAT.RowMap(),1);
00467   Epetra_MultiVector Z(MAT.RowMap(),1);
00468   X.SetSeed(24601); X.Random();
00469   MAT.Apply(X,Y);
00470   Y.Norm2(&norm0);
00471   
00472   Cheby.ApplyInverse(Y,Z);
00473   X.Update(1.0,Z,-1.0);
00474   X.Norm2(&norm2);
00475   return norm2 / norm0;
00476 }
00477 
00478 
00479 //============================================= 
00480 bool TestBlockChebyshev(const Epetra_Comm & Comm){
00481   const int NUM_ROWS=100;
00482   
00483   bool TestPassed=true;
00484   Epetra_CrsMatrix *MAT;
00485   int NumBlocks, *Blockstart_,*Blockids_;
00486   double norm2;
00487 
00488   // Test #1 - Local, Contiguous matrix w/ diagonal precond
00489   Build_Local_Contiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00490   delete [] Blockstart_; delete [] Blockids_;
00491   Build_DiagonalStructure(MAT->RowMap(),NumBlocks,Blockstart_,Blockids_,true);
00492   norm2=Test_Cheby(*MAT,NumBlocks,Blockstart_,Blockids_,100,true);
00493   if(norm2 > 1e-12) TestPassed=false;
00494   if(!Comm.MyPID()) cout<<"Cheby LC-D   nrm-red = "<<norm2<<endl;
00495   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00496 
00497   // Test #2 - Local, Non-Contiguous matrix w/ diagonal precond
00498   Build_Local_NonContiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00499   delete [] Blockstart_; delete [] Blockids_;
00500   Build_DiagonalStructure(MAT->RowMap(),NumBlocks,Blockstart_,Blockids_,true);
00501   norm2=Test_Cheby(*MAT,NumBlocks,Blockstart_,Blockids_,100,true);
00502   if(norm2 > 1e-12) TestPassed=false;
00503   if(!Comm.MyPID()) cout<<"Cheby LNC-D  nrm-red = "<<norm2<<endl;
00504   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00505 
00506   // TEST #3 - Non-Local matrix w/ diagonal precond
00507   Build_NonLocal_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00508   delete [] Blockstart_; delete [] Blockids_;
00509   Build_DiagonalStructure(MAT->RowMap(),NumBlocks,Blockstart_,Blockids_,false);
00510   norm2=Test_Cheby(*MAT,NumBlocks,Blockstart_,Blockids_,100,false);
00511   if(norm2 > 1e-12) TestPassed=false;
00512   if(!Comm.MyPID()) cout<<"Cheby NL-D   nrm-red = "<<norm2<<endl;
00513   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00514 
00515   // Test #4 - Local, Contiguous matrix w/ exact precond
00516   Build_Local_Contiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00517   norm2=Test_Cheby(*MAT,NumBlocks,Blockstart_,Blockids_,1,true);
00518   if(norm2 > 1e-12) TestPassed=false;
00519   if(!Comm.MyPID()) cout<<"Cheby LC-E   nrm-red = "<<norm2<<endl;
00520   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00521 
00522   // Test #5 - Local, Non-Contiguous matrix w/ exact precond
00523   Build_Local_NonContiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00524   norm2=Test_Cheby(*MAT,NumBlocks,Blockstart_,Blockids_,1,true);
00525   if(norm2 > 1e-12) TestPassed=false;  
00526   if(!Comm.MyPID()) cout<<"Cheby LNC-E  nrm-red = "<<norm2<<endl;
00527   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00528   
00529   // TEST #6 - Non-Local matrix w/ exact precond
00530   Build_NonLocal_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00531   norm2=Test_Cheby(*MAT,NumBlocks,Blockstart_,Blockids_,1,false);
00532   if(norm2 > 1e-12) TestPassed=false;
00533   if(!Comm.MyPID()) cout<<"Cheby NL-E   nrm-red = "<<norm2<<endl;
00534   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00535   
00536   // Test #7 - Local, Contiguous matrix w/ diagonal precond (contiguous mode)
00537   Build_Local_Contiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00538   norm2=Test_Cheby_C(*MAT,1,100);
00539   if(norm2 > 1e-12) TestPassed=false;
00540   if(!Comm.MyPID()) cout<<"Cheby LC-Dc  nrm-red = "<<norm2<<endl;
00541   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00542 
00543   // Test #8 - Local, Contiguous matrix w/ exact precond (contiguous mode)
00544   Build_Local_Contiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00545   norm2=Test_Cheby_C(*MAT,2,1);
00546   if(norm2 > 1e-12) TestPassed=false;
00547   if(!Comm.MyPID()) cout<<"Cheby LC-Ec  nrm-red = "<<norm2<<endl;
00548   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00549 
00550   return TestPassed;
00551 }
00552 
00553 //============================================= 
00554 //============================================= 
00555 //============================================= 
00556 int main(int argc, char *argv[])
00557 {
00558 
00559 #ifdef HAVE_MPI
00560   MPI_Init(&argc,&argv);
00561   Epetra_MpiComm Comm( MPI_COMM_WORLD );
00562 #else
00563   Epetra_SerialComm Comm;
00564 #endif
00565   bool TestPassed=true;
00566 
00567   // BlockDiagMatrix test
00568   TestPassed=TestPassed && TestBlockDiagMatrix(Comm);
00569 
00570   // PointToBlockDiagPermute tests
00571   TestPassed=TestPassed && TestPointToBlockDiagPermute(Comm);
00572 
00573   // Block Chebyshev Tests
00574   TestPassed=TestPassed && TestBlockChebyshev(Comm);
00575   
00576   // ============ //
00577   // final output //
00578   // ============ //
00579 
00580   if (!TestPassed) {
00581     if(!Comm.MyPID()) cout << "Test `BlockCheby.exe' FAILED!" << endl;
00582     exit(EXIT_FAILURE);
00583   }
00584   
00585 #ifdef HAVE_MPI
00586   MPI_Finalize(); 
00587 #endif
00588   if(!Comm.MyPID()) cout << "Test `BlockCheby.exe' passed!" << endl;
00589   exit(EXIT_SUCCESS);
00590 }
00591 
00592 #else
00593 
00594 #ifdef HAVE_MPI
00595 #include "Epetra_MpiComm.h"
00596 #else
00597 #include "Epetra_SerialComm.h"
00598 #endif
00599 
00600 int main(int argc, char *argv[])
00601 {
00602 
00603 #ifdef HAVE_MPI
00604   MPI_Init(&argc,&argv);
00605   Epetra_MpiComm Comm( MPI_COMM_WORLD );
00606 #else
00607   Epetra_SerialComm Comm;
00608 #endif
00609 
00610   puts("please configure IFPACK with --eanble-aztecoo --enable-teuchos --enable-epetraext");
00611   puts("to run this test");
00612 
00613 #ifdef HAVE_MPI
00614   MPI_Finalize() ;
00615 #endif
00616   return(EXIT_SUCCESS);
00617 }
00618 
00619 #endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines