petsc_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_PETSC_NONLINEAR_SOLVER_H
21 #define LIBMESH_PETSC_NONLINEAR_SOLVER_H
22 
23 #include "libmesh/libmesh_config.h"
24 
25 // Petsc include files.
26 #ifdef LIBMESH_HAVE_PETSC
27 
28 // Local includes
30 #include "libmesh/petsc_macro.h"
31 
32 // PETSc includes
33 # include <petscsnes.h>
34 
35 namespace libMesh
36 {
37 class ResidualContext;
38 
39 // Allow users access to these functions in case they want to reuse them. Users shouldn't
40 // need access to these most of the time as they are used internally by this object.
41 extern "C"
42 {
43  PetscErrorCode __libmesh_petsc_snes_monitor (SNES, PetscInt its, PetscReal fnorm, void *);
44  PetscErrorCode __libmesh_petsc_snes_residual (SNES, Vec x, Vec r, void * ctx);
45  PetscErrorCode __libmesh_petsc_snes_fd_residual (SNES, Vec x, Vec r, void * ctx);
46  PetscErrorCode __libmesh_petsc_snes_mffd_interface (void * ctx, Vec x, Vec r);
47  PetscErrorCode __libmesh_petsc_snes_mffd_residual (SNES, Vec x, Vec r, void * ctx);
48 #if PETSC_RELEASE_LESS_THAN(3,5,0)
49  PetscErrorCode __libmesh_petsc_snes_jacobian (SNES, Vec x, Mat * jac, Mat * pc, MatStructure * msflag, void * ctx);
50 #else
51  PetscErrorCode __libmesh_petsc_snes_jacobian (SNES, Vec x, Mat jac, Mat pc, void * ctx);
52 #endif
53 
54  PetscErrorCode __libmesh_petsc_snes_postcheck(
55 #if PETSC_VERSION_LESS_THAN(3,3,0)
56  SNES, Vec x, Vec y, Vec w, void * context, PetscBool * changed_y, PetscBool * changed_w
57 #else
58  SNESLineSearch, Vec x, Vec y, Vec w, PetscBool * changed_y, PetscBool * changed_w, void * context
59 #endif
60  );
61 }
62 
71 template <typename T>
73 {
74 public:
79 
83  explicit
84  PetscNonlinearSolver (sys_type & system);
85 
90 
94  virtual void clear () libmesh_override;
95 
100  virtual void init (const char * name = libmesh_nullptr) libmesh_override;
101 
105  SNES snes() { this->init(); return _snes; }
106 
111  virtual std::pair<unsigned int, Real>
112  solve (SparseMatrix<T> &, // System Jacobian Matrix
113  NumericVector<T> &, // Solution vector
114  NumericVector<T> &, // Residual vector
115  const double, // Stopping tolerance
116  const unsigned int) libmesh_override; // N. Iterations
117 
122  virtual void print_converged_reason() libmesh_override;
123 
131  SNESConvergedReason get_converged_reason();
132 
136  virtual int get_total_linear_iterations() libmesh_override;
137 
143  virtual unsigned get_current_nonlinear_iteration_number() const libmesh_override
145 
149  void set_residual_zero_out(bool state) { _zero_out_residual = state; }
150 
154  void set_jacobian_zero_out(bool state) { _zero_out_jacobian = state; }
155 
159  void use_default_monitor(bool state) { _default_monitor = state; }
160 
161 protected:
162 
166  SNES _snes;
167 
178  SNESConvergedReason _reason;
179 
184 
189 
194 
199 
204 
205 #if !PETSC_VERSION_LESS_THAN(3,3,0)
207  void (*)(std::vector<NumericVector<Number> *> &, sys_type &),
208  MatNullSpace *);
209 #endif
210 private:
211  friend ResidualContext libmesh_petsc_snes_residual_helper (SNES snes, Vec x, void * ctx);
212  friend PetscErrorCode __libmesh_petsc_snes_residual (SNES snes, Vec x, Vec r, void * ctx);
213  friend PetscErrorCode __libmesh_petsc_snes_fd_residual (SNES snes, Vec x, Vec r, void * ctx);
214  friend PetscErrorCode __libmesh_petsc_snes_mffd_interface (void * ctx, Vec x, Vec r);
215  friend PetscErrorCode __libmesh_petsc_snes_mffd_residual (SNES snes, Vec x, Vec r, void * ctx);
216 #if PETSC_RELEASE_LESS_THAN(3,5,0)
217  friend PetscErrorCode __libmesh_petsc_snes_jacobian (SNES snes, Vec x, Mat * jac, Mat * pc, MatStructure * msflag, void * ctx);
218 #else
219  friend PetscErrorCode __libmesh_petsc_snes_jacobian (SNES snes, Vec x, Mat jac, Mat pc, void * ctx);
220 #endif
221 };
222 
223 
224 
225 } // namespace libMesh
226 
227 
228 #endif // #ifdef LIBMESH_HAVE_PETSC
229 #endif // LIBMESH_PETSC_NONLINEAR_SOLVER_H
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
virtual void clear() libmesh_override
PetscErrorCode __libmesh_petsc_snes_mffd_residual(SNES snes, Vec x, Vec r, void *ctx)
friend PetscErrorCode __libmesh_petsc_snes_mffd_interface(void *ctx, Vec x, Vec r)
PetscErrorCode __libmesh_petsc_snes_residual(SNES snes, Vec x, Vec r, void *ctx)
PetscErrorCode __libmesh_petsc_snes_monitor(SNES, PetscInt its, PetscReal fnorm, void *)
const class libmesh_nullptr_t libmesh_nullptr
friend PetscErrorCode __libmesh_petsc_snes_mffd_residual(SNES snes, Vec x, Vec r, void *ctx)
PetscErrorCode __libmesh_petsc_snes_postcheck(#if PETSC_VERSION_LESS_THAN(3, 3, 0) SNES, Vec x, Vec y, Vec w, void *context, PetscBool *changed_y, PetscBool *changed_w#else SNESLineSearch, Vec x, Vec y, Vec w, PetscBool *changed_y, PetscBool *changed_w, void *context#endif)
PetscErrorCode __libmesh_petsc_snes_mffd_interface(void *ctx, Vec x, Vec r)
friend ResidualContext libmesh_petsc_snes_residual_helper(SNES snes, Vec x, void *ctx)
PetscErrorCode __libmesh_petsc_snes_jacobian(#if PETSC_RELEASE_LESS_THAN(3, 5, 0) SNES snes, Vec x, Mat *jac, Mat *pc, MatStructure *msflag, void *ctx#else SNES snes, Vec x, Mat jac, Mat pc, void *ctx#endif)
Used for solving nonlinear implicit systems of equations.
virtual void init(const char *name=libmesh_nullptr) libmesh_override
friend PetscErrorCode __libmesh_petsc_snes_residual(SNES snes, Vec x, Vec r, void *ctx)
void build_mat_null_space(NonlinearImplicitSystem::ComputeVectorSubspace *computeSubspaceObject, void(*)(std::vector< NumericVector< Number > * > &, sys_type &), MatNullSpace *)
SNESConvergedReason get_converged_reason()
friend PetscErrorCode __libmesh_petsc_snes_jacobian(SNES snes, Vec x, Mat *jac, Mat *pc, MatStructure *msflag, void *ctx)
PetscTruth PetscBool
Definition: petsc_macro.h:64
const sys_type & system() const
virtual void print_converged_reason() libmesh_override
virtual std::pair< unsigned int, Real > solve(SparseMatrix< T > &, NumericVector< T > &, NumericVector< T > &, const double, const unsigned int) libmesh_override
friend PetscErrorCode __libmesh_petsc_snes_fd_residual(SNES snes, Vec x, Vec r, void *ctx)
virtual unsigned get_current_nonlinear_iteration_number() const libmesh_override
virtual int get_total_linear_iterations() libmesh_override
NonlinearImplicitSystem sys_type
PetscErrorCode __libmesh_petsc_snes_fd_residual(SNES snes, Vec x, Vec r, void *ctx)