trilinos_aztec_linear_solver.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #include "libmesh/libmesh_common.h"
21 
22 #ifdef LIBMESH_HAVE_AZTECOO
23 
24 
25 // C++ includes
26 
27 // Local Includes
29 #include "libmesh/string_to_enum.h"
33 
34 namespace libMesh
35 {
36 
37 
38 /*----------------------- functions ----------------------------------*/
39 template <typename T>
41 {
42  if (this->initialized())
43  {
44  this->_is_initialized = false;
45 
46  // Mimic PETSc default solver and preconditioner
47  this->_solver_type = GMRES;
48 
49  if (this->n_processors() == 1)
50  this->_preconditioner_type = ILU_PRECOND;
51  else
52  this->_preconditioner_type = BLOCK_JACOBI_PRECOND;
53  }
54 }
55 
56 
57 
58 template <typename T>
60 {
61  // Initialize the data structures if not done so already.
62  if (!this->initialized())
63  {
64  this->_is_initialized = true;
65 
66  _linear_solver = new AztecOO();
67 
68  set_solver_type();
69 
70  switch(this->_preconditioner_type)
71  {
72  case ILU_PRECOND:
73  _linear_solver->SetAztecOption(AZ_precond,AZ_dom_decomp);
74  _linear_solver->SetAztecOption(AZ_subdomain_solve,AZ_ilu);
75  break;
76 
78  _linear_solver->SetAztecOption(AZ_precond,AZ_Jacobi);
79  break;
80 
81  case ICC_PRECOND:
82  _linear_solver->SetAztecOption(AZ_precond,AZ_dom_decomp);
83  _linear_solver->SetAztecOption(AZ_subdomain_solve,AZ_icc);
84  break;
85 
86  case LU_PRECOND:
87  _linear_solver->SetAztecOption(AZ_precond,AZ_dom_decomp);
88  _linear_solver->SetAztecOption(AZ_subdomain_solve,AZ_lu);
89  break;
90 
91  default:
92  _linear_solver->SetAztecOption(AZ_precond,AZ_dom_decomp);
93  _linear_solver->SetAztecOption(AZ_subdomain_solve,AZ_ilu);
94  }
95  }
96 }
97 
98 
99 
100 
101 template <typename T>
102 std::pair<unsigned int, Real>
104  SparseMatrix<T>& precond_in,
105  NumericVector<T>& solution_in,
106  NumericVector<T>& rhs_in,
107  const double tol,
108  const unsigned int m_its)
109 {
110  START_LOG("solve()", "AztecLinearSolver");
111 
112  // Make sure the data passed in are really of Epetra types
113  EpetraMatrix<T>* matrix = libmesh_cast_ptr<EpetraMatrix<T>*>(&matrix_in);
114  EpetraMatrix<T>* precond = libmesh_cast_ptr<EpetraMatrix<T>*>(&precond_in);
115  EpetraVector<T>* solution = libmesh_cast_ptr<EpetraVector<T>*>(&solution_in);
116  EpetraVector<T>* rhs = libmesh_cast_ptr<EpetraVector<T>*>(&rhs_in);
117 
118  this->init();
119 
120  // Close the matrices and vectors in case this wasn't already done.
121  matrix->close ();
122  precond->close ();
123  solution->close ();
124  rhs->close ();
125 
126  _linear_solver->SetAztecOption(AZ_max_iter,m_its);
127  _linear_solver->SetAztecParam(AZ_tol,tol);
128 
129  Epetra_FECrsMatrix * emat = matrix->mat();
130  Epetra_Vector * esol = solution->vec();
131  Epetra_Vector * erhs = rhs->vec();
132 
133  _linear_solver->Iterate(emat, esol, erhs, m_its, tol);
134 
135  STOP_LOG("solve()", "AztecLinearSolver");
136 
137  // return the # of its. and the final residual norm.
138  return std::make_pair(_linear_solver->NumIters(), _linear_solver->TrueResidual());
139 }
140 
141 
142 
143 template <typename T>
144 std::pair<unsigned int, Real>
148  const double,
149  const unsigned int)
150 //AztecLinearSolver<T>::solve (const ShellMatrix<T>& shell_matrix,
151 // NumericVector<T>& solution_in,
152 // NumericVector<T>& rhs_in,
153 // const double tol,
154 // const unsigned int m_its)
155 {
156  libmesh_not_implemented();
157 }
158 
159 
160 
161 template <typename T>
162 std::pair<unsigned int, Real>
164  const SparseMatrix<T>&,
167  const double,
168  const unsigned int)
169 //AztecLinearSolver<T>::solve (const ShellMatrix<T>& shell_matrix,
170 // const SparseMatrix<T>& precond_matrix,
171 // NumericVector<T> &solution_in,
172 // NumericVector<T> &rhs_in,
173 // const double tol,
174 // const unsigned int m_its)
175 {
176  libmesh_not_implemented();
177 }
178 
179 
180 
181 template <typename T>
182 void AztecLinearSolver<T>::get_residual_history(std::vector<double>& /* hist */)
183 {
184  libmesh_not_implemented();
185 }
186 
187 
188 
189 
190 template <typename T>
192 {
193  return _linear_solver->TrueResidual();
194 }
195 
196 
197 
198 template <typename T>
200 {
201  libmesh_not_implemented();
202 }
203 
204 template <typename T>
206 {
207  switch (this->_solver_type)
208  {
209  case CG:
210  _linear_solver->SetAztecOption(AZ_solver, AZ_cg); return;
211 
212  case CGS:
213  _linear_solver->SetAztecOption(AZ_solver, AZ_cgs); return;
214 
215  case TFQMR:
216  _linear_solver->SetAztecOption(AZ_solver, AZ_tfqmr); return;
217 
218  case BICGSTAB:
219  _linear_solver->SetAztecOption(AZ_solver, AZ_bicgstab); return;
220 
221  case GMRES:
222  _linear_solver->SetAztecOption(AZ_solver, AZ_gmres); return;
223 
224  default:
225  libMesh::err << "ERROR: Unsupported AztecOO Solver: "
226  << Utility::enum_to_string(this->_solver_type) << std::endl
227  << "Continuing with AztecOO defaults" << std::endl;
228  }
229 }
230 
231 //------------------------------------------------------------------
232 // Explicit instantiations
233 template class AztecLinearSolver<Number>;
234 
235 } // namespace libMesh
236 
237 
238 
239 #endif // #ifdef LIBMESH_HAVE_AZTECOO

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:07 UTC

Hosted By:
SourceForge.net Logo