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   int MyPID=Comm.MyPID();
00159   
00160   NumBlocks=MAT->NumMyRows()/2;
00161   Blockstart_=new int [NumBlocks+1];
00162   Blockids_=new int [MAT->NumMyRows()];
00163   Blockstart_[0]=0;
00164   int curr_idx=0,curr_block=0;
00165 
00166   for(int i=0;i<MAT->NumMyRows();i++){
00167     // local contiguous blocks of constant size 2
00168     int row_in_block=i%2;
00169     if(row_in_block==0){
00170       Blockstart_[curr_block]=i;      
00171       indices[0]=i; indices[1]=i+1;
00172       Blockids_[i]=i;
00173       MAT->InsertGlobalValues(Map.GID(i),2,&values[1],&indices[0]);
00174     }
00175     else if(row_in_block==1){
00176       indices[0]=i-1; indices[1]=i;      
00177       MAT->InsertGlobalValues(Map.GID(i),2,&values[0],&indices[0]);
00178       Blockids_[i]=i;
00179       curr_block++;
00180     }
00181   }
00182   Blockstart_[NumBlocks]=MAT->NumMyRows();
00183   
00184   MAT->FillComplete();
00185 }
00186 
00187 //============================================= 
00188 void Build_Local_NonContiguous_Size2_BlockMatrix(const Epetra_Comm & Comm, int NUM_ROWS,int &NumBlocks,
00189                                                  int *&Blockstart_, int *&Blockids_,Epetra_CrsMatrix *&MAT){
00190   double values[3]={-1.0,2.0,-1.0};
00191   int indices[3];
00192 
00193   // Build me a CrsMatrix
00194   Epetra_Map Map(-1,NUM_ROWS,0,Comm);
00195   MAT=new Epetra_CrsMatrix(Copy,Map,2);
00196   assert(MAT->NumMyRows()%2==0);
00197 
00198   int MyPID=Comm.MyPID();
00199   
00200   NumBlocks=MAT->NumMyRows()/2;
00201   Blockstart_=new int [NumBlocks+1];
00202   Blockids_=new int [MAT->NumMyRows()];
00203   Blockstart_[0]=0;
00204   int curr_idx=0,curr_block=0;
00205 
00206   for(int i=0;i<MAT->NumMyRows();i++){
00207     int row_in_block=(i%2)?1:0;
00208     if(row_in_block==0){
00209       int idx=i/2;
00210       Blockstart_[curr_block]=i;      
00211       indices[0]=Map.GID(idx); indices[1]=Map.GID(idx+NumBlocks);
00212       Blockids_[i]=idx;
00213       MAT->InsertGlobalValues(Map.GID(idx),2,&values[1],&indices[0]);
00214 
00215     }
00216     else if(row_in_block==1){
00217       int idx=(i-1)/2+NumBlocks;
00218       indices[0]=Map.GID(idx-NumBlocks); indices[1]=Map.GID(idx);      
00219       MAT->InsertGlobalValues(Map.GID(idx),2,&values[0],&indices[0]);
00220       Blockids_[i]=idx;
00221       curr_block++;
00222     }
00223   }
00224   Blockstart_[NumBlocks]=MAT->NumMyRows();
00225   
00226   MAT->FillComplete();
00227 }
00228 //============================================= 
00229 void Build_NonLocal_BlockMatrix(const Epetra_Comm & Comm, int NUM_ROWS,int &NumBlocks,
00230                                 int *&Blockstart_, int *&Blockids_,Epetra_CrsMatrix *&MAT){
00231   double values[3]={-1.0,2.0,-1.0};
00232   int indices[3];
00233   int NumProcs=Comm.NumProc();
00234   
00235   // Build me a CrsMatrix
00236   Epetra_Map Map(-1,NUM_ROWS,0,Comm);
00237   MAT=new Epetra_CrsMatrix(Copy,Map,2); 
00238   int MyPID=Comm.MyPID();
00239 
00240   for(int i=0;i<MAT->NumMyRows();i++){
00241     int GID=Map.GID(i);
00242     indices[0]=GID;     
00243     if(i==0) values[1]=2.0;
00244     else values[1]=NumProcs+1;
00245     
00246     MAT->InsertGlobalValues(GID,1,&values[1],&indices[0]);    
00247 
00248     // A little off-diagonal for good measure
00249     if(NumProcs>1 && MyPID==0 && i>0){
00250       indices[1]=GID;
00251       for(int j=1;j<NumProcs;j++){
00252         indices[1]+=NUM_ROWS;//HAQ
00253         MAT->InsertGlobalValues(GID,1,&values[0],&indices[1]);
00254       }
00255     }
00256     else if(NumProcs>1 && MyPID!=0 && i>0){
00257       indices[1]=GID%NUM_ROWS;//HAQ
00258       MAT->InsertGlobalValues(GID,1,&values[0],&indices[1]);
00259     }    
00260   }
00261   MAT->FillComplete();
00262 
00263   // Build me some block structure
00264   // The first row on each proc is a solo block.  All others get blocked to
00265   // PID=0's second block  
00266   if(MyPID==0) NumBlocks=NUM_ROWS;
00267   else NumBlocks=1;  
00268   Blockstart_=new int [NumBlocks+1];
00269   Blockstart_[0]=0;
00270   int curr_idx,curr_block;
00271   
00272   if(MyPID){
00273     // PID > 0
00274     Blockids_=new int[1];
00275     Blockstart_[0]=0; Blockstart_[1]=1;
00276     Blockids_[0]=Map.GID(0);    
00277   }
00278   else{
00279     // PID 0
00280     int nnz=NumProcs*NumBlocks;
00281     Blockids_=new int[nnz+1];
00282     Blockstart_[0]=0;
00283     Blockids_[0]=Map.GID(0);
00284     curr_idx=1; curr_block=1;
00285     for(int j=1;j<NUM_ROWS;j++){
00286       Blockstart_[curr_block]=curr_idx;
00287       for(int i=0;i<NumProcs;i++){
00288         Blockids_[curr_idx]=NUM_ROWS*i+j;//FIX: THIS IS A HACK
00289         curr_idx++;
00290       }
00291       curr_block++;
00292     }
00293     Blockstart_[curr_block]=curr_idx;
00294   }
00295 }
00296 
00297 //============================================= 
00298 void Build_DiagonalStructure(const Epetra_Map &Map,int &NumBlocks,int *&Blockstart_, int *&Blockids_,bool local_ids){
00299   NumBlocks=Map.NumMyElements();
00300   Blockstart_=new int[NumBlocks+1];
00301   Blockids_=new int[NumBlocks];
00302   for(int i=0;i<NumBlocks;i++){
00303     Blockstart_[i]=i;
00304     if(local_ids) Blockids_[i]=i;
00305     else Blockids_[i]=Map.GID(i);
00306   }
00307   Blockstart_[NumBlocks]=NumBlocks;
00308 }
00309 
00310   
00311 
00312 //============================================= 
00313 //============================================= 
00314 //============================================= 
00315 double Test_PTBDP(const Epetra_CrsMatrix& MAT, int NumBlocks,int* Blockstart_,int* Blockids_,bool is_lid){
00316   // Build the block lists
00317   Teuchos::ParameterList List,Sublist;
00318   List.set("number of local blocks",NumBlocks);
00319   List.set("block start index",Blockstart_);
00320   if(is_lid) List.set("block entry lids",Blockids_);
00321   else List.set("block entry gids",Blockids_);
00322 
00323   Sublist.set("apply mode","invert");
00324   //Sublist.set("apply mode","multiply");
00325   List.set("blockdiagmatrix: list",Sublist);
00326   
00327   EpetraExt_PointToBlockDiagPermute Perm(MAT);
00328   Perm.SetParameters(List);
00329 
00330   Perm.Compute();
00331   Epetra_MultiVector X(MAT.RowMap(),1);
00332   Epetra_MultiVector Y(MAT.RowMap(),1);
00333   Epetra_MultiVector Z(MAT.RowMap(),1);
00334   X.SetSeed(24601); X.Random();  
00335 
00336   double norm2;
00337   Perm.ApplyInverse(X,Y);
00338   MAT.Apply(Y,Z);
00339   X.Update(1.0,Z,-1.0);
00340   X.Norm2(&norm2);
00341   return norm2;
00342 }
00343 
00344 //============================================= 
00345 bool TestPointToBlockDiagPermute(const Epetra_Comm & Comm){
00346   const int NUM_ROWS=64;
00347   
00348   bool TestPassed=true;
00349   Epetra_CrsMatrix *MAT;
00350   int NumBlocks, *Blockstart_,*Blockids_;
00351   double norm2;
00352 
00353   // TEST #1 - Local, Contiguous
00354   Build_Local_Contiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00355   norm2=Test_PTBDP(*MAT,NumBlocks,Blockstart_,Blockids_,true);
00356   if(norm2 > 1e-12) TestPassed=false;
00357   if(!Comm.MyPID()) cout<<"P2BDP LCMAT    Error = "<<norm2<<endl;
00358   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00359 
00360   // TEST #2 - Local, Non-Contiguous
00361   Build_Local_NonContiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00362   norm2=Test_PTBDP(*MAT,NumBlocks,Blockstart_,Blockids_,true);
00363   if(norm2 > 1e-12) TestPassed=false;
00364   if(!Comm.MyPID()) cout<<"P2BDP LNCMat   Error = "<<norm2<<endl;
00365   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00366   
00367   // TEST #3 - Non-Local
00368   Build_NonLocal_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00369   norm2=Test_PTBDP(*MAT,NumBlocks,Blockstart_,Blockids_,false);
00370   if(norm2 > 1e-12) TestPassed=false;
00371   if(!Comm.MyPID()) cout<<"P2BDP NLMat    Error = "<<norm2<<endl;
00372   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00373 
00374   return TestPassed;
00375 }
00376 
00377 
00378 //============================================= 
00379 double Test_Cheby(const Epetra_CrsMatrix& MAT, int NumBlocks,int* Blockstart_,int* Blockids_,int maxits,bool is_lid){
00380   double norm2,norm0;
00381   // Build the block lists  
00382   Teuchos::ParameterList ChebyList,List,Sublist;
00383   List.set("number of local blocks",NumBlocks);
00384   List.set("block start index",Blockstart_);
00385   if(is_lid) List.set("block entry lids",Blockids_);
00386   else List.set("block entry gids",Blockids_);
00387 
00388   Sublist.set("apply mode","invert");
00389   List.set("blockdiagmatrix: list",Sublist);
00390 
00391   ChebyList.set("chebyshev: use block mode",true);
00392   ChebyList.set("chebyshev: block list",List);
00393   ChebyList.set("chebyshev: eigenvalue autocompute ratio",30.0);//HAQ
00394   ChebyList.set("chebyshev: degree",maxits);
00395 
00396   // Build a Chebyshev
00397   Ifpack_Chebyshev Cheby(&MAT);
00398   Cheby.SetParameters(ChebyList);
00399   Cheby.Compute();
00400 
00401   Epetra_MultiVector X(MAT.RowMap(),1);
00402   Epetra_MultiVector Y(MAT.RowMap(),1);
00403   Epetra_MultiVector Z(MAT.RowMap(),1);
00404   X.SetSeed(24601); X.Random();
00405   MAT.Apply(X,Y);
00406   Y.Norm2(&norm0);
00407   
00408   Cheby.ApplyInverse(Y,Z);
00409   X.Update(1.0,Z,-1.0);
00410   X.Norm2(&norm2);
00411   return norm2 / norm0;
00412 }
00413 
00414 
00415 //============================================= 
00416 bool TestBlockChebyshev(const Epetra_Comm & Comm){
00417   const int NUM_ROWS=100;
00418   
00419   bool TestPassed=true;
00420   Epetra_CrsMatrix *MAT;
00421   int NumBlocks, *Blockstart_,*Blockids_;
00422   double norm2;
00423 
00424   // Test #1 - Local, Contiguous matrix w/ diagonal precond
00425   Build_Local_Contiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00426   delete [] Blockstart_; delete [] Blockids_;
00427   Build_DiagonalStructure(MAT->RowMap(),NumBlocks,Blockstart_,Blockids_,true);
00428   norm2=Test_Cheby(*MAT,NumBlocks,Blockstart_,Blockids_,100,true);
00429   if(norm2 > 1e-12) TestPassed=false;
00430   if(!Comm.MyPID()) cout<<"Cheby LC-D   nrm-red = "<<norm2<<endl;
00431   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00432 
00433   // Test #2 - Local, Non-Contiguous matrix w/ diagonal precond
00434   Build_Local_NonContiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00435   delete [] Blockstart_; delete [] Blockids_;
00436   Build_DiagonalStructure(MAT->RowMap(),NumBlocks,Blockstart_,Blockids_,true);
00437   norm2=Test_Cheby(*MAT,NumBlocks,Blockstart_,Blockids_,100,true);
00438   if(norm2 > 1e-12) TestPassed=false;
00439   if(!Comm.MyPID()) cout<<"Cheby LNC-D  nrm-red = "<<norm2<<endl;
00440   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00441 
00442   // TEST #3 - Non-Local matrix w/ diagonal precond
00443   Build_NonLocal_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00444   delete [] Blockstart_; delete [] Blockids_;
00445   Build_DiagonalStructure(MAT->RowMap(),NumBlocks,Blockstart_,Blockids_,false);
00446   norm2=Test_Cheby(*MAT,NumBlocks,Blockstart_,Blockids_,100,false);
00447   if(norm2 > 1e-12) TestPassed=false;
00448   if(!Comm.MyPID()) cout<<"Cheby NL-D   nrm-red = "<<norm2<<endl;
00449   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00450 
00451   // Test #4 - Local, Contiguous matrix w/ exact precond
00452   Build_Local_Contiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00453   norm2=Test_Cheby(*MAT,NumBlocks,Blockstart_,Blockids_,1,true);
00454   if(norm2 > 1e-12) TestPassed=false;
00455   if(!Comm.MyPID()) cout<<"Cheby LC-E   nrm-red = "<<norm2<<endl;
00456   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00457 
00458   // Test #5 - Local, Non-Contiguous matrix w/ exact precond
00459   Build_Local_NonContiguous_Size2_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00460   norm2=Test_Cheby(*MAT,NumBlocks,Blockstart_,Blockids_,1,true);
00461   if(norm2 > 1e-12) TestPassed=false;  
00462   if(!Comm.MyPID()) cout<<"Cheby LNC-E  nrm-red = "<<norm2<<endl;
00463   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00464   
00465   // TEST #6 - Non-Local matrix w/ exact precond
00466   Build_NonLocal_BlockMatrix(Comm,NUM_ROWS,NumBlocks,Blockstart_,Blockids_,MAT);
00467   norm2=Test_Cheby(*MAT,NumBlocks,Blockstart_,Blockids_,1,false);
00468   if(norm2 > 1e-12) TestPassed=false;
00469   if(!Comm.MyPID()) cout<<"Cheby NL-E   nrm-red = "<<norm2<<endl;
00470   delete MAT; delete [] Blockstart_; delete [] Blockids_;
00471   
00472   return TestPassed;
00473 }
00474 
00475 //============================================= 
00476 //============================================= 
00477 //============================================= 
00478 int main(int argc, char *argv[])
00479 {
00480 
00481 #ifdef HAVE_MPI
00482   MPI_Init(&argc,&argv);
00483   Epetra_MpiComm Comm( MPI_COMM_WORLD );
00484 #else
00485   Epetra_SerialComm Comm;
00486 #endif
00487   bool TestPassed=true;
00488 
00489   // BlockDiagMatrix test
00490   TestPassed=TestPassed && TestBlockDiagMatrix(Comm);
00491 
00492   // PointToBlockDiagPermute tests
00493   TestPassed=TestPassed && TestPointToBlockDiagPermute(Comm);
00494 
00495   // Block Chebyshev Tests
00496   TestPassed=TestPassed && TestBlockChebyshev(Comm);
00497   
00498   // ============ //
00499   // final output //
00500   // ============ //
00501 
00502   if (!TestPassed) {
00503     if(!Comm.MyPID()) cout << "Test `BlockCheby.exe' FAILED!" << endl;
00504     exit(EXIT_FAILURE);
00505   }
00506   
00507 #ifdef HAVE_MPI
00508   MPI_Finalize(); 
00509 #endif
00510   if(!Comm.MyPID()) cout << "Test `BlockCheby.exe' passed!" << endl;
00511   exit(EXIT_SUCCESS);
00512 }
00513 
00514 #else
00515 
00516 #ifdef HAVE_MPI
00517 #include "Epetra_MpiComm.h"
00518 #else
00519 #include "Epetra_SerialComm.h"
00520 #endif
00521 
00522 int main(int argc, char *argv[])
00523 {
00524 
00525 #ifdef HAVE_MPI
00526   MPI_Init(&argc,&argv);
00527   Epetra_MpiComm Comm( MPI_COMM_WORLD );
00528 #else
00529   Epetra_SerialComm Comm;
00530 #endif
00531 
00532   puts("please configure IFPACK with --eanble-aztecoo --enable-teuchos --enable-epetraext");
00533   puts("to run this test");
00534 
00535 #ifdef HAVE_MPI
00536   MPI_Finalize() ;
00537 #endif
00538   return(EXIT_SUCCESS);
00539 }
00540 
00541 #endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
Generated on Wed Apr 13 10:05:33 2011 for Ifpack Package Browser (Single Doxygen Collection) by  doxygen 1.6.3