Ifpack Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
test/TestAll_LL/cxx_main.cpp
Go to the documentation of this file.
1 /*@HEADER
2 // ***********************************************************************
3 //
4 // Ifpack: Object-Oriented Algebraic Preconditioner Package
5 // Copyright (2002) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ***********************************************************************
40 //@HEADER
41 */
42 
43 #include "Ifpack_ConfigDefs.h"
44 
45 #ifdef HAVE_MPI
46 #include "Epetra_MpiComm.h"
47 #else
48 #include "Epetra_SerialComm.h"
49 #endif
50 #include "Epetra_CrsMatrix.h"
51 #include "Epetra_Vector.h"
52 #include "Epetra_LinearProblem.h"
53 #include "Galeri_Maps.h"
54 #include "Galeri_CrsMatrices.h"
55 #include "Teuchos_ParameterList.hpp"
56 #include "Teuchos_RefCountPtr.hpp"
57 #include "Ifpack_AdditiveSchwarz.h"
58 #include "AztecOO.h"
59 #include "Ifpack_PointRelaxation.h"
60 #include "Ifpack_BlockRelaxation.h"
61 #include "Ifpack_DenseContainer.h"
62 #include "Ifpack_SparseContainer.h"
63 #include "Ifpack_Amesos.h"
64 #include "Ifpack_Utils.h"
65 #include "Ifpack_Chebyshev.h"
66 #include "Ifpack_Polynomial.h"
67 #include "Ifpack_Krylov.h"
68 
69 template <class T>
70 bool Test(const Teuchos::RefCountPtr<Epetra_RowMatrix>& Matrix, Teuchos::ParameterList& List)
71 {
72  using std::cout;
73  using std::endl;
74 
75  int NumVectors = 1;
76  bool UseTranspose = false;
77 
78  Epetra_MultiVector LHS(Matrix->OperatorDomainMap(),NumVectors);
79  Epetra_MultiVector RHS(Matrix->OperatorRangeMap(),NumVectors);
80  Epetra_MultiVector LHSexact(Matrix->OperatorDomainMap(),NumVectors);
81 
82  LHS.PutScalar(0.0);
83  LHSexact.Random();
84  Matrix->Multiply(UseTranspose,LHSexact,RHS);
85 
86  Epetra_LinearProblem Problem(&*Matrix,&LHS,&RHS);
87 
88  Teuchos::RefCountPtr<T> Prec;
89 
90  Prec = Teuchos::rcp( new T(&*Matrix) );
91  assert(Prec != Teuchos::null);
92 
93  IFPACK_CHK_ERR(Prec->SetParameters(List));
94  IFPACK_CHK_ERR(Prec->Initialize());
95  IFPACK_CHK_ERR(Prec->Compute());
96 
97  // create the AztecOO solver
98  AztecOO AztecOOSolver(Problem);
99 
100  // specify solver
101  AztecOOSolver.SetAztecOption(AZ_solver,AZ_gmres);
102  AztecOOSolver.SetAztecOption(AZ_output,32);
103 
104  AztecOOSolver.SetPrecOperator(&*Prec);
105 
106  // solver. The solver should converge in one iteration,
107  // or maximum two (numerical errors)
108  AztecOOSolver.Iterate(1550,1e-8);
109 
110  cout << *Prec;
111 
112  std::vector<double> Norm(NumVectors);
113  LHS.Update(1.0,LHSexact,-1.0);
114  LHS.Norm2(&Norm[0]);
115  for (int i = 0 ; i < NumVectors ; ++i) {
116  cout << "Norm[" << i << "] = " << Norm[i] << endl;
117  if (Norm[i] > 1e-3)
118  return(false);
119  }
120  return(true);
121 
122 }
123 
124 int main(int argc, char *argv[])
125 {
126 
127 #ifdef HAVE_MPI
128  MPI_Init(&argc,&argv);
129  Epetra_MpiComm Comm(MPI_COMM_WORLD);
130 #else
131  Epetra_SerialComm Comm;
132 #endif
133 
134  bool verbose = (Comm.MyPID() == 0);
135 
136  Teuchos::ParameterList GaleriList;
137  const int nx = 30;
138  GaleriList.set("n", nx * nx);
139  GaleriList.set("nx", nx);
140  GaleriList.set("ny", nx);
141  Teuchos::RefCountPtr<Epetra_Map> Map = Teuchos::rcp( Galeri::CreateMap64("Linear", Comm, GaleriList) );
142  Teuchos::RefCountPtr<Epetra_RowMatrix> Matrix = Teuchos::rcp( Galeri::CreateCrsMatrix("Laplace2D", &*Map, GaleriList) );
143 
144  Teuchos::ParameterList List, DefaultList;
145 
146  // test the preconditioner
147  int TestPassed = true;
148 
149  if (!Test<Ifpack_Chebyshev>(Matrix,List))
150  {
151  TestPassed = false;
152  }
153 
154  List.set("polynomial: degree",3);
155  if (!Test<Ifpack_Polynomial>(Matrix,List))
156  {
157  TestPassed = false;
158  }
159 
160  List = DefaultList;
161  List.set("krylov: tolerance", 1e-14);
162  List.set("krylov: iterations", 100);
163  List.set("krylov: preconditioner", 2);
164  List.set("krylov: block size", 9);
165  List.set("krylov: number of sweeps", 2);
166  if (!Test<Ifpack_Krylov>(Matrix,List))
167  {
168  TestPassed = false;
169  }
170 
171  if (!Test< Ifpack_AdditiveSchwarz<Ifpack_Krylov> >(Matrix,List)) {
172  TestPassed = false;
173  }
174 
175  if (!Test<Ifpack_Amesos>(Matrix,List))
176  {
177  TestPassed = false;
178  }
179 
180  // FIXME
181 #if 0
183  TestPassed = false;
184  }
185 #endif
186 
187  // this is ok as long as just one sweep is applied
188  List = DefaultList;
189  List.set("relaxation: type", "Gauss-Seidel");
190  if (!Test<Ifpack_PointRelaxation>(Matrix,List)) {
191  TestPassed = false;
192  }
193 
194  // this is ok as long as just one sweep is applied
195  List = DefaultList;
196  List.set("relaxation: type", "symmetric Gauss-Seidel");
197  List.set("relaxation: sweeps", 5);
198  List.set("partitioner: local parts", 128);
199  List.set("partitioner: type", "linear");
201  TestPassed = false;
202  }
203 
204  // this is ok as long as just one sweep is applied
205  List = DefaultList;
206  List.set("relaxation: type", "symmetric Gauss-Seidel");
207  List.set("partitioner: local parts", 128);
208  List.set("partitioner: type", "linear");
210  TestPassed = false;
211  }
212 
213  // this is ok as long as just one sweep is applied
214  List = DefaultList;
215  List.set("relaxation: type", "symmetric Gauss-Seidel");
216  List.set("partitioner: local parts", 128);
217  List.set("partitioner: type", "linear");
219  TestPassed = false;
220  }
221  if (!TestPassed) {
222  cerr << "Test `TestAll.exe' FAILED!" << endl;
223  exit(EXIT_FAILURE);
224  }
225 
226 #ifdef HAVE_MPI
227  MPI_Finalize();
228 #endif
229  if (verbose)
230  cout << "Test `TestAll.exe' passed!" << endl;
231 
232  return(EXIT_SUCCESS);
233 }
bool Test(const Teuchos::RefCountPtr< Epetra_RowMatrix > &Matrix, Teuchos::ParameterList &List)
Ifpack_BlockRelaxation: a class to define block relaxation preconditioners of Epetra_RowMatrix&#39;s.
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
const int NumVectors
Definition: performance.cpp:71
int MyPID() const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Ifpack_AdditiveSchwarz: a class to define Additive Schwarz preconditioners of Epetra_RowMatrix&#39;s.
Ifpack_SparseContainer: a class for storing and solving linear systems using sparse matrices...
int main(int argc, char *argv[])
#define IFPACK_CHK_ERR(ifpack_err)
#define RHS(a)
Definition: MatGenFD.c:60