nonlinear_solver.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2018 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 #include "libmesh/auto_ptr.h" // deprecated
31 
32 // C++ includes
33 #include <cstddef>
34 #include <memory>
35 
36 namespace libMesh
37 {
38 
39 // forward declarations
40 template <typename T> class SparseMatrix;
41 template <typename T> class NumericVector;
42 template <typename T> class Preconditioner;
43 class SolverConfiguration;
44 
52 template <typename T>
53 class NonlinearSolver : public ReferenceCountedObject<NonlinearSolver<T>>,
54  public ParallelObject
55 {
56 public:
61 
65  explicit
66  NonlinearSolver (sys_type & s);
67 
71  virtual ~NonlinearSolver ();
72 
77  static std::unique_ptr<NonlinearSolver<T>> build(sys_type & s,
78  const SolverPackage solver_package = libMesh::default_solver_package());
79 
84  bool initialized () const { return _is_initialized; }
85 
89  virtual void clear () {}
90 
95  virtual void init (const char * name = libmesh_nullptr) = 0;
96 
100  virtual std::pair<unsigned int, Real> solve (SparseMatrix<T> &, // System Jacobian Matrix
101  NumericVector<T> &, // Solution vector
102  NumericVector<T> &, // Residual vector
103  const double, // Stopping tolerance
104  const unsigned int) = 0; // N. Iterations
105 
110  virtual void print_converged_reason() { libmesh_not_implemented(); }
111 
115  virtual int get_total_linear_iterations() = 0;
116 
124  virtual unsigned get_current_nonlinear_iteration_number() const = 0;
125 
130  void (* residual) (const NumericVector<Number> & X,
132  sys_type & S);
133 
139 
145 
152 
157  void (* jacobian) (const NumericVector<Number> & X,
159  sys_type & S);
160 
166 
174  void (* matvec) (const NumericVector<Number> & X,
177  sys_type & S);
178 
187 
193  sys_type & S);
198 
205  void (* nullspace) (std::vector<NumericVector<Number> *> & sp, sys_type & S);
206 
214 
220  void (* transpose_nullspace) (std::vector<NumericVector<Number> *> & sp, sys_type & S);
221 
228 
234  void (* nearnullspace) (std::vector<NumericVector<Number> *> & sp, sys_type & S);
235 
242 
247  void (* user_presolve)(sys_type & S);
248 
255  void (* postcheck) (const NumericVector<Number> & old_soln,
256  NumericVector<Number> & search_direction,
257  NumericVector<Number> & new_soln,
258  bool & changed_search_direction,
259  bool & changed_new_soln,
260  sys_type & S);
261 
268 
272  const sys_type & system () const { return _system; }
273 
277  sys_type & system () { return _system; }
278 
282  void attach_preconditioner(Preconditioner<T> * preconditioner);
283 
288 
293 
306 
320 
325  unsigned int max_linear_iterations;
326 
332 
337 
342  bool converged;
343 
347  void set_solver_configuration(SolverConfiguration & solver_configuration);
348 
349 protected:
353  sys_type & _system;
354 
359 
364 
370 };
371 
372 
373 
374 
375 /*----------------------- inline functions ----------------------------------*/
376 template <typename T>
377 inline
379  ParallelObject (s),
408  converged(false),
409  _system(s),
410  _is_initialized (false),
413 {
414 }
415 
416 
417 
418 template <typename T>
419 inline
421 {
422  this->clear ();
423 }
424 
425 
426 } // namespace libMesh
427 
428 
429 #endif // LIBMESH_NONLINEAR_SOLVER_H
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
void(* residual)(const NumericVector< Number > &X, NumericVector< Number > &R, sys_type &S)
NonlinearImplicitSystem::ComputeResidualandJacobian * residual_and_jacobian_object
void(* transpose_nullspace)(std::vector< NumericVector< Number > * > &sp, sys_type &S)
unsigned int max_function_evaluations
Preconditioner< T > * _preconditioner
virtual void init(const char *name=libmesh_nullptr)=0
void(* user_presolve)(sys_type &S)
void(* nullspace)(std::vector< NumericVector< Number > * > &sp, sys_type &S)
void attach_preconditioner(Preconditioner< T > *preconditioner)
const class libmesh_nullptr_t libmesh_nullptr
virtual unsigned get_current_nonlinear_iteration_number() const =0
unsigned int max_linear_iterations
SolverConfiguration * _solver_configuration
virtual std::pair< unsigned int, Real > solve(SparseMatrix< T > &, NumericVector< T > &, NumericVector< T > &, const double, const unsigned int)=0
NonlinearImplicitSystem::ComputeResidual * residual_object
void(* jacobian)(const NumericVector< Number > &X, SparseMatrix< Number > &J, sys_type &S)
SolverPackage default_solver_package()
Definition: libmesh.C:994
NonlinearImplicitSystem::ComputeJacobian * jacobian_object
void set_solver_configuration(SolverConfiguration &solver_configuration)
virtual int get_total_linear_iterations()=0
NonlinearImplicitSystem::ComputeBounds * bounds_object
Used for solving nonlinear implicit systems of equations.
const sys_type & system() const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void(* nearnullspace)(std::vector< NumericVector< Number > * > &sp, sys_type &S)
NonlinearImplicitSystem::ComputePostCheck * postcheck_object
virtual void print_converged_reason()
static std::unique_ptr< NonlinearSolver< T > > build(sys_type &s, const SolverPackage solver_package=libMesh::default_solver_package())
void(* bounds)(NumericVector< Number > &XL, NumericVector< Number > &XU, sys_type &S)
NonlinearImplicitSystem::ComputeResidual * mffd_residual_object
void(* postcheck)(const NumericVector< Number > &old_soln, NumericVector< Number > &search_direction, NumericVector< Number > &new_soln, bool &changed_search_direction, bool &changed_new_soln, sys_type &S)
NonlinearImplicitSystem::ComputeVectorSubspace * nearnullspace_object
NonlinearImplicitSystem sys_type
unsigned int max_nonlinear_iterations
NonlinearImplicitSystem::ComputeVectorSubspace * transpose_nullspace_object
void(* matvec)(const NumericVector< Number > &X, NumericVector< Number > *R, SparseMatrix< Number > *J, sys_type &S)
NonlinearImplicitSystem::ComputeResidual * fd_residual_object
NonlinearImplicitSystem::ComputeVectorSubspace * nullspace_object