test/TestAll/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 
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 "Galeri_Maps.h"
00040 #include "Galeri_CrsMatrices.h"
00041 #include "Teuchos_ParameterList.hpp"
00042 #include "Ifpack_AdditiveSchwarz.h"
00043 #include "AztecOO.h"
00044 #include "Ifpack_PointRelaxation.h"
00045 #include "Ifpack_BlockRelaxation.h"
00046 #include "Ifpack_DenseContainer.h"
00047 #include "Ifpack_SparseContainer.h"
00048 #include "Ifpack_Amesos.h"
00049 #include "Ifpack_Utils.h"
00050 #include "Ifpack_Chebyshev.h"
00051 
00052 template <class T>
00053 bool Test(Epetra_RowMatrix* Matrix, Teuchos::ParameterList& List)
00054 {
00055 
00056   int NumVectors = 1;
00057   bool UseTranspose = false;
00058 
00059   Epetra_MultiVector LHS(Matrix->OperatorDomainMap(),NumVectors);
00060   Epetra_MultiVector RHS(Matrix->OperatorRangeMap(),NumVectors);
00061   Epetra_MultiVector LHSexact(Matrix->OperatorDomainMap(),NumVectors);
00062 
00063   LHS.PutScalar(0.0);
00064   LHSexact.Random();
00065   Matrix->Multiply(UseTranspose,LHSexact,RHS);
00066 
00067   Epetra_LinearProblem Problem(Matrix,&LHS,&RHS);
00068 
00069   T* Prec;
00070   
00071   Prec = new T(Matrix);
00072   assert(Prec != 0);
00073 
00074   IFPACK_CHK_ERR(Prec->SetParameters(List));
00075   IFPACK_CHK_ERR(Prec->Initialize());
00076   IFPACK_CHK_ERR(Prec->Compute());
00077 
00078   // create the AztecOO solver
00079   AztecOO AztecOOSolver(Problem);
00080 
00081   // specify solver
00082   AztecOOSolver.SetAztecOption(AZ_solver,AZ_gmres);
00083   AztecOOSolver.SetAztecOption(AZ_output,32);
00084 
00085   AztecOOSolver.SetPrecOperator(Prec);
00086 
00087   // solver. The solver should converge in one iteration,
00088   // or maximum two (numerical errors)
00089   AztecOOSolver.Iterate(1550,1e-8);
00090 
00091   cout << *Prec;
00092   delete Prec;
00093   
00094   vector<double> Norm(NumVectors);
00095   LHS.Update(1.0,LHSexact,-1.0);
00096   LHS.Norm2(&Norm[0]);
00097   for (int i = 0 ; i < NumVectors ; ++i) {
00098     cout << "Norm[" << i << "] = " << Norm[i] << endl;
00099     if (Norm[i] > 1e-3)
00100       return(false);
00101   }
00102   return(true);
00103 
00104 }
00105 
00106 int main(int argc, char *argv[])
00107 {
00108 
00109 #ifdef HAVE_MPI
00110   MPI_Init(&argc,&argv);
00111   Epetra_MpiComm Comm(MPI_COMM_WORLD);
00112 #else
00113   Epetra_SerialComm Comm;
00114 #endif
00115 
00116   bool verbose = (Comm.MyPID() == 0);
00117 
00118   Teuchos::ParameterList GaleriList;
00119   const int nx = 30;
00120   GaleriList.set("n", nx * nx);
00121   GaleriList.set("nx", nx);
00122   GaleriList.set("ny", nx);
00123   Epetra_Map* Map = Galeri::CreateMap("Linear", Comm, GaleriList);
00124   Epetra_RowMatrix* Matrix = Galeri::CreateCrsMatrix("Laplace2D", Map, GaleriList);
00125 
00126   Teuchos::ParameterList List, DefaultList;
00127 
00128   // test the preconditioner
00129   int TestPassed = true;
00130 
00131   if (!Test<Ifpack_Chebyshev>(Matrix,List)) 
00132   {
00133     TestPassed = false;
00134   }
00135 
00136   if (!Test<Ifpack_Amesos>(Matrix,List)) 
00137   {
00138     TestPassed = false;
00139   }
00140 
00141   // FIXME
00142 #if 0
00143   if (!Test<Ifpack_AdditiveSchwarz<Ifpack_BlockRelaxation<Ifpack_DenseContainer> > >(Matrix,List)) {
00144     TestPassed = false;
00145   }
00146 #endif
00147 
00148   // this is ok as long as just one sweep is applied
00149   List = DefaultList;
00150   List.set("relaxation: type", "Gauss-Seidel");
00151   if (!Test<Ifpack_PointRelaxation>(Matrix,List)) {
00152     TestPassed = false;
00153   }
00154 
00155   // this is ok as long as just one sweep is applied
00156   List = DefaultList;
00157   List.set("relaxation: type", "symmetric Gauss-Seidel");
00158   List.set("relaxation: sweeps", 5);
00159   List.set("partitioner: local parts", 128);
00160   List.set("partitioner: type", "linear");
00161   if (!Test<Ifpack_BlockRelaxation<Ifpack_DenseContainer> >(Matrix,List)) {
00162     TestPassed = false;
00163   }
00164  
00165   // this is ok as long as just one sweep is applied
00166   List = DefaultList;
00167   List.set("relaxation: type", "symmetric Gauss-Seidel");
00168   List.set("partitioner: local parts", 128);
00169   List.set("partitioner: type", "linear");
00170   if (!Test<Ifpack_BlockRelaxation<Ifpack_SparseContainer<Ifpack_Amesos> > >(Matrix,List)) {
00171     TestPassed = false;
00172   }
00173 
00174   // this is ok as long as just one sweep is applied
00175   List = DefaultList;
00176   List.set("relaxation: type", "symmetric Gauss-Seidel");
00177   List.set("partitioner: local parts", 128);
00178   List.set("partitioner: type", "linear");
00179   if (!Test<Ifpack_AdditiveSchwarz<Ifpack_BlockRelaxation<Ifpack_SparseContainer<Ifpack_Amesos> > > >(Matrix,List)) {
00180     TestPassed = false;
00181   }
00182   if (!TestPassed) {
00183     cerr << "Test `TestAll.exe' FAILED!" << endl;
00184     exit(EXIT_FAILURE);
00185   }
00186 
00187   delete Matrix;
00188   delete Map;
00189 
00190 #ifdef HAVE_MPI
00191   MPI_Finalize(); 
00192 #endif
00193   if (verbose)
00194     cout << "Test `TestAll.exe' passed!" << endl;
00195 
00196   return(EXIT_SUCCESS);
00197 }

Generated on Thu Sep 18 12:37:21 2008 for Ifpack Package Browser (Single Doxygen Collection) by doxygen 1.3.9.1