linear_solver.h
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 #ifndef LIBMESH_LINEAR_SOLVER_H
21 #define LIBMESH_LINEAR_SOLVER_H
22 
23 // Local includes
24 #include "libmesh/libmesh_common.h"
30 #include "libmesh/libmesh.h"
32 
33 // C++ includes
34 #include <cstddef>
35 #include <vector>
36 
37 namespace libMesh
38 {
39 
40 // forward declarations
41 template <typename T> class AutoPtr;
42 template <typename T> class SparseMatrix;
43 template <typename T> class NumericVector;
44 template <typename T> class ShellMatrix;
45 template <typename T> class Preconditioner;
46 
47 
56 template <typename T>
57 class LinearSolver : public ReferenceCountedObject<LinearSolver<T> >,
58  public ParallelObject
59 {
60 public:
61 
66  LIBMESH_CAN_DEFAULT_TO_COMMWORLD);
67 
71  virtual ~LinearSolver ();
72 
77  static AutoPtr<LinearSolver<T> > build(const
79  const SolverPackage solver_package = libMesh::default_solver_package());
80 
85  bool initialized () const { return _is_initialized; }
86 
90  virtual void clear () {}
91 
95  virtual void init () = 0;
96 
100  SolverType solver_type () const { return _solver_type; }
101 
105  void set_solver_type (const SolverType st)
106  { _solver_type = st; }
107 
112 
117 
121  void attach_preconditioner(Preconditioner<T> * preconditioner);
122 
123  virtual void reuse_preconditioner(bool );
124 
126 
134  virtual void restrict_solve_to (const std::vector<unsigned int>* const dofs,
135  const SubsetSolveMode subset_solve_mode=SUBSET_ZERO);
136 
143  virtual std::pair<unsigned int, Real> solve (SparseMatrix<T>&, // System Matrix
144  NumericVector<T>&, // Solution vector
145  NumericVector<T>&, // RHS vector
146  const double, // Stopping tolerance
147  const unsigned int) = 0; // N. Iterations
148 
154  virtual std::pair<unsigned int, Real> adjoint_solve (SparseMatrix<T>&, // System Matrix
155  NumericVector<T>&, // Solution vector
156  NumericVector<T>&, // RHS vector
157  const double, // Stopping tolerance
158  const unsigned int); // N. Iterations
159 
165  virtual std::pair<unsigned int, Real> solve (SparseMatrix<T>&, // System Matrix
166  SparseMatrix<T>&, // Preconditioning Matrix
167  NumericVector<T>&, // Solution vector
168  NumericVector<T>&, // RHS vector
169  const double, // Stopping tolerance
170  const unsigned int) = 0; // N. Iterations
171 
178  std::pair<unsigned int, Real> solve (SparseMatrix<T>& matrix,
179  SparseMatrix<T>* precond_matrix,
180  NumericVector<T>&, // Solution vector
181  NumericVector<T>&, // RHS vector
182  const double, // Stopping tolerance
183  const unsigned int); // N. Iterations
184 
185 
186 
190  virtual std::pair<unsigned int, Real> solve (const ShellMatrix<T>& shell_matrix,
191  NumericVector<T>&, // Solution vector
192  NumericVector<T>&, // RHS vector
193  const double, // Stopping tolerance
194  const unsigned int) = 0; // N. Iterations
195 
196 
197 
203  virtual std::pair<unsigned int, Real> solve (const ShellMatrix<T>& shell_matrix,
204  const SparseMatrix<T>& precond_matrix,
205  NumericVector<T>&, // Solution vector
206  NumericVector<T>&, // RHS vector
207  const double, // Stopping tolerance
208  const unsigned int) = 0; // N. Iterations
209 
210 
215  std::pair<unsigned int, Real> solve (const ShellMatrix<T>& matrix,
216  const SparseMatrix<T>* precond_matrix,
217  NumericVector<T>&, // Solution vector
218  NumericVector<T>&, // RHS vector
219  const double, // Stopping tolerance
220  const unsigned int); // N. Iterations
221 
222 
227  virtual void print_converged_reason() = 0;
228 
229 
230 protected:
231 
232 
237 
242 
247 
252 
260 
261 };
262 
263 
264 
265 
266 /*----------------------- inline functions ----------------------------------*/
267 template <typename T>
268 inline
270  ParallelObject (comm_in),
271  _solver_type (GMRES),
272  _preconditioner_type (ILU_PRECOND),
273  _is_initialized (false),
274  _preconditioner (NULL),
275  same_preconditioner (false)
276 {
277 }
278 
279 
280 
281 template <typename T>
282 inline
284 {
285  this->clear ();
286 }
287 
288 template <typename T>
289 inline
291 {
292  return same_preconditioner;
293 }
294 
295 template <typename T>
296 inline
297 std::pair<unsigned int, Real>
299  SparseMatrix<T>* pc_mat,
300  NumericVector<T>& sol,
301  NumericVector<T>& rhs,
302  const double tol,
303  const unsigned int n_iter)
304 {
305  if (pc_mat)
306  return this->solve(mat, *pc_mat, sol, rhs, tol, n_iter);
307  else
308  return this->solve(mat, sol, rhs, tol, n_iter);
309 }
310 
311 
312 template <typename T>
313 inline
314 std::pair<unsigned int, Real>
316  const SparseMatrix<T>* pc_mat,
317  NumericVector<T>& sol,
318  NumericVector<T>& rhs,
319  const double tol,
320  const unsigned int n_iter)
321 {
322  if (pc_mat)
323  return this->solve(mat, *pc_mat, sol, rhs, tol, n_iter);
324  else
325  return this->solve(mat, sol, rhs, tol, n_iter);
326 }
327 
328 } // namespace libMesh
329 
330 
331 #endif // LIBMESH_LINEAR_SOLVER_H

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

Hosted By:
SourceForge.net Logo