nonlinear_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_NONLINEAR_SOLVER_H
21 #define LIBMESH_NONLINEAR_SOLVER_H
22 
23 // Local includes
24 #include "libmesh/libmesh_common.h"
28 #include "libmesh/libmesh.h"
30 
31 // C++ includes
32 #include <cstddef>
33 
34 namespace libMesh
35 {
36 
37 // forward declarations
38 template <typename T> class AutoPtr;
39 template <typename T> class SparseMatrix;
40 template <typename T> class NumericVector;
41 template <typename T> class Preconditioner;
42 
43 
44 
45 
54 template <typename T>
55 class NonlinearSolver : public ReferenceCountedObject<NonlinearSolver<T> >,
56  public ParallelObject
57 {
58 public:
63 
67  explicit
69 
73  virtual ~NonlinearSolver ();
74 
80  const SolverPackage solver_package = libMesh::default_solver_package());
81 
86  bool initialized () const { return _is_initialized; }
87 
91  virtual void clear () {}
92 
96  virtual void init () = 0;
97 
101  virtual std::pair<unsigned int, Real> solve (SparseMatrix<T>&, // System Jacobian Matrix
102  NumericVector<T>&, // Solution vector
103  NumericVector<T>&, // Residual vector
104  const double, // Stopping tolerance
105  const unsigned int) = 0; // N. Iterations
106 
111  virtual void print_converged_reason() { libmesh_not_implemented(); }
112 
116  virtual int get_total_linear_iterations() = 0;
117 
123  virtual unsigned get_current_nonlinear_iteration_number() const = 0;
124 
129  void (* residual) (const NumericVector<Number>& X,
131  sys_type& S);
132 
138 
143  void (* jacobian) (const NumericVector<Number>& X,
145  sys_type& S);
146 
152 
159  void (* matvec) (const NumericVector<Number>& X,
162  sys_type& S);
163 
171 
177  sys_type& S);
182 
189  void (* nullspace) (std::vector<NumericVector<Number>*>& sp, sys_type& S);
190 
198 
204  void (* nearnullspace) (std::vector<NumericVector<Number>*>& sp, sys_type& S);
205 
212 
213 
214  void (* user_presolve)(sys_type& S);
215 
219  const sys_type & system () const { return _system; }
220 
224  sys_type & system () { return _system; }
225 
229  void attach_preconditioner(Preconditioner<T> * preconditioner);
230 
235 
240 
253 
267 
272  unsigned int max_linear_iterations;
273 
279 
284 
289  bool converged;
290 
291 protected:
296 
301 
306 };
307 
308 
309 
310 
311 /*----------------------- inline functions ----------------------------------*/
312 template <typename T>
313 inline
315  ParallelObject (s),
316  residual (NULL),
317  residual_object (NULL),
318  jacobian (NULL),
319  jacobian_object (NULL),
320  matvec (NULL),
321  residual_and_jacobian_object (NULL),
322  bounds (NULL),
323  bounds_object (NULL),
324  nullspace (NULL),
325  nullspace_object (NULL),
326  nearnullspace (NULL),
327  nearnullspace_object (NULL),
328  user_presolve (NULL),
329  max_nonlinear_iterations(0),
330  max_function_evaluations(0),
331  absolute_residual_tolerance(0),
332  relative_residual_tolerance(0),
333  absolute_step_tolerance(0),
334  relative_step_tolerance(0),
335  max_linear_iterations(0),
336  initial_linear_tolerance(0),
337  minimum_linear_tolerance(0),
338  _system(s),
339  _is_initialized (false),
340  _preconditioner (NULL)
341 {
342 }
343 
344 
345 
346 template <typename T>
347 inline
349 {
350  this->clear ();
351 }
352 
353 
354 } // namespace libMesh
355 
356 
357 #endif // LIBMESH_NONLINEAR_SOLVER_H

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

Hosted By:
SourceForge.net Logo