Ifpack Package Browser (Single Doxygen Collection) Development
test/IHSS-SORa/cxx_main.cpp
Go to the documentation of this file.
00001 
00002 // @HEADER
00003 // ***********************************************************************
00004 // 
00005 //            Trilinos: An Object-Oriented Solver Framework
00006 //                 Copyright (2001) Sandia Corporation
00007 // 
00008 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00009 // license for use of this work by or on behalf of the U.S. Government.
00010 // 
00011 // This library is free software; you can redistribute it and/or modify
00012 // it under the terms of the GNU Lesser General Public License as
00013 // published by the Free Software Foundation; either version 2.1 of the
00014 // License, or (at your option) any later version.
00015 //  
00016 // This library is distributed in the hope that it will be useful, but
00017 // WITHOUT ANY WARRANTY; without even the implied warranty of
00018 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00019 // Lesser General Public License for more details.
00020 //  
00021 // You should have received a copy of the GNU Lesser General Public
00022 // License along with this library; if not, write to the Free Software
00023 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00024 // USA
00025 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00026 // 
00027 // ***********************************************************************
00028 // @HEADER
00029 
00030 #include "Ifpack_ConfigDefs.h"
00031 
00032 #include "Epetra_ConfigDefs.h"
00033 #ifdef HAVE_MPI
00034 #include "mpi.h"
00035 #include "Epetra_MpiComm.h"
00036 #else
00037 #include "Epetra_SerialComm.h"
00038 #endif
00039 #include "Epetra_Comm.h"
00040 #include "Epetra_Map.h"
00041 #include "Epetra_Time.h"
00042 #include "Epetra_BlockMap.h"
00043 #include "Epetra_MultiVector.h"
00044 #include "Epetra_Vector.h"
00045 #include "Epetra_Export.h"
00046 #include "AztecOO.h"
00047 #include "Galeri_Maps.h"
00048 #include "Galeri_CrsMatrices.h"
00049 #include "Ifpack_IHSS.h"
00050 #include "Ifpack_SORa.h"
00051 #include "Ifpack.h"
00052 #include "Teuchos_RefCountPtr.hpp"
00053 
00054 // function for fancy output
00055 
00056 string toString(const int& x) {
00057   char s[100];
00058   sprintf(s, "%d", x);
00059   return string(s);
00060 }
00061 
00062 string toString(const double& x) {
00063   char s[100];
00064   sprintf(s, "%g", x);
00065   return string(s);
00066 }
00067 
00068 // main driver
00069 
00070 int main(int argc, char *argv[]) {
00071 
00072 #ifdef HAVE_MPI
00073   MPI_Init(&argc,&argv);
00074   Epetra_MpiComm Comm (MPI_COMM_WORLD);
00075 #else
00076   Epetra_SerialComm Comm;
00077 #endif
00078 
00079   int MyPID = Comm.MyPID();
00080   bool verbose = false; 
00081   if (MyPID==0) verbose = true;
00082 
00083   Teuchos::ParameterList GaleriList;
00084   int nx = 30; 
00085 
00086   GaleriList.set("nx", nx);
00087   //  GaleriList.set("ny", nx * Comm.NumProc());
00088   GaleriList.set("ny", nx);
00089   GaleriList.set("mx", 1);
00090   GaleriList.set("my", Comm.NumProc());
00091   GaleriList.set("alpha", .0);
00092   GaleriList.set("diff", 1.0);
00093   GaleriList.set("conv", 100.0);
00094 
00095   Teuchos::RefCountPtr<Epetra_Map> Map = Teuchos::rcp( Galeri::CreateMap("Cartesian2D", Comm, GaleriList) );
00096   Teuchos::RefCountPtr<Epetra_CrsMatrix> A = Teuchos::rcp( Galeri::CreateCrsMatrix("UniFlow2D", &*Map, GaleriList) );
00097   Teuchos::RefCountPtr<Epetra_MultiVector> LHS = Teuchos::rcp( new Epetra_MultiVector(*Map, 1) );
00098   Teuchos::RefCountPtr<Epetra_MultiVector> RHS = Teuchos::rcp( new Epetra_MultiVector(*Map, 1) );
00099   LHS->PutScalar(0.0); RHS->Random();
00100   Ifpack Factory;  
00101   int Niters = 100;
00102 
00103   // ============================= //
00104   // Construct IHSS preconditioner //
00105   // ============================= //
00106   Teuchos::RefCountPtr<Ifpack_Preconditioner> Prec = Teuchos::rcp( Factory.Create("IHSS", &*A,0) );
00107   Teuchos::ParameterList List;
00108   List.set("ihss: hermetian type","ILU");
00109   List.set("ihss: skew hermetian type","ILU");
00110   List.set("ihss: ratio eigenvalue",100.0);
00111   // Could set sublist values here to better control the ILU, but this isn't needed for this example.
00112   IFPACK_CHK_ERR(Prec->SetParameters(List));
00113   IFPACK_CHK_ERR(Prec->Compute());
00114 
00115   // ============================= //
00116   // Create solver Object          //
00117   // ============================= //
00118 
00119   AztecOO solver;
00120   solver.SetUserMatrix(&*A);
00121   solver.SetLHS(&*LHS);
00122   solver.SetRHS(&*RHS);
00123   solver.SetAztecOption(AZ_solver,AZ_gmres);
00124   solver.SetPrecOperator(&*Prec);
00125   solver.SetAztecOption(AZ_output, 1); 
00126   solver.Iterate(Niters, 1e-8);
00127 
00128   // ============================= //
00129   // Construct SORa preconditioner //
00130   // ============================= //
00131   Teuchos::RefCountPtr<Ifpack_Preconditioner> Prec2 = Teuchos::rcp( Factory.Create("SORa", &*A,0) );
00132   Teuchos::ParameterList List2;
00133   List2.set("sora: sweeps",1);
00134   // Could set sublist values here to better control the ILU, but this isn't needed for this example.
00135   IFPACK_CHK_ERR(Prec2->SetParameters(List2));
00136   IFPACK_CHK_ERR(Prec2->Compute());
00137 
00138   // ============================= //
00139   // Create solver Object          //
00140   // ============================= //
00141   AztecOO solver2;
00142   LHS->PutScalar(0.0);
00143   solver2.SetUserMatrix(&*A);
00144   solver2.SetLHS(&*LHS);
00145   solver2.SetRHS(&*RHS);
00146   solver2.SetAztecOption(AZ_solver,AZ_gmres);
00147   solver2.SetPrecOperator(&*Prec2);
00148   solver2.SetAztecOption(AZ_output, 1); 
00149   solver2.Iterate(Niters, 1e-8);
00150 
00151 #ifdef HAVE_MPI
00152   MPI_Finalize() ;
00153 #endif
00154 
00155   return(EXIT_SUCCESS);
00156 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines