Actual source code: petscsnes.h
1: /*
2: User interface for the nonlinear solvers package.
3: */
6: #include petscksp.h
9: /*S
10: SNES - Abstract PETSc object that manages all nonlinear solves
12: Level: beginner
14: Concepts: nonlinear solvers
16: .seealso: SNESCreate(), SNESSetType(), SNESType, TS, KSP, KSP, PC
17: S*/
18: typedef struct _p_SNES* SNES;
20: /*E
21: SNESType - String with the name of a PETSc SNES method or the creation function
22: with an optional dynamic library name, for example
23: http://www.mcs.anl.gov/petsc/lib.a:mysnescreate()
25: Level: beginner
27: .seealso: SNESSetType(), SNES
28: E*/
29: #define SNESType char*
30: #define SNESLS "ls"
31: #define SNESTR "tr"
32: #define SNESPYTHON "python"
33: #define SNESTEST "test"
34: #define SNESPICARD "picard"
35: #define SNESKSPONLY "ksponly"
36: #define SNESVI "vi"
37: #define SNESNGMRES "ngmres"
38: #define SNESSORQN "sorqn"
40: /* Logging support */
70: /*MC
71: SNESRegisterDynamic - Adds a method to the nonlinear solver package.
73: Synopsis:
74: PetscErrorCode SNESRegisterDynamic(const char *name_solver,const char *path,const char *name_create,PetscErrorCode (*routine_create)(SNES))
76: Not collective
78: Input Parameters:
79: + name_solver - name of a new user-defined solver
80: . path - path (either absolute or relative) the library containing this solver
81: . name_create - name of routine to create method context
82: - routine_create - routine to create method context
84: Notes:
85: SNESRegisterDynamic() may be called multiple times to add several user-defined solvers.
87: If dynamic libraries are used, then the fourth input argument (routine_create)
88: is ignored.
90: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR},
91: and others of the form ${any_environmental_variable} occuring in pathname will be
92: replaced with appropriate values.
94: Sample usage:
95: .vb
96: SNESRegisterDynamic("my_solver",/home/username/my_lib/lib/libg/solaris/mylib.a,
97: "MySolverCreate",MySolverCreate);
98: .ve
100: Then, your solver can be chosen with the procedural interface via
101: $ SNESSetType(snes,"my_solver")
102: or at runtime via the option
103: $ -snes_type my_solver
105: Level: advanced
107: Note: If your function is not being put into a shared library then use SNESRegister() instead
109: .keywords: SNES, nonlinear, register
111: .seealso: SNESRegisterAll(), SNESRegisterDestroy()
112: M*/
113: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
114: #define SNESRegisterDynamic(a,b,c,d) SNESRegister(a,b,c,0)
115: #else
116: #define SNESRegisterDynamic(a,b,c,d) SNESRegister(a,b,c,d)
117: #endif
187: /*E
188: SNESConvergedReason - reason a SNES method was said to
189: have converged or diverged
191: Level: beginner
193: The two most common reasons for divergence are
194: $ 1) an incorrectly coded or computed Jacobian or
195: $ 2) failure or lack of convergence in the linear system (in this case we recommend
196: $ testing with -pc_type lu to eliminate the linear solver as the cause of the problem).
198: Diverged Reasons:
199: . SNES_DIVERGED_LOCAL_MIN - this can only occur when using the line-search variant of SNES.
200: The line search wants to minimize Q(alpha) = 1/2 || F(x + alpha s) ||^2_2 this occurs
201: at Q'(alpha) = s^T F'(x+alpha s)^T F(x+alpha s) = 0. If s is the Newton direction - F'(x)^(-1)F(x) then
202: you get Q'(alpha) = -F(x)^T F'(x)^(-1)^T F'(x+alpha s)F(x+alpha s); when alpha = 0
203: Q'(0) = - ||F(x)||^2_2 which is always NEGATIVE if F'(x) is invertible. This means the Newton
204: direction is a descent direction and the line search should succeed if alpha is small enough.
206: If F'(x) is NOT invertible AND F'(x)^T F(x) = 0 then Q'(0) = 0 and the Newton direction
207: is NOT a descent direction so the line search will fail. All one can do at this point
208: is change the initial guess and try again.
210: An alternative explanation: Newton's method can be regarded as replacing the function with
211: its linear approximation and minimizing the 2-norm of that. That is F(x+s) approx F(x) + F'(x)s
212: so we minimize || F(x) + F'(x) s ||^2_2; do this using Least Squares. If F'(x) is invertible then
213: s = - F'(x)^(-1)F(x) otherwise F'(x)^T F'(x) s = -F'(x)^T F(x). If F'(x)^T F(x) is NOT zero then there
214: exists a nontrival (that is F'(x)s != 0) solution to the equation and this direction is
215: s = - [F'(x)^T F'(x)]^(-1) F'(x)^T F(x) so Q'(0) = - F(x)^T F'(x) [F'(x)^T F'(x)]^(-T) F'(x)^T F(x)
216: = - (F'(x)^T F(x)) [F'(x)^T F'(x)]^(-T) (F'(x)^T F(x)). Since we are assuming (F'(x)^T F(x)) != 0
217: and F'(x)^T F'(x) has no negative eigenvalues Q'(0) < 0 so s is a descent direction and the line
218: search should succeed for small enough alpha.
220: Note that this RARELY happens in practice. Far more likely the linear system is not being solved
221: (well enough?) or the Jacobian is wrong.
222:
223: SNES_DIVERGED_MAX_IT means that the solver reached the maximum number of iterations without satisfying any
224: convergence criteria. SNES_CONVERGED_ITS means that SNESSkipConverged() was chosen as the convergence test;
225: thus the usual convergence criteria have not been checked and may or may not be satisfied.
227: Developer Notes: this must match finclude/petscsnes.h
229: The string versions of these are in SNESConvergedReason, if you change any value here you must
230: also adjust that array.
232: Each reason has its own manual page.
234: .seealso: SNESSolve(), SNESGetConvergedReason(), KSPConvergedReason, SNESSetConvergenceTest()
235: E*/
236: typedef enum {/* converged */
237: SNES_CONVERGED_FNORM_ABS = 2, /* ||F|| < atol */
238: SNES_CONVERGED_FNORM_RELATIVE = 3, /* ||F|| < rtol*||F_initial|| */
239: SNES_CONVERGED_PNORM_RELATIVE = 4, /* Newton computed step size small; || delta x || < stol */
240: SNES_CONVERGED_ITS = 5, /* maximum iterations reached */
241: SNES_CONVERGED_TR_DELTA = 7,
242: /* diverged */
243: SNES_DIVERGED_FUNCTION_DOMAIN = -1, /* the new x location passed the function is not in the domain of F */
244: SNES_DIVERGED_FUNCTION_COUNT = -2,
245: SNES_DIVERGED_LINEAR_SOLVE = -3, /* the linear solve failed */
246: SNES_DIVERGED_FNORM_NAN = -4,
247: SNES_DIVERGED_MAX_IT = -5,
248: SNES_DIVERGED_LINE_SEARCH = -6, /* the line search failed */
249: SNES_DIVERGED_LOCAL_MIN = -8, /* || J^T b || is small, implies converged to local minimum of F() */
250: SNES_CONVERGED_ITERATING = 0} SNESConvergedReason;
253: /*MC
254: SNES_CONVERGED_FNORM_ABS - 2-norm(F) <= abstol
256: Level: beginner
258: .seealso: SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()
260: M*/
262: /*MC
263: SNES_CONVERGED_FNORM_RELATIVE - 2-norm(F) <= rtol*2-norm(F(x_0)) where x_0 is the initial guess
265: Level: beginner
267: .seealso: SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()
269: M*/
271: /*MC
272: SNES_CONVERGED_PNORM_RELATIVE - The 2-norm of the last step <= stol * 2-norm(x) where x is the current
273: solution and stol is the 4th argument to SNESSetTolerances()
275: Level: beginner
277: .seealso: SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()
279: M*/
281: /*MC
282: SNES_DIVERGED_FUNCTION_COUNT - The user provided function has been called more times then the final
283: argument to SNESSetTolerances()
285: Level: beginner
287: .seealso: SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()
289: M*/
291: /*MC
292: SNES_DIVERGED_FNORM_NAN - the 2-norm of the current function evaluation is not-a-number (NaN), this
293: is usually caused by a division of 0 by 0.
295: Level: beginner
297: .seealso: SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()
299: M*/
301: /*MC
302: SNES_DIVERGED_MAX_IT - SNESSolve() has reached the maximum number of iterations requested
304: Level: beginner
306: .seealso: SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()
308: M*/
310: /*MC
311: SNES_DIVERGED_LINE_SEARCH - The line search has failed. This only occurs for a SNESType of SNESLS
313: Level: beginner
315: .seealso: SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()
317: M*/
319: /*MC
320: SNES_DIVERGED_LOCAL_MIN - the algorithm seems to have stagnated at a local minimum that is not zero.
321: See the manual page for SNESConvergedReason for more details
323: Level: beginner
325: .seealso: SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()
327: M*/
329: /*MC
330: SNES_CONERGED_ITERATING - this only occurs if SNESGetConvergedReason() is called during the SNESSolve()
332: Level: beginner
334: .seealso: SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()
336: M*/
351: /* --------- Solving systems of nonlinear equations --------------- */
352: typedef PetscErrorCode (*SNESFunction)(SNES,Vec,Vec,void*);
353: typedef PetscErrorCode (*SNESJacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
363: /* --------- Routines specifically for line search methods --------------- */
376: /* Routines for VI solver */
381: #define SNES_VI_INF 1.0e20
382: #define SNES_VI_NINF -1.0e20
386: /* Should this routine be private? */
394: /* Routines for Multiblock solver */
395: PetscErrorCode SNESMultiblockSetFields(SNES, const char [], PetscInt, const PetscInt *);
396: PetscErrorCode SNESMultiblockSetIS(SNES, const char [], IS);
397: PetscErrorCode SNESMultiblockSetBlockSize(SNES, PetscInt);
398: PetscErrorCode SNESMultiblockSetType(SNES, PCCompositeType);
401: #endif