Actual source code: snes.c
1: #include <petsc/private/snesimpl.h>
2: #include <petsc/private/linesearchimpl.h>
3: #include <petscdmshell.h>
4: #include <petscdraw.h>
5: #include <petscds.h>
6: #include <petscdmadaptor.h>
7: #include <petscconvest.h>
8: #include <petscpreload.h>
10: PetscBool SNESRegisterAllCalled = PETSC_FALSE;
11: PetscFunctionList SNESList = NULL;
13: /* Logging support */
14: PetscClassId SNES_CLASSID, DMSNES_CLASSID;
15: PetscLogEvent SNES_Solve, SNES_SetUp, SNES_FunctionEval, SNES_JacobianEval, SNES_NGSEval, SNES_NGSFuncEval, SNES_NewtonALEval, SNES_NPCSolve, SNES_ObjectiveEval;
17: /*@
18: SNESSetErrorIfNotConverged - Causes `SNESSolve()` to generate an error immediately if the solver has not converged.
20: Logically Collective
22: Input Parameters:
23: + snes - iterative context obtained from `SNESCreate()`
24: - flg - `PETSC_TRUE` indicates you want the error generated
26: Options Database Key:
27: . -snes_error_if_not_converged (true|false) - cause an immediate error condition and stop the program if the solver does not converge
29: Level: intermediate
31: Note:
32: Normally PETSc continues if a solver fails to converge, you can call `SNESGetConvergedReason()` after a `SNESSolve()`
33: to determine if it has converged. Otherwise the solution may be inaccurate or wrong
35: .seealso: [](ch_snes), `SNES`, `SNESGetErrorIfNotConverged()`, `KSPGetErrorIfNotConverged()`, `KSPSetErrorIfNotConverged()`
36: @*/
37: PetscErrorCode SNESSetErrorIfNotConverged(SNES snes, PetscBool flg)
38: {
39: PetscFunctionBegin;
42: snes->errorifnotconverged = flg;
43: PetscFunctionReturn(PETSC_SUCCESS);
44: }
46: /*@
47: SNESGetErrorIfNotConverged - Indicates if `SNESSolve()` will generate an error if the solver does not converge?
49: Not Collective
51: Input Parameter:
52: . snes - iterative context obtained from `SNESCreate()`
54: Output Parameter:
55: . flag - `PETSC_TRUE` if it will generate an error, else `PETSC_FALSE`
57: Level: intermediate
59: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESSetErrorIfNotConverged()`, `KSPGetErrorIfNotConverged()`, `KSPSetErrorIfNotConverged()`
60: @*/
61: PetscErrorCode SNESGetErrorIfNotConverged(SNES snes, PetscBool *flag)
62: {
63: PetscFunctionBegin;
65: PetscAssertPointer(flag, 2);
66: *flag = snes->errorifnotconverged;
67: PetscFunctionReturn(PETSC_SUCCESS);
68: }
70: /*@
71: SNESSetAlwaysComputesFinalResidual - tells the `SNES` to always compute the residual (nonlinear function value) at the final solution
73: Logically Collective
75: Input Parameters:
76: + snes - the shell `SNES`
77: - flg - `PETSC_TRUE` to always compute the residual
79: Level: advanced
81: Note:
82: Some solvers (such as smoothers in a `SNESFAS`) do not need the residual computed at the final solution so skip computing it
83: to save time.
85: .seealso: [](ch_snes), `SNES`, `SNESFAS`, `SNESSolve()`, `SNESGetAlwaysComputesFinalResidual()`
86: @*/
87: PetscErrorCode SNESSetAlwaysComputesFinalResidual(SNES snes, PetscBool flg)
88: {
89: PetscFunctionBegin;
91: snes->alwayscomputesfinalresidual = flg;
92: PetscFunctionReturn(PETSC_SUCCESS);
93: }
95: /*@
96: SNESGetAlwaysComputesFinalResidual - checks if the `SNES` always computes the residual at the final solution
98: Logically Collective
100: Input Parameter:
101: . snes - the `SNES` context
103: Output Parameter:
104: . flg - `PETSC_TRUE` if the residual is computed
106: Level: advanced
108: .seealso: [](ch_snes), `SNES`, `SNESFAS`, `SNESSolve()`, `SNESSetAlwaysComputesFinalResidual()`
109: @*/
110: PetscErrorCode SNESGetAlwaysComputesFinalResidual(SNES snes, PetscBool *flg)
111: {
112: PetscFunctionBegin;
114: *flg = snes->alwayscomputesfinalresidual;
115: PetscFunctionReturn(PETSC_SUCCESS);
116: }
118: /*@
119: SNESSetFunctionDomainError - tells `SNES` that the input vector, a proposed new solution, to your function you provided to `SNESSetFunction()` is not
120: in the function's domain. For example, a step with negative pressure.
122: Not Collective
124: Input Parameter:
125: . snes - the `SNES` context
127: Level: advanced
129: Notes:
130: This does not need to be called by all processes in the `SNES` MPI communicator.
132: A few solvers will try to cut the step size to avoid the domain error but for other solvers `SNESSolve()` stops iterating and
133: returns with a `SNESConvergedReason` of `SNES_DIVERGED_FUNCTION_DOMAIN`
135: You can direct `SNES` to avoid certain steps by using `SNESVISetVariableBounds()`, `SNESVISetComputeVariableBounds()` or
136: `SNESLineSearchSetPreCheck()`, `SNESLineSearchSetPostCheck()`
138: You should always call `SNESGetConvergedReason()` after each `SNESSolve()` and verify if the iteration converged (positive result) or diverged (negative result).
140: You can call `SNESSetJacobianDomainError()` during a Jacobian computation to indicate the proposed solution is not in the domain.
142: Developer Note:
143: This value is used by `SNESCheckFunctionDomainError()` to determine if the `SNESConvergedReason` is set to `SNES_DIVERGED_FUNCTION_DOMAIN`
145: .seealso: [](ch_snes), `SNESCreate()`, `SNESSetFunction()`, `SNESFunctionFn`, `SNESSetJacobianDomainError()`, `SNESVISetVariableBounds()`,
146: `SNESVISetComputeVariableBounds()`, `SNESLineSearchSetPreCheck()`, `SNESLineSearchSetPostCheck()`, `SNESConvergedReason`, `SNESGetConvergedReason()`,
147: `SNES_DIVERGED_FUNCTION_DOMAIN`, `SNESSetObjectiveDomainError()`, `SNES_DIVERGED_OBJECTIVE_DOMAIN`
148: @*/
149: PetscErrorCode SNESSetFunctionDomainError(SNES snes)
150: {
151: PetscFunctionBegin;
153: snes->functiondomainerror = PETSC_TRUE;
154: PetscFunctionReturn(PETSC_SUCCESS);
155: }
157: /*@
158: SNESSetObjectiveDomainError - tells `SNES` that the input vector, a proposed new solution, to your function you provided to `SNESSetObjective()` is not
159: in the function's domain. For example, a step with negative pressure.
161: Not Collective
163: Input Parameter:
164: . snes - the `SNES` context
166: Level: advanced
168: Notes:
169: This does not need to be called by all processes in the `SNES` MPI communicator.
171: A few solvers will try to cut the step size to avoid the domain error but for other solvers `SNESSolve()` stops iterating and
172: returns with a `SNESConvergedReason` of `SNES_DIVERGED_OBJECTIVE_DOMAIN`
174: You can direct `SNES` to avoid certain steps by using `SNESVISetVariableBounds()`, `SNESVISetComputeVariableBounds()` or
175: `SNESLineSearchSetPreCheck()`, `SNESLineSearchSetPostCheck()`
177: You should always call `SNESGetConvergedReason()` after each `SNESSolve()` and verify if the iteration converged (positive result) or diverged (negative result).
179: You can call `SNESSetJacobianDomainError()` during a Jacobian computation to indicate the proposed solution is not in the domain.
181: Developer Note:
182: This value is used by `SNESCheckObjectiveDomainError()` to determine if the `SNESConvergedReason` is set to `SNES_DIVERGED_OBJECTIVE_DOMAIN`
184: .seealso: [](ch_snes), `SNESCreate()`, `SNESSetFunction()`, `SNESFunctionFn`, `SNESSetJacobianDomainError()`, `SNESVISetVariableBounds()`,
185: `SNESVISetComputeVariableBounds()`, `SNESLineSearchSetPreCheck()`, `SNESLineSearchSetPostCheck()`, `SNESConvergedReason`, `SNESGetConvergedReason()`,
186: `SNES_DIVERGED_OBJECTIVE_DOMAIN`, `SNESSetFunctionDomainError()`, `SNES_DIVERGED_FUNCTION_DOMAIN`
187: @*/
188: PetscErrorCode SNESSetObjectiveDomainError(SNES snes)
189: {
190: PetscFunctionBegin;
192: snes->objectivedomainerror = PETSC_TRUE;
193: PetscFunctionReturn(PETSC_SUCCESS);
194: }
196: /*@
197: SNESSetJacobianDomainError - tells `SNES` that the function you provided to `SNESSetJacobian()` at the proposed step. For example there is a negative element transformation.
199: Logically Collective
201: Input Parameter:
202: . snes - the `SNES` context
204: Level: advanced
206: Notes:
207: If this is called the `SNESSolve()` stops iterating and returns with a `SNESConvergedReason` of `SNES_DIVERGED_JACOBIAN_DOMAIN`
209: You should always call `SNESGetConvergedReason()` after each `SNESSolve()` and verify if the iteration converged (positive result) or diverged (negative result).
211: You can direct `SNES` to avoid certain steps by using `SNESVISetVariableBounds()`, `SNESVISetComputeVariableBounds()` or
212: `SNESLineSearchSetPreCheck()`, `SNESLineSearchSetPostCheck()`
214: .seealso: [](ch_snes), `SNESCreate()`, `SNESSetFunction()`, `SNESFunctionFn`, `SNESSetFunctionDomainError()`, `SNESVISetVariableBounds()`,
215: `SNESVISetComputeVariableBounds()`, `SNESLineSearchSetPreCheck()`, `SNESLineSearchSetPostCheck()`, `SNESConvergedReason`, `SNESGetConvergedReason()`
216: @*/
217: PetscErrorCode SNESSetJacobianDomainError(SNES snes)
218: {
219: PetscFunctionBegin;
221: PetscCheck(!snes->errorifnotconverged, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "User code indicates computeJacobian does not make sense");
222: snes->jacobiandomainerror = PETSC_TRUE;
223: PetscFunctionReturn(PETSC_SUCCESS);
224: }
226: /*@
227: SNESSetCheckJacobianDomainError - tells `SNESSolve()` whether to check if the user called `SNESSetJacobianDomainError()` to indicate a Jacobian domain error after
228: each Jacobian evaluation.
230: Logically Collective
232: Input Parameters:
233: + snes - the `SNES` context
234: - flg - indicates if or not to check Jacobian domain error after each Jacobian evaluation
236: Level: advanced
238: Notes:
239: By default, it checks for the Jacobian domain error in the debug mode, and does not check it in the optimized mode.
241: Checks require one extra parallel synchronization for each Jacobian evaluation
243: .seealso: [](ch_snes), `SNES`, `SNESConvergedReason`, `SNESCreate()`, `SNESSetFunction()`, `SNESFunctionFn`, `SNESSetFunctionDomainError()`, `SNESGetCheckJacobianDomainError()`
244: @*/
245: PetscErrorCode SNESSetCheckJacobianDomainError(SNES snes, PetscBool flg)
246: {
247: PetscFunctionBegin;
249: snes->checkjacdomainerror = flg;
250: PetscFunctionReturn(PETSC_SUCCESS);
251: }
253: /*@
254: SNESGetCheckJacobianDomainError - Get an indicator whether or not `SNES` is checking Jacobian domain errors after each Jacobian evaluation.
256: Logically Collective
258: Input Parameter:
259: . snes - the `SNES` context
261: Output Parameter:
262: . flg - `PETSC_FALSE` indicates that it is not checking Jacobian domain errors after each Jacobian evaluation
264: Level: advanced
266: .seealso: [](ch_snes), `SNES`, `SNESCreate()`, `SNESSetFunction()`, `SNESFunctionFn`, `SNESSetFunctionDomainError()`, `SNESSetCheckJacobianDomainError()`
267: @*/
268: PetscErrorCode SNESGetCheckJacobianDomainError(SNES snes, PetscBool *flg)
269: {
270: PetscFunctionBegin;
272: PetscAssertPointer(flg, 2);
273: *flg = snes->checkjacdomainerror;
274: PetscFunctionReturn(PETSC_SUCCESS);
275: }
277: /*@
278: SNESLoad - Loads a `SNES` that has been stored in `PETSCVIEWERBINARY` with `SNESView()`.
280: Collective
282: Input Parameters:
283: + snes - the newly loaded `SNES`, this needs to have been created with `SNESCreate()` or
284: some related function before a call to `SNESLoad()`.
285: - viewer - binary file viewer, obtained from `PetscViewerBinaryOpen()`
287: Level: intermediate
289: Note:
290: The `SNESType` is determined by the data in the file, any type set into the `SNES` before this call is ignored.
292: .seealso: [](ch_snes), `SNES`, `PetscViewer`, `SNESCreate()`, `SNESType`, `PetscViewerBinaryOpen()`, `SNESView()`, `MatLoad()`, `VecLoad()`
293: @*/
294: PetscErrorCode SNESLoad(SNES snes, PetscViewer viewer)
295: {
296: PetscBool isbinary;
297: PetscInt classid;
298: char type[256];
299: KSP ksp;
300: DM dm;
301: DMSNES dmsnes;
303: PetscFunctionBegin;
306: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
307: PetscCheck(isbinary, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid viewer; open viewer with PetscViewerBinaryOpen()");
309: PetscCall(PetscViewerBinaryRead(viewer, &classid, 1, NULL, PETSC_INT));
310: PetscCheck(classid == SNES_FILE_CLASSID, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONG, "Not SNES next in file");
311: PetscCall(PetscViewerBinaryRead(viewer, type, 256, NULL, PETSC_CHAR));
312: PetscCall(SNESSetType(snes, type));
313: PetscTryTypeMethod(snes, load, viewer);
314: PetscCall(SNESGetDM(snes, &dm));
315: PetscCall(DMGetDMSNES(dm, &dmsnes));
316: PetscCall(DMSNESLoad(dmsnes, viewer));
317: PetscCall(SNESGetKSP(snes, &ksp));
318: PetscCall(KSPLoad(ksp, viewer));
319: PetscFunctionReturn(PETSC_SUCCESS);
320: }
322: #include <petscdraw.h>
323: #if PetscDefined(HAVE_SAWS)
324: #include <petscviewersaws.h>
325: #endif
327: /*@
328: SNESViewFromOptions - View a `SNES` based on values in the options database
330: Collective
332: Input Parameters:
333: + A - the `SNES` context
334: . obj - optional object that provides the options prefix for the checks, pass `NULL` to use the options prefix of `A`
335: - name - command line option
337: Options Database Key:
338: . -name viewer_specification - See `PetscOptionsCreateViewer()` for the values of `viewer_specification`
340: Level: intermediate
342: Note:
343: This checks the options database, creates the viewer on-the-fly, uses it and then destroys it. Hence it should not be called in heavily used routines,
344: rather `PetscOptionsCreateViewer()` should be used to construct the viewer once which can then be utilized in the heavily used routine.
346: .seealso: [](ch_snes), `SNES`, `SNESView()`, `PetscObjectViewFromOptions()`, `SNESCreate()`, `PetscOptionsCreateViewer()`
347: @*/
348: PetscErrorCode SNESViewFromOptions(SNES A, PetscObject obj, const char name[])
349: {
350: PetscFunctionBegin;
352: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
353: PetscFunctionReturn(PETSC_SUCCESS);
354: }
356: PETSC_EXTERN PetscErrorCode SNESComputeJacobian_DMDA(SNES, Vec, Mat, Mat, void *);
358: /*@
359: SNESView - Prints or visualizes the `SNES` data structure.
361: Collective
363: Input Parameters:
364: + snes - the `SNES` context
365: - viewer - the `PetscViewer`
367: Options Database Key:
368: . -snes_view viewer_specification - Calls `SNESView()` at end of `SNESSolve()`, see `PetscOptionsCreateViewer()` for the format of `viewer_specification`
370: Level: beginner
372: Notes:
373: The available visualization contexts include
374: + `PETSC_VIEWER_STDOUT_SELF` - standard output (default)
375: - `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard
376: output where only the first processor opens
377: the file. All other processors send their
378: data to the first processor to print.
380: The available formats include
381: + `PETSC_VIEWER_DEFAULT` - standard output (default)
382: - `PETSC_VIEWER_ASCII_INFO_DETAIL` - more verbose output for `SNESNASM`
384: The user can open an alternative visualization context with
385: `PetscViewerASCIIOpen()` - output to a specified file.
387: In the debugger you can do "call `SNESView`(snes,0)" to display the `SNES` solver. (The same holds for any PETSc object viewer).
389: .seealso: [](ch_snes), `SNES`, `SNESLoad()`, `SNESCreate()`, `PetscViewerASCIIOpen()`, `SNESViewFromOptions()`, `PetscOptionsCreateViewer()`
390: @*/
391: PetscErrorCode SNESView(SNES snes, PetscViewer viewer)
392: {
393: SNESKSPEW *kctx;
394: KSP ksp;
395: Vec u;
396: SNESLineSearch linesearch;
397: PetscBool isascii, isstring, isbinary, isdraw;
398: DMSNES dmsnes;
399: #if PetscDefined(HAVE_SAWS)
400: PetscBool issaws;
401: #endif
403: PetscFunctionBegin;
405: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)snes), &viewer));
407: PetscCheckSameComm(snes, 1, viewer, 2);
409: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
410: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSTRING, &isstring));
411: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
412: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERDRAW, &isdraw));
413: #if PetscDefined(HAVE_SAWS)
414: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSAWS, &issaws));
415: #endif
416: if (isascii) {
417: SNESNormSchedule normschedule;
418: DM dm;
419: SNESJacobianFn *cJ;
420: void *ctx;
421: const char *pre = "";
423: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)snes, viewer));
424: if (!snes->setupcalled) PetscCall(PetscViewerASCIIPrintf(viewer, " SNES has not been set up so information may be incomplete\n"));
425: if (snes->ops->view) {
426: PetscCall(PetscViewerASCIIPushTab(viewer));
427: PetscUseTypeMethod(snes, view, viewer);
428: PetscCall(PetscViewerASCIIPopTab(viewer));
429: }
430: if (snes->max_funcs == PETSC_UNLIMITED) {
431: PetscCall(PetscViewerASCIIPrintf(viewer, " maximum iterations=%" PetscInt_FMT ", maximum function evaluations=unlimited\n", snes->max_its));
432: } else {
433: PetscCall(PetscViewerASCIIPrintf(viewer, " maximum iterations=%" PetscInt_FMT ", maximum function evaluations=%" PetscInt_FMT "\n", snes->max_its, snes->max_funcs));
434: }
435: PetscCall(PetscViewerASCIIPrintf(viewer, " tolerances: relative=%g, absolute=%g, solution=%g\n", (double)snes->rtol, (double)snes->abstol, (double)snes->stol));
436: if (snes->usesksp) PetscCall(PetscViewerASCIIPrintf(viewer, " total number of linear solver iterations=%" PetscInt_FMT "\n", snes->linear_its));
437: PetscCall(PetscViewerASCIIPrintf(viewer, " total number of function evaluations=%" PetscInt_FMT "\n", snes->nfuncs));
438: PetscCall(SNESGetNormSchedule(snes, &normschedule));
439: if (normschedule > 0) PetscCall(PetscViewerASCIIPrintf(viewer, " norm schedule %s\n", SNESNormSchedules[normschedule]));
440: if (snes->gridsequence) PetscCall(PetscViewerASCIIPrintf(viewer, " total number of grid sequence refinements=%" PetscInt_FMT "\n", snes->gridsequence));
441: if (snes->ksp_ewconv) {
442: kctx = (SNESKSPEW *)snes->kspconvctx;
443: if (kctx) {
444: PetscCall(PetscViewerASCIIPrintf(viewer, " Eisenstat-Walker computation of KSP relative tolerance (version %" PetscInt_FMT ")\n", kctx->version));
445: PetscCall(PetscViewerASCIIPrintf(viewer, " rtol_0=%g, rtol_max=%g, threshold=%g\n", (double)kctx->rtol_0, (double)kctx->rtol_max, (double)kctx->threshold));
446: PetscCall(PetscViewerASCIIPrintf(viewer, " gamma=%g, alpha=%g, alpha2=%g\n", (double)kctx->gamma, (double)kctx->alpha, (double)kctx->alpha2));
447: }
448: }
449: if (snes->lagpreconditioner == -1) {
450: PetscCall(PetscViewerASCIIPrintf(viewer, " Preconditioned is never rebuilt\n"));
451: } else if (snes->lagpreconditioner > 1) {
452: PetscCall(PetscViewerASCIIPrintf(viewer, " Preconditioned is rebuilt every %" PetscInt_FMT " new Jacobians\n", snes->lagpreconditioner));
453: }
454: if (snes->lagjacobian == -1) {
455: PetscCall(PetscViewerASCIIPrintf(viewer, " Jacobian is never rebuilt\n"));
456: } else if (snes->lagjacobian > 1) {
457: PetscCall(PetscViewerASCIIPrintf(viewer, " Jacobian is rebuilt every %" PetscInt_FMT " SNES iterations\n", snes->lagjacobian));
458: }
459: PetscCall(SNESGetDM(snes, &dm));
460: PetscCall(DMSNESGetJacobian(dm, &cJ, &ctx));
461: if (snes->mf_operator) {
462: PetscCall(PetscViewerASCIIPrintf(viewer, " Jacobian is applied matrix-free with differencing\n"));
463: pre = "Preconditioning ";
464: }
465: if (cJ == SNESComputeJacobianDefault) {
466: PetscCall(PetscViewerASCIIPrintf(viewer, " %sJacobian is built using finite differences one column at a time\n", pre));
467: } else if (cJ == SNESComputeJacobianDefaultColor) {
468: PetscCall(PetscViewerASCIIPrintf(viewer, " %sJacobian is built using finite differences with coloring\n", pre));
469: /* it slightly breaks data encapsulation for access the DMDA information directly */
470: } else if (cJ == SNESComputeJacobian_DMDA) {
471: MatFDColoring fdcoloring;
472: PetscCall(PetscObjectQuery((PetscObject)dm, "DMDASNES_FDCOLORING", (PetscObject *)&fdcoloring));
473: if (fdcoloring) {
474: PetscCall(PetscViewerASCIIPrintf(viewer, " %sJacobian is built using colored finite differences on a DMDA\n", pre));
475: } else {
476: PetscCall(PetscViewerASCIIPrintf(viewer, " %sJacobian is built using a DMDA local Jacobian\n", pre));
477: }
478: } else if (snes->mf && !snes->mf_operator) {
479: PetscCall(PetscViewerASCIIPrintf(viewer, " Jacobian is applied matrix-free with differencing, no explicit Jacobian\n"));
480: }
481: } else if (isstring) {
482: const char *type;
483: PetscCall(SNESGetType(snes, &type));
484: PetscCall(PetscViewerStringSPrintf(viewer, " SNESType: %-7.7s", type));
485: PetscTryTypeMethod(snes, view, viewer);
486: } else if (isbinary) {
487: PetscInt classid = SNES_FILE_CLASSID;
488: MPI_Comm comm;
489: PetscMPIInt rank;
490: char type[256];
492: PetscCall(PetscObjectGetComm((PetscObject)snes, &comm));
493: PetscCallMPI(MPI_Comm_rank(comm, &rank));
494: if (rank == 0) {
495: PetscCall(PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT));
496: PetscCall(PetscStrncpy(type, ((PetscObject)snes)->type_name, sizeof(type)));
497: PetscCall(PetscViewerBinaryWrite(viewer, type, sizeof(type), PETSC_CHAR));
498: }
499: PetscTryTypeMethod(snes, view, viewer);
500: } else if (isdraw) {
501: PetscDraw draw;
502: char str[36];
503: PetscReal x, y, bottom, h;
505: PetscCall(PetscViewerDrawGetDraw(viewer, 0, &draw));
506: PetscCall(PetscDrawGetCurrentPoint(draw, &x, &y));
507: PetscCall(PetscStrncpy(str, "SNES: ", sizeof(str)));
508: PetscCall(PetscStrlcat(str, ((PetscObject)snes)->type_name, sizeof(str)));
509: PetscCall(PetscDrawStringBoxed(draw, x, y, PETSC_DRAW_BLUE, PETSC_DRAW_BLACK, str, NULL, &h));
510: bottom = y - h;
511: PetscCall(PetscDrawPushCurrentPoint(draw, x, bottom));
512: PetscTryTypeMethod(snes, view, viewer);
513: #if PetscDefined(HAVE_SAWS)
514: } else if (issaws) {
515: PetscMPIInt rank;
516: const char *name;
518: PetscCall(PetscObjectGetName((PetscObject)snes, &name));
519: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
520: if (!((PetscObject)snes)->amsmem && rank == 0) {
521: char dir[1024];
523: PetscCall(PetscObjectViewSAWs((PetscObject)snes, viewer));
524: PetscCall(PetscSNPrintf(dir, 1024, "/PETSc/Objects/%s/its", name));
525: PetscCallSAWs(SAWs_Register, dir, &snes->iter, 1, SAWs_READ, SAWs_INT);
526: if (!snes->conv_hist) PetscCall(SNESSetConvergenceHistory(snes, NULL, NULL, PETSC_DECIDE, PETSC_TRUE));
527: PetscCall(PetscSNPrintf(dir, 1024, "/PETSc/Objects/%s/conv_hist", name));
528: PetscCallSAWs(SAWs_Register, dir, snes->conv_hist, 10, SAWs_READ, SAWs_DOUBLE);
529: }
530: #endif
531: }
532: if (snes->linesearch) {
533: PetscCall(SNESGetLineSearch(snes, &linesearch));
534: PetscCall(PetscViewerASCIIPushTab(viewer));
535: PetscCall(SNESLineSearchView(linesearch, viewer));
536: PetscCall(PetscViewerASCIIPopTab(viewer));
537: }
538: if (snes->npc && snes->usesnpc) {
539: PetscCall(PetscViewerASCIIPushTab(viewer));
540: PetscCall(SNESView(snes->npc, viewer));
541: PetscCall(PetscViewerASCIIPopTab(viewer));
542: }
543: PetscCall(PetscViewerASCIIPushTab(viewer));
544: PetscCall(DMGetDMSNES(snes->dm, &dmsnes));
545: PetscCall(DMSNESView(dmsnes, viewer));
546: PetscCall(PetscViewerASCIIPopTab(viewer));
547: if (snes->usesksp) {
548: PetscCall(SNESGetKSP(snes, &ksp));
549: PetscCall(PetscViewerASCIIPushTab(viewer));
550: PetscCall(KSPView(ksp, viewer));
551: PetscCall(PetscViewerASCIIPopTab(viewer));
552: } else {
553: PetscViewerFormat format;
555: PetscCall(SNESGetSolution(snes, &u));
556: PetscCall(PetscViewerGetFormat(viewer, &format));
557: if (u && isascii) {
558: if (format != PETSC_VIEWER_ASCII_INFO_DETAIL) PetscCall(PetscViewerPushFormat(viewer, PETSC_VIEWER_ASCII_INFO));
559: PetscCall(PetscViewerASCIIPrintf(viewer, "solution vector:\n"));
560: PetscCall(PetscViewerASCIIPushTab(viewer));
561: PetscCall(VecView(u, viewer));
562: PetscCall(PetscViewerASCIIPopTab(viewer));
563: if (format != PETSC_VIEWER_ASCII_INFO_DETAIL) PetscCall(PetscViewerPopFormat(viewer));
564: }
565: }
566: if (isdraw) {
567: PetscDraw draw;
568: PetscCall(PetscViewerDrawGetDraw(viewer, 0, &draw));
569: PetscCall(PetscDrawPopCurrentPoint(draw));
570: }
571: PetscFunctionReturn(PETSC_SUCCESS);
572: }
574: /*
575: We retain a list of functions that also take SNES command
576: line options. These are called at the end SNESSetFromOptions()
577: */
578: #define MAXSETFROMOPTIONS 5
579: static PetscInt numberofsetfromoptions;
580: static PetscErrorCode (*othersetfromoptions[MAXSETFROMOPTIONS])(SNES);
582: /*@
583: SNESAddOptionsChecker - Adds an additional function to check for `SNES` options.
585: Not Collective
587: Input Parameter:
588: . snescheck - function that checks for options
590: Calling sequence of `snescheck`:
591: . snes - the `SNES` object for which it is checking options
593: Level: developer
595: .seealso: [](ch_snes), `SNES`, `SNESSetFromOptions()`
596: @*/
597: PetscErrorCode SNESAddOptionsChecker(PetscErrorCode (*snescheck)(SNES snes))
598: {
599: PetscFunctionBegin;
600: PetscCheck(numberofsetfromoptions < MAXSETFROMOPTIONS, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Too many options checkers, only %d allowed", MAXSETFROMOPTIONS);
601: othersetfromoptions[numberofsetfromoptions++] = snescheck;
602: PetscFunctionReturn(PETSC_SUCCESS);
603: }
605: static PetscErrorCode SNESSetUpMatrixFree_Private(SNES snes, PetscBool hasOperator, PetscInt version)
606: {
607: Mat J;
608: MatNullSpace nullsp;
610: PetscFunctionBegin;
613: if (!snes->vec_func && (snes->jacobian || snes->jacobian_pre)) {
614: Mat A = snes->jacobian, B = snes->jacobian_pre;
615: PetscCall(MatCreateVecs(A ? A : B, NULL, &snes->vec_func));
616: }
618: PetscCheck(version == 1 || version == 2, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "matrix-free operator routines, only version 1 and 2");
619: if (version == 1) {
620: PetscCall(MatCreateSNESMF(snes, &J));
621: PetscCall(MatMFFDSetOptionsPrefix(J, ((PetscObject)snes)->prefix));
622: PetscCall(MatSetFromOptions(J));
623: /* TODO: the version 2 code should be merged into the MatCreateSNESMF() and MatCreateMFFD() infrastructure and then removed */
624: } else /* if (version == 2) */ {
625: PetscCheck(snes->vec_func, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "SNESSetFunction() must be called first");
626: #if !PetscDefined(USE_COMPLEX) && !PetscDefined(USE_REAL_SINGLE) && !PetscDefined(USE_REAL___FLOAT128) && !PetscDefined(USE_REAL___FP16)
627: PetscCall(MatCreateSNESMFMore(snes, snes->vec_func, &J));
628: #else
629: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "matrix-free operator routines (version 2)");
630: #endif
631: }
633: /* attach any user provided null space that was on Amat to the newly created matrix-free matrix */
634: if (snes->jacobian) {
635: PetscCall(MatGetNullSpace(snes->jacobian, &nullsp));
636: if (nullsp) PetscCall(MatSetNullSpace(J, nullsp));
637: }
639: PetscCall(PetscInfo(snes, "Setting default matrix-free operator routines (version %" PetscInt_FMT ")\n", version));
640: if (hasOperator) {
641: /* This version replaces the user provided Jacobian matrix with a
642: matrix-free version but still employs the user-provided matrix used for computing the preconditioner. */
643: PetscCall(SNESSetJacobian(snes, J, NULL, NULL, NULL));
644: } else {
645: /* This version replaces both the user-provided Jacobian and the user-
646: provided preconditioner Jacobian with the default matrix-free version. */
647: if (snes->npcside == PC_LEFT && snes->npc) {
648: if (!snes->jacobian) PetscCall(SNESSetJacobian(snes, J, NULL, NULL, NULL));
649: } else PetscCall(SNESSetJacobian(snes, J, J, MatMFFDComputeJacobian, NULL));
650: }
651: PetscCall(MatDestroy(&J));
652: PetscFunctionReturn(PETSC_SUCCESS);
653: }
655: static PetscErrorCode DMRestrictHook_SNESVecSol(DM dmfine, Mat Restrict, Vec Rscale, Mat Inject, DM dmcoarse, PetscCtx ctx)
656: {
657: SNES snes = (SNES)ctx;
658: Vec Xfine, Xfine_named = NULL, Xcoarse;
660: PetscFunctionBegin;
661: if (PetscLogPrintInfo) {
662: PetscInt finelevel, coarselevel, fineclevel, coarseclevel;
663: PetscCall(DMGetRefineLevel(dmfine, &finelevel));
664: PetscCall(DMGetCoarsenLevel(dmfine, &fineclevel));
665: PetscCall(DMGetRefineLevel(dmcoarse, &coarselevel));
666: PetscCall(DMGetCoarsenLevel(dmcoarse, &coarseclevel));
667: PetscCall(PetscInfo(dmfine, "Restricting SNES solution vector from level %" PetscInt_FMT "-%" PetscInt_FMT " to level %" PetscInt_FMT "-%" PetscInt_FMT "\n", finelevel, fineclevel, coarselevel, coarseclevel));
668: }
669: if (dmfine == snes->dm) Xfine = snes->vec_sol;
670: else {
671: PetscCall(DMGetNamedGlobalVector(dmfine, "SNESVecSol", &Xfine_named));
672: Xfine = Xfine_named;
673: }
674: PetscCall(DMGetNamedGlobalVector(dmcoarse, "SNESVecSol", &Xcoarse));
675: if (Inject) {
676: PetscCall(MatRestrict(Inject, Xfine, Xcoarse));
677: } else {
678: PetscCall(MatRestrict(Restrict, Xfine, Xcoarse));
679: PetscCall(VecPointwiseMult(Xcoarse, Xcoarse, Rscale));
680: }
681: PetscCall(DMRestoreNamedGlobalVector(dmcoarse, "SNESVecSol", &Xcoarse));
682: if (Xfine_named) PetscCall(DMRestoreNamedGlobalVector(dmfine, "SNESVecSol", &Xfine_named));
683: PetscFunctionReturn(PETSC_SUCCESS);
684: }
686: static PetscErrorCode DMCoarsenHook_SNESVecSol(DM dm, DM dmc, PetscCtx ctx)
687: {
688: PetscFunctionBegin;
689: PetscCall(DMCoarsenHookAdd(dmc, DMCoarsenHook_SNESVecSol, DMRestrictHook_SNESVecSol, ctx));
690: PetscFunctionReturn(PETSC_SUCCESS);
691: }
693: /* This may be called to rediscretize the operator on levels of linear multigrid. The DM shuffle is so the user can
694: * safely call SNESGetDM() in their residual evaluation routine. */
695: static PetscErrorCode KSPComputeOperators_SNES(KSP ksp, Mat A, Mat B, PetscCtx ctx)
696: {
697: SNES snes = (SNES)ctx;
698: DMSNES sdm;
699: Vec X, Xnamed = NULL;
700: DM dmsave;
701: void *ctxsave;
702: SNESJacobianFn *jac = NULL;
704: PetscFunctionBegin;
705: dmsave = snes->dm;
706: PetscCall(KSPGetDM(ksp, &snes->dm));
707: if (dmsave == snes->dm) X = snes->vec_sol; /* We are on the finest level */
708: else {
709: PetscBool has;
711: /* We are on a coarser level, this vec was initialized using a DM restrict hook */
712: PetscCall(DMHasNamedGlobalVector(snes->dm, "SNESVecSol", &has));
713: PetscCheck(has, PetscObjectComm((PetscObject)snes->dm), PETSC_ERR_PLIB, "Missing SNESVecSol");
714: PetscCall(DMGetNamedGlobalVector(snes->dm, "SNESVecSol", &Xnamed));
715: X = Xnamed;
716: PetscCall(SNESGetJacobian(snes, NULL, NULL, &jac, &ctxsave));
717: /* If the DM's don't match up, the MatFDColoring context needed for the jacobian won't match up either -- fixit. */
718: if (jac == SNESComputeJacobianDefaultColor) PetscCall(SNESSetJacobian(snes, NULL, NULL, SNESComputeJacobianDefaultColor, NULL));
719: }
721: /* Compute the operators */
722: PetscCall(DMGetDMSNES(snes->dm, &sdm));
723: if (Xnamed && sdm->ops->computefunction) {
724: /* The SNES contract with the user is that ComputeFunction is always called before ComputeJacobian.
725: We make sure of this here. Disable affine shift since it is for the finest level */
726: Vec F, saverhs = snes->vec_rhs;
728: snes->vec_rhs = NULL;
729: PetscCall(DMGetGlobalVector(snes->dm, &F));
730: PetscCall(SNESComputeFunction(snes, X, F));
731: PetscCall(DMRestoreGlobalVector(snes->dm, &F));
732: snes->vec_rhs = saverhs;
733: snes->nfuncs--; /* Do not log coarser level evaluations */
734: }
735: /* Make sure KSP DM has the Jacobian computation routine */
736: if (!sdm->ops->computejacobian) PetscCall(DMCopyDMSNES(dmsave, snes->dm));
737: PetscCall(SNESComputeJacobian(snes, X, A, B)); /* cannot handle previous SNESSetJacobianDomainError() calls */
739: /* Put the previous context back */
740: if (snes->dm != dmsave && jac == SNESComputeJacobianDefaultColor) PetscCall(SNESSetJacobian(snes, NULL, NULL, jac, ctxsave));
742: if (Xnamed) PetscCall(DMRestoreNamedGlobalVector(snes->dm, "SNESVecSol", &Xnamed));
743: snes->dm = dmsave;
744: PetscFunctionReturn(PETSC_SUCCESS);
745: }
747: /*@
748: SNESSetUpMatrices - ensures that matrices are available for `SNES` Newton-like methods, this is called by `SNESSetUp_XXX()`
750: Collective
752: Input Parameter:
753: . snes - `SNES` object to configure
755: Level: developer
757: Note:
758: If the matrices do not yet exist it attempts to create them based on options previously set for the `SNES` such as `-snes_mf`
760: Developer Note:
761: The functionality of this routine overlaps in a confusing way with the functionality of `SNESSetUpMatrixFree_Private()` which is called by
762: `SNESSetUp()` but sometimes `SNESSetUpMatrices()` is called without `SNESSetUp()` being called. A refactorization to simplify the
763: logic that handles the matrix-free case is desirable.
765: .seealso: [](ch_snes), `SNES`, `SNESSetUp()`
766: @*/
767: PetscErrorCode SNESSetUpMatrices(SNES snes)
768: {
769: DM dm;
770: DMSNES sdm;
772: PetscFunctionBegin;
773: PetscCall(SNESGetDM(snes, &dm));
774: PetscCall(DMGetDMSNES(dm, &sdm));
775: if (!snes->jacobian && snes->mf && !snes->mf_operator && !snes->jacobian_pre) {
776: Mat J;
777: void *functx;
778: PetscCall(MatCreateSNESMF(snes, &J));
779: PetscCall(MatMFFDSetOptionsPrefix(J, ((PetscObject)snes)->prefix));
780: PetscCall(MatSetFromOptions(J));
781: PetscCall(SNESGetFunction(snes, NULL, NULL, &functx));
782: PetscCall(SNESSetJacobian(snes, J, J, NULL, NULL));
783: PetscCall(MatDestroy(&J));
784: } else if (snes->mf_operator && !snes->jacobian_pre && !snes->jacobian) {
785: Mat J, B;
786: PetscCall(MatCreateSNESMF(snes, &J));
787: PetscCall(MatMFFDSetOptionsPrefix(J, ((PetscObject)snes)->prefix));
788: PetscCall(MatSetFromOptions(J));
789: PetscCall(DMCreateMatrix(snes->dm, &B));
790: /* sdm->computejacobian was already set to reach here */
791: PetscCall(SNESSetJacobian(snes, J, B, NULL, NULL));
792: PetscCall(MatDestroy(&J));
793: PetscCall(MatDestroy(&B));
794: } else if (!snes->jacobian_pre) {
795: PetscDS prob;
796: Mat J, B;
797: PetscBool hasPrec = PETSC_FALSE;
799: J = snes->jacobian;
800: PetscCall(DMGetDS(dm, &prob));
801: if (prob) PetscCall(PetscDSHasJacobianPreconditioner(prob, &hasPrec));
802: if (!J && hasPrec) PetscCall(DMCreateMatrix(snes->dm, &J));
803: else PetscCall(PetscObjectReference((PetscObject)J));
804: PetscCall(DMCreateMatrix(snes->dm, &B));
805: PetscCall(SNESSetJacobian(snes, J ? J : B, B, NULL, NULL));
806: PetscCall(MatDestroy(&J));
807: PetscCall(MatDestroy(&B));
808: }
809: {
810: KSP ksp;
811: PetscCall(SNESGetKSP(snes, &ksp));
812: PetscCall(KSPSetComputeOperators(ksp, KSPComputeOperators_SNES, snes));
813: PetscCall(DMCoarsenHookAdd(snes->dm, DMCoarsenHook_SNESVecSol, DMRestrictHook_SNESVecSol, snes));
814: }
815: PetscFunctionReturn(PETSC_SUCCESS);
816: }
818: PETSC_SINGLE_LIBRARY_INTERN PetscErrorCode PetscMonitorPauseFinal_Internal(PetscInt, PetscCtx);
820: static PetscErrorCode SNESMonitorPauseFinal_Internal(SNES snes)
821: {
822: PetscFunctionBegin;
823: if (!snes->pauseFinal) PetscFunctionReturn(PETSC_SUCCESS);
824: PetscCall(PetscMonitorPauseFinal_Internal(snes->numbermonitors, snes->monitorcontext));
825: PetscFunctionReturn(PETSC_SUCCESS);
826: }
828: /*@
829: SNESMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type indicated by the user
831: Collective
833: Input Parameters:
834: + snes - `SNES` object you wish to monitor
835: . name - the monitor type one is seeking
836: . help - message indicating what monitoring is done
837: . manual - manual page for the monitor
838: . monitor - the monitor function, this must use a `PetscViewerFormat` as its context
839: - monitorsetup - a function that is called once ONLY if the user selected this monitor that may set additional features of the `SNES` or `PetscViewer` objects
841: Calling sequence of `monitor`:
842: + snes - the nonlinear solver context
843: . it - the current iteration
844: . r - the current function norm
845: - vf - a `PetscViewerAndFormat` struct that contains the `PetscViewer` and `PetscViewerFormat` to use
847: Calling sequence of `monitorsetup`:
848: + snes - the nonlinear solver context
849: - vf - a `PetscViewerAndFormat` struct that contains the `PetscViewer` and `PetscViewerFormat` to use
851: Options Database Key:
852: . -name - trigger the use of this monitor in `SNESSetFromOptions()`
854: Level: advanced
856: .seealso: [](ch_snes), `PetscOptionsCreateViewer()`, `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`,
857: `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`,
858: `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`,
859: `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`,
860: `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`,
861: `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`,
862: `PetscOptionsFList()`, `PetscOptionsEList()`
863: @*/
864: PetscErrorCode SNESMonitorSetFromOptions(SNES snes, const char name[], const char help[], const char manual[], PetscErrorCode (*monitor)(SNES snes, PetscInt it, PetscReal r, PetscViewerAndFormat *vf), PetscErrorCode (*monitorsetup)(SNES snes, PetscViewerAndFormat *vf))
865: {
866: PetscViewer viewer;
867: PetscViewerFormat format;
868: PetscBool flg;
870: PetscFunctionBegin;
871: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, name, &viewer, &format, &flg));
872: if (flg) {
873: PetscViewerAndFormat *vf;
874: PetscCall(PetscViewerAndFormatCreate(viewer, format, &vf));
875: PetscCall(PetscViewerDestroy(&viewer));
876: if (monitorsetup) PetscCall((*monitorsetup)(snes, vf));
877: PetscCall(SNESMonitorSet(snes, (PetscErrorCode (*)(SNES, PetscInt, PetscReal, PetscCtx))monitor, vf, (PetscCtxDestroyFn *)PetscViewerAndFormatDestroy));
878: }
879: PetscFunctionReturn(PETSC_SUCCESS);
880: }
882: PetscErrorCode SNESEWSetFromOptions_Private(SNESKSPEW *kctx, PetscBool print_api, MPI_Comm comm, const char *prefix)
883: {
884: const char *api = print_api ? "SNESKSPSetParametersEW" : NULL;
886: PetscFunctionBegin;
887: PetscOptionsBegin(comm, prefix, "Eisenstat and Walker type forcing options", "KSP");
888: PetscCall(PetscOptionsInt("-ksp_ew_version", "Version 1, 2 or 3", api, kctx->version, &kctx->version, NULL));
889: PetscCall(PetscOptionsReal("-ksp_ew_rtol0", "0 <= rtol0 < 1", api, kctx->rtol_0, &kctx->rtol_0, NULL));
890: kctx->rtol_max = PetscMax(kctx->rtol_0, kctx->rtol_max);
891: PetscCall(PetscOptionsReal("-ksp_ew_rtolmax", "0 <= rtolmax < 1", api, kctx->rtol_max, &kctx->rtol_max, NULL));
892: PetscCall(PetscOptionsReal("-ksp_ew_gamma", "0 <= gamma <= 1", api, kctx->gamma, &kctx->gamma, NULL));
893: PetscCall(PetscOptionsReal("-ksp_ew_alpha", "1 < alpha <= 2", api, kctx->alpha, &kctx->alpha, NULL));
894: PetscCall(PetscOptionsReal("-ksp_ew_alpha2", "alpha2", NULL, kctx->alpha2, &kctx->alpha2, NULL));
895: PetscCall(PetscOptionsReal("-ksp_ew_threshold", "0 < threshold < 1", api, kctx->threshold, &kctx->threshold, NULL));
896: PetscCall(PetscOptionsReal("-ksp_ew_v4_p1", "p1", NULL, kctx->v4_p1, &kctx->v4_p1, NULL));
897: PetscCall(PetscOptionsReal("-ksp_ew_v4_p2", "p2", NULL, kctx->v4_p2, &kctx->v4_p2, NULL));
898: PetscCall(PetscOptionsReal("-ksp_ew_v4_p3", "p3", NULL, kctx->v4_p3, &kctx->v4_p3, NULL));
899: PetscCall(PetscOptionsReal("-ksp_ew_v4_m1", "Scaling when rk-1 in [p2,p3)", NULL, kctx->v4_m1, &kctx->v4_m1, NULL));
900: PetscCall(PetscOptionsReal("-ksp_ew_v4_m2", "Scaling when rk-1 in [p3,+infty)", NULL, kctx->v4_m2, &kctx->v4_m2, NULL));
901: PetscCall(PetscOptionsReal("-ksp_ew_v4_m3", "Threshold for successive rtol (0.1 in Eq.7)", NULL, kctx->v4_m3, &kctx->v4_m3, NULL));
902: PetscCall(PetscOptionsReal("-ksp_ew_v4_m4", "Adaptation scaling (0.5 in Eq.7)", NULL, kctx->v4_m4, &kctx->v4_m4, NULL));
903: PetscOptionsEnd();
904: PetscFunctionReturn(PETSC_SUCCESS);
905: }
907: /*@
908: SNESSetFromOptions - Sets various `SNES` and `KSP` parameters from user options.
910: Collective
912: Input Parameter:
913: . snes - the `SNES` context
915: Options Database Keys:
916: + -snes_view_pre viewer_specification - view information about the nonlinear solver before the solve
917: . -snes_view viewer_specification - view information about the nonlinear solver after the solve
918: . -snes_type type - newtonls, newtontr, ngmres, ncg, nrichardson, qn, vi, fas, `SNESType` for complete list
919: . -snes_rtol rtol - relative decrease in tolerance norm from initial
920: . -snes_atol abstol - absolute tolerance of residual norm
921: . -snes_stol stol - convergence tolerance in terms of the norm of the change in the solution between steps
922: . -snes_divergence_tolerance divtol - if the residual goes above divtol*rnorm0, exit with divergence
923: . -snes_max_it max_it - maximum number of iterations
924: . -snes_max_funcs max_funcs - maximum number of function evaluations
925: . -snes_force_iteration force - force `SNESSolve()` to take at least one iteration
926: . -snes_max_fail max_fail - maximum number of line search failures allowed before stopping, default is none
927: . -snes_max_linear_solve_fail num - number of linear solver failures allowed before `SNESSolve()` is stopped
928: . -snes_lag_preconditioner lag - how often preconditioner is rebuilt (use -1 to never rebuild)
929: . -snes_lag_preconditioner_persists (true|false) - retains the `-snes_lag_preconditioner` information across multiple `SNESSolve()`
930: . -snes_lag_jacobian lag - how often Jacobian is rebuilt (use -1 to never rebuild)
931: . -snes_lag_jacobian_persists (true|false) - retains the `-snes_lag_jacobian` information across multiple `SNESSolve()`
932: . -snes_convergence_test (default|skip|correct_pressure) - convergence test in nonlinear solver. default `SNESConvergedDefault()`. `skip` uses `SNESConvergedSkip()` and means continue
933: iterating until `max_it` or some other criterion is reached. `correct_pressure` uses
934: `SNESConvergedCorrectPressure()` and has special handling of a pressure null space.
935: . -snes_monitor viewer_specification - displays residual norm at each iteration. if no filename given prints to stdout
936: . -snes_monitor_solution viewer_specification - displays solution at each iteration
937: . -snes_monitor_residual viewer_specification - displays residual (not its norm) at each iteration
938: . -snes_monitor_solution_update viewer_specification - displays update to solution at each iteration
939: . -snes_monitor_lg_range (true|false) - plots function range at each iteration
940: . -snes_monitor_pause_final (true|false) - Pauses all monitor drawing after the solver ends
941: . -snes_fd (true|false) - use finite differences to compute Jacobian; very slow, only for testing
942: . -snes_fd_color (true|false) - use finite differences with coloring to compute Jacobian
943: . -snes_converged_reason viewer_specification - print the reason for convergence/divergence after each solve
944: . -npc_snes_type type - the `SNES` type to use as a nonlinear preconditioner
945: . -snes_test_jacobian [threshold] - compare the user provided Jacobian with one computed via finite differences to check for errors.
946: If a threshold is given, display only those entries whose difference is greater than the threshold.
947: - -snes_test_jacobian_view viewer_specification - display the user provided Jacobian, the finite difference Jacobian and the difference between them
948: to help users detect the location of errors in the user provided Jacobian.
950: Options Database Keys for Eisenstat-Walker method:
951: + -snes_ksp_ew - use Eisenstat-Walker method for determining linear system convergence
952: . -snes_ksp_ew_version ver - version of Eisenstat-Walker method
953: . -snes_ksp_ew_rtol0 rtol0 - Sets rtol0
954: . -snes_ksp_ew_rtolmax rtolmax - Sets rtolmax
955: . -snes_ksp_ew_gamma gamma - Sets gamma
956: . -snes_ksp_ew_alpha alpha - Sets alpha
957: . -snes_ksp_ew_alpha2 alpha2 - Sets alpha2
958: - -snes_ksp_ew_threshold threshold - Sets threshold
960: Level: beginner
962: Notes:
963: See `PetscOptionsCreateViewer()` for the format of `viewer_specification`
965: See `KSPSetFromOptions()` for the `KSP` options database keys
967: To see all options, run your program with the -help option or consult the users manual
969: See `SNESLineSearchSetFromOptions()` for all the line search options available
971: `SNES` supports three approaches for computing (approximate) Jacobians: user provided via `SNESSetJacobian()`, matrix-free using `MatCreateSNESMF()`,
972: and computing explicitly with
973: finite differences and coloring using `MatFDColoring`. It is also possible to use automatic differentiation and the `MatFDColoring` object.
975: .seealso: [](ch_snes), `SNESType`, `SNESSetOptionsPrefix()`, `SNESResetFromOptions()`, `SNES`, `SNESCreate()`, `MatCreateSNESMF()`, `MatFDColoring`, `SNESLineSearchSetFromOptions()`
976: @*/
977: PetscErrorCode SNESSetFromOptions(SNES snes)
978: {
979: PetscBool flg, pcset, persist, set;
980: PetscInt i, indx, lag, grids, max_its, max_funcs;
981: const char *deft = SNESNEWTONLS;
982: const char *convtests[] = {"default", "skip", "correct_pressure"};
983: SNESKSPEW *kctx = NULL;
984: char type[256], monfilename[PETSC_MAX_PATH_LEN], ewprefix[256];
985: PCSide pcside;
986: const char *optionsprefix;
987: PetscReal rtol, abstol, stol;
989: PetscFunctionBegin;
991: PetscCall(SNESRegisterAll());
992: PetscObjectOptionsBegin((PetscObject)snes);
993: if (((PetscObject)snes)->type_name) deft = ((PetscObject)snes)->type_name;
994: PetscCall(PetscOptionsFList("-snes_type", "Nonlinear solver method", "SNESSetType", SNESList, deft, type, sizeof(type), &flg));
995: if (flg) PetscCall(SNESSetType(snes, type));
996: else if (!((PetscObject)snes)->type_name) PetscCall(SNESSetType(snes, deft));
998: abstol = snes->abstol;
999: rtol = snes->rtol;
1000: stol = snes->stol;
1001: max_its = snes->max_its;
1002: max_funcs = snes->max_funcs;
1003: PetscCall(PetscOptionsReal("-snes_rtol", "Stop if decrease in function norm less than", "SNESSetTolerances", snes->rtol, &rtol, NULL));
1004: PetscCall(PetscOptionsReal("-snes_atol", "Stop if function norm less than", "SNESSetTolerances", snes->abstol, &abstol, NULL));
1005: PetscCall(PetscOptionsReal("-snes_stol", "Stop if step length less than", "SNESSetTolerances", snes->stol, &stol, NULL));
1006: PetscCall(PetscOptionsInt("-snes_max_it", "Maximum iterations", "SNESSetTolerances", snes->max_its, &max_its, NULL));
1007: PetscCall(PetscOptionsInt("-snes_max_funcs", "Maximum function evaluations", "SNESSetTolerances", snes->max_funcs, &max_funcs, NULL));
1008: PetscCall(SNESSetTolerances(snes, abstol, rtol, stol, max_its, max_funcs));
1010: PetscCall(PetscOptionsReal("-snes_divergence_tolerance", "Stop if residual norm increases by this factor", "SNESSetDivergenceTolerance", snes->divtol, &snes->divtol, &flg));
1011: if (flg) PetscCall(SNESSetDivergenceTolerance(snes, snes->divtol));
1013: PetscCall(PetscOptionsInt("-snes_max_fail", "Maximum nonlinear step failures", "SNESSetMaxNonlinearStepFailures", snes->maxFailures, &snes->maxFailures, &flg));
1014: if (flg) PetscCall(SNESSetMaxNonlinearStepFailures(snes, snes->maxFailures));
1016: PetscCall(PetscOptionsInt("-snes_max_linear_solve_fail", "Maximum failures in linear solves allowed", "SNESSetMaxLinearSolveFailures", snes->maxLinearSolveFailures, &snes->maxLinearSolveFailures, &flg));
1017: if (flg) PetscCall(SNESSetMaxLinearSolveFailures(snes, snes->maxLinearSolveFailures));
1019: PetscCall(PetscOptionsBool("-snes_error_if_not_converged", "Generate error if solver does not converge", "SNESSetErrorIfNotConverged", snes->errorifnotconverged, &snes->errorifnotconverged, NULL));
1020: PetscCall(PetscOptionsBool("-snes_force_iteration", "Force SNESSolve() to take at least one iteration", "SNESSetForceIteration", snes->forceiteration, &snes->forceiteration, NULL));
1021: PetscCall(PetscOptionsBool("-snes_check_jacobian_domain_error", "Check Jacobian domain error after Jacobian evaluation", "SNESCheckJacobianDomainError", snes->checkjacdomainerror, &snes->checkjacdomainerror, NULL));
1023: PetscCall(PetscOptionsInt("-snes_lag_preconditioner", "How often to rebuild preconditioner", "SNESSetLagPreconditioner", snes->lagpreconditioner, &lag, &flg));
1024: if (flg) {
1025: PetscCheck(lag != -1, PetscObjectComm((PetscObject)snes), PETSC_ERR_USER, "Cannot set the lag to -1 from the command line since the preconditioner must be built as least once, perhaps you mean -2");
1026: PetscCall(SNESSetLagPreconditioner(snes, lag));
1027: }
1028: PetscCall(PetscOptionsBool("-snes_lag_preconditioner_persists", "Preconditioner lagging through multiple SNES solves", "SNESSetLagPreconditionerPersists", snes->lagjac_persist, &persist, &flg));
1029: if (flg) PetscCall(SNESSetLagPreconditionerPersists(snes, persist));
1030: PetscCall(PetscOptionsInt("-snes_lag_jacobian", "How often to rebuild Jacobian", "SNESSetLagJacobian", snes->lagjacobian, &lag, &flg));
1031: if (flg) {
1032: PetscCheck(lag != -1, PetscObjectComm((PetscObject)snes), PETSC_ERR_USER, "Cannot set the lag to -1 from the command line since the Jacobian must be built as least once, perhaps you mean -2");
1033: PetscCall(SNESSetLagJacobian(snes, lag));
1034: }
1035: PetscCall(PetscOptionsBool("-snes_lag_jacobian_persists", "Jacobian lagging through multiple SNES solves", "SNESSetLagJacobianPersists", snes->lagjac_persist, &persist, &flg));
1036: if (flg) PetscCall(SNESSetLagJacobianPersists(snes, persist));
1038: PetscCall(PetscOptionsInt("-snes_grid_sequence", "Use grid sequencing to generate initial guess", "SNESSetGridSequence", snes->gridsequence, &grids, &flg));
1039: if (flg) PetscCall(SNESSetGridSequence(snes, grids));
1041: PetscCall(PetscOptionsEList("-snes_convergence_test", "Convergence test", "SNESSetConvergenceTest", convtests, PETSC_STATIC_ARRAY_LENGTH(convtests), "default", &indx, &flg));
1042: if (flg) {
1043: switch (indx) {
1044: case 0:
1045: PetscCall(SNESSetConvergenceTest(snes, SNESConvergedDefault, NULL, NULL));
1046: break;
1047: case 1:
1048: PetscCall(SNESSetConvergenceTest(snes, SNESConvergedSkip, NULL, NULL));
1049: break;
1050: case 2:
1051: PetscCall(SNESSetConvergenceTest(snes, SNESConvergedCorrectPressure, NULL, NULL));
1052: break;
1053: }
1054: }
1056: PetscCall(PetscOptionsEList("-snes_norm_schedule", "SNES Norm schedule", "SNESSetNormSchedule", SNESNormSchedules, 5, "function", &indx, &flg));
1057: if (flg) PetscCall(SNESSetNormSchedule(snes, (SNESNormSchedule)indx));
1059: PetscCall(PetscOptionsEList("-snes_function_type", "SNES Norm schedule", "SNESSetFunctionType", SNESFunctionTypes, 2, "unpreconditioned", &indx, &flg));
1060: if (flg) PetscCall(SNESSetFunctionType(snes, (SNESFunctionType)indx));
1062: kctx = (SNESKSPEW *)snes->kspconvctx;
1064: PetscCall(PetscOptionsBool("-snes_ksp_ew", "Use Eisentat-Walker linear system convergence test", "SNESKSPSetUseEW", snes->ksp_ewconv, &snes->ksp_ewconv, NULL));
1066: PetscCall(SNESGetOptionsPrefix(snes, &optionsprefix));
1067: PetscCall(PetscSNPrintf(ewprefix, sizeof(ewprefix), "%s%s", optionsprefix ? optionsprefix : "", "snes_"));
1068: PetscCall(SNESEWSetFromOptions_Private(kctx, PETSC_TRUE, PetscObjectComm((PetscObject)snes), ewprefix));
1070: flg = PETSC_FALSE;
1071: PetscCall(PetscOptionsBool("-snes_monitor_cancel", "Remove all monitors", "SNESMonitorCancel", flg, &flg, &set));
1072: if (set && flg) PetscCall(SNESMonitorCancel(snes));
1074: PetscCall(PetscOptionsDeprecated("-snes_monitor_short", "-snes_monitor", "3.26", NULL));
1075: PetscCall(SNESMonitorSetFromOptions(snes, "-snes_monitor", "Monitor norm of function", "SNESMonitorDefault", SNESMonitorDefault, SNESMonitorDefaultSetUp));
1076: PetscCall(SNESMonitorSetFromOptions(snes, "-snes_monitor_range", "Monitor range of elements of function", "SNESMonitorRange", SNESMonitorRange, NULL));
1078: PetscCall(SNESMonitorSetFromOptions(snes, "-snes_monitor_ratio", "Monitor ratios of the norm of function for consecutive steps", "SNESMonitorRatio", SNESMonitorRatio, SNESMonitorRatioSetUp));
1079: PetscCall(SNESMonitorSetFromOptions(snes, "-snes_monitor_field", "Monitor norm of function (split into fields)", "SNESMonitorDefaultField", SNESMonitorDefaultField, NULL));
1080: PetscCall(SNESMonitorSetFromOptions(snes, "-snes_monitor_solution", "View solution at each iteration", "SNESMonitorSolution", SNESMonitorSolution, NULL));
1081: PetscCall(SNESMonitorSetFromOptions(snes, "-snes_monitor_solution_update", "View correction at each iteration", "SNESMonitorSolutionUpdate", SNESMonitorSolutionUpdate, NULL));
1082: PetscCall(SNESMonitorSetFromOptions(snes, "-snes_monitor_residual", "View residual at each iteration", "SNESMonitorResidual", SNESMonitorResidual, NULL));
1083: PetscCall(SNESMonitorSetFromOptions(snes, "-snes_monitor_jacupdate_spectrum", "Print the change in the spectrum of the Jacobian", "SNESMonitorJacUpdateSpectrum", SNESMonitorJacUpdateSpectrum, NULL));
1084: PetscCall(SNESMonitorSetFromOptions(snes, "-snes_monitor_fields", "Monitor norm of function per field", "SNESMonitorSet", SNESMonitorFields, NULL));
1085: PetscCall(PetscOptionsBool("-snes_monitor_pause_final", "Pauses all draw monitors at the final iterate", "SNESMonitorPauseFinal_Internal", PETSC_FALSE, &snes->pauseFinal, NULL));
1087: PetscCall(PetscOptionsString("-snes_monitor_python", "Use Python function", "SNESMonitorSet", NULL, monfilename, sizeof(monfilename), &flg));
1088: if (flg) PetscCall(PetscPythonMonitorSet((PetscObject)snes, monfilename));
1090: flg = PETSC_FALSE;
1091: PetscCall(PetscOptionsBool("-snes_monitor_lg_range", "Plot function range at each iteration", "SNESMonitorLGRange", flg, &flg, NULL));
1092: if (flg) {
1093: PetscViewer ctx;
1095: PetscCall(PetscViewerDrawOpen(PetscObjectComm((PetscObject)snes), NULL, NULL, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &ctx));
1096: PetscCall(SNESMonitorSet(snes, SNESMonitorLGRange, ctx, (PetscCtxDestroyFn *)PetscViewerDestroy));
1097: }
1099: PetscCall(PetscViewerDestroy(&snes->convergedreasonviewer));
1100: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_converged_reason", &snes->convergedreasonviewer, &snes->convergedreasonformat, NULL));
1101: flg = PETSC_FALSE;
1102: PetscCall(PetscOptionsBool("-snes_converged_reason_view_cancel", "Remove all converged reason viewers", "SNESConvergedReasonViewCancel", flg, &flg, &set));
1103: if (set && flg) PetscCall(SNESConvergedReasonViewCancel(snes));
1105: flg = PETSC_FALSE;
1106: PetscCall(PetscOptionsBool("-snes_fd", "Use finite differences (slow) to compute Jacobian", "SNESComputeJacobianDefault", flg, &flg, NULL));
1107: if (flg) {
1108: void *functx;
1109: DM dm;
1110: PetscCall(SNESGetDM(snes, &dm));
1111: PetscCall(DMSNESUnsetJacobianContext_Internal(dm));
1112: PetscCall(SNESGetFunction(snes, NULL, NULL, &functx));
1113: PetscCall(SNESSetJacobian(snes, snes->jacobian, snes->jacobian_pre, SNESComputeJacobianDefault, functx));
1114: PetscCall(PetscInfo(snes, "Setting default finite difference Jacobian matrix\n"));
1115: }
1117: flg = PETSC_FALSE;
1118: PetscCall(PetscOptionsBool("-snes_fd_function", "Use finite differences (slow) to compute function from user objective", "SNESObjectiveComputeFunctionDefaultFD", flg, &flg, NULL));
1119: if (flg) PetscCall(SNESSetFunction(snes, NULL, SNESObjectiveComputeFunctionDefaultFD, NULL));
1121: flg = PETSC_FALSE;
1122: PetscCall(PetscOptionsBool("-snes_fd_color", "Use finite differences with coloring to compute Jacobian", "SNESComputeJacobianDefaultColor", flg, &flg, NULL));
1123: if (flg) {
1124: DM dm;
1125: PetscCall(SNESGetDM(snes, &dm));
1126: PetscCall(DMSNESUnsetJacobianContext_Internal(dm));
1127: PetscCall(SNESSetJacobian(snes, snes->jacobian, snes->jacobian_pre, SNESComputeJacobianDefaultColor, NULL));
1128: PetscCall(PetscInfo(snes, "Setting default finite difference coloring Jacobian matrix\n"));
1129: }
1131: flg = PETSC_FALSE;
1132: PetscCall(PetscOptionsBool("-snes_mf_operator", "Use a Matrix-Free Jacobian with user-provided matrix for computing the preconditioner", "SNESSetUseMatrixFree", PETSC_FALSE, &snes->mf_operator, &flg));
1133: if (flg && snes->mf_operator) {
1134: snes->mf_operator = PETSC_TRUE;
1135: snes->mf = PETSC_TRUE;
1136: }
1137: flg = PETSC_FALSE;
1138: PetscCall(PetscOptionsBool("-snes_mf", "Use a Matrix-Free Jacobian with no preconditioner by default", "SNESSetUseMatrixFree", PETSC_FALSE, &snes->mf, &flg));
1139: if (!flg && snes->mf_operator) snes->mf = PETSC_TRUE;
1140: PetscCall(PetscOptionsInt("-snes_mf_version", "Matrix-Free routines version 1 or 2", "None", snes->mf_version, &snes->mf_version, NULL));
1142: PetscCall(PetscOptionsName("-snes_test_function", "Compare hand-coded and finite difference functions", "None", &snes->testFunc));
1143: PetscCall(PetscOptionsName("-snes_test_jacobian", "Compare hand-coded and finite difference Jacobians", "None", &snes->testJac));
1145: flg = PETSC_FALSE;
1146: PetscCall(SNESGetNPCSide(snes, &pcside));
1147: PetscCall(PetscOptionsEnum("-snes_npc_side", "SNES nonlinear preconditioner side", "SNESSetNPCSide", PCSides, (PetscEnum)pcside, (PetscEnum *)&pcside, &flg));
1148: if (flg) PetscCall(SNESSetNPCSide(snes, pcside));
1150: #if PetscDefined(HAVE_SAWS)
1151: /*
1152: Publish convergence information using SAWs
1153: */
1154: flg = PETSC_FALSE;
1155: PetscCall(PetscOptionsBool("-snes_monitor_saws", "Publish SNES progress using SAWs", "SNESMonitorSet", flg, &flg, NULL));
1156: if (flg) {
1157: PetscCtx ctx;
1158: PetscCall(SNESMonitorSAWsCreate(snes, &ctx));
1159: PetscCall(SNESMonitorSet(snes, SNESMonitorSAWs, ctx, SNESMonitorSAWsDestroy));
1160: }
1161: #endif
1162: #if PetscDefined(HAVE_SAWS)
1163: {
1164: PetscBool set;
1165: flg = PETSC_FALSE;
1166: PetscCall(PetscOptionsBool("-snes_saws_block", "Block for SAWs at end of SNESSolve", "PetscObjectSAWsBlock", ((PetscObject)snes)->amspublishblock, &flg, &set));
1167: if (set) PetscCall(PetscObjectSAWsSetBlock((PetscObject)snes, flg));
1168: }
1169: #endif
1171: for (i = 0; i < numberofsetfromoptions; i++) PetscCall((*othersetfromoptions[i])(snes));
1173: PetscTryTypeMethod(snes, setfromoptions, PetscOptionsObject);
1175: /* process any options handlers added with PetscObjectAddOptionsHandler() */
1176: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)snes, PetscOptionsObject));
1177: PetscOptionsEnd();
1179: if (snes->linesearch) {
1180: PetscCall(SNESGetLineSearch(snes, &snes->linesearch));
1181: PetscCall(SNESLineSearchSetFromOptions(snes->linesearch));
1182: }
1184: /* if user has set the SNES NPC type via options database, create it. */
1185: PetscCall(SNESGetOptionsPrefix(snes, &optionsprefix));
1186: PetscCall(PetscOptionsHasName(((PetscObject)snes)->options, optionsprefix, "-npc_snes_type", &pcset));
1187: if (pcset && !snes->npc) PetscCall(SNESGetNPC(snes, &snes->npc));
1189: if (snes->usesksp) {
1190: PC pc;
1191: PCType pctype;
1193: if (!snes->ksp) PetscCall(SNESGetKSP(snes, &snes->ksp));
1194: PetscCall(KSPSetOperators(snes->ksp, snes->jacobian, snes->jacobian_pre));
1195: PetscCall(KSPGetPC(snes->ksp, &pc));
1196: PetscCall(PCGetType(pc, &pctype));
1197: /* If the first two conditions in the following conditional are true, we know a matrix-free Mat
1198: will be used eventually with the PC, but we cannot provide the matrix-free Mat to the PC here
1199: since we do not have enough information to construct it here (it is constructed after the
1200: start of SNESSetUp()). If we do not set the PCNONE here, then the PCSetFromOptions() called
1201: from KSPSetFromOptions() below will use PCGetDefaultType_Private() to set a PCType
1202: appropriate for the current pc->pmat that will likely not work for the matrix-free Mat, thus
1203: producing a later confusing error message. A significant refactoring of how SNES handles
1204: matrix-free Mat would be needed to eliminate the next line of code. Note that if the PC type
1205: has already been set (third condition), we do not override it. The fourth condition exempts
1206: left-side nonlinear preconditioners, which require a real PC */
1207: if (snes->mf && !snes->mf_operator && !pctype && !(snes->npcside == PC_LEFT && snes->npc)) {
1208: PetscCall(PetscInfo(snes, "Setting PCNONE since no PC type was set and the Jacobian will be matrix-free\n"));
1209: PetscCall(PCSetType(pc, PCNONE));
1210: }
1211: PetscCall(KSPSetFromOptions(snes->ksp));
1212: }
1214: if (snes->npc) PetscCall(SNESSetFromOptions(snes->npc));
1215: snes->setfromoptionscalled++;
1216: PetscFunctionReturn(PETSC_SUCCESS);
1217: }
1219: /*@
1220: SNESResetFromOptions - Sets various `SNES` and `KSP` parameters from user options ONLY if the `SNESSetFromOptions()` was previously called
1222: Collective
1224: Input Parameter:
1225: . snes - the `SNES` context
1227: Level: advanced
1229: .seealso: [](ch_snes), `SNES`, `SNESSetFromOptions()`, `SNESSetOptionsPrefix()`
1230: @*/
1231: PetscErrorCode SNESResetFromOptions(SNES snes)
1232: {
1233: PetscFunctionBegin;
1234: if (snes->setfromoptionscalled) PetscCall(SNESSetFromOptions(snes));
1235: PetscFunctionReturn(PETSC_SUCCESS);
1236: }
1238: /*@
1239: SNESSetComputeApplicationContext - Sets an optional function to compute a user-defined context for
1240: the nonlinear solvers.
1242: Logically Collective; No Fortran Support
1244: Input Parameters:
1245: + snes - the `SNES` context
1246: . compute - function to compute the context
1247: - destroy - function to destroy the context, see `PetscCtxDestroyFn` for the calling sequence
1249: Calling sequence of `compute`:
1250: + snes - the `SNES` context
1251: - ctx - context to be computed
1253: Level: intermediate
1255: Note:
1256: This routine is useful if you are performing grid sequencing or using `SNESFAS` and need the appropriate context generated for each level.
1258: Use `SNESSetApplicationContext()` to see the context immediately
1260: .seealso: [](ch_snes), `SNESGetApplicationContext()`, `SNESSetApplicationContext()`, `PetscCtxDestroyFn`
1261: @*/
1262: PetscErrorCode SNESSetComputeApplicationContext(SNES snes, PetscErrorCode (*compute)(SNES snes, PetscCtxRt ctx), PetscCtxDestroyFn *destroy)
1263: {
1264: PetscFunctionBegin;
1266: snes->ops->ctxcompute = compute;
1267: snes->ops->ctxdestroy = destroy;
1268: PetscFunctionReturn(PETSC_SUCCESS);
1269: }
1271: /*@
1272: SNESSetApplicationContext - Sets the optional user-defined context for the nonlinear solvers.
1274: Logically Collective
1276: Input Parameters:
1277: + snes - the `SNES` context
1278: - ctx - the application context
1280: Level: intermediate
1282: Notes:
1283: Users can provide a context when constructing the `SNES` options and then access it inside their function, Jacobian computation, or other evaluation function
1284: with `SNESGetApplicationContext()`
1286: To provide a function that computes the context for you use `SNESSetComputeApplicationContext()`
1288: Fortran Note:
1289: This only works when `ctx` is a Fortran derived type (it cannot be a `PetscObject`), we recommend writing a Fortran interface definition for this
1290: function that tells the Fortran compiler the derived data type that is passed in as the `ctx` argument. See `SNESGetApplicationContext()` for
1291: an example.
1293: .seealso: [](ch_snes), `SNES`, `SNESSetComputeApplicationContext()`, `SNESGetApplicationContext()`
1294: @*/
1295: PetscErrorCode SNESSetApplicationContext(SNES snes, PetscCtx ctx)
1296: {
1297: KSP ksp;
1299: PetscFunctionBegin;
1301: PetscCall(SNESGetKSP(snes, &ksp));
1302: PetscCall(KSPSetApplicationContext(ksp, ctx));
1303: snes->ctx = ctx;
1304: PetscFunctionReturn(PETSC_SUCCESS);
1305: }
1307: /*@
1308: SNESGetApplicationContext - Gets the user-defined context for the
1309: nonlinear solvers set with `SNESGetApplicationContext()` or `SNESSetComputeApplicationContext()`
1311: Not Collective
1313: Input Parameter:
1314: . snes - `SNES` context
1316: Output Parameter:
1317: . ctx - the application context
1319: Level: intermediate
1321: Fortran Notes:
1322: This only works when the context is a Fortran derived type or a `PetscObject`. Declare `ctx` with
1323: .vb
1324: type(tUsertype), pointer :: ctx
1325: .ve
1327: .seealso: [](ch_snes), `SNESSetApplicationContext()`, `SNESSetComputeApplicationContext()`
1328: @*/
1329: PetscErrorCode SNESGetApplicationContext(SNES snes, PetscCtxRt ctx)
1330: {
1331: PetscFunctionBegin;
1333: *(void **)ctx = snes->ctx;
1334: PetscFunctionReturn(PETSC_SUCCESS);
1335: }
1337: /*@
1338: SNESSetUseMatrixFree - indicates that `SNES` should use matrix-free finite difference matrix-vector products to apply the Jacobian.
1340: Logically Collective
1342: Input Parameters:
1343: + snes - `SNES` context
1344: . mf_operator - use matrix-free only for the Amat used by `SNESSetJacobian()`, this means the user provided Pmat will continue to be used
1345: - mf - use matrix-free for both the Amat and Pmat used by `SNESSetJacobian()`, both the Amat and Pmat set in `SNESSetJacobian()` will be ignored. With
1346: this option no matrix-element based preconditioners can be used in the linear solve since the matrix won't be explicitly available
1348: Options Database Keys:
1349: + -snes_mf_operator - use matrix-free only for the mat operator
1350: . -snes_mf - use matrix-free for both the mat and pmat operator
1351: . -snes_fd_color - compute the Jacobian via coloring and finite differences.
1352: - -snes_fd - compute the Jacobian via finite differences (slow)
1354: Level: intermediate
1356: Notes:
1357: `SNES` supports three approaches for computing (approximate) Jacobians: user provided via `SNESSetJacobian()`, matrix-free using `MatCreateSNESMF()`,
1358: and computing explicitly with
1359: finite differences and coloring using `MatFDColoring`. It is also possible to use automatic differentiation and the `MatFDColoring` object.
1361: When `mf` is used, `SNESSetFromOptions()` sets the `KSP`'s `PC` to `PCNONE` unless a `PC` type has already been selected.
1363: .seealso: [](ch_snes), `SNES`, `SNESGetUseMatrixFree()`, `MatCreateSNESMF()`, `SNESComputeJacobianDefaultColor()`, `MatFDColoring`
1364: @*/
1365: PetscErrorCode SNESSetUseMatrixFree(SNES snes, PetscBool mf_operator, PetscBool mf)
1366: {
1367: PetscFunctionBegin;
1371: snes->mf = mf_operator ? PETSC_TRUE : mf;
1372: snes->mf_operator = mf_operator;
1373: PetscFunctionReturn(PETSC_SUCCESS);
1374: }
1376: /*@
1377: SNESGetUseMatrixFree - indicates if the `SNES` uses matrix-free finite difference matrix vector products to apply the Jacobian.
1379: Not Collective, but the resulting flags will be the same on all MPI processes
1381: Input Parameter:
1382: . snes - `SNES` context
1384: Output Parameters:
1385: + mf_operator - use matrix-free only for the Amat used by `SNESSetJacobian()`, this means the user provided Pmat will continue to be used
1386: - mf - use matrix-free for both the Amat and Pmat used by `SNESSetJacobian()`, both the Amat and Pmat set in `SNESSetJacobian()` will be ignored
1388: Level: intermediate
1390: .seealso: [](ch_snes), `SNES`, `SNESSetUseMatrixFree()`, `MatCreateSNESMF()`
1391: @*/
1392: PetscErrorCode SNESGetUseMatrixFree(SNES snes, PetscBool *mf_operator, PetscBool *mf)
1393: {
1394: PetscFunctionBegin;
1396: if (mf) *mf = snes->mf;
1397: if (mf_operator) *mf_operator = snes->mf_operator;
1398: PetscFunctionReturn(PETSC_SUCCESS);
1399: }
1401: /*@
1402: SNESGetIterationNumber - Gets the number of nonlinear iterations completed in the current or most recent `SNESSolve()`
1404: Not Collective
1406: Input Parameter:
1407: . snes - `SNES` context
1409: Output Parameter:
1410: . iter - iteration number
1412: Level: intermediate
1414: Notes:
1415: For example, during the computation of iteration 2 this would return 1.
1417: This is useful for using lagged Jacobians (where one does not recompute the
1418: Jacobian at each `SNES` iteration). For example, the code
1419: .vb
1420: ierr = SNESGetIterationNumber(snes,&it);
1421: if (!(it % 2)) {
1422: [compute Jacobian here]
1423: }
1424: .ve
1425: can be used in your function that computes the Jacobian to cause the Jacobian to be
1426: recomputed every second `SNES` iteration. See also `SNESSetLagJacobian()`
1428: After the `SNES` solve is complete this will return the number of nonlinear iterations used.
1430: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESSetLagJacobian()`, `SNESGetLinearSolveIterations()`, `SNESSetMonitor()`
1431: @*/
1432: PetscErrorCode SNESGetIterationNumber(SNES snes, PetscInt *iter)
1433: {
1434: PetscFunctionBegin;
1436: PetscAssertPointer(iter, 2);
1437: *iter = snes->iter;
1438: PetscFunctionReturn(PETSC_SUCCESS);
1439: }
1441: /*@
1442: SNESSetIterationNumber - Sets the current iteration number.
1444: Not Collective
1446: Input Parameters:
1447: + snes - `SNES` context
1448: - iter - iteration number
1450: Level: developer
1452: Note:
1453: This should only be called inside a `SNES` nonlinear solver.
1455: .seealso: [](ch_snes), `SNESGetLinearSolveIterations()`
1456: @*/
1457: PetscErrorCode SNESSetIterationNumber(SNES snes, PetscInt iter)
1458: {
1459: PetscFunctionBegin;
1461: PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
1462: snes->iter = iter;
1463: PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
1464: PetscFunctionReturn(PETSC_SUCCESS);
1465: }
1467: /*@
1468: SNESGetNonlinearStepFailures - Gets the number of unsuccessful steps
1469: taken by the nonlinear solver in the current or most recent `SNESSolve()` .
1471: Not Collective
1473: Input Parameter:
1474: . snes - `SNES` context
1476: Output Parameter:
1477: . nfails - number of unsuccessful steps attempted
1479: Level: intermediate
1481: Notes:
1482: A failed step is a step that was generated and taken but did not satisfy the requested step criteria. For example,
1483: the `SNESLineSearchApply()` could not generate a sufficient decrease in the function norm (in fact it may have produced an increase).
1485: Taken steps that produce a infinity or NaN in the function evaluation or generate a `SNESSetFunctionDomainError()`
1486: will always immediately terminate the `SNESSolve()` regardless of the value of `maxFails`.
1488: `SNESSetMaxNonlinearStepFailures()` determines how many unsuccessful steps are allowed before the `SNESSolve()` terminates
1490: This counter is reset to zero for each successive call to `SNESSolve()`.
1492: .seealso: [](ch_snes), `SNES`, `SNESGetMaxLinearSolveFailures()`, `SNESGetLinearSolveIterations()`, `SNESSetMaxLinearSolveFailures()`, `SNESGetLinearSolveFailures()`,
1493: `SNESSetMaxNonlinearStepFailures()`, `SNESGetMaxNonlinearStepFailures()`
1494: @*/
1495: PetscErrorCode SNESGetNonlinearStepFailures(SNES snes, PetscInt *nfails)
1496: {
1497: PetscFunctionBegin;
1499: PetscAssertPointer(nfails, 2);
1500: *nfails = snes->numFailures;
1501: PetscFunctionReturn(PETSC_SUCCESS);
1502: }
1504: /*@
1505: SNESSetMaxNonlinearStepFailures - Sets the maximum number of unsuccessful steps
1506: attempted by the nonlinear solver before it gives up and returns unconverged or generates an error
1508: Not Collective
1510: Input Parameters:
1511: + snes - `SNES` context
1512: - maxFails - maximum of unsuccessful steps allowed, use `PETSC_UNLIMITED` to have no limit on the number of failures
1514: Options Database Key:
1515: . -snes_max_fail n - maximum number of unsuccessful steps allowed
1517: Level: intermediate
1519: Note:
1520: A failed step is a step that was generated and taken but did not satisfy the requested criteria. For example,
1521: the `SNESLineSearchApply()` could not generate a sufficient decrease in the function norm (in fact it may have produced an increase).
1523: Taken steps that produce a infinity or NaN in the function evaluation or generate a `SNESSetFunctionDomainError()`
1524: will always immediately terminate the `SNESSolve()` regardless of the value of `maxFails`.
1526: Developer Note:
1527: The options database key is wrong for this function name
1529: .seealso: [](ch_snes), `SNESSetErrorIfNotConverged()`, `SNESGetMaxLinearSolveFailures()`, `SNESGetLinearSolveIterations()`, `SNESSetMaxLinearSolveFailures()`,
1530: `SNESGetLinearSolveFailures()`, `SNESGetMaxNonlinearStepFailures()`, `SNESGetNonlinearStepFailures()`, `SNESCheckLineSearchFailure()`
1531: @*/
1532: PetscErrorCode SNESSetMaxNonlinearStepFailures(SNES snes, PetscInt maxFails)
1533: {
1534: PetscFunctionBegin;
1537: if (maxFails == PETSC_UNLIMITED) {
1538: snes->maxFailures = PETSC_INT_MAX;
1539: } else {
1540: PetscCheck(maxFails >= 0, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_OUTOFRANGE, "Cannot have a negative maximum number of failures");
1541: snes->maxFailures = maxFails;
1542: }
1543: PetscFunctionReturn(PETSC_SUCCESS);
1544: }
1546: /*@
1547: SNESGetMaxNonlinearStepFailures - Gets the maximum number of unsuccessful steps
1548: attempted by the nonlinear solver before it gives up and returns unconverged or generates an error
1550: Not Collective
1552: Input Parameter:
1553: . snes - `SNES` context
1555: Output Parameter:
1556: . maxFails - maximum of unsuccessful steps
1558: Level: intermediate
1560: .seealso: [](ch_snes), `SNESSetErrorIfNotConverged()`, `SNESGetMaxLinearSolveFailures()`, `SNESGetLinearSolveIterations()`, `SNESSetMaxLinearSolveFailures()`, `SNESGetLinearSolveFailures()`,
1561: `SNESSetMaxNonlinearStepFailures()`, `SNESGetNonlinearStepFailures()`
1562: @*/
1563: PetscErrorCode SNESGetMaxNonlinearStepFailures(SNES snes, PetscInt *maxFails)
1564: {
1565: PetscFunctionBegin;
1567: PetscAssertPointer(maxFails, 2);
1568: *maxFails = snes->maxFailures;
1569: PetscFunctionReturn(PETSC_SUCCESS);
1570: }
1572: /*@
1573: SNESGetNumberFunctionEvals - Gets the number of user provided function evaluations
1574: done by the `SNES` object in the current or most recent `SNESSolve()`
1576: Not Collective
1578: Input Parameter:
1579: . snes - `SNES` context
1581: Output Parameter:
1582: . nfuncs - number of evaluations
1584: Level: intermediate
1586: Note:
1587: Reset every time `SNESSolve()` is called unless `SNESSetCountersReset()` is used.
1589: .seealso: [](ch_snes), `SNES`, `SNESGetMaxLinearSolveFailures()`, `SNESGetLinearSolveIterations()`, `SNESSetMaxLinearSolveFailures()`, `SNESGetLinearSolveFailures()`, `SNESSetCountersReset()`
1590: @*/
1591: PetscErrorCode SNESGetNumberFunctionEvals(SNES snes, PetscInt *nfuncs)
1592: {
1593: PetscFunctionBegin;
1595: PetscAssertPointer(nfuncs, 2);
1596: *nfuncs = snes->nfuncs;
1597: PetscFunctionReturn(PETSC_SUCCESS);
1598: }
1600: /*@
1601: SNESGetLinearSolveFailures - Gets the number of failed (non-converged)
1602: linear solvers in the current or most recent `SNESSolve()`
1604: Not Collective
1606: Input Parameter:
1607: . snes - `SNES` context
1609: Output Parameter:
1610: . nfails - number of failed solves
1612: Options Database Key:
1613: . -snes_max_linear_solve_fail num - The number of failures before the solve is terminated
1615: Level: intermediate
1617: Note:
1618: This counter is reset to zero for each successive call to `SNESSolve()`.
1620: .seealso: [](ch_snes), `SNESGetMaxLinearSolveFailures()`, `SNESGetLinearSolveIterations()`, `SNESSetMaxLinearSolveFailures()`
1621: @*/
1622: PetscErrorCode SNESGetLinearSolveFailures(SNES snes, PetscInt *nfails)
1623: {
1624: PetscFunctionBegin;
1626: PetscAssertPointer(nfails, 2);
1627: *nfails = snes->numLinearSolveFailures;
1628: PetscFunctionReturn(PETSC_SUCCESS);
1629: }
1631: /*@
1632: SNESSetMaxLinearSolveFailures - the number of failed linear solve attempts
1633: allowed before `SNES` returns with a diverged reason of `SNES_DIVERGED_LINEAR_SOLVE`
1635: Logically Collective
1637: Input Parameters:
1638: + snes - `SNES` context
1639: - maxFails - maximum allowed linear solve failures, use `PETSC_UNLIMITED` to have no limit on the number of failures
1641: Options Database Key:
1642: . -snes_max_linear_solve_fail num - The number of failures before the solve is terminated
1644: Level: intermediate
1646: Note:
1647: By default this is 0; that is `SNES` returns on the first failed linear solve
1649: Developer Note:
1650: The options database key is wrong for this function name
1652: .seealso: [](ch_snes), `SNESSetErrorIfNotConverged()`, `SNESGetLinearSolveFailures()`, `SNESGetMaxLinearSolveFailures()`, `SNESGetLinearSolveIterations()`
1653: @*/
1654: PetscErrorCode SNESSetMaxLinearSolveFailures(SNES snes, PetscInt maxFails)
1655: {
1656: PetscFunctionBegin;
1660: if (maxFails == PETSC_UNLIMITED) {
1661: snes->maxLinearSolveFailures = PETSC_INT_MAX;
1662: } else {
1663: PetscCheck(maxFails >= 0, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_OUTOFRANGE, "Cannot have a negative maximum number of failures");
1664: snes->maxLinearSolveFailures = maxFails;
1665: }
1666: PetscFunctionReturn(PETSC_SUCCESS);
1667: }
1669: /*@
1670: SNESGetMaxLinearSolveFailures - gets the maximum number of linear solve failures that
1671: are allowed before `SNES` returns as unsuccessful
1673: Not Collective
1675: Input Parameter:
1676: . snes - `SNES` context
1678: Output Parameter:
1679: . maxFails - maximum of unsuccessful solves allowed
1681: Level: intermediate
1683: Note:
1684: By default this is 1; that is `SNES` returns on the first failed linear solve
1686: .seealso: [](ch_snes), `SNESSetErrorIfNotConverged()`, `SNESGetLinearSolveFailures()`, `SNESGetLinearSolveIterations()`, `SNESSetMaxLinearSolveFailures()`
1687: @*/
1688: PetscErrorCode SNESGetMaxLinearSolveFailures(SNES snes, PetscInt *maxFails)
1689: {
1690: PetscFunctionBegin;
1692: PetscAssertPointer(maxFails, 2);
1693: *maxFails = snes->maxLinearSolveFailures;
1694: PetscFunctionReturn(PETSC_SUCCESS);
1695: }
1697: /*@
1698: SNESGetLinearSolveIterations - Gets the total number of linear iterations
1699: used by the nonlinear solver in the most recent `SNESSolve()`
1701: Not Collective
1703: Input Parameter:
1704: . snes - `SNES` context
1706: Output Parameter:
1707: . lits - number of linear iterations
1709: Level: intermediate
1711: Notes:
1712: This counter is reset to zero for each successive call to `SNESSolve()` unless `SNESSetCountersReset()` is used.
1714: If the linear solver fails inside the `SNESSolve()` the iterations for that call to the linear solver are not included. If you wish to count them
1715: then call `KSPGetIterationNumber()` after the failed solve.
1717: .seealso: [](ch_snes), `SNES`, `SNESGetIterationNumber()`, `SNESGetLinearSolveFailures()`, `SNESGetMaxLinearSolveFailures()`, `SNESSetCountersReset()`
1718: @*/
1719: PetscErrorCode SNESGetLinearSolveIterations(SNES snes, PetscInt *lits)
1720: {
1721: PetscFunctionBegin;
1723: PetscAssertPointer(lits, 2);
1724: *lits = snes->linear_its;
1725: PetscFunctionReturn(PETSC_SUCCESS);
1726: }
1728: /*@
1729: SNESSetCountersReset - Sets whether or not the counters for linear iterations and function evaluations
1730: are reset every time `SNESSolve()` is called.
1732: Logically Collective
1734: Input Parameters:
1735: + snes - `SNES` context
1736: - reset - whether to reset the counters or not, defaults to `PETSC_TRUE`
1738: Level: developer
1740: .seealso: [](ch_snes), `SNESGetNumberFunctionEvals()`, `SNESGetLinearSolveIterations()`, `SNESGetNPC()`
1741: @*/
1742: PetscErrorCode SNESSetCountersReset(SNES snes, PetscBool reset)
1743: {
1744: PetscFunctionBegin;
1747: snes->counters_reset = reset;
1748: PetscFunctionReturn(PETSC_SUCCESS);
1749: }
1751: /*@
1752: SNESResetCounters - Reset counters for linear iterations and function evaluations.
1754: Logically Collective
1756: Input Parameters:
1757: . snes - `SNES` context
1759: Level: developer
1761: Note:
1762: It honors the flag set with `SNESSetCountersReset()`
1764: .seealso: [](ch_snes), `SNESGetNumberFunctionEvals()`, `SNESGetLinearSolveIterations()`, `SNESGetNPC()`
1765: @*/
1766: PetscErrorCode SNESResetCounters(SNES snes)
1767: {
1768: PetscFunctionBegin;
1770: if (snes->counters_reset) {
1771: snes->nfuncs = 0;
1772: snes->linear_its = 0;
1773: snes->numFailures = 0;
1774: }
1775: PetscFunctionReturn(PETSC_SUCCESS);
1776: }
1778: /*@
1779: SNESSetKSP - Sets a `KSP` context for the `SNES` object to use
1781: Not Collective, but the `SNES` and `KSP` objects must live on the same `MPI_Comm`
1783: Input Parameters:
1784: + snes - the `SNES` context
1785: - ksp - the `KSP` context
1787: Level: developer
1789: Notes:
1790: The `SNES` object already has its `KSP` object, you can obtain with `SNESGetKSP()`
1791: so this routine is rarely needed.
1793: The `KSP` object that is already in the `SNES` object has its reference count
1794: decreased by one when this is called.
1796: .seealso: [](ch_snes), `SNES`, `KSP`, `KSPGetPC()`, `SNESCreate()`, `KSPCreate()`
1797: @*/
1798: PetscErrorCode SNESSetKSP(SNES snes, KSP ksp)
1799: {
1800: PetscFunctionBegin;
1803: PetscCheckSameComm(snes, 1, ksp, 2);
1804: PetscCall(PetscObjectReference((PetscObject)ksp));
1805: PetscCall(PetscObjectDereference((PetscObject)snes->ksp));
1806: snes->ksp = ksp;
1807: PetscFunctionReturn(PETSC_SUCCESS);
1808: }
1810: /*@
1811: SNESParametersInitialize - Sets the base defaults for parameters in `snes`, updating a parameter's current value when it matches its previously recorded default.
1813: Logically collective
1815: Input Parameter:
1816: . snes - the `SNES` object
1818: Level: developer
1820: Notes:
1822: The base defaults are the non-type-specific values established when the `SNES` is created. A `SNESType` constructor may subsequently replace them with type-specific defaults.
1824: Developer Notes:
1826: `SNESCreate()` calls this routine to establish the base defaults. `SNESSetType()` calls it before constructing a new `SNESType`, so the recorded defaults associated with the previous type are replaced before the new type installs its own defaults.
1828: Default tracking is based on value equality, not on whether a setter was called. Consequently, an explicitly assigned value that equals the recorded default may be updated when the type changes.
1830: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESDestroy()`, `SNESSetLagPreconditioner()`, `SNESSetLagJacobian()`,
1831: `PetscObjectParameterSetDefault()`
1832: @*/
1833: PetscErrorCode SNESParametersInitialize(SNES snes)
1834: {
1835: PetscObjectParameterSetDefault(snes, max_its, 50);
1836: PetscObjectParameterSetDefault(snes, max_funcs, 10000);
1837: PetscObjectParameterSetDefault(snes, rtol, PetscDefined(USE_REAL_SINGLE) ? 1.e-5 : 1.e-8);
1838: PetscObjectParameterSetDefault(snes, abstol, PetscDefined(USE_REAL_SINGLE) ? 1.e-25 : 1.e-50);
1839: PetscObjectParameterSetDefault(snes, stol, PetscDefined(USE_REAL_SINGLE) ? 1.e-5 : 1.e-8);
1840: PetscObjectParameterSetDefault(snes, divtol, 1.e4);
1841: return PETSC_SUCCESS;
1842: }
1844: /*@
1845: SNESCreate - Creates a nonlinear solver context used to manage a set of nonlinear solves
1847: Collective
1849: Input Parameter:
1850: . comm - MPI communicator
1852: Output Parameter:
1853: . outsnes - the new `SNES` context
1855: Options Database Keys:
1856: + -snes_mf - Activates default matrix-free Jacobian-vector products, with no preconditioner by default
1857: . -snes_mf_operator - Activates default matrix-free Jacobian-vector products, and a user-provided matrix as set by `SNESSetJacobian()`
1858: . -snes_fd_coloring - uses a relative fast computation of the Jacobian using finite differences and a graph coloring
1859: - -snes_fd - Uses (slow!) finite differences to compute Jacobian
1861: Level: beginner
1863: Developer Notes:
1864: `SNES` always creates a `KSP` object even though many `SNES` methods do not use it. This is
1865: unfortunate and should be fixed at some point. The flag snes->usesksp indicates if the
1866: particular method does use `KSP` and regulates if the information about the `KSP` is printed
1867: in `SNESView()`.
1869: `TSSetFromOptions()` does call `SNESSetFromOptions()` which can lead to users being confused
1870: by help messages about meaningless `SNES` options.
1872: `SNES` always creates the `snes->kspconvctx` even though it is used by only one type. This should be fixed.
1874: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESDestroy()`, `SNESSetLagPreconditioner()`, `SNESSetLagJacobian()`
1875: @*/
1876: PetscErrorCode SNESCreate(MPI_Comm comm, SNES *outsnes)
1877: {
1878: SNES snes;
1879: SNESKSPEW *kctx;
1881: PetscFunctionBegin;
1882: PetscAssertPointer(outsnes, 2);
1883: PetscCall(SNESInitializePackage());
1885: PetscCall(PetscHeaderCreate(snes, SNES_CLASSID, "SNES", "Nonlinear solver", "SNES", comm, SNESDestroy, SNESView));
1886: snes->ops->converged = SNESConvergedDefault;
1887: snes->usesksp = PETSC_TRUE;
1888: snes->norm = 0.0;
1889: snes->xnorm = 0.0;
1890: snes->ynorm = 0.0;
1891: snes->normschedule = SNES_NORM_ALWAYS;
1892: snes->functype = SNES_FUNCTION_DEFAULT;
1893: snes->ttol = 0.0;
1895: snes->rnorm0 = 0;
1896: snes->nfuncs = 0;
1897: snes->numFailures = 0;
1898: snes->maxFailures = 1;
1899: snes->linear_its = 0;
1900: snes->lagjacobian = 1;
1901: snes->jac_iter = 0;
1902: snes->lagjac_persist = PETSC_FALSE;
1903: snes->lagpreconditioner = 1;
1904: snes->pre_iter = 0;
1905: snes->lagpre_persist = PETSC_FALSE;
1906: snes->numbermonitors = 0;
1907: snes->numberreasonviews = 0;
1908: snes->data = NULL;
1909: snes->setupcalled = PETSC_FALSE;
1910: snes->ksp_ewconv = PETSC_FALSE;
1911: snes->nwork = 0;
1912: snes->work = NULL;
1913: snes->nvwork = 0;
1914: snes->vwork = NULL;
1915: snes->conv_hist_len = 0;
1916: snes->conv_hist_max = 0;
1917: snes->conv_hist = NULL;
1918: snes->conv_hist_its = NULL;
1919: snes->conv_hist_reset = PETSC_TRUE;
1920: snes->counters_reset = PETSC_TRUE;
1921: snes->vec_func_init_set = PETSC_FALSE;
1922: snes->reason = SNES_CONVERGED_ITERATING;
1923: snes->npcside = PC_RIGHT;
1924: snes->setfromoptionscalled = 0;
1926: snes->mf = PETSC_FALSE;
1927: snes->mf_operator = PETSC_FALSE;
1928: snes->mf_version = 1;
1930: snes->numLinearSolveFailures = 0;
1931: snes->maxLinearSolveFailures = 1;
1933: snes->vizerotolerance = 1.e-8;
1934: snes->checkjacdomainerror = PetscDefined(USE_DEBUG) ? PETSC_TRUE : PETSC_FALSE;
1936: /* Set this to true if the implementation of SNESSolve_XXX does compute the residual at the final solution. */
1937: snes->alwayscomputesfinalresidual = PETSC_FALSE;
1939: /* Create context to compute Eisenstat-Walker relative tolerance for KSP */
1940: PetscCall(PetscNew(&kctx));
1942: snes->kspconvctx = kctx;
1943: kctx->version = 2;
1944: kctx->rtol_0 = 0.3; /* Eisenstat and Walker suggest rtol_0=.5, but
1945: this was too large for some test cases */
1946: kctx->rtol_last = 0.0;
1947: kctx->rtol_max = 0.9;
1948: kctx->gamma = 1.0;
1949: kctx->alpha = 0.5 * (1.0 + PetscSqrtReal(5.0));
1950: kctx->alpha2 = kctx->alpha;
1951: kctx->threshold = 0.1;
1952: kctx->lresid_last = 0.0;
1953: kctx->norm_last = 0.0;
1955: kctx->rk_last = 0.0;
1956: kctx->rk_last_2 = 0.0;
1957: kctx->rtol_last_2 = 0.0;
1958: kctx->v4_p1 = 0.1;
1959: kctx->v4_p2 = 0.4;
1960: kctx->v4_p3 = 0.7;
1961: kctx->v4_m1 = 0.8;
1962: kctx->v4_m2 = 0.5;
1963: kctx->v4_m3 = 0.1;
1964: kctx->v4_m4 = 0.5;
1966: PetscCall(SNESParametersInitialize(snes));
1967: *outsnes = snes;
1968: PetscFunctionReturn(PETSC_SUCCESS);
1969: }
1971: /*@
1972: SNESSetFunction - Sets the function evaluation routine and function
1973: vector for use by the `SNES` routines in solving systems of nonlinear
1974: equations.
1976: Logically Collective
1978: Input Parameters:
1979: + snes - the `SNES` context
1980: . r - vector to store function values, may be `NULL`
1981: . f - function evaluation routine; for calling sequence see `SNESFunctionFn`
1982: - ctx - [optional] user-defined context for private data for the
1983: function evaluation routine (may be `NULL`)
1985: Level: beginner
1987: .seealso: [](ch_snes), `SNES`, `SNESGetFunction()`, `SNESComputeFunction()`, `SNESSetJacobian()`, `SNESSetPicard()`, `SNESFunctionFn`
1988: @*/
1989: PetscErrorCode SNESSetFunction(SNES snes, Vec r, SNESFunctionFn *f, PetscCtx ctx)
1990: {
1991: DM dm;
1993: PetscFunctionBegin;
1995: if (r) {
1997: PetscCheckSameComm(snes, 1, r, 2);
1998: PetscCall(PetscObjectReference((PetscObject)r));
1999: PetscCall(VecDestroy(&snes->vec_func));
2000: snes->vec_func = r;
2001: }
2002: /* update DMSNES
2003: We support incremental information; so update the function context only if r is not specified
2004: (which allows to disable the callbacks when both f and ctx are NULL),
2005: or, if r is specified, when at least one of f and ctx is not NULL */
2006: PetscCall(SNESGetDM(snes, &dm));
2007: if (!r || f || ctx) PetscCall(DMSNESSetFunction(dm, f, ctx));
2008: if (f == SNESPicardComputeFunction) PetscCall(DMSNESSetMFFunction(dm, SNESPicardComputeMFFunction, ctx));
2009: PetscFunctionReturn(PETSC_SUCCESS);
2010: }
2012: /*@
2013: SNESSetInitialFunction - Set an already computed function evaluation at the initial guess to be reused by `SNESSolve()`.
2015: Logically Collective
2017: Input Parameters:
2018: + snes - the `SNES` context
2019: - f - vector to store function value
2021: Level: developer
2023: Notes:
2024: This should not be modified during the solution procedure.
2026: This is used extensively in the `SNESFAS` hierarchy and in nonlinear preconditioning.
2028: .seealso: [](ch_snes), `SNES`, `SNESFAS`, `SNESSetFunction()`, `SNESComputeFunction()`, `SNESSetInitialFunctionNorm()`
2029: @*/
2030: PetscErrorCode SNESSetInitialFunction(SNES snes, Vec f)
2031: {
2032: Vec vec_func;
2034: PetscFunctionBegin;
2037: PetscCheckSameComm(snes, 1, f, 2);
2038: if (snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_PRECONDITIONED) {
2039: snes->vec_func_init_set = PETSC_FALSE;
2040: PetscFunctionReturn(PETSC_SUCCESS);
2041: }
2042: PetscCall(SNESGetFunction(snes, &vec_func, NULL, NULL));
2043: PetscCall(VecCopy(f, vec_func));
2045: snes->vec_func_init_set = PETSC_TRUE;
2046: PetscFunctionReturn(PETSC_SUCCESS);
2047: }
2049: /*@
2050: SNESSetNormSchedule - Sets the `SNESNormSchedule` used in convergence and monitoring
2051: of the `SNES` method, when norms are computed in the solving process
2053: Logically Collective
2055: Input Parameters:
2056: + snes - the `SNES` context
2057: - normschedule - the frequency of norm computation
2059: Options Database Key:
2060: . -snes_norm_schedule (none|always|initialonly|finalonly|initialfinalonly) - set the schedule
2062: Level: advanced
2064: Notes:
2065: Only certain `SNES` methods support certain `SNESNormSchedules`. Most require evaluation
2066: of the nonlinear function and the taking of its norm at every iteration to
2067: even ensure convergence at all. However, methods such as custom Gauss-Seidel methods
2068: `SNESNGS` and the like do not require the norm of the function to be computed, and therefore
2069: may either be monitored for convergence or not. As these are often used as nonlinear
2070: preconditioners, monitoring the norm of their error is not a useful enterprise within
2071: their solution.
2073: .seealso: [](ch_snes), `SNESNormSchedule`, `SNESGetNormSchedule()`, `SNESComputeFunction()`, `VecNorm()`, `SNESSetFunction()`, `SNESSetInitialFunction()`
2074: @*/
2075: PetscErrorCode SNESSetNormSchedule(SNES snes, SNESNormSchedule normschedule)
2076: {
2077: PetscFunctionBegin;
2080: snes->normschedule = normschedule;
2081: PetscFunctionReturn(PETSC_SUCCESS);
2082: }
2084: /*@
2085: SNESGetNormSchedule - Gets the `SNESNormSchedule` used in convergence and monitoring
2086: of the `SNES` method.
2088: Logically Collective
2090: Input Parameters:
2091: + snes - the `SNES` context
2092: - normschedule - the type of the norm used
2094: Level: advanced
2096: .seealso: [](ch_snes), `SNES`, `SNESSetNormSchedule()`, `SNESComputeFunction()`, `VecNorm()`, `SNESSetFunction()`, `SNESSetInitialFunction()`, `SNESNormSchedule`
2097: @*/
2098: PetscErrorCode SNESGetNormSchedule(SNES snes, SNESNormSchedule *normschedule)
2099: {
2100: PetscFunctionBegin;
2102: *normschedule = snes->normschedule;
2103: PetscFunctionReturn(PETSC_SUCCESS);
2104: }
2106: /*@
2107: SNESSetFunctionNorm - Sets the last computed residual norm.
2109: Logically Collective
2111: Input Parameters:
2112: + snes - the `SNES` context
2113: - norm - the value of the norm
2115: Level: developer
2117: .seealso: [](ch_snes), `SNES`, `SNESGetNormSchedule()`, `SNESComputeFunction()`, `VecNorm()`, `SNESSetFunction()`, `SNESSetInitialFunction()`, `SNESNormSchedule`
2118: @*/
2119: PetscErrorCode SNESSetFunctionNorm(SNES snes, PetscReal norm)
2120: {
2121: PetscFunctionBegin;
2123: snes->norm = norm;
2124: PetscFunctionReturn(PETSC_SUCCESS);
2125: }
2127: /*@
2128: SNESGetFunctionNorm - Gets the last computed norm of the residual
2130: Not Collective
2132: Input Parameter:
2133: . snes - the `SNES` context
2135: Output Parameter:
2136: . norm - the last computed residual norm
2138: Level: developer
2140: .seealso: [](ch_snes), `SNES`, `SNESSetNormSchedule()`, `SNESComputeFunction()`, `VecNorm()`, `SNESSetFunction()`, `SNESSetInitialFunction()`, `SNESNormSchedule`
2141: @*/
2142: PetscErrorCode SNESGetFunctionNorm(SNES snes, PetscReal *norm)
2143: {
2144: PetscFunctionBegin;
2146: PetscAssertPointer(norm, 2);
2147: *norm = snes->norm;
2148: PetscFunctionReturn(PETSC_SUCCESS);
2149: }
2151: /*@
2152: SNESGetUpdateNorm - Gets the last computed norm of the solution update
2154: Not Collective
2156: Input Parameter:
2157: . snes - the `SNES` context
2159: Output Parameter:
2160: . ynorm - the last computed update norm
2162: Level: developer
2164: Note:
2165: The new solution is the current solution plus the update, so this norm is an indication of the size of the update
2167: .seealso: [](ch_snes), `SNES`, `SNESSetNormSchedule()`, `SNESComputeFunction()`, `SNESGetFunctionNorm()`
2168: @*/
2169: PetscErrorCode SNESGetUpdateNorm(SNES snes, PetscReal *ynorm)
2170: {
2171: PetscFunctionBegin;
2173: PetscAssertPointer(ynorm, 2);
2174: *ynorm = snes->ynorm;
2175: PetscFunctionReturn(PETSC_SUCCESS);
2176: }
2178: /*@
2179: SNESGetSolutionNorm - Gets the last computed norm of the solution
2181: Not Collective
2183: Input Parameter:
2184: . snes - the `SNES` context
2186: Output Parameter:
2187: . xnorm - the last computed solution norm
2189: Level: developer
2191: .seealso: [](ch_snes), `SNES`, `SNESSetNormSchedule()`, `SNESComputeFunction()`, `SNESGetFunctionNorm()`, `SNESGetUpdateNorm()`
2192: @*/
2193: PetscErrorCode SNESGetSolutionNorm(SNES snes, PetscReal *xnorm)
2194: {
2195: PetscFunctionBegin;
2197: PetscAssertPointer(xnorm, 2);
2198: *xnorm = snes->xnorm;
2199: PetscFunctionReturn(PETSC_SUCCESS);
2200: }
2202: /*@
2203: SNESSetFunctionType - Sets the `SNESFunctionType`
2204: of the `SNES` method.
2206: Logically Collective
2208: Input Parameters:
2209: + snes - the `SNES` context
2210: - type - the function type
2212: Level: developer
2214: Values of the function type\:
2215: + `SNES_FUNCTION_DEFAULT` - the default for the given `SNESType`
2216: . `SNES_FUNCTION_UNPRECONDITIONED` - an unpreconditioned function evaluation (this is the function provided with `SNESSetFunction()`
2217: - `SNES_FUNCTION_PRECONDITIONED` - a transformation of the function provided with `SNESSetFunction()`
2219: Note:
2220: Different `SNESType`s use this value in different ways
2222: .seealso: [](ch_snes), `SNES`, `SNESFunctionType`, `SNESGetNormSchedule()`, `SNESComputeFunction()`, `VecNorm()`, `SNESSetFunction()`, `SNESSetInitialFunction()`, `SNESNormSchedule`
2223: @*/
2224: PetscErrorCode SNESSetFunctionType(SNES snes, SNESFunctionType type)
2225: {
2226: PetscFunctionBegin;
2228: snes->functype = type;
2229: PetscFunctionReturn(PETSC_SUCCESS);
2230: }
2232: /*@
2233: SNESGetFunctionType - Gets the `SNESFunctionType` used in convergence and monitoring set with `SNESSetFunctionType()`
2234: of the SNES method.
2236: Logically Collective
2238: Input Parameters:
2239: + snes - the `SNES` context
2240: - type - the type of the function evaluation, see `SNESSetFunctionType()`
2242: Level: advanced
2244: .seealso: [](ch_snes), `SNESSetFunctionType()`, `SNESFunctionType`, `SNESSetNormSchedule()`, `SNESComputeFunction()`, `VecNorm()`, `SNESSetFunction()`, `SNESSetInitialFunction()`, `SNESNormSchedule`
2245: @*/
2246: PetscErrorCode SNESGetFunctionType(SNES snes, SNESFunctionType *type)
2247: {
2248: PetscFunctionBegin;
2250: *type = snes->functype;
2251: PetscFunctionReturn(PETSC_SUCCESS);
2252: }
2254: /*@
2255: SNESSetNGS - Sets the user nonlinear Gauss-Seidel routine for
2256: use with composed nonlinear solvers.
2258: Input Parameters:
2259: + snes - the `SNES` context, usually of the `SNESType` `SNESNGS`
2260: . f - function evaluation routine to apply Gauss-Seidel, see `SNESNGSFn` for calling sequence
2261: - ctx - [optional] user-defined context for private data for the smoother evaluation routine (may be `NULL`)
2263: Level: intermediate
2265: Note:
2266: The `SNESNGS` routines are used by the composed nonlinear solver to generate
2267: a problem appropriate update to the solution, particularly `SNESFAS`.
2269: .seealso: [](ch_snes), `SNESNGS`, `SNESGetNGS()`, `SNESNCG`, `SNESGetFunction()`, `SNESComputeNGS()`, `SNESNGSFn`
2270: @*/
2271: PetscErrorCode SNESSetNGS(SNES snes, SNESNGSFn *f, PetscCtx ctx)
2272: {
2273: DM dm;
2275: PetscFunctionBegin;
2277: PetscCall(SNESGetDM(snes, &dm));
2278: PetscCall(DMSNESSetNGS(dm, f, ctx));
2279: PetscFunctionReturn(PETSC_SUCCESS);
2280: }
2282: /*@
2283: SNESPicardComputeMFFunction - Matrix-free residual $A(x) x - b(x)$ used by `SNESSetPicard()` when the operator is applied through `-snes_mf_operator`
2285: Collective
2287: Input Parameters:
2288: + snes - the `SNES` context
2289: . x - the current iterate
2290: - ctx - unused application context; the Picard callbacks are retrieved from the attached `DMSNES`
2292: Output Parameter:
2293: . f - the residual vector
2295: Level: developer
2297: Note:
2298: Uses a duplicate of `snes->jacobian_pre` because `snes->jacobian_pre` cannot be changed during the `KSPSolve()`.
2300: .seealso: [](ch_snes), `SNES`, `SNESSetPicard()`, `SNESPicardComputeFunction()`, `SNESPicardComputeJacobian()`
2301: @*/
2302: PetscErrorCode SNESPicardComputeMFFunction(SNES snes, Vec x, Vec f, PetscCtx ctx)
2303: {
2304: DM dm;
2305: DMSNES sdm;
2307: PetscFunctionBegin;
2308: PetscCall(SNESGetDM(snes, &dm));
2309: PetscCall(DMGetDMSNES(dm, &sdm));
2310: /* A(x)*x - b(x) */
2311: if (sdm->ops->computepfunction) {
2312: PetscCallBack("SNES Picard callback function", (*sdm->ops->computepfunction)(snes, x, f, sdm->pctx));
2313: PetscCall(VecScale(f, -1.0));
2314: /* Cannot share nonzero pattern because of the possible use of SNESComputeJacobianDefault() */
2315: if (!snes->picard) PetscCall(MatDuplicate(snes->jacobian_pre, MAT_DO_NOT_COPY_VALUES, &snes->picard));
2316: PetscCallBack("SNES Picard callback Jacobian", (*sdm->ops->computepjacobian)(snes, x, snes->picard, snes->picard, sdm->pctx));
2317: PetscCall(MatMultAdd(snes->picard, x, f, f));
2318: } else {
2319: PetscCallBack("SNES Picard callback Jacobian", (*sdm->ops->computepjacobian)(snes, x, snes->picard, snes->picard, sdm->pctx));
2320: PetscCall(MatMult(snes->picard, x, f));
2321: }
2322: PetscFunctionReturn(PETSC_SUCCESS);
2323: }
2325: /*@
2326: SNESPicardComputeFunction - Compute the residual $A(x) x - b(x)$ using the callbacks registered by `SNESSetPicard()`
2328: Collective
2330: Input Parameters:
2331: + snes - the `SNES` context
2332: . x - the current iterate
2333: - ctx - unused application context; the Picard callbacks are retrieved from the attached `DMSNES`
2335: Output Parameter:
2336: . f - the residual vector
2338: Level: developer
2340: .seealso: [](ch_snes), `SNES`, `SNESSetPicard()`, `SNESPicardComputeMFFunction()`, `SNESPicardComputeJacobian()`
2341: @*/
2342: PetscErrorCode SNESPicardComputeFunction(SNES snes, Vec x, Vec f, PetscCtx ctx)
2343: {
2344: DM dm;
2345: DMSNES sdm;
2347: PetscFunctionBegin;
2348: PetscCall(SNESGetDM(snes, &dm));
2349: PetscCall(DMGetDMSNES(dm, &sdm));
2350: /* A(x)*x - b(x) */
2351: if (sdm->ops->computepfunction) {
2352: PetscCallBack("SNES Picard callback function", (*sdm->ops->computepfunction)(snes, x, f, sdm->pctx));
2353: PetscCall(VecScale(f, -1.0));
2354: PetscCallBack("SNES Picard callback Jacobian", (*sdm->ops->computepjacobian)(snes, x, snes->jacobian, snes->jacobian_pre, sdm->pctx));
2355: PetscCall(MatMultAdd(snes->jacobian_pre, x, f, f));
2356: } else {
2357: PetscCallBack("SNES Picard callback Jacobian", (*sdm->ops->computepjacobian)(snes, x, snes->jacobian, snes->jacobian_pre, sdm->pctx));
2358: PetscCall(MatMult(snes->jacobian_pre, x, f));
2359: }
2360: PetscFunctionReturn(PETSC_SUCCESS);
2361: }
2363: /*@
2364: SNESPicardComputeJacobian - Trivial Jacobian assembly callback used by `SNESSetPicard()`; the Picard operator is filled in by `SNESPicardComputeFunction()`
2366: Collective
2368: Input Parameters:
2369: + snes - the `SNES` context
2370: . x1 - the current iterate (unused)
2371: . J - the Jacobian matrix to assemble
2372: . B - the preconditioning matrix (unused)
2373: - ctx - unused application context
2375: Level: developer
2377: Note:
2378: Only calls `MatAssemblyBegin()`/`MatAssemblyEnd()` on `J`, because the Picard iteration reuses the operator already assembled by `SNESPicardComputeFunction()`.
2380: .seealso: [](ch_snes), `SNES`, `SNESSetPicard()`, `SNESPicardComputeFunction()`, `SNESPicardComputeMFFunction()`
2381: @*/
2382: PetscErrorCode SNESPicardComputeJacobian(SNES snes, Vec x1, Mat J, Mat B, PetscCtx ctx)
2383: {
2384: PetscFunctionBegin;
2385: /* the jacobian matrix should be pre-filled in SNESPicardComputeFunction */
2386: /* must assembly if matrix-free to get the last SNES solution */
2387: PetscCall(MatAssemblyBegin(J, MAT_FINAL_ASSEMBLY));
2388: PetscCall(MatAssemblyEnd(J, MAT_FINAL_ASSEMBLY));
2389: PetscFunctionReturn(PETSC_SUCCESS);
2390: }
2392: /*@
2393: SNESSetPicard - Use `SNES` to solve the system $A(x) x = bp(x) + b$ via a Picard type iteration (Picard linearization)
2395: Logically Collective
2397: Input Parameters:
2398: + snes - the `SNES` context
2399: . r - vector to store function values, may be `NULL`
2400: . bp - function evaluation routine, may be `NULL`, for the calling sequence see `SNESFunctionFn`
2401: . Amat - matrix with which $A(x) x - bp(x) - b$ is to be computed
2402: . Pmat - matrix from which preconditioner is computed (usually the same as `Amat`)
2403: . J - function to compute matrix values, for the calling sequence see `SNESJacobianFn`
2404: - ctx - [optional] user-defined context for private data for the function evaluation routine (may be `NULL`)
2406: Level: intermediate
2408: Notes:
2409: It is often better to provide the nonlinear function $F()$ and some approximation to its Jacobian directly and use
2410: an approximate Newton solver. This interface is provided to allow porting/testing a previous Picard based code in PETSc before converting it to approximate Newton.
2412: One can call `SNESSetPicard()` or `SNESSetFunction()` (and possibly `SNESSetJacobian()`) but cannot call both
2414: Solves the equation $A(x) x = bp(x) + b$ via the defect correction algorithm $A(x^{n}) (x^{n+1} - x^{n}) = bp(x^{n}) + b - A(x^{n})x^{n}$.
2415: When an exact solver is used this corresponds to the "classic" Picard $A(x^{n}) x^{n+1} = bp(x^{n}) + b$ iteration.
2417: Run with `-snes_mf_operator` to solve the system with Newton's method using $A(x^{n})$ to construct the preconditioner.
2419: We implement the defect correction form of the Picard iteration because it converges much more generally when inexact linear solvers are used then
2420: the direct Picard iteration $A(x^n) x^{n+1} = bp(x^n) + b$
2422: There is some controversity over the definition of a Picard iteration for nonlinear systems but almost everyone agrees that it involves a linear solve and some
2423: believe it is the iteration $A(x^{n}) x^{n+1} = b(x^{n})$ hence we use the name Picard. If anyone has an authoritative reference that defines the Picard iteration
2424: different please contact us at petsc-dev@mcs.anl.gov and we'll have an entirely new argument \:-).
2426: When used with `-snes_mf_operator` this will run matrix-free Newton's method where the matrix-vector product is of the true Jacobian of $A(x)x - bp(x) - b$ and
2427: $A(x^{n})$ is used to build the preconditioner
2429: When used with `-snes_fd` this will compute the true Jacobian (very slowly one column at a time) and thus represent Newton's method.
2431: When used with `-snes_fd_coloring` this will compute the Jacobian via coloring and thus represent a faster implementation of Newton's method. But the
2432: the nonzero structure of the Jacobian is, in general larger than that of the Picard matrix $A$ so you must provide in $A$ the needed nonzero structure for the correct
2433: coloring. When using `DMDA` this may mean creating the matrix $A$ with `DMCreateMatrix()` using a wider stencil than strictly needed for $A$ or with a `DMDA_STENCIL_BOX`.
2434: See the comment in src/snes/tutorials/ex15.c.
2436: .seealso: [](ch_snes), `SNES`, `SNESGetFunction()`, `SNESSetFunction()`, `SNESComputeFunction()`, `SNESSetJacobian()`, `SNESGetPicard()`, `SNESLineSearchPreCheckPicard()`,
2437: `SNESFunctionFn`, `SNESJacobianFn`
2438: @*/
2439: PetscErrorCode SNESSetPicard(SNES snes, Vec r, SNESFunctionFn *bp, Mat Amat, Mat Pmat, SNESJacobianFn *J, PetscCtx ctx)
2440: {
2441: DM dm;
2443: PetscFunctionBegin;
2445: PetscCall(SNESGetDM(snes, &dm));
2446: PetscCall(DMSNESSetPicard(dm, bp, J, ctx));
2447: PetscCall(DMSNESSetMFFunction(dm, SNESPicardComputeMFFunction, ctx));
2448: PetscCall(SNESSetFunction(snes, r, SNESPicardComputeFunction, ctx));
2449: PetscCall(SNESSetJacobian(snes, Amat, Pmat, SNESPicardComputeJacobian, ctx));
2450: PetscFunctionReturn(PETSC_SUCCESS);
2451: }
2453: /*@
2454: SNESGetPicard - Returns the context for the Picard iteration
2456: Not Collective, but `Vec` is parallel if `SNES` is parallel. Collective if `Vec` is requested, but has not been created yet.
2458: Input Parameter:
2459: . snes - the `SNES` context
2461: Output Parameters:
2462: + r - the function (or `NULL`)
2463: . f - the function (or `NULL`); for calling sequence see `SNESFunctionFn`
2464: . Amat - the matrix used to defined the operation A(x) x - b(x) (or `NULL`)
2465: . Pmat - the matrix from which the preconditioner will be constructed (or `NULL`)
2466: . J - the function for matrix evaluation (or `NULL`); for calling sequence see `SNESJacobianFn`
2467: - ctx - the function context (or `NULL`)
2469: Level: advanced
2471: .seealso: [](ch_snes), `SNESSetFunction()`, `SNESSetPicard()`, `SNESGetFunction()`, `SNESGetJacobian()`, `SNESGetDM()`, `SNESFunctionFn`, `SNESJacobianFn`
2472: @*/
2473: PetscErrorCode SNESGetPicard(SNES snes, Vec *r, SNESFunctionFn **f, Mat *Amat, Mat *Pmat, SNESJacobianFn **J, PetscCtxRt ctx)
2474: {
2475: DM dm;
2477: PetscFunctionBegin;
2479: PetscCall(SNESGetFunction(snes, r, NULL, NULL));
2480: PetscCall(SNESGetJacobian(snes, Amat, Pmat, NULL, NULL));
2481: PetscCall(SNESGetDM(snes, &dm));
2482: PetscCall(DMSNESGetPicard(dm, f, J, ctx));
2483: PetscFunctionReturn(PETSC_SUCCESS);
2484: }
2486: /*@
2487: SNESSetComputeInitialGuess - Sets a routine used to compute an initial guess for the nonlinear problem
2489: Logically Collective
2491: Input Parameters:
2492: + snes - the `SNES` context
2493: . func - function evaluation routine, see `SNESInitialGuessFn` for the calling sequence
2494: - ctx - [optional] user-defined context for private data for the
2495: function evaluation routine (may be `NULL`)
2497: Level: intermediate
2499: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESSetFunction()`, `SNESGetFunction()`, `SNESComputeFunction()`, `SNESSetJacobian()`, `SNESInitialGuessFn`
2500: @*/
2501: PetscErrorCode SNESSetComputeInitialGuess(SNES snes, SNESInitialGuessFn *func, PetscCtx ctx)
2502: {
2503: PetscFunctionBegin;
2505: if (func) snes->ops->computeinitialguess = func;
2506: if (ctx) snes->initialguessP = ctx;
2507: PetscFunctionReturn(PETSC_SUCCESS);
2508: }
2510: /*@
2511: SNESGetRhs - Gets the vector for solving F(x) = `rhs`. If `rhs` is not set
2512: it assumes a zero right-hand side.
2514: Logically Collective
2516: Input Parameter:
2517: . snes - the `SNES` context
2519: Output Parameter:
2520: . rhs - the right-hand side vector or `NULL` if there is no right-hand side vector
2522: Level: intermediate
2524: .seealso: [](ch_snes), `SNES`, `SNESGetSolution()`, `SNESGetFunction()`, `SNESComputeFunction()`, `SNESSetJacobian()`, `SNESSetFunction()`
2525: @*/
2526: PetscErrorCode SNESGetRhs(SNES snes, Vec *rhs)
2527: {
2528: PetscFunctionBegin;
2530: PetscAssertPointer(rhs, 2);
2531: *rhs = snes->vec_rhs;
2532: PetscFunctionReturn(PETSC_SUCCESS);
2533: }
2535: /*@
2536: SNESComputeFunction - Calls the function that has been set with `SNESSetFunction()`.
2538: Collective
2540: Input Parameters:
2541: + snes - the `SNES` context
2542: - x - input vector
2544: Output Parameter:
2545: . f - function vector, as set by `SNESSetFunction()`
2547: Level: developer
2549: Notes:
2550: `SNESComputeFunction()` is typically used within nonlinear solvers
2551: implementations, so users would not generally call this routine themselves.
2553: When solving for $F(x) = b$, this routine computes $f = F(x) - b$.
2555: This function usually appears in the pattern.
2556: .vb
2557: SNESComputeFunction(snes, x, f);
2558: VecNorm(f, &fnorm);
2559: SNESCheckFunctionDomainError(snes, fnorm); or SNESLineSearchCheckFunctionDomainError(ls, fnorm);
2560: .ve
2561: to collectively handle the use of `SNESSetFunctionDomainError()` in the provided callback function.
2563: .seealso: [](ch_snes), `SNES`, `SNESSetFunction()`, `SNESGetFunction()`, `SNESComputeMFFunction()`, `SNESSetFunctionDomainError()`
2564: @*/
2565: PetscErrorCode SNESComputeFunction(SNES snes, Vec x, Vec f)
2566: {
2567: DM dm;
2568: DMSNES sdm;
2570: PetscFunctionBegin;
2574: PetscCheckSameComm(snes, 1, x, 2);
2575: PetscCheckSameComm(snes, 1, f, 3);
2576: PetscCall(VecValidValues_Internal(x, 2, PETSC_TRUE));
2578: PetscCall(SNESGetDM(snes, &dm));
2579: PetscCall(DMGetDMSNES(dm, &sdm));
2580: PetscCheck(sdm->ops->computefunction || snes->vec_rhs, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetFunction() or SNESSetDM() before SNESComputeFunction(), likely called from SNESSolve().");
2581: if (sdm->ops->computefunction) {
2582: if (sdm->ops->computefunction != SNESObjectiveComputeFunctionDefaultFD) PetscCall(PetscLogEventBegin(SNES_FunctionEval, snes, x, f, 0));
2583: PetscCall(VecLockReadPush(x));
2584: /* ensure domainerror is false prior to computefunction evaluation (may not have been reset) */
2585: snes->functiondomainerror = PETSC_FALSE;
2586: {
2587: void *ctx;
2588: SNESFunctionFn *computefunction;
2589: PetscCall(DMSNESGetFunction(dm, &computefunction, &ctx));
2590: PetscCallBack("SNES callback function", (*computefunction)(snes, x, f, ctx));
2591: }
2592: PetscCall(VecLockReadPop(x));
2593: if (sdm->ops->computefunction != SNESObjectiveComputeFunctionDefaultFD) PetscCall(PetscLogEventEnd(SNES_FunctionEval, snes, x, f, 0));
2594: } else /* if (snes->vec_rhs) */ {
2595: PetscCall(MatMult(snes->jacobian, x, f));
2596: }
2597: if (snes->vec_rhs) PetscCall(VecAXPY(f, -1.0, snes->vec_rhs));
2598: snes->nfuncs++;
2599: /*
2600: domainerror might not be set on all processes; so we tag vector locally with infinity and the next inner product or norm will
2601: propagate the value to all processes
2602: */
2603: PetscCall(VecFlag(f, snes->functiondomainerror));
2604: PetscFunctionReturn(PETSC_SUCCESS);
2605: }
2607: /*@
2608: SNESComputeMFFunction - Calls the function that has been set with `DMSNESSetMFFunction()`.
2610: Collective
2612: Input Parameters:
2613: + snes - the `SNES` context
2614: - x - input vector
2616: Output Parameter:
2617: . y - output vector
2619: Level: developer
2621: Notes:
2622: `SNESComputeMFFunction()` is used within the matrix-vector products called by the matrix created with `MatCreateSNESMF()`
2623: so users would not generally call this routine themselves.
2625: Since this function is intended for use with finite differencing it does not subtract the right-hand side vector provided with `SNESSolve()`
2626: while `SNESComputeFunction()` does. As such, this routine cannot be used with `MatMFFDSetBase()` with a provided F function value even if it applies the
2627: same function as `SNESComputeFunction()` if a `SNESSolve()` right-hand side vector is use because the two functions difference would include this right hand side function.
2629: .seealso: [](ch_snes), `SNES`, `SNESSetFunction()`, `SNESGetFunction()`, `SNESComputeFunction()`, `MatCreateSNESMF()`, `DMSNESSetMFFunction()`
2630: @*/
2631: PetscErrorCode SNESComputeMFFunction(SNES snes, Vec x, Vec y)
2632: {
2633: DM dm;
2634: DMSNES sdm;
2636: PetscFunctionBegin;
2640: PetscCheckSameComm(snes, 1, x, 2);
2641: PetscCheckSameComm(snes, 1, y, 3);
2642: PetscCall(VecValidValues_Internal(x, 2, PETSC_TRUE));
2644: PetscCall(SNESGetDM(snes, &dm));
2645: PetscCall(DMGetDMSNES(dm, &sdm));
2646: PetscCall(PetscLogEventBegin(SNES_FunctionEval, snes, x, y, 0));
2647: PetscCall(VecLockReadPush(x));
2648: /* ensure domainerror is false prior to computefunction evaluation (may not have been reset) */
2649: snes->functiondomainerror = PETSC_FALSE;
2650: PetscCallBack("SNES callback function", (*sdm->ops->computemffunction)(snes, x, y, sdm->mffunctionctx));
2651: PetscCall(VecLockReadPop(x));
2652: PetscCall(PetscLogEventEnd(SNES_FunctionEval, snes, x, y, 0));
2653: snes->nfuncs++;
2654: /*
2655: domainerror might not be set on all processes; so we tag vector locally with infinity and the next inner product or norm will
2656: propagate the value to all processes
2657: */
2658: PetscCall(VecFlag(y, snes->functiondomainerror));
2659: PetscFunctionReturn(PETSC_SUCCESS);
2660: }
2662: /*@
2663: SNESComputeNGS - Calls the Gauss-Seidel function that has been set with `SNESSetNGS()`.
2665: Collective
2667: Input Parameters:
2668: + snes - the `SNES` context
2669: . x - input vector
2670: - b - rhs vector
2672: Output Parameter:
2673: . x - new solution vector
2675: Level: developer
2677: Note:
2678: `SNESComputeNGS()` is typically used within composed nonlinear solver
2679: implementations, so most users would not generally call this routine
2680: themselves.
2682: .seealso: [](ch_snes), `SNESNGSFn`, `SNESSetNGS()`, `SNESComputeFunction()`, `SNESNGS`
2683: @*/
2684: PetscErrorCode SNESComputeNGS(SNES snes, Vec b, Vec x)
2685: {
2686: DM dm;
2687: DMSNES sdm;
2689: PetscFunctionBegin;
2693: PetscCheckSameComm(snes, 1, x, 3);
2694: if (b) PetscCheckSameComm(snes, 1, b, 2);
2695: if (b) PetscCall(VecValidValues_Internal(b, 2, PETSC_TRUE));
2696: PetscCall(PetscLogEventBegin(SNES_NGSEval, snes, x, b, 0));
2697: PetscCall(SNESGetDM(snes, &dm));
2698: PetscCall(DMGetDMSNES(dm, &sdm));
2699: PetscCheck(sdm->ops->computegs, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetNGS() before SNESComputeNGS(), likely called from SNESSolve().");
2700: if (b) PetscCall(VecLockReadPush(b));
2701: PetscCallBack("SNES callback NGS", (*sdm->ops->computegs)(snes, x, b, sdm->gsctx));
2702: if (b) PetscCall(VecLockReadPop(b));
2703: PetscCall(PetscLogEventEnd(SNES_NGSEval, snes, x, b, 0));
2704: PetscFunctionReturn(PETSC_SUCCESS);
2705: }
2707: static PetscErrorCode SNESComputeFunction_FD(SNES snes, Vec Xin, Vec G)
2708: {
2709: Vec X;
2710: PetscScalar *g;
2711: PetscReal f, f2;
2712: PetscInt low, high, N, i;
2713: PetscBool flg;
2714: PetscReal h = .5 * PETSC_SQRT_MACHINE_EPSILON;
2716: PetscFunctionBegin;
2717: PetscCall(PetscOptionsGetReal(((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_fd_delta", &h, &flg));
2718: PetscCall(VecDuplicate(Xin, &X));
2719: PetscCall(VecCopy(Xin, X));
2720: PetscCall(VecGetSize(X, &N));
2721: PetscCall(VecGetOwnershipRange(X, &low, &high));
2722: PetscCall(VecSetOption(X, VEC_IGNORE_OFF_PROC_ENTRIES, PETSC_TRUE));
2723: PetscCall(VecGetArray(G, &g));
2724: for (i = 0; i < N; i++) {
2725: PetscCall(VecSetValue(X, i, -h, ADD_VALUES));
2726: PetscCall(VecAssemblyBegin(X));
2727: PetscCall(VecAssemblyEnd(X));
2728: PetscCall(SNESComputeObjective(snes, X, &f));
2729: PetscCall(VecSetValue(X, i, 2.0 * h, ADD_VALUES));
2730: PetscCall(VecAssemblyBegin(X));
2731: PetscCall(VecAssemblyEnd(X));
2732: PetscCall(SNESComputeObjective(snes, X, &f2));
2733: PetscCall(VecSetValue(X, i, -h, ADD_VALUES));
2734: PetscCall(VecAssemblyBegin(X));
2735: PetscCall(VecAssemblyEnd(X));
2736: if (i >= low && i < high) g[i - low] = (f2 - f) / (2.0 * h);
2737: }
2738: PetscCall(VecRestoreArray(G, &g));
2739: PetscCall(VecDestroy(&X));
2740: PetscFunctionReturn(PETSC_SUCCESS);
2741: }
2743: /*@
2744: SNESTestFunction - Computes the difference between the computed and finite-difference functions
2746: Collective
2748: Input Parameter:
2749: . snes - the `SNES` context
2751: Options Database Keys:
2752: + -snes_test_function - compare the user provided function with one compute via finite differences to check for errors.
2753: - -snes_test_function_view - display the user provided function, the finite difference function and the difference
2755: Level: developer
2757: .seealso: [](ch_snes), `SNESTestJacobian()`, `SNESSetFunction()`, `SNESComputeFunction()`
2758: @*/
2759: PetscErrorCode SNESTestFunction(SNES snes)
2760: {
2761: Vec x, g1, g2, g3;
2762: PetscBool complete_print = PETSC_FALSE;
2763: PetscReal hcnorm, fdnorm, hcmax, fdmax, diffmax, diffnorm;
2764: PetscScalar dot;
2765: MPI_Comm comm;
2766: PetscViewer viewer, mviewer;
2767: PetscViewerFormat format;
2768: PetscInt tabs;
2769: static PetscBool directionsprinted = PETSC_FALSE;
2770: SNESObjectiveFn *objective;
2772: PetscFunctionBegin;
2773: PetscCall(SNESGetObjective(snes, &objective, NULL));
2774: if (!objective) PetscFunctionReturn(PETSC_SUCCESS);
2776: PetscObjectOptionsBegin((PetscObject)snes);
2777: PetscCall(PetscOptionsViewer("-snes_test_function_view", "View difference between hand-coded and finite difference function element entries", "None", &mviewer, &format, &complete_print));
2778: PetscOptionsEnd();
2780: PetscCall(PetscObjectGetComm((PetscObject)snes, &comm));
2781: PetscCall(PetscViewerASCIIGetStdout(comm, &viewer));
2782: PetscCall(PetscViewerASCIIGetTab(viewer, &tabs));
2783: PetscCall(PetscViewerASCIISetTab(viewer, ((PetscObject)snes)->tablevel));
2784: PetscCall(PetscViewerASCIIPrintf(viewer, " ---------- Testing Function -------------\n"));
2785: if (!complete_print && !directionsprinted) {
2786: PetscCall(PetscViewerASCIIPrintf(viewer, " Run with -snes_test_function_view and optionally -snes_test_function <threshold> to show difference\n"));
2787: PetscCall(PetscViewerASCIIPrintf(viewer, " of hand-coded and finite difference function entries greater than <threshold>.\n"));
2788: }
2789: if (!directionsprinted) {
2790: PetscCall(PetscViewerASCIIPrintf(viewer, " Testing hand-coded Function, if (for double precision runs) ||F - Ffd||/||F|| is\n"));
2791: PetscCall(PetscViewerASCIIPrintf(viewer, " O(1.e-8), the hand-coded Function is probably correct.\n"));
2792: directionsprinted = PETSC_TRUE;
2793: }
2794: if (complete_print) PetscCall(PetscViewerPushFormat(mviewer, format));
2796: PetscCall(SNESGetSolution(snes, &x));
2797: PetscCall(VecDuplicate(x, &g1));
2798: PetscCall(VecDuplicate(x, &g2));
2799: PetscCall(VecDuplicate(x, &g3));
2800: PetscCall(SNESComputeFunction(snes, x, g1)); /* does not handle use of SNESSetFunctionDomainError() correctly */
2801: PetscCall(SNESComputeFunction_FD(snes, x, g2));
2803: PetscCall(VecNorm(g2, NORM_2, &fdnorm));
2804: PetscCall(VecNorm(g1, NORM_2, &hcnorm));
2805: PetscCall(VecNorm(g2, NORM_INFINITY, &fdmax));
2806: PetscCall(VecNorm(g1, NORM_INFINITY, &hcmax));
2807: PetscCall(VecDot(g1, g2, &dot));
2808: PetscCall(VecCopy(g1, g3));
2809: PetscCall(VecAXPY(g3, -1.0, g2));
2810: PetscCall(VecNorm(g3, NORM_2, &diffnorm));
2811: PetscCall(VecNorm(g3, NORM_INFINITY, &diffmax));
2812: PetscCall(PetscViewerASCIIPrintf(viewer, " ||Ffd|| %g, ||F|| = %g, angle cosine = (Ffd'F)/||Ffd||||F|| = %g\n", (double)fdnorm, (double)hcnorm, (double)(PetscRealPart(dot) / (fdnorm * hcnorm))));
2813: PetscCall(PetscViewerASCIIPrintf(viewer, " 2-norm ||F - Ffd||/||F|| = %g, ||F - Ffd|| = %g\n", (double)(diffnorm / PetscMax(hcnorm, fdnorm)), (double)diffnorm));
2814: PetscCall(PetscViewerASCIIPrintf(viewer, " max-norm ||F - Ffd||/||F|| = %g, ||F - Ffd|| = %g\n", (double)(diffmax / PetscMax(hcmax, fdmax)), (double)diffmax));
2816: if (complete_print) {
2817: PetscCall(PetscViewerASCIIPrintf(viewer, " Hand-coded function ----------\n"));
2818: PetscCall(VecView(g1, mviewer));
2819: PetscCall(PetscViewerASCIIPrintf(viewer, " Finite difference function ----------\n"));
2820: PetscCall(VecView(g2, mviewer));
2821: PetscCall(PetscViewerASCIIPrintf(viewer, " Hand-coded minus finite-difference function ----------\n"));
2822: PetscCall(VecView(g3, mviewer));
2823: }
2824: PetscCall(VecDestroy(&g1));
2825: PetscCall(VecDestroy(&g2));
2826: PetscCall(VecDestroy(&g3));
2828: if (complete_print) {
2829: PetscCall(PetscViewerPopFormat(mviewer));
2830: PetscCall(PetscViewerDestroy(&mviewer));
2831: }
2832: PetscCall(PetscViewerASCIISetTab(viewer, tabs));
2833: PetscFunctionReturn(PETSC_SUCCESS);
2834: }
2836: /*@
2837: SNESTestJacobian - Computes the difference between the computed and finite-difference Jacobians
2839: Collective
2841: Input Parameter:
2842: . snes - the `SNES` context
2844: Output Parameters:
2845: + Jnorm - the Frobenius norm of the computed Jacobian, or `NULL`
2846: - diffNorm - the Frobenius norm of the difference of the computed and finite-difference Jacobians, or `NULL`
2848: Options Database Keys:
2849: + -snes_test_jacobian [threshold] - compare the user provided Jacobian with one compute via finite differences to check for errors.
2850: If a threshold is given, display only those entries whose difference is greater than the threshold.
2851: - -snes_test_jacobian_view viewer_specification - display the user provided Jacobian, the finite difference Jacobian and the difference, see `PetscOptionsCreateViewer()` for the
2852: format of `viewer_specification`
2854: Level: developer
2856: Note:
2857: Directions and norms are printed to stdout if `diffNorm` is `NULL`.
2859: .seealso: [](ch_snes), `SNESTestFunction()`, `SNESSetJacobian()`, `SNESComputeJacobian()`
2860: @*/
2861: PetscErrorCode SNESTestJacobian(SNES snes, PetscReal *Jnorm, PetscReal *diffNorm)
2862: {
2863: Mat A, B, C, D, jacobian;
2864: Vec x = snes->vec_sol, f;
2865: PetscReal nrm, gnorm;
2866: PetscReal threshold = 1.e-5;
2867: void *functx;
2868: PetscBool complete_print = PETSC_FALSE, threshold_print = PETSC_FALSE, flg, istranspose;
2869: PetscBool silent = diffNorm != PETSC_NULLPTR ? PETSC_TRUE : PETSC_FALSE;
2870: PetscViewer viewer, mviewer;
2871: MPI_Comm comm;
2872: PetscInt tabs;
2873: static PetscBool directionsprinted = PETSC_FALSE;
2874: PetscViewerFormat format;
2876: PetscFunctionBegin;
2877: PetscObjectOptionsBegin((PetscObject)snes);
2878: PetscCall(PetscOptionsReal("-snes_test_jacobian", "Threshold for element difference between hand-coded and finite difference being meaningful", "None", threshold, &threshold, NULL));
2879: PetscCall(PetscOptionsDeprecated("-snes_test_jacobian_display", "-snes_test_jacobian_view", "3.13", NULL));
2880: PetscCall(PetscOptionsViewer("-snes_test_jacobian_view", "View difference between hand-coded and finite difference Jacobians element entries", "None", &mviewer, &format, &complete_print));
2881: PetscCall(PetscOptionsDeprecated("-snes_test_jacobian_display_threshold", "-snes_test_jacobian", "3.13", "-snes_test_jacobian accepts an optional threshold (since v3.10)"));
2882: PetscCall(PetscOptionsReal("-snes_test_jacobian_display_threshold", "Display difference between hand-coded and finite difference Jacobians which exceed input threshold", "None", threshold, &threshold, &threshold_print));
2883: PetscOptionsEnd();
2885: PetscCall(PetscObjectGetComm((PetscObject)snes, &comm));
2886: PetscCall(PetscViewerASCIIGetStdout(comm, &viewer));
2887: PetscCall(PetscViewerASCIIGetTab(viewer, &tabs));
2888: PetscCall(PetscViewerASCIISetTab(viewer, ((PetscObject)snes)->tablevel));
2889: if (!silent) PetscCall(PetscViewerASCIIPrintf(viewer, " ---------- Testing Jacobian -------------\n"));
2890: if (!complete_print && !silent && !directionsprinted) {
2891: PetscCall(PetscViewerASCIIPrintf(viewer, " Run with -snes_test_jacobian_view and optionally -snes_test_jacobian <threshold> to show difference\n"));
2892: PetscCall(PetscViewerASCIIPrintf(viewer, " of hand-coded and finite difference Jacobian entries greater than <threshold>.\n"));
2893: }
2894: if (!directionsprinted && !silent) {
2895: PetscCall(PetscViewerASCIIPrintf(viewer, " Testing hand-coded Jacobian, if (for double precision runs) ||J - Jfd||_F/||J||_F is\n"));
2896: PetscCall(PetscViewerASCIIPrintf(viewer, " O(1.e-8), the hand-coded Jacobian is probably correct.\n"));
2897: directionsprinted = PETSC_TRUE;
2898: }
2899: if (complete_print) PetscCall(PetscViewerPushFormat(mviewer, format));
2901: PetscCall(PetscObjectTypeCompare((PetscObject)snes->jacobian, MATMFFD, &flg));
2902: if (!flg) jacobian = snes->jacobian;
2903: else jacobian = snes->jacobian_pre;
2905: if (!x) PetscCall(MatCreateVecs(jacobian, &x, NULL));
2906: else PetscCall(PetscObjectReference((PetscObject)x));
2907: PetscCall(VecDuplicate(x, &f));
2909: /* evaluate the function at this point because SNESComputeJacobianDefault() assumes that the function has been evaluated and put into snes->vec_func */
2910: PetscCall(SNESComputeFunction(snes, x, f));
2911: PetscCall(VecDestroy(&f));
2912: PetscCall(PetscObjectTypeCompare((PetscObject)snes, SNESKSPTRANSPOSEONLY, &istranspose));
2913: while (jacobian) {
2914: Mat JT = NULL, Jsave = NULL;
2916: if (istranspose) {
2917: PetscCall(MatCreateTranspose(jacobian, &JT));
2918: Jsave = jacobian;
2919: jacobian = JT;
2920: }
2921: PetscCall(PetscObjectBaseTypeCompareAny((PetscObject)jacobian, &flg, MATSEQAIJ, MATMPIAIJ, MATSEQDENSE, MATMPIDENSE, MATSEQBAIJ, MATMPIBAIJ, MATSEQSBAIJ, MATMPISBAIJ, ""));
2922: if (flg) {
2923: A = jacobian;
2924: PetscCall(PetscObjectReference((PetscObject)A));
2925: } else {
2926: PetscCall(MatComputeOperator(jacobian, MATAIJ, &A));
2927: }
2929: PetscCall(MatDuplicate(A, MAT_DO_NOT_COPY_VALUES, &B));
2930: PetscCall(MatSetOption(B, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE));
2932: PetscCall(SNESGetFunction(snes, NULL, NULL, &functx));
2933: PetscCall(SNESComputeJacobianDefault(snes, x, B, B, functx));
2935: PetscCall(MatDuplicate(B, MAT_COPY_VALUES, &D));
2936: PetscCall(MatAYPX(D, -1.0, A, DIFFERENT_NONZERO_PATTERN));
2937: PetscCall(MatNorm(D, NORM_FROBENIUS, &nrm));
2938: PetscCall(MatNorm(A, NORM_FROBENIUS, &gnorm));
2939: PetscCall(MatDestroy(&D));
2940: if (!gnorm) gnorm = 1; /* just in case */
2941: if (!silent) PetscCall(PetscViewerASCIIPrintf(viewer, " ||J - Jfd||_F/||J||_F = %g, ||J - Jfd||_F = %g\n", (double)(nrm / gnorm), (double)nrm));
2942: if (complete_print) {
2943: PetscCall(PetscViewerASCIIPrintf(viewer, " Hand-coded Jacobian ----------\n"));
2944: PetscCall(MatView(A, mviewer));
2945: PetscCall(PetscViewerASCIIPrintf(viewer, " Finite difference Jacobian ----------\n"));
2946: PetscCall(MatView(B, mviewer));
2947: }
2949: if (threshold_print || complete_print) {
2950: PetscInt Istart, Iend, *ccols, bncols, cncols, j, row;
2951: PetscScalar *cvals;
2952: const PetscInt *bcols;
2953: const PetscScalar *bvals;
2955: PetscCall(MatDuplicate(A, MAT_DO_NOT_COPY_VALUES, &C));
2956: PetscCall(MatSetOption(C, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE));
2958: PetscCall(MatAYPX(B, -1.0, A, DIFFERENT_NONZERO_PATTERN));
2959: PetscCall(MatGetOwnershipRange(B, &Istart, &Iend));
2961: for (row = Istart; row < Iend; row++) {
2962: PetscCall(MatGetRow(B, row, &bncols, &bcols, &bvals));
2963: PetscCall(PetscMalloc2(bncols, &ccols, bncols, &cvals));
2964: for (j = 0, cncols = 0; j < bncols; j++) {
2965: if (PetscAbsScalar(bvals[j]) > threshold) {
2966: ccols[cncols] = bcols[j];
2967: cvals[cncols] = bvals[j];
2968: cncols += 1;
2969: }
2970: }
2971: if (cncols) PetscCall(MatSetValues(C, 1, &row, cncols, ccols, cvals, INSERT_VALUES));
2972: PetscCall(MatRestoreRow(B, row, &bncols, &bcols, &bvals));
2973: PetscCall(PetscFree2(ccols, cvals));
2974: }
2975: PetscCall(MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY));
2976: PetscCall(MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY));
2977: PetscCall(PetscViewerASCIIPrintf(viewer, " Hand-coded minus finite-difference Jacobian with tolerance %g ----------\n", (double)threshold));
2978: PetscCall(MatView(C, complete_print ? mviewer : viewer));
2979: PetscCall(MatDestroy(&C));
2980: }
2981: PetscCall(MatDestroy(&A));
2982: PetscCall(MatDestroy(&B));
2983: PetscCall(MatDestroy(&JT));
2984: if (Jsave) jacobian = Jsave;
2985: if (jacobian != snes->jacobian_pre) {
2986: jacobian = snes->jacobian_pre;
2987: if (!silent) PetscCall(PetscViewerASCIIPrintf(viewer, " ---------- Testing Jacobian for preconditioner -------------\n"));
2988: } else jacobian = NULL;
2989: }
2990: PetscCall(VecDestroy(&x));
2991: if (complete_print) PetscCall(PetscViewerPopFormat(mviewer));
2992: PetscCall(PetscViewerDestroy(&mviewer));
2993: PetscCall(PetscViewerASCIISetTab(viewer, tabs));
2995: if (Jnorm) *Jnorm = gnorm;
2996: if (diffNorm) *diffNorm = nrm;
2997: PetscFunctionReturn(PETSC_SUCCESS);
2998: }
3000: /*@
3001: SNESComputeJacobian - Computes the Jacobian matrix that has been set with `SNESSetJacobian()`.
3003: Collective
3005: Input Parameters:
3006: + snes - the `SNES` context
3007: - X - input vector
3009: Output Parameters:
3010: + A - Jacobian matrix
3011: - B - optional matrix for building the preconditioner, usually the same as `A`
3013: Options Database Keys:
3014: + -snes_lag_preconditioner lag - how often to rebuild preconditioner
3015: . -snes_lag_jacobian lag - how often to rebuild Jacobian
3016: . -snes_test_jacobian [threshold] - compare the user provided Jacobian with one compute via finite differences to check for errors.
3017: If a threshold is given, display only those entries whose difference is greater than the threshold.
3018: . -snes_test_jacobian_view viewer_specification - display the user provided Jacobian, the finite difference Jacobian and the difference between them to help users detect the location of errors in the user provided Jacobian.
3019: See `PetscOptionsCreateViewer()` for the format of `viewer_specification`
3020: . -snes_compare_explicit - compare the computed Jacobian to the finite difference Jacobian and output the differences
3021: . -snes_compare_explicit_draw - compare the computed Jacobian to the finite difference Jacobian and draw the result
3022: . -snes_compare_explicit_draw_contour - compare the computed Jacobian to the finite difference Jacobian and draw a contour plot with the result
3023: . -snes_compare_operator - make the comparison options above use the operator instead of the matrix used to construct the preconditioner
3024: . -snes_compare_coloring - compute the finite difference Jacobian using coloring and display norms of difference
3025: . -snes_compare_coloring_display - compute the finite difference Jacobian using coloring and display verbose differences
3026: . -snes_compare_coloring_threshold - display only those matrix entries that differ by more than a given threshold
3027: . -snes_compare_coloring_threshold_atol - absolute tolerance for difference in matrix entries to be displayed by `-snes_compare_coloring_threshold`
3028: . -snes_compare_coloring_threshold_rtol - relative tolerance for difference in matrix entries to be displayed by `-snes_compare_coloring_threshold`
3029: . -snes_compare_coloring_draw - compute the finite difference Jacobian using coloring and draw differences
3030: - -snes_compare_coloring_draw_contour - compute the finite difference Jacobian using coloring and show contours of matrices and differences
3032: Level: developer
3034: Note:
3035: Most users should not need to explicitly call this routine, as it
3036: is used internally within the nonlinear solvers.
3038: Developer Note:
3039: This has duplicative ways of checking the accuracy of the user provided Jacobian (see the options above). This is for historical reasons.
3041: .seealso: [](ch_snes), `SNESSetJacobian()`, `KSPSetOperators()`, `MatStructure`, `SNESSetLagPreconditioner()`, `SNESSetLagJacobian()`,
3042: `SNESSetJacobianDomainError()`, `SNESCheckJacobianDomainError()`, `SNESSetCheckJacobianDomainError()`
3043: @*/
3044: PetscErrorCode SNESComputeJacobian(SNES snes, Vec X, Mat A, Mat B)
3045: {
3046: PetscBool flag;
3047: DM dm;
3048: DMSNES sdm;
3049: KSP ksp;
3051: PetscFunctionBegin;
3054: PetscCheckSameComm(snes, 1, X, 2);
3055: PetscCall(VecValidValues_Internal(X, 2, PETSC_TRUE));
3056: PetscCall(SNESGetDM(snes, &dm));
3057: PetscCall(DMGetDMSNES(dm, &sdm));
3059: /* make sure that MatAssemblyBegin/End() is called on A matrix if it is matrix-free */
3060: if (snes->lagjacobian == -2) {
3061: snes->lagjacobian = -1;
3063: PetscCall(PetscInfo(snes, "Recomputing Jacobian/preconditioner because lag is -2 (means compute Jacobian, but then never again) \n"));
3064: } else if (snes->lagjacobian == -1) {
3065: PetscCall(PetscInfo(snes, "Reusing Jacobian/preconditioner because lag is -1\n"));
3066: PetscCall(PetscObjectTypeCompare((PetscObject)A, MATMFFD, &flag));
3067: if (flag) {
3068: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
3069: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
3070: }
3071: PetscFunctionReturn(PETSC_SUCCESS);
3072: } else if (snes->lagjacobian > 1 && (snes->iter + snes->jac_iter) % snes->lagjacobian) {
3073: PetscCall(PetscInfo(snes, "Reusing Jacobian/preconditioner because lag is %" PetscInt_FMT " and SNES iteration is %" PetscInt_FMT "\n", snes->lagjacobian, snes->iter));
3074: PetscCall(PetscObjectTypeCompare((PetscObject)A, MATMFFD, &flag));
3075: if (flag) {
3076: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
3077: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
3078: }
3079: PetscFunctionReturn(PETSC_SUCCESS);
3080: }
3081: if (snes->npc && snes->npcside == PC_LEFT) {
3082: /* SNESASPIN uses SNESNASM as the nonlinear preconditioner. When SNESNASM
3083: is done solving the sub-systems it calls the user-provided Jacobian function
3084: (corresponding to the unpreconditioned residual) retrieved through the DM.
3085: Consequently it would be redundant to call the Jacobian function here. In
3086: the future we may move the outer Jacobian function call out of SNESNASM
3087: in which case no special casing will be required here. */
3088: PetscCall(PetscObjectTypeCompare((PetscObject)snes, SNESASPIN, &flag));
3089: if (flag) {
3090: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
3091: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
3092: PetscFunctionReturn(PETSC_SUCCESS);
3093: }
3094: }
3096: PetscCall(PetscLogEventBegin(SNES_JacobianEval, snes, X, A, B));
3097: PetscCall(VecLockReadPush(X));
3098: {
3099: void *ctx;
3100: SNESJacobianFn *J;
3101: PetscCall(DMSNESGetJacobian(dm, &J, &ctx));
3102: PetscCallBack("SNES callback Jacobian", (*J)(snes, X, A, B, ctx));
3103: }
3104: PetscCall(VecLockReadPop(X));
3105: PetscCall(PetscLogEventEnd(SNES_JacobianEval, snes, X, A, B));
3107: /* attach latest linearization point to the matrix used to construct the preconditioner */
3108: PetscCall(PetscObjectCompose((PetscObject)B, "__SNES_latest_X", (PetscObject)X));
3110: /* the next line ensures that snes->ksp exists */
3111: PetscCall(SNESGetKSP(snes, &ksp));
3112: if (snes->lagpreconditioner == -2) {
3113: PetscCall(PetscInfo(snes, "Rebuilding preconditioner exactly once since lag is -2\n"));
3114: PetscCall(KSPSetReusePreconditioner(snes->ksp, PETSC_FALSE));
3115: snes->lagpreconditioner = -1;
3116: } else if (snes->lagpreconditioner == -1) {
3117: PetscCall(PetscInfo(snes, "Reusing preconditioner because lag is -1\n"));
3118: PetscCall(KSPSetReusePreconditioner(snes->ksp, PETSC_TRUE));
3119: } else if (snes->lagpreconditioner > 1 && (snes->iter + snes->pre_iter) % snes->lagpreconditioner) {
3120: PetscCall(PetscInfo(snes, "Reusing preconditioner because lag is %" PetscInt_FMT " and SNES iteration is %" PetscInt_FMT "\n", snes->lagpreconditioner, snes->iter));
3121: PetscCall(KSPSetReusePreconditioner(snes->ksp, PETSC_TRUE));
3122: } else {
3123: PetscCall(PetscInfo(snes, "Rebuilding preconditioner\n"));
3124: PetscCall(KSPSetReusePreconditioner(snes->ksp, PETSC_FALSE));
3125: }
3127: /* monkey business to allow testing Jacobians in multilevel solvers.
3128: This is needed because the SNESTestXXX interface does not accept vectors and matrices */
3129: {
3130: Vec xsave = snes->vec_sol;
3131: Mat jacobiansave = snes->jacobian;
3132: Mat jacobian_presave = snes->jacobian_pre;
3134: snes->vec_sol = X;
3135: snes->jacobian = A;
3136: snes->jacobian_pre = B;
3137: if (snes->testFunc) PetscCall(SNESTestFunction(snes));
3138: if (snes->testJac) PetscCall(SNESTestJacobian(snes, NULL, NULL));
3140: snes->vec_sol = xsave;
3141: snes->jacobian = jacobiansave;
3142: snes->jacobian_pre = jacobian_presave;
3143: }
3145: {
3146: PetscBool flag = PETSC_FALSE, flag_draw = PETSC_FALSE, flag_contour = PETSC_FALSE, flag_operator = PETSC_FALSE;
3147: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_compare_explicit", NULL, NULL, &flag));
3148: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_compare_explicit_draw", NULL, NULL, &flag_draw));
3149: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_compare_explicit_draw_contour", NULL, NULL, &flag_contour));
3150: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_compare_operator", NULL, NULL, &flag_operator));
3151: if (flag || flag_draw || flag_contour) {
3152: Mat Bexp_mine = NULL, Bexp, FDexp;
3153: PetscViewer vdraw, vstdout;
3154: PetscBool flg;
3155: if (flag_operator) {
3156: PetscCall(MatComputeOperator(A, MATAIJ, &Bexp_mine));
3157: Bexp = Bexp_mine;
3158: } else {
3159: /* See if the matrix used to construct the preconditioner can be viewed and added directly */
3160: PetscCall(PetscObjectBaseTypeCompareAny((PetscObject)B, &flg, MATSEQAIJ, MATMPIAIJ, MATSEQDENSE, MATMPIDENSE, MATSEQBAIJ, MATMPIBAIJ, MATSEQSBAIJ, MATMPISBAIJ, ""));
3161: if (flg) Bexp = B;
3162: else {
3163: /* If the "preconditioning" matrix is itself MATSHELL or some other type without direct support */
3164: PetscCall(MatComputeOperator(B, MATAIJ, &Bexp_mine));
3165: Bexp = Bexp_mine;
3166: }
3167: }
3168: PetscCall(MatConvert(Bexp, MATSAME, MAT_INITIAL_MATRIX, &FDexp));
3169: PetscCall(SNESComputeJacobianDefault(snes, X, FDexp, FDexp, NULL));
3170: PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)snes), &vstdout));
3171: if (flag_draw || flag_contour) {
3172: PetscCall(PetscViewerDrawOpen(PetscObjectComm((PetscObject)snes), NULL, "Explicit Jacobians", PETSC_DECIDE, PETSC_DECIDE, 300, 300, &vdraw));
3173: if (flag_contour) PetscCall(PetscViewerPushFormat(vdraw, PETSC_VIEWER_DRAW_CONTOUR));
3174: } else vdraw = NULL;
3175: PetscCall(PetscViewerASCIIPrintf(vstdout, "Explicit %s\n", flag_operator ? "Jacobian" : "preconditioning Jacobian"));
3176: if (flag) PetscCall(MatView(Bexp, vstdout));
3177: if (vdraw) PetscCall(MatView(Bexp, vdraw));
3178: PetscCall(PetscViewerASCIIPrintf(vstdout, "Finite difference Jacobian\n"));
3179: if (flag) PetscCall(MatView(FDexp, vstdout));
3180: if (vdraw) PetscCall(MatView(FDexp, vdraw));
3181: PetscCall(MatAYPX(FDexp, -1.0, Bexp, SAME_NONZERO_PATTERN));
3182: PetscCall(PetscViewerASCIIPrintf(vstdout, "User-provided matrix minus finite difference Jacobian\n"));
3183: if (flag) PetscCall(MatView(FDexp, vstdout));
3184: if (vdraw) { /* Always use contour for the difference */
3185: PetscCall(PetscViewerPushFormat(vdraw, PETSC_VIEWER_DRAW_CONTOUR));
3186: PetscCall(MatView(FDexp, vdraw));
3187: PetscCall(PetscViewerPopFormat(vdraw));
3188: }
3189: if (flag_contour) PetscCall(PetscViewerPopFormat(vdraw));
3190: PetscCall(PetscViewerDestroy(&vdraw));
3191: PetscCall(MatDestroy(&Bexp_mine));
3192: PetscCall(MatDestroy(&FDexp));
3193: }
3194: }
3195: {
3196: PetscBool flag = PETSC_FALSE, flag_display = PETSC_FALSE, flag_draw = PETSC_FALSE, flag_contour = PETSC_FALSE, flag_threshold = PETSC_FALSE;
3197: PetscReal threshold_atol = PETSC_SQRT_MACHINE_EPSILON, threshold_rtol = 10 * PETSC_SQRT_MACHINE_EPSILON;
3198: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_compare_coloring", NULL, NULL, &flag));
3199: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_compare_coloring_display", NULL, NULL, &flag_display));
3200: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_compare_coloring_draw", NULL, NULL, &flag_draw));
3201: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_compare_coloring_draw_contour", NULL, NULL, &flag_contour));
3202: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_compare_coloring_threshold", NULL, NULL, &flag_threshold));
3203: if (flag_threshold) {
3204: PetscCall(PetscOptionsGetReal(((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_compare_coloring_threshold_rtol", &threshold_rtol, NULL));
3205: PetscCall(PetscOptionsGetReal(((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_compare_coloring_threshold_atol", &threshold_atol, NULL));
3206: }
3207: if (flag || flag_display || flag_draw || flag_contour || flag_threshold) {
3208: Mat Bfd;
3209: PetscViewer vdraw, vstdout;
3210: MatColoring coloring;
3211: ISColoring iscoloring;
3212: MatFDColoring matfdcoloring;
3213: SNESFunctionFn *func;
3214: void *funcctx;
3215: PetscReal norm1, norm2, normmax;
3217: PetscCall(MatDuplicate(B, MAT_DO_NOT_COPY_VALUES, &Bfd));
3218: PetscCall(MatColoringCreate(Bfd, &coloring));
3219: PetscCall(MatColoringSetType(coloring, MATCOLORINGSL));
3220: PetscCall(MatColoringSetFromOptions(coloring));
3221: PetscCall(MatColoringApply(coloring, &iscoloring));
3222: PetscCall(MatColoringDestroy(&coloring));
3223: PetscCall(MatFDColoringCreate(Bfd, iscoloring, &matfdcoloring));
3224: PetscCall(MatFDColoringSetFromOptions(matfdcoloring));
3225: PetscCall(MatFDColoringSetUp(Bfd, iscoloring, matfdcoloring));
3226: PetscCall(ISColoringDestroy(&iscoloring));
3228: /* This method of getting the function is currently unreliable since it doesn't work for DM local functions. */
3229: PetscCall(SNESGetFunction(snes, NULL, &func, &funcctx));
3230: PetscCall(MatFDColoringSetFunction(matfdcoloring, (MatFDColoringFn *)func, funcctx));
3231: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)matfdcoloring, ((PetscObject)snes)->prefix));
3232: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)matfdcoloring, "coloring_"));
3233: PetscCall(MatFDColoringSetFromOptions(matfdcoloring));
3234: PetscCall(MatFDColoringApply(Bfd, matfdcoloring, X, snes));
3235: PetscCall(MatFDColoringDestroy(&matfdcoloring));
3237: PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)snes), &vstdout));
3238: if (flag_draw || flag_contour) {
3239: PetscCall(PetscViewerDrawOpen(PetscObjectComm((PetscObject)snes), NULL, "Colored Jacobians", PETSC_DECIDE, PETSC_DECIDE, 300, 300, &vdraw));
3240: if (flag_contour) PetscCall(PetscViewerPushFormat(vdraw, PETSC_VIEWER_DRAW_CONTOUR));
3241: } else vdraw = NULL;
3242: PetscCall(PetscViewerASCIIPrintf(vstdout, "Explicit preconditioning Jacobian\n"));
3243: if (flag_display) PetscCall(MatView(B, vstdout));
3244: if (vdraw) PetscCall(MatView(B, vdraw));
3245: PetscCall(PetscViewerASCIIPrintf(vstdout, "Colored Finite difference Jacobian\n"));
3246: if (flag_display) PetscCall(MatView(Bfd, vstdout));
3247: if (vdraw) PetscCall(MatView(Bfd, vdraw));
3248: PetscCall(MatAYPX(Bfd, -1.0, B, SAME_NONZERO_PATTERN));
3249: PetscCall(MatNorm(Bfd, NORM_1, &norm1));
3250: PetscCall(MatNorm(Bfd, NORM_FROBENIUS, &norm2));
3251: PetscCall(MatNorm(Bfd, NORM_MAX, &normmax));
3252: PetscCall(PetscViewerASCIIPrintf(vstdout, "User-provided matrix minus finite difference Jacobian, norm1=%g normFrob=%g normmax=%g\n", (double)norm1, (double)norm2, (double)normmax));
3253: if (flag_display) PetscCall(MatView(Bfd, vstdout));
3254: if (vdraw) { /* Always use contour for the difference */
3255: PetscCall(PetscViewerPushFormat(vdraw, PETSC_VIEWER_DRAW_CONTOUR));
3256: PetscCall(MatView(Bfd, vdraw));
3257: PetscCall(PetscViewerPopFormat(vdraw));
3258: }
3259: if (flag_contour) PetscCall(PetscViewerPopFormat(vdraw));
3261: if (flag_threshold) {
3262: PetscInt bs, rstart, rend, i;
3263: PetscCall(MatGetBlockSize(B, &bs));
3264: PetscCall(MatGetOwnershipRange(B, &rstart, &rend));
3265: for (i = rstart; i < rend; i++) {
3266: const PetscScalar *ba, *ca;
3267: const PetscInt *bj, *cj;
3268: PetscInt bn, cn, j, maxentrycol = -1, maxdiffcol = -1, maxrdiffcol = -1;
3269: PetscReal maxentry = 0, maxdiff = 0, maxrdiff = 0;
3270: PetscCall(MatGetRow(B, i, &bn, &bj, &ba));
3271: PetscCall(MatGetRow(Bfd, i, &cn, &cj, &ca));
3272: PetscCheck(bn == cn, ((PetscObject)A)->comm, PETSC_ERR_PLIB, "Unexpected different nonzero pattern in -snes_compare_coloring_threshold");
3273: for (j = 0; j < bn; j++) {
3274: PetscReal rdiff = PetscAbsScalar(ca[j]) / (threshold_atol + threshold_rtol * PetscAbsScalar(ba[j]));
3275: if (PetscAbsScalar(ba[j]) > PetscAbs(maxentry)) {
3276: maxentrycol = bj[j];
3277: maxentry = PetscRealPart(ba[j]);
3278: }
3279: if (PetscAbsScalar(ca[j]) > PetscAbs(maxdiff)) {
3280: maxdiffcol = bj[j];
3281: maxdiff = PetscRealPart(ca[j]);
3282: }
3283: if (rdiff > maxrdiff) {
3284: maxrdiffcol = bj[j];
3285: maxrdiff = rdiff;
3286: }
3287: }
3288: if (maxrdiff > 1) {
3289: PetscCall(PetscViewerASCIIPrintf(vstdout, "row %" PetscInt_FMT " (maxentry=%g at %" PetscInt_FMT ", maxdiff=%g at %" PetscInt_FMT ", maxrdiff=%g at %" PetscInt_FMT "):", i, (double)maxentry, maxentrycol, (double)maxdiff, maxdiffcol, (double)maxrdiff, maxrdiffcol));
3290: for (j = 0; j < bn; j++) {
3291: PetscReal rdiff;
3292: rdiff = PetscAbsScalar(ca[j]) / (threshold_atol + threshold_rtol * PetscAbsScalar(ba[j]));
3293: if (rdiff > 1) PetscCall(PetscViewerASCIIPrintf(vstdout, " (%" PetscInt_FMT ",%g:%g)", bj[j], (double)PetscRealPart(ba[j]), (double)PetscRealPart(ca[j])));
3294: }
3295: PetscCall(PetscViewerASCIIPrintf(vstdout, "\n"));
3296: }
3297: PetscCall(MatRestoreRow(B, i, &bn, &bj, &ba));
3298: PetscCall(MatRestoreRow(Bfd, i, &cn, &cj, &ca));
3299: }
3300: }
3301: PetscCall(PetscViewerDestroy(&vdraw));
3302: PetscCall(MatDestroy(&Bfd));
3303: }
3304: }
3305: PetscFunctionReturn(PETSC_SUCCESS);
3306: }
3308: /*@
3309: SNESSetJacobian - Sets the function to compute Jacobian as well as the
3310: location to store the matrix.
3312: Logically Collective
3314: Input Parameters:
3315: + snes - the `SNES` context
3316: . Amat - the matrix that defines the (approximate) Jacobian
3317: . Pmat - the matrix to be used in constructing the preconditioner, usually the same as `Amat`.
3318: . J - Jacobian evaluation routine (if `NULL` then `SNES` retains any previously set value), see `SNESJacobianFn` for details
3319: - ctx - [optional] user-defined context for private data for the
3320: Jacobian evaluation routine (may be `NULL`) (if `NULL` then `SNES` retains any previously set value)
3322: Level: beginner
3324: Notes:
3325: If the `Amat` matrix and `Pmat` matrix are different you must call `MatAssemblyBegin()`/`MatAssemblyEnd()` on
3326: each matrix.
3328: If you know the operator `Amat` has a null space you can use `MatSetNullSpace()` and `MatSetTransposeNullSpace()` to supply the null
3329: space to `Amat` and the `KSP` solvers will automatically use that null space as needed during the solution process.
3331: If using `SNESComputeJacobianDefaultColor()` to assemble a Jacobian, the `ctx` argument
3332: must be a `MatFDColoring`.
3334: Other defect-correction schemes can be used by computing a different matrix in place of the Jacobian. One common
3335: example is to use the "Picard linearization" which only differentiates through the highest order parts of each term using `SNESSetPicard()`
3337: .seealso: [](ch_snes), `SNES`, `KSPSetOperators()`, `SNESSetFunction()`, `MatMFFDComputeJacobian()`, `SNESComputeJacobianDefaultColor()`, `MatStructure`,
3338: `SNESSetPicard()`, `SNESJacobianFn`, `SNESFunctionFn`
3339: @*/
3340: PetscErrorCode SNESSetJacobian(SNES snes, Mat Amat, Mat Pmat, SNESJacobianFn *J, PetscCtx ctx)
3341: {
3342: DM dm;
3344: PetscFunctionBegin;
3348: if (Amat) PetscCheckSameComm(snes, 1, Amat, 2);
3349: if (Pmat) PetscCheckSameComm(snes, 1, Pmat, 3);
3350: /* update DMSNES
3351: We support incremental information; so update the function context only if both Amat and Pmat are not specified
3352: (which allows to disable the callbacks when both J and ctx are NULL),
3353: or, if any of the mats is specified, when at least one of J and ctx is not NULL */
3354: PetscCall(SNESGetDM(snes, &dm));
3355: if ((!Amat && !Pmat) || J || ctx) PetscCall(DMSNESSetJacobian(dm, J, ctx));
3356: if (Amat) {
3357: PetscCall(PetscObjectReference((PetscObject)Amat));
3358: PetscCall(MatDestroy(&snes->jacobian));
3360: snes->jacobian = Amat;
3361: }
3362: if (Pmat) {
3363: PetscCall(PetscObjectReference((PetscObject)Pmat));
3364: PetscCall(MatDestroy(&snes->jacobian_pre));
3366: snes->jacobian_pre = Pmat;
3367: }
3368: PetscFunctionReturn(PETSC_SUCCESS);
3369: }
3371: /*@
3372: SNESGetJacobian - Returns the Jacobian matrix and optionally the user
3373: provided context for evaluating the Jacobian.
3375: Not Collective, but `Mat` object will be parallel if `SNES` is
3377: Input Parameter:
3378: . snes - the nonlinear solver context
3380: Output Parameters:
3381: + Amat - location to stash (approximate) Jacobian matrix (or `NULL`)
3382: . Pmat - location to stash matrix used to compute the preconditioner (or `NULL`)
3383: . J - location to put Jacobian function (or `NULL`), for calling sequence see `SNESJacobianFn`
3384: - ctx - location to stash Jacobian ctx (or `NULL`)
3386: Level: advanced
3388: .seealso: [](ch_snes), `SNES`, `Mat`, `SNESSetJacobian()`, `SNESComputeJacobian()`, `SNESJacobianFn`, `SNESGetFunction()`
3389: @*/
3390: PetscErrorCode SNESGetJacobian(SNES snes, Mat *Amat, Mat *Pmat, SNESJacobianFn **J, PetscCtxRt ctx)
3391: {
3392: DM dm;
3394: PetscFunctionBegin;
3396: if (Amat) *Amat = snes->jacobian;
3397: if (Pmat) *Pmat = snes->jacobian_pre;
3398: PetscCall(SNESGetDM(snes, &dm));
3399: PetscCall(DMSNESGetJacobian(dm, J, ctx));
3400: PetscFunctionReturn(PETSC_SUCCESS);
3401: }
3403: static PetscErrorCode SNESSetDefaultComputeJacobian(SNES snes)
3404: {
3405: DM dm;
3406: DMSNES sdm;
3408: PetscFunctionBegin;
3409: PetscCall(SNESGetDM(snes, &dm));
3410: PetscCall(DMGetDMSNES(dm, &sdm));
3411: if (!sdm->ops->computejacobian && snes->jacobian_pre) {
3412: DM dm;
3413: PetscBool isdense, ismf;
3415: PetscCall(SNESGetDM(snes, &dm));
3416: PetscCall(PetscObjectTypeCompareAny((PetscObject)snes->jacobian_pre, &isdense, MATSEQDENSE, MATMPIDENSE, MATDENSE, NULL));
3417: PetscCall(PetscObjectTypeCompareAny((PetscObject)snes->jacobian_pre, &ismf, MATMFFD, MATSHELL, NULL));
3418: if (isdense) PetscCall(DMSNESSetJacobian(dm, SNESComputeJacobianDefault, NULL));
3419: else if (!ismf) PetscCall(DMSNESSetJacobian(dm, SNESComputeJacobianDefaultColor, NULL));
3420: }
3421: PetscFunctionReturn(PETSC_SUCCESS);
3422: }
3424: /*@
3425: SNESSetUp - Sets up the internal data structures for the later use
3426: of a nonlinear solver `SNESSolve()`.
3428: Collective
3430: Input Parameter:
3431: . snes - the `SNES` context
3433: Level: advanced
3435: Note:
3436: For basic use of the `SNES` solvers the user does not need to explicitly call
3437: `SNESSetUp()`, since these actions will automatically occur during
3438: the call to `SNESSolve()`. However, if one wishes to control this
3439: phase separately, `SNESSetUp()` should be called after `SNESCreate()`
3440: and optional routines of the form SNESSetXXX(), but before `SNESSolve()`.
3442: .seealso: [](ch_snes), `SNES`, `SNESCreate()`, `SNESSolve()`, `SNESDestroy()`, `SNESSetFromOptions()`
3443: @*/
3444: PetscErrorCode SNESSetUp(SNES snes)
3445: {
3446: DM dm;
3447: DMSNES sdm;
3448: SNESLineSearch linesearch, pclinesearch;
3449: void *lsprectx, *lspostctx;
3450: PetscBool mf_operator, mf;
3451: Vec f, fpc;
3452: void *funcctx;
3453: void *jacctx, *appctx;
3454: Mat j, jpre;
3455: PetscErrorCode (*precheck)(SNESLineSearch, Vec, Vec, PetscBool *, PetscCtx);
3456: PetscErrorCode (*postcheck)(SNESLineSearch, Vec, Vec, Vec, PetscBool *, PetscBool *, PetscCtx);
3457: SNESFunctionFn *func;
3458: SNESJacobianFn *jac;
3460: PetscFunctionBegin;
3462: if (snes->setupcalled) PetscFunctionReturn(PETSC_SUCCESS);
3463: PetscCall(PetscLogEventBegin(SNES_SetUp, snes, 0, 0, 0));
3465: if (!((PetscObject)snes)->type_name) PetscCall(SNESSetType(snes, SNESNEWTONLS));
3467: PetscCall(SNESGetFunction(snes, &snes->vec_func, NULL, NULL));
3469: PetscCall(SNESGetDM(snes, &dm));
3470: PetscCall(DMGetDMSNES(dm, &sdm));
3471: PetscCall(SNESSetDefaultComputeJacobian(snes));
3473: if (!snes->vec_func) PetscCall(DMCreateGlobalVector(dm, &snes->vec_func));
3475: if (snes->usesksp && !snes->ksp) PetscCall(SNESGetKSP(snes, &snes->ksp));
3477: if (snes->linesearch) {
3478: PetscCall(SNESGetLineSearch(snes, &snes->linesearch));
3479: PetscCall(SNESLineSearchSetFunction(snes->linesearch, SNESComputeFunction));
3480: }
3482: PetscCall(SNESGetUseMatrixFree(snes, &mf_operator, &mf));
3483: if (snes->npc && snes->npcside == PC_LEFT) {
3484: snes->mf = PETSC_TRUE;
3485: snes->mf_operator = PETSC_FALSE;
3486: }
3487: if (snes->ops->ctxcompute && !snes->ctx) PetscCallBack("SNES callback compute application context", (*snes->ops->ctxcompute)(snes, &snes->ctx));
3488: if (snes->mf) PetscCall(SNESSetUpMatrixFree_Private(snes, snes->mf_operator, snes->mf_version));
3490: if (snes->npc) {
3491: SNESNormSchedule npc_norm_schedule;
3493: /* copy the DM over and the functions if NPC DM is not present */
3494: if (!snes->npc->dm) {
3495: PetscCall(SNESGetDM(snes, &dm));
3496: PetscCall(SNESSetDM(snes->npc, dm));
3498: PetscCall(SNESGetFunction(snes, &f, &func, &funcctx));
3499: PetscCall(VecDuplicate(f, &fpc));
3500: PetscCall(SNESSetFunction(snes->npc, fpc, func, funcctx));
3501: PetscCall(SNESGetJacobian(snes, &j, &jpre, &jac, &jacctx));
3502: PetscCall(SNESSetJacobian(snes->npc, j, jpre, jac, jacctx));
3503: PetscCall(SNESSetUseMatrixFree(snes->npc, mf_operator, mf));
3504: PetscCall(VecDestroy(&fpc));
3506: /* copy the function pointers over */
3507: PetscCall(PetscObjectCopyFortranFunctionPointers((PetscObject)snes, (PetscObject)snes->npc));
3509: /* Propagate app context if not present */
3510: PetscCall(SNESGetApplicationContext(snes->npc, &appctx));
3511: if (!appctx && !snes->npc->ops->ctxcompute) {
3512: if (snes->ops->ctxcompute) {
3513: PetscCall(SNESSetComputeApplicationContext(snes->npc, snes->ops->ctxcompute, snes->ops->ctxdestroy));
3514: } else {
3515: PetscCall(SNESGetApplicationContext(snes, &appctx));
3516: PetscCall(SNESSetApplicationContext(snes->npc, appctx));
3517: }
3518: }
3519: }
3521: /* Set default norm schedule for NPC if not yet set */
3522: PetscCall(SNESGetNormSchedule(snes->npc, &npc_norm_schedule));
3523: if (npc_norm_schedule == SNES_NORM_DEFAULT) {
3524: if (snes->npcside == PC_RIGHT) PetscCall(SNESSetNormSchedule(snes->npc, SNES_NORM_FINAL_ONLY));
3525: else if (snes->npcside == PC_LEFT) PetscCall(SNESSetNormSchedule(snes->npc, SNES_NORM_NONE));
3526: }
3528: PetscCall(SNESSetFromOptions(snes->npc));
3530: /* copy the line search context over */
3531: if (snes->dm == snes->npc->dm && snes->linesearch && snes->npc->linesearch) {
3532: PetscCall(SNESGetLineSearch(snes, &linesearch));
3533: PetscCall(SNESGetLineSearch(snes->npc, &pclinesearch));
3534: PetscCall(SNESLineSearchGetPreCheck(linesearch, &precheck, &lsprectx));
3535: PetscCall(SNESLineSearchGetPostCheck(linesearch, &postcheck, &lspostctx));
3536: PetscCall(SNESLineSearchSetPreCheck(pclinesearch, precheck, lsprectx));
3537: PetscCall(SNESLineSearchSetPostCheck(pclinesearch, postcheck, lspostctx));
3538: PetscCall(PetscObjectCopyFortranFunctionPointers((PetscObject)linesearch, (PetscObject)pclinesearch));
3539: }
3540: }
3542: snes->jac_iter = 0;
3543: snes->pre_iter = 0;
3545: PetscTryTypeMethod(snes, setup);
3547: PetscCall(SNESSetDefaultComputeJacobian(snes));
3549: if (snes->npc && snes->npcside == PC_LEFT) {
3550: if (snes->functype == SNES_FUNCTION_PRECONDITIONED) {
3551: if (snes->linesearch) {
3552: PetscCall(SNESGetLineSearch(snes, &linesearch));
3553: PetscCall(SNESLineSearchSetFunction(linesearch, SNESComputeFunctionDefaultNPC));
3554: }
3555: }
3556: }
3557: PetscCall(PetscLogEventEnd(SNES_SetUp, snes, 0, 0, 0));
3558: snes->setupcalled = PETSC_TRUE;
3559: PetscFunctionReturn(PETSC_SUCCESS);
3560: }
3562: /*@
3563: SNESReset - Resets a `SNES` context to the state it was in before `SNESSetUp()` was called and removes any allocated `Vec` and `Mat` from its data structures
3565: Collective
3567: Input Parameter:
3568: . snes - the nonlinear iterative solver context obtained from `SNESCreate()`
3570: Level: intermediate
3572: Notes:
3573: Any options set on the `SNES` object, including those set with `SNESSetFromOptions()` remain.
3575: Call this if you wish to reuse a `SNES` but with different size vectors
3577: Also calls the application context destroy routine set with `SNESSetComputeApplicationContext()`
3579: .seealso: [](ch_snes), `SNES`, `SNESDestroy()`, `SNESCreate()`, `SNESSetUp()`, `SNESSolve()`
3580: @*/
3581: PetscErrorCode SNESReset(SNES snes)
3582: {
3583: PetscFunctionBegin;
3585: if (snes->ops->ctxdestroy && snes->ctx) {
3586: PetscCallBack("SNES callback destroy application context", (*snes->ops->ctxdestroy)(&snes->ctx));
3587: snes->ctx = NULL;
3588: }
3589: if (snes->npc) PetscCall(SNESReset(snes->npc));
3591: PetscTryTypeMethod(snes, reset);
3592: if (snes->ksp) PetscCall(KSPReset(snes->ksp));
3594: if (snes->linesearch) PetscCall(SNESLineSearchReset(snes->linesearch));
3596: PetscCall(VecDestroy(&snes->vec_rhs));
3597: PetscCall(VecDestroy(&snes->vec_sol));
3598: PetscCall(VecDestroy(&snes->vec_sol_update));
3599: PetscCall(VecDestroy(&snes->vec_func));
3600: PetscCall(MatDestroy(&snes->jacobian));
3601: PetscCall(MatDestroy(&snes->jacobian_pre));
3602: PetscCall(MatDestroy(&snes->picard));
3603: PetscCall(VecDestroyVecs(snes->nwork, &snes->work));
3604: PetscCall(VecDestroyVecs(snes->nvwork, &snes->vwork));
3606: snes->alwayscomputesfinalresidual = PETSC_FALSE;
3608: snes->nwork = snes->nvwork = 0;
3609: snes->setupcalled = PETSC_FALSE;
3610: PetscFunctionReturn(PETSC_SUCCESS);
3611: }
3613: /*@
3614: SNESConvergedReasonViewCancel - Clears all the reason view functions for a `SNES` object provided with `SNESConvergedReasonViewSet()` also
3615: removes the default viewer.
3617: Collective
3619: Input Parameter:
3620: . snes - the nonlinear iterative solver context obtained from `SNESCreate()`
3622: Level: intermediate
3624: .seealso: [](ch_snes), `SNES`, `SNESCreate()`, `SNESDestroy()`, `SNESReset()`, `SNESConvergedReasonViewSet()`
3625: @*/
3626: PetscErrorCode SNESConvergedReasonViewCancel(SNES snes)
3627: {
3628: PetscFunctionBegin;
3630: for (PetscInt i = 0; i < snes->numberreasonviews; i++) {
3631: if (snes->reasonviewdestroy[i]) PetscCall((*snes->reasonviewdestroy[i])(&snes->reasonviewcontext[i]));
3632: }
3633: snes->numberreasonviews = 0;
3634: PetscCall(PetscViewerDestroy(&snes->convergedreasonviewer));
3635: PetscFunctionReturn(PETSC_SUCCESS);
3636: }
3638: /*@
3639: SNESDestroy - Destroys the nonlinear solver context that was created
3640: with `SNESCreate()`.
3642: Collective
3644: Input Parameter:
3645: . snes - the `SNES` context
3647: Level: beginner
3649: .seealso: [](ch_snes), `SNES`, `SNESCreate()`, `SNESSolve()`
3650: @*/
3651: PetscErrorCode SNESDestroy(SNES *snes)
3652: {
3653: DM dm;
3655: PetscFunctionBegin;
3656: if (!*snes) PetscFunctionReturn(PETSC_SUCCESS);
3658: if (--((PetscObject)*snes)->refct > 0) {
3659: *snes = NULL;
3660: PetscFunctionReturn(PETSC_SUCCESS);
3661: }
3663: PetscCall(SNESReset(*snes));
3664: PetscCall(SNESDestroy(&(*snes)->npc));
3666: /* if memory was published with SAWs then destroy it */
3667: PetscCall(PetscObjectSAWsViewOff((PetscObject)*snes));
3668: PetscTryTypeMethod(*snes, destroy);
3670: dm = (*snes)->dm;
3671: while (dm) {
3672: PetscCall(DMCoarsenHookRemove(dm, DMCoarsenHook_SNESVecSol, DMRestrictHook_SNESVecSol, *snes));
3673: PetscCall(DMGetCoarseDM(dm, &dm));
3674: }
3676: PetscCall(DMDestroy(&(*snes)->dm));
3677: PetscCall(KSPDestroy(&(*snes)->ksp));
3678: PetscCall(SNESLineSearchDestroy(&(*snes)->linesearch));
3680: PetscCall(PetscFree((*snes)->kspconvctx));
3681: if ((*snes)->ops->convergeddestroy) PetscCall((*(*snes)->ops->convergeddestroy)(&(*snes)->cnvP));
3682: if ((*snes)->conv_hist_alloc) PetscCall(PetscFree2((*snes)->conv_hist, (*snes)->conv_hist_its));
3683: PetscCall(SNESMonitorCancel(*snes));
3684: PetscCall(SNESConvergedReasonViewCancel(*snes));
3685: PetscCall(PetscHeaderDestroy(snes));
3686: PetscFunctionReturn(PETSC_SUCCESS);
3687: }
3689: /* ----------- Routines to set solver parameters ---------- */
3691: /*@
3692: SNESSetLagPreconditioner - Sets when the preconditioner is rebuilt in the nonlinear solve `SNESSolve()`.
3694: Logically Collective
3696: Input Parameters:
3697: + snes - the `SNES` context
3698: - lag - 1 means rebuild every time the Jacobian is computed within a single nonlinear solve, 2 means every second time
3699: the Jacobian is built etc. -2 indicates rebuild preconditioner at next chance but then never rebuild after that
3701: Options Database Keys:
3702: + -snes_lag_jacobian_persists (true|false) - sets the persistence through multiple `SNESSolve()`
3703: . -snes_lag_jacobian (-2|1|2|...) - sets the lag
3704: . -snes_lag_preconditioner_persists (true|false) - sets the persistence through multiple `SNESSolve()`
3705: - -snes_lag_preconditioner (-2|1|2|...) - sets the lag
3707: Level: intermediate
3709: Notes:
3710: The default is 1
3712: The preconditioner is ALWAYS built in the first iteration of a nonlinear solve unless lag is -1 or `SNESSetLagPreconditionerPersists()` was called
3714: `SNESSetLagPreconditionerPersists()` allows using the same uniform lagging (for example every second linear solve) across multiple nonlinear solves.
3716: .seealso: [](ch_snes), `SNESGetLagPreconditioner()`, `SNESSetLagJacobian()`, `SNESGetLagJacobian()`, `SNESSetLagPreconditionerPersists()`,
3717: `SNESSetLagJacobianPersists()`, `SNES`, `SNESSolve()`
3718: @*/
3719: PetscErrorCode SNESSetLagPreconditioner(SNES snes, PetscInt lag)
3720: {
3721: PetscFunctionBegin;
3724: PetscCheck(lag >= -2, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Lag must be -2, -1, 1 or greater");
3725: PetscCheck(lag, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Lag cannot be 0");
3726: snes->lagpreconditioner = lag;
3727: PetscFunctionReturn(PETSC_SUCCESS);
3728: }
3730: /*@
3731: SNESSetGridSequence - sets the number of steps of grid sequencing that `SNES` will do
3733: Logically Collective
3735: Input Parameters:
3736: + snes - the `SNES` context
3737: - steps - the number of refinements to do, defaults to 0
3739: Options Database Key:
3740: . -snes_grid_sequence steps - Use grid sequencing to generate initial guess
3742: Level: intermediate
3744: Notes:
3745: Once grid sequencing is turned on `SNESSolve()` will automatically perform the solve on each grid refinement.
3747: Use `SNESGetSolution()` to extract the fine grid solution after grid sequencing.
3749: .seealso: [](ch_snes), `SNES`, `SNESGetLagPreconditioner()`, `SNESSetLagJacobian()`, `SNESGetLagJacobian()`, `SNESGetGridSequence()`,
3750: `SNESSetDM()`, `SNESSolve()`
3751: @*/
3752: PetscErrorCode SNESSetGridSequence(SNES snes, PetscInt steps)
3753: {
3754: PetscFunctionBegin;
3757: snes->gridsequence = steps;
3758: PetscFunctionReturn(PETSC_SUCCESS);
3759: }
3761: /*@
3762: SNESGetGridSequence - gets the number of steps of grid sequencing that `SNES` will do
3764: Logically Collective
3766: Input Parameter:
3767: . snes - the `SNES` context
3769: Output Parameter:
3770: . steps - the number of refinements to do, defaults to 0
3772: Level: intermediate
3774: .seealso: [](ch_snes), `SNESGetLagPreconditioner()`, `SNESSetLagJacobian()`, `SNESGetLagJacobian()`, `SNESSetGridSequence()`
3775: @*/
3776: PetscErrorCode SNESGetGridSequence(SNES snes, PetscInt *steps)
3777: {
3778: PetscFunctionBegin;
3780: *steps = snes->gridsequence;
3781: PetscFunctionReturn(PETSC_SUCCESS);
3782: }
3784: /*@
3785: SNESGetLagPreconditioner - Return how often the preconditioner is rebuilt
3787: Not Collective
3789: Input Parameter:
3790: . snes - the `SNES` context
3792: Output Parameter:
3793: . lag - -1 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is computed within a single nonlinear solve, 2 means every second time
3794: the Jacobian is built etc. -2 indicates rebuild preconditioner at next chance but then never rebuild after that
3796: Level: intermediate
3798: Notes:
3799: The default is 1
3801: The preconditioner is ALWAYS built in the first iteration of a nonlinear solve unless lag is -1
3803: .seealso: [](ch_snes), `SNES`, `SNESSetLagPreconditioner()`, `SNESSetLagJacobianPersists()`, `SNESSetLagPreconditionerPersists()`
3804: @*/
3805: PetscErrorCode SNESGetLagPreconditioner(SNES snes, PetscInt *lag)
3806: {
3807: PetscFunctionBegin;
3809: *lag = snes->lagpreconditioner;
3810: PetscFunctionReturn(PETSC_SUCCESS);
3811: }
3813: /*@
3814: SNESSetLagJacobian - Set when the Jacobian is rebuilt in the nonlinear solve. See `SNESSetLagPreconditioner()` for determining how
3815: often the preconditioner is rebuilt.
3817: Logically Collective
3819: Input Parameters:
3820: + snes - the `SNES` context
3821: - lag - -1 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is computed within a single nonlinear solve, 2 means every second time
3822: the Jacobian is built etc. -2 means rebuild at next chance but then never again
3824: Options Database Keys:
3825: + -snes_lag_jacobian_persists (true|false) - sets the persistence through multiple SNES solves
3826: . -snes_lag_jacobian (-2|1|2|...) - sets the lag
3827: . -snes_lag_preconditioner_persists (true|false) - sets the persistence through multiple SNES solves
3828: - -snes_lag_preconditioner (-2|1|2|...) - sets the lag.
3830: Level: intermediate
3832: Notes:
3833: The default is 1
3835: The Jacobian is ALWAYS built in the first iteration of a nonlinear solve unless lag is -1
3837: If -1 is used before the very first nonlinear solve the CODE WILL FAIL! because no Jacobian is used, use -2 to indicate you want it recomputed
3838: at the next Newton step but never again (unless it is reset to another value)
3840: .seealso: [](ch_snes), `SNES`, `SNESGetLagPreconditioner()`, `SNESSetLagPreconditioner()`, `SNESGetLagJacobianPersists()`, `SNESSetLagPreconditionerPersists()`
3841: @*/
3842: PetscErrorCode SNESSetLagJacobian(SNES snes, PetscInt lag)
3843: {
3844: PetscFunctionBegin;
3846: PetscCheck(lag >= -2, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Lag must be -2, -1, 1 or greater");
3847: PetscCheck(lag, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Lag cannot be 0");
3849: snes->lagjacobian = lag;
3850: PetscFunctionReturn(PETSC_SUCCESS);
3851: }
3853: /*@
3854: SNESGetLagJacobian - Get how often the Jacobian is rebuilt. See `SNESGetLagPreconditioner()` to determine when the preconditioner is rebuilt
3856: Not Collective
3858: Input Parameter:
3859: . snes - the `SNES` context
3861: Output Parameter:
3862: . lag - -1 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is computed within a single nonlinear solve, 2 means every second time
3863: the Jacobian is built etc.
3865: Level: intermediate
3867: Notes:
3868: The default is 1
3870: The jacobian is ALWAYS built in the first iteration of a nonlinear solve unless lag is -1 or `SNESSetLagJacobianPersists()` was called.
3872: .seealso: [](ch_snes), `SNES`, `SNESSetLagJacobian()`, `SNESSetLagPreconditioner()`, `SNESGetLagPreconditioner()`, `SNESSetLagJacobianPersists()`, `SNESSetLagPreconditionerPersists()`
3873: @*/
3874: PetscErrorCode SNESGetLagJacobian(SNES snes, PetscInt *lag)
3875: {
3876: PetscFunctionBegin;
3878: *lag = snes->lagjacobian;
3879: PetscFunctionReturn(PETSC_SUCCESS);
3880: }
3882: /*@
3883: SNESSetLagJacobianPersists - Set whether or not the Jacobian lagging persists through multiple nonlinear solves
3885: Logically collective
3887: Input Parameters:
3888: + snes - the `SNES` context
3889: - flg - jacobian lagging persists if true
3891: Options Database Keys:
3892: + -snes_lag_jacobian_persists (true|false) - sets the persistence through multiple SNES solves
3893: . -snes_lag_jacobian (-2|1|2|...) - sets the lag
3894: . -snes_lag_preconditioner_persists (true|false) - sets the persistence through multiple SNES solves
3895: - -snes_lag_preconditioner (-2|1|2|...) - sets the lag
3897: Level: advanced
3899: Notes:
3900: Normally when `SNESSetLagJacobian()` is used, the Jacobian is always rebuilt at the beginning of each new nonlinear solve, this removes that behavior
3902: This is useful both for nonlinear preconditioning, where it's appropriate to have the Jacobian be stale by
3903: several solves, and for implicit time-stepping, where Jacobian lagging in the inner nonlinear solve over several
3904: timesteps may present huge efficiency gains.
3906: .seealso: [](ch_snes), `SNES`, `SNESSetLagPreconditionerPersists()`, `SNESSetLagJacobian()`, `SNESGetLagJacobian()`, `SNESGetNPC()`
3907: @*/
3908: PetscErrorCode SNESSetLagJacobianPersists(SNES snes, PetscBool flg)
3909: {
3910: PetscFunctionBegin;
3913: snes->lagjac_persist = flg;
3914: PetscFunctionReturn(PETSC_SUCCESS);
3915: }
3917: /*@
3918: SNESSetLagPreconditionerPersists - Set whether or not the preconditioner lagging persists through multiple nonlinear solves
3920: Logically Collective
3922: Input Parameters:
3923: + snes - the `SNES` context
3924: - flg - preconditioner lagging persists if true
3926: Options Database Keys:
3927: + -snes_lag_jacobian_persists (true|false) - sets the persistence through multiple SNES solves
3928: . -snes_lag_jacobian (-2|1|2|...) - sets the lag
3929: . -snes_lag_preconditioner_persists (true|false) - sets the persistence through multiple SNES solves
3930: - -snes_lag_preconditioner (-2|1|2|...) - sets the lag
3932: Level: developer
3934: Notes:
3935: Normally when `SNESSetLagPreconditioner()` is used, the preconditioner is always rebuilt at the beginning of each new nonlinear solve, this removes that behavior
3937: This is useful both for nonlinear preconditioning, where it's appropriate to have the preconditioner be stale
3938: by several solves, and for implicit time-stepping, where preconditioner lagging in the inner nonlinear solve over
3939: several timesteps may present huge efficiency gains.
3941: .seealso: [](ch_snes), `SNES`, `SNESSetLagJacobianPersists()`, `SNESSetLagJacobian()`, `SNESGetLagJacobian()`, `SNESGetNPC()`, `SNESSetLagPreconditioner()`
3942: @*/
3943: PetscErrorCode SNESSetLagPreconditionerPersists(SNES snes, PetscBool flg)
3944: {
3945: PetscFunctionBegin;
3948: snes->lagpre_persist = flg;
3949: PetscFunctionReturn(PETSC_SUCCESS);
3950: }
3952: /*@
3953: SNESSetForceIteration - force `SNESSolve()` to take at least one iteration regardless of the initial residual norm
3955: Logically Collective
3957: Input Parameters:
3958: + snes - the `SNES` context
3959: - force - `PETSC_TRUE` require at least one iteration
3961: Options Database Key:
3962: . -snes_force_iteration force - Sets forcing an iteration
3964: Level: intermediate
3966: Note:
3967: This is used sometimes with `TS` to prevent `TS` from detecting a false steady state solution
3969: .seealso: [](ch_snes), `SNES`, `TS`, `SNESSetDivergenceTolerance()`
3970: @*/
3971: PetscErrorCode SNESSetForceIteration(SNES snes, PetscBool force)
3972: {
3973: PetscFunctionBegin;
3975: snes->forceiteration = force;
3976: PetscFunctionReturn(PETSC_SUCCESS);
3977: }
3979: /*@
3980: SNESGetForceIteration - Check whether or not `SNESSolve()` take at least one iteration regardless of the initial residual norm
3982: Logically Collective
3984: Input Parameter:
3985: . snes - the `SNES` context
3987: Output Parameter:
3988: . force - `PETSC_TRUE` requires at least one iteration.
3990: Level: intermediate
3992: .seealso: [](ch_snes), `SNES`, `SNESSetForceIteration()`, `SNESSetDivergenceTolerance()`
3993: @*/
3994: PetscErrorCode SNESGetForceIteration(SNES snes, PetscBool *force)
3995: {
3996: PetscFunctionBegin;
3998: *force = snes->forceiteration;
3999: PetscFunctionReturn(PETSC_SUCCESS);
4000: }
4002: /*@
4003: SNESSetTolerances - Sets various parameters used in `SNES` convergence tests.
4005: Logically Collective
4007: Input Parameters:
4008: + snes - the `SNES` context
4009: . abstol - the absolute convergence tolerance, $F(x^n) \le abstol$
4010: . rtol - the relative convergence tolerance, $F(x^n) \le reltol \, F(x^0)$
4011: . stol - convergence tolerance in terms of the norm of the change in the solution between steps, || delta x || < stol*|| x ||
4012: . maxit - the maximum number of iterations allowed in the solver, default 50.
4013: - maxf - the maximum number of function evaluations allowed in the solver (use `PETSC_UNLIMITED` indicates no limit), default 10,000
4015: Options Database Keys:
4016: + -snes_atol abstol - Sets `abstol`
4017: . -snes_rtol rtol - Sets `rtol`
4018: . -snes_stol stol - Sets `stol`
4019: . -snes_max_it maxit - Sets `maxit`
4020: - -snes_max_funcs maxf - Sets `maxf` (use `unlimited` to have no maximum)
4022: Level: intermediate
4024: Note:
4025: All parameters must be non-negative
4027: Use `PETSC_CURRENT` to retain the current value of any parameter and `PETSC_DETERMINE` to use the default value for the given `SNES`.
4028: The default value is the value in the object when its type is set.
4030: Use `PETSC_UNLIMITED` on `maxit` or `maxf` to indicate there is no bound on the number of iterations or number of function evaluations.
4032: Fortran Note:
4033: Use `PETSC_CURRENT_INTEGER`, `PETSC_CURRENT_REAL`, `PETSC_UNLIMITED_INTEGER`, `PETSC_DETERMINE_INTEGER`, or `PETSC_DETERMINE_REAL`
4035: .seealso: [](ch_snes), `SNESSolve()`, `SNES`, `SNESSetDivergenceTolerance()`, `SNESSetForceIteration()`
4036: @*/
4037: PetscErrorCode SNESSetTolerances(SNES snes, PetscReal abstol, PetscReal rtol, PetscReal stol, PetscInt maxit, PetscInt maxf)
4038: {
4039: PetscFunctionBegin;
4047: if (abstol == (PetscReal)PETSC_DETERMINE) {
4048: snes->abstol = snes->default_abstol;
4049: } else if (abstol != (PetscReal)PETSC_CURRENT) {
4050: PetscCheck(abstol >= 0.0, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_OUTOFRANGE, "Absolute tolerance %g must be non-negative", (double)abstol);
4051: snes->abstol = abstol;
4052: }
4054: if (rtol == (PetscReal)PETSC_DETERMINE) {
4055: snes->rtol = snes->default_rtol;
4056: } else if (rtol != (PetscReal)PETSC_CURRENT) {
4057: PetscCheck(rtol >= 0.0 && 1.0 > rtol, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_OUTOFRANGE, "Relative tolerance %g must be non-negative and less than 1.0", (double)rtol);
4058: snes->rtol = rtol;
4059: }
4061: if (stol == (PetscReal)PETSC_DETERMINE) {
4062: snes->stol = snes->default_stol;
4063: } else if (stol != (PetscReal)PETSC_CURRENT) {
4064: PetscCheck(stol >= 0.0, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_OUTOFRANGE, "Step tolerance %g must be non-negative", (double)stol);
4065: snes->stol = stol;
4066: }
4068: if (maxit == PETSC_DETERMINE) {
4069: snes->max_its = snes->default_max_its;
4070: } else if (maxit == PETSC_UNLIMITED) {
4071: snes->max_its = PETSC_INT_MAX;
4072: } else if (maxit != PETSC_CURRENT) {
4073: PetscCheck(maxit >= 0, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_OUTOFRANGE, "Maximum number of iterations %" PetscInt_FMT " must be non-negative", maxit);
4074: snes->max_its = maxit;
4075: }
4077: if (maxf == PETSC_DETERMINE) {
4078: snes->max_funcs = snes->default_max_funcs;
4079: } else if (maxf == PETSC_UNLIMITED || maxf == -1) {
4080: snes->max_funcs = PETSC_UNLIMITED;
4081: } else if (maxf != PETSC_CURRENT) {
4082: PetscCheck(maxf >= 0, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_OUTOFRANGE, "Maximum number of function evaluations %" PetscInt_FMT " must be nonnegative", maxf);
4083: snes->max_funcs = maxf;
4084: }
4085: PetscFunctionReturn(PETSC_SUCCESS);
4086: }
4088: /*@
4089: SNESSetDivergenceTolerance - Sets the divergence tolerance used for the `SNES` divergence test.
4091: Logically Collective
4093: Input Parameters:
4094: + snes - the `SNES` context
4095: - divtol - the divergence tolerance. Use `PETSC_UNLIMITED` to deactivate the test. If the residual norm $F(x^n) \ge divtol \, F(x^0)$ the solver
4096: is stopped due to divergence.
4098: Options Database Key:
4099: . -snes_divergence_tolerance divtol - Sets `divtol`
4101: Level: intermediate
4103: Notes:
4104: Use `PETSC_DETERMINE` to use the default value from when the object's type was set.
4106: Fortran Note:
4107: Use ``PETSC_DETERMINE_REAL` or `PETSC_UNLIMITED_REAL`
4109: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESSetTolerances()`, `SNESGetDivergenceTolerance()`
4110: @*/
4111: PetscErrorCode SNESSetDivergenceTolerance(SNES snes, PetscReal divtol)
4112: {
4113: PetscFunctionBegin;
4117: if (divtol == (PetscReal)PETSC_DETERMINE) {
4118: snes->divtol = snes->default_divtol;
4119: } else if (divtol == (PetscReal)PETSC_UNLIMITED || divtol == -1) {
4120: snes->divtol = PETSC_UNLIMITED;
4121: } else if (divtol != (PetscReal)PETSC_CURRENT) {
4122: PetscCheck(divtol >= 1.0, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_OUTOFRANGE, "Divergence tolerance %g must be greater than 1.0", (double)divtol);
4123: snes->divtol = divtol;
4124: }
4125: PetscFunctionReturn(PETSC_SUCCESS);
4126: }
4128: /*@
4129: SNESGetTolerances - Gets various parameters used in `SNES` convergence tests.
4131: Not Collective
4133: Input Parameter:
4134: . snes - the `SNES` context
4136: Output Parameters:
4137: + atol - the absolute convergence tolerance
4138: . rtol - the relative convergence tolerance
4139: . stol - convergence tolerance in terms of the norm of the change in the solution between steps
4140: . maxit - the maximum number of iterations allowed
4141: - maxf - the maximum number of function evaluations allowed, `PETSC_UNLIMITED` indicates no bound
4143: Level: intermediate
4145: Notes:
4146: See `SNESSetTolerances()` for details on the parameters.
4148: The user can specify `NULL` for any parameter that is not needed.
4150: .seealso: [](ch_snes), `SNES`, `SNESSetTolerances()`
4151: @*/
4152: PetscErrorCode SNESGetTolerances(SNES snes, PetscReal *atol, PetscReal *rtol, PetscReal *stol, PetscInt *maxit, PetscInt *maxf)
4153: {
4154: PetscFunctionBegin;
4156: if (atol) *atol = snes->abstol;
4157: if (rtol) *rtol = snes->rtol;
4158: if (stol) *stol = snes->stol;
4159: if (maxit) *maxit = snes->max_its;
4160: if (maxf) *maxf = snes->max_funcs;
4161: PetscFunctionReturn(PETSC_SUCCESS);
4162: }
4164: /*@
4165: SNESGetDivergenceTolerance - Gets divergence tolerance used in divergence test.
4167: Not Collective
4169: Input Parameters:
4170: + snes - the `SNES` context
4171: - divtol - divergence tolerance
4173: Level: intermediate
4175: .seealso: [](ch_snes), `SNES`, `SNESSetDivergenceTolerance()`
4176: @*/
4177: PetscErrorCode SNESGetDivergenceTolerance(SNES snes, PetscReal *divtol)
4178: {
4179: PetscFunctionBegin;
4181: if (divtol) *divtol = snes->divtol;
4182: PetscFunctionReturn(PETSC_SUCCESS);
4183: }
4185: PETSC_INTERN PetscErrorCode SNESMonitorRange_Private(SNES, PetscInt, PetscReal *);
4187: /*@
4188: SNESMonitorLGRange - Line-graph monitor that plots the residual norm together with residual-range statistics for a `SNESSolve()`
4190: Collective
4192: Input Parameters:
4193: + snes - the `SNES` context
4194: . n - the iteration number
4195: . rnorm - the 2-norm of the residual
4196: - monctx - a `PetscViewer` of type `PETSCVIEWERDRAW` set up with `PetscViewerMonitorLGSetUp()`
4198: Level: intermediate
4200: Note:
4201: Plots four line graphs in the viewer: the residual norm (log scale), the fraction of residual entries larger than 20% of the maximum entry, the relative decrease `(prev - rnorm)/prev`, and their product.
4203: .seealso: [](ch_snes), `SNES`, `SNESMonitorSet()`, `SNESMonitorDefault()`, `PetscViewerDrawGetDrawLG()`, `PetscDrawLG`
4204: @*/
4205: PetscErrorCode SNESMonitorLGRange(SNES snes, PetscInt n, PetscReal rnorm, PetscCtx monctx)
4206: {
4207: PetscDrawLG lg;
4208: PetscReal x, y, per;
4209: PetscViewer v = (PetscViewer)monctx;
4210: static PetscReal prev; /* should be in the context */
4211: PetscDraw draw;
4213: PetscFunctionBegin;
4215: PetscCall(PetscViewerDrawGetDrawLG(v, 0, &lg));
4216: if (!n) PetscCall(PetscDrawLGReset(lg));
4217: PetscCall(PetscDrawLGGetDraw(lg, &draw));
4218: PetscCall(PetscDrawSetTitle(draw, "Residual norm"));
4219: x = (PetscReal)n;
4220: if (rnorm > 0.0) y = PetscLog10Real(rnorm);
4221: else y = -15.0;
4222: PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
4223: if (n < 20 || !(n % 5) || snes->reason) {
4224: PetscCall(PetscDrawLGDraw(lg));
4225: PetscCall(PetscDrawLGSave(lg));
4226: }
4228: PetscCall(PetscViewerDrawGetDrawLG(v, 1, &lg));
4229: if (!n) PetscCall(PetscDrawLGReset(lg));
4230: PetscCall(PetscDrawLGGetDraw(lg, &draw));
4231: PetscCall(PetscDrawSetTitle(draw, "% elements > .2*max element"));
4232: PetscCall(SNESMonitorRange_Private(snes, n, &per));
4233: x = (PetscReal)n;
4234: y = 100.0 * per;
4235: PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
4236: if (n < 20 || !(n % 5) || snes->reason) {
4237: PetscCall(PetscDrawLGDraw(lg));
4238: PetscCall(PetscDrawLGSave(lg));
4239: }
4241: PetscCall(PetscViewerDrawGetDrawLG(v, 2, &lg));
4242: if (!n) {
4243: prev = rnorm;
4244: PetscCall(PetscDrawLGReset(lg));
4245: }
4246: PetscCall(PetscDrawLGGetDraw(lg, &draw));
4247: PetscCall(PetscDrawSetTitle(draw, "(norm -oldnorm)/oldnorm"));
4248: x = (PetscReal)n;
4249: y = (prev - rnorm) / prev;
4250: PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
4251: if (n < 20 || !(n % 5) || snes->reason) {
4252: PetscCall(PetscDrawLGDraw(lg));
4253: PetscCall(PetscDrawLGSave(lg));
4254: }
4256: PetscCall(PetscViewerDrawGetDrawLG(v, 3, &lg));
4257: if (!n) PetscCall(PetscDrawLGReset(lg));
4258: PetscCall(PetscDrawLGGetDraw(lg, &draw));
4259: PetscCall(PetscDrawSetTitle(draw, "(norm -oldnorm)/oldnorm*(% > .2 max)"));
4260: x = (PetscReal)n;
4261: y = (prev - rnorm) / (prev * per);
4262: if (n > 2) { /*skip initial crazy value */
4263: PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
4264: }
4265: if (n < 20 || !(n % 5) || snes->reason) {
4266: PetscCall(PetscDrawLGDraw(lg));
4267: PetscCall(PetscDrawLGSave(lg));
4268: }
4269: prev = rnorm;
4270: PetscFunctionReturn(PETSC_SUCCESS);
4271: }
4273: /*@
4274: SNESConverged - Run the convergence test and update the `SNESConvergedReason`.
4276: Collective
4278: Input Parameters:
4279: + snes - the `SNES` context
4280: . it - current iteration
4281: . xnorm - 2-norm of current iterate
4282: . snorm - 2-norm of current step
4283: - fnorm - 2-norm of function
4285: Level: developer
4287: Note:
4288: This routine is called by the `SNESSolve()` implementations.
4289: It does not typically need to be called by the user.
4291: .seealso: [](ch_snes), `SNES`, `SNESSolve`, `SNESSetConvergenceTest()`
4292: @*/
4293: PetscErrorCode SNESConverged(SNES snes, PetscInt it, PetscReal xnorm, PetscReal snorm, PetscReal fnorm)
4294: {
4295: PetscFunctionBegin;
4296: if (!snes->reason) {
4297: if (snes->normschedule == SNES_NORM_ALWAYS) PetscUseTypeMethod(snes, converged, it, xnorm, snorm, fnorm, &snes->reason, snes->cnvP);
4298: if (it == snes->max_its && !snes->reason) {
4299: if (snes->normschedule == SNES_NORM_ALWAYS) {
4300: PetscCall(PetscInfo(snes, "Maximum number of iterations has been reached: %" PetscInt_FMT "\n", snes->max_its));
4301: snes->reason = SNES_DIVERGED_MAX_IT;
4302: } else snes->reason = SNES_CONVERGED_ITS;
4303: }
4304: }
4305: PetscFunctionReturn(PETSC_SUCCESS);
4306: }
4308: /*@
4309: SNESMonitor - runs any `SNES` monitor routines provided with `SNESMonitor()` or the options database
4311: Collective
4313: Input Parameters:
4314: + snes - nonlinear solver context obtained from `SNESCreate()`
4315: . iter - current iteration number
4316: - rnorm - current relative norm of the residual
4318: Level: developer
4320: Note:
4321: This routine is called by the `SNESSolve()` implementations.
4322: It does not typically need to be called by the user.
4324: .seealso: [](ch_snes), `SNES`, `SNESMonitorSet()`
4325: @*/
4326: PetscErrorCode SNESMonitor(SNES snes, PetscInt iter, PetscReal rnorm)
4327: {
4328: PetscInt i, n = snes->numbermonitors;
4330: PetscFunctionBegin;
4331: PetscCall(VecLockReadPush(snes->vec_sol));
4332: for (i = 0; i < n; i++) PetscCall((*snes->monitor[i])(snes, iter, rnorm, snes->monitorcontext[i]));
4333: PetscCall(VecLockReadPop(snes->vec_sol));
4334: PetscFunctionReturn(PETSC_SUCCESS);
4335: }
4337: /* ------------ Routines to set performance monitoring options ----------- */
4339: /*MC
4340: SNESMonitorFunction - functional form passed to `SNESMonitorSet()` to monitor convergence of nonlinear solver
4342: Synopsis:
4343: #include <petscsnes.h>
4344: PetscErrorCode SNESMonitorFunction(SNES snes, PetscInt its, PetscReal norm, PetscCtx mctx)
4346: Collective
4348: Input Parameters:
4349: + snes - the `SNES` context
4350: . its - iteration number
4351: . norm - 2-norm function value (may be estimated)
4352: - mctx - [optional] monitoring context
4354: Level: advanced
4356: .seealso: [](ch_snes), `SNESMonitorSet()`, `PetscCtx`
4357: M*/
4359: /*@
4360: SNESMonitorSet - Sets an ADDITIONAL function that is to be used at every
4361: iteration of the `SNES` nonlinear solver to display the iteration's
4362: progress.
4364: Logically Collective
4366: Input Parameters:
4367: + snes - the `SNES` context
4368: . f - the monitor function, for the calling sequence see `SNESMonitorFunction`
4369: . mctx - [optional] user-defined context for private data for the monitor routine (use `NULL` if no context is desired)
4370: - monitordestroy - [optional] routine that frees monitor context (may be `NULL`), see `PetscCtxDestroyFn` for the calling sequence
4372: Calling sequence of f:
4373: + snes - the `SNES` object
4374: . it - the current iteration
4375: . rnorm - norm of the residual
4376: - mctx - the optional monitor context
4378: Options Database Keys:
4379: + -snes_monitor - sets `SNESMonitorDefault()`
4380: . -snes_monitor draw::draw_lg - sets line graph monitor
4381: - -snes_monitor_cancel - cancels all monitors that have been hardwired into a code by calls to `SNESMonitorSet()`, but does not cancel those set via
4382: the options database.
4384: Level: intermediate
4386: Note:
4387: Several different monitoring routines may be set by calling
4388: `SNESMonitorSet()` multiple times; all will be called in the
4389: order in which they were set.
4391: Fortran Note:
4392: Only a single monitor function can be set for each `SNES` object
4394: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESMonitorDefault()`, `SNESMonitorCancel()`, `SNESMonitorFunction`, `PetscCtxDestroyFn`
4395: @*/
4396: PetscErrorCode SNESMonitorSet(SNES snes, PetscErrorCode (*f)(SNES snes, PetscInt it, PetscReal rnorm, PetscCtx mctx), PetscCtx mctx, PetscCtxDestroyFn *monitordestroy)
4397: {
4398: PetscFunctionBegin;
4400: for (PetscInt i = 0; i < snes->numbermonitors; i++) {
4401: PetscBool identical;
4403: PetscCall(PetscMonitorCompare((PetscErrorCode (*)(void))(PetscVoidFn *)f, mctx, monitordestroy, (PetscErrorCode (*)(void))(PetscVoidFn *)snes->monitor[i], snes->monitorcontext[i], snes->monitordestroy[i], &identical));
4404: if (identical) PetscFunctionReturn(PETSC_SUCCESS);
4405: }
4406: PetscCheck(snes->numbermonitors < MAXSNESMONITORS, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Too many monitors set");
4407: snes->monitor[snes->numbermonitors] = f;
4408: snes->monitordestroy[snes->numbermonitors] = monitordestroy;
4409: snes->monitorcontext[snes->numbermonitors++] = mctx;
4410: PetscFunctionReturn(PETSC_SUCCESS);
4411: }
4413: /*@
4414: SNESMonitorCancel - Clears all the monitor functions for a `SNES` object.
4416: Logically Collective
4418: Input Parameter:
4419: . snes - the `SNES` context
4421: Options Database Key:
4422: . -snes_monitor_cancel - cancels all monitors that have been hardwired
4423: into a code by calls to `SNESMonitorSet()`, but does not cancel those
4424: set via the options database
4426: Level: intermediate
4428: Note:
4429: There is no way to clear one specific monitor from a `SNES` object.
4431: .seealso: [](ch_snes), `SNES`, `SNESMonitorDefault()`, `SNESMonitorSet()`
4432: @*/
4433: PetscErrorCode SNESMonitorCancel(SNES snes)
4434: {
4435: PetscFunctionBegin;
4437: for (PetscInt i = 0; i < snes->numbermonitors; i++) {
4438: if (snes->monitordestroy[i]) PetscCall((*snes->monitordestroy[i])(&snes->monitorcontext[i]));
4439: }
4440: snes->numbermonitors = 0;
4441: PetscFunctionReturn(PETSC_SUCCESS);
4442: }
4444: /*@
4445: SNESSetConvergenceTest - Sets the function that is to be used
4446: to test for convergence of the nonlinear iterative solution.
4448: Logically Collective
4450: Input Parameters:
4451: + snes - the `SNES` context
4452: . func - routine to test for convergence
4453: . ctx - [optional] context for private data for the convergence routine (may be `NULL`)
4454: - destroy - [optional] destructor for the context (may be `NULL`; `PETSC_NULL_FUNCTION` in Fortran)
4456: Calling sequence of func:
4457: + snes - the `SNES` context
4458: . it - the current iteration number
4459: . xnorm - the norm of the new solution
4460: . snorm - the norm of the step
4461: . fnorm - the norm of the function value
4462: . reason - output, the reason convergence or divergence as declared
4463: - ctx - the optional convergence test context
4465: Level: advanced
4467: .seealso: [](ch_snes), `SNES`, `SNESConvergedDefault()`, `SNESConvergedSkip()`
4468: @*/
4469: PetscErrorCode SNESSetConvergenceTest(SNES snes, PetscErrorCode (*func)(SNES snes, PetscInt it, PetscReal xnorm, PetscReal snorm, PetscReal fnorm, SNESConvergedReason *reason, PetscCtx ctx), PetscCtx ctx, PetscCtxDestroyFn *destroy)
4470: {
4471: PetscFunctionBegin;
4473: if (!func) func = SNESConvergedSkip;
4474: if (snes->ops->convergeddestroy) PetscCall((*snes->ops->convergeddestroy)(&snes->cnvP));
4475: snes->ops->converged = func;
4476: snes->ops->convergeddestroy = destroy;
4477: snes->cnvP = ctx;
4478: PetscFunctionReturn(PETSC_SUCCESS);
4479: }
4481: /*@
4482: SNESGetConvergedReason - Gets the reason the `SNES` iteration was stopped, which may be due to convergence, divergence, or stagnation
4484: Not Collective
4486: Input Parameter:
4487: . snes - the `SNES` context
4489: Output Parameter:
4490: . reason - negative value indicates diverged, positive value converged, see `SNESConvergedReason` for the individual convergence tests for complete lists
4492: Options Database Key:
4493: . -snes_converged_reason - prints the reason to standard out
4495: Level: intermediate
4497: Note:
4498: Should only be called after the call the `SNESSolve()` is complete, if it is called earlier it returns the value `SNES__CONVERGED_ITERATING`.
4500: .seealso: [](ch_snes), `SNESSolve()`, `SNESSetConvergenceTest()`, `SNESSetConvergedReason()`, `SNESConvergedReason`, `SNESGetConvergedReasonString()`
4501: @*/
4502: PetscErrorCode SNESGetConvergedReason(SNES snes, SNESConvergedReason *reason)
4503: {
4504: PetscFunctionBegin;
4506: PetscAssertPointer(reason, 2);
4507: *reason = snes->reason;
4508: PetscFunctionReturn(PETSC_SUCCESS);
4509: }
4511: /*@
4512: SNESGetConvergedReasonString - Return a human readable string for `SNESConvergedReason`
4514: Not Collective
4516: Input Parameter:
4517: . snes - the `SNES` context
4519: Output Parameter:
4520: . strreason - a human readable string that describes `SNES` converged reason
4522: Level: beginner
4524: .seealso: [](ch_snes), `SNES`, `SNESGetConvergedReason()`
4525: @*/
4526: PetscErrorCode SNESGetConvergedReasonString(SNES snes, const char *strreason[])
4527: {
4528: PetscFunctionBegin;
4530: PetscAssertPointer(strreason, 2);
4531: *strreason = SNESConvergedReasons[snes->reason];
4532: PetscFunctionReturn(PETSC_SUCCESS);
4533: }
4535: /*@
4536: SNESSetConvergedReason - Sets the reason the `SNES` iteration was stopped.
4538: Logically Collective
4540: Input Parameters:
4541: + snes - the `SNES` context
4542: - reason - negative value indicates diverged, positive value converged, see `SNESConvergedReason` or the
4543: manual pages for the individual convergence tests for complete lists
4545: Level: developer
4547: Developer Note:
4548: Called inside the various `SNESSolve()` implementations
4550: .seealso: [](ch_snes), `SNESGetConvergedReason()`, `SNESSetConvergenceTest()`, `SNESConvergedReason`
4551: @*/
4552: PetscErrorCode SNESSetConvergedReason(SNES snes, SNESConvergedReason reason)
4553: {
4554: PetscFunctionBegin;
4557: snes->reason = reason;
4558: PetscFunctionReturn(PETSC_SUCCESS);
4559: }
4561: /*@
4562: SNESSetConvergenceHistory - Sets the arrays used to hold the convergence history.
4564: Logically Collective
4566: Input Parameters:
4567: + snes - iterative context obtained from `SNESCreate()`
4568: . a - array to hold history, this array will contain the function norms computed at each step
4569: . its - integer array holds the number of linear iterations for each solve.
4570: . na - size of `a` and `its`
4571: - reset - `PETSC_TRUE` indicates each new nonlinear solve resets the history counter to zero,
4572: else it continues storing new values for new nonlinear solves after the old ones
4574: Level: intermediate
4576: Notes:
4577: If 'a' and 'its' are `NULL` then space is allocated for the history. If 'na' is `PETSC_DECIDE` (or, deprecated, `PETSC_DEFAULT`) then a
4578: default array of length 1,000 is allocated.
4580: This routine is useful, e.g., when running a code for purposes
4581: of accurate performance monitoring, when no I/O should be done
4582: during the section of code that is being timed.
4584: If the arrays run out of space after a number of iterations then the later values are not saved in the history
4586: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergenceHistory()`
4587: @*/
4588: PetscErrorCode SNESSetConvergenceHistory(SNES snes, PetscReal a[], PetscInt its[], PetscInt na, PetscBool reset)
4589: {
4590: PetscFunctionBegin;
4592: if (a) PetscAssertPointer(a, 2);
4593: if (its) PetscAssertPointer(its, 3);
4594: if (!a) {
4595: if (na == PETSC_DECIDE) na = 1000;
4596: PetscCall(PetscCalloc2(na, &a, na, &its));
4597: snes->conv_hist_alloc = PETSC_TRUE;
4598: }
4599: snes->conv_hist = a;
4600: snes->conv_hist_its = its;
4601: snes->conv_hist_max = (size_t)na;
4602: snes->conv_hist_len = 0;
4603: snes->conv_hist_reset = reset;
4604: PetscFunctionReturn(PETSC_SUCCESS);
4605: }
4607: #if PetscDefined(HAVE_MATLAB)
4608: #include <engine.h> /* MATLAB include file */
4609: #include <mex.h> /* MATLAB include file */
4611: PETSC_EXTERN mxArray *SNESGetConvergenceHistoryMatlab(SNES snes)
4612: {
4613: mxArray *mat;
4614: PetscReal *ar;
4616: mat = mxCreateDoubleMatrix(snes->conv_hist_len, 1, mxREAL);
4617: ar = (PetscReal *)mxGetData(mat);
4618: for (PetscInt i = 0; i < snes->conv_hist_len; i++) ar[i] = snes->conv_hist[i];
4619: return mat;
4620: }
4621: #endif
4623: /*@
4624: SNESGetConvergenceHistory - Gets the arrays used to hold the convergence history.
4626: Not Collective
4628: Input Parameter:
4629: . snes - iterative context obtained from `SNESCreate()`
4631: Output Parameters:
4632: + a - array to hold history, usually was set with `SNESSetConvergenceHistory()`
4633: . its - integer array holds the number of linear iterations (or
4634: negative if not converged) for each solve.
4635: - na - size of `a` and `its`
4637: Level: intermediate
4639: Note:
4640: This routine is useful, e.g., when running a code for purposes
4641: of accurate performance monitoring, when no I/O should be done
4642: during the section of code that is being timed.
4644: Fortran Notes:
4645: Return the arrays with ``SNESRestoreConvergenceHistory()`
4647: Use the arguments
4648: .vb
4649: PetscReal, pointer :: a(:)
4650: PetscInt, pointer :: its(:)
4651: .ve
4653: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESSetConvergenceHistory()`
4654: @*/
4655: PetscErrorCode SNESGetConvergenceHistory(SNES snes, PetscReal *a[], PetscInt *its[], PetscInt *na)
4656: {
4657: PetscFunctionBegin;
4659: if (a) *a = snes->conv_hist;
4660: if (its) *its = snes->conv_hist_its;
4661: if (na) *na = (PetscInt)snes->conv_hist_len;
4662: PetscFunctionReturn(PETSC_SUCCESS);
4663: }
4665: /*@
4666: SNESSetUpdate - Sets the general-purpose update function called
4667: at the beginning of every iteration of the nonlinear solve. Specifically
4668: it is called just before the Jacobian is "evaluated" and after the function
4669: evaluation.
4671: Logically Collective
4673: Input Parameters:
4674: + snes - The nonlinear solver context
4675: - func - The update function; for calling sequence see `SNESUpdateFn`
4677: Level: advanced
4679: Notes:
4680: This is NOT what one uses to update the ghost points before a function evaluation, that should be done at the beginning of your function provided
4681: to `SNESSetFunction()`, or `SNESSetPicard()`
4682: This is not used by most users, and it is intended to provide a general hook that is run
4683: right before the direction step is computed.
4685: Users are free to modify the current residual vector,
4686: the current linearization point, or any other vector associated to the specific solver used.
4687: If such modifications take place, it is the user responsibility to update all the relevant
4688: vectors. For example, if one is adjusting the model parameters at each Newton step their code may look like
4689: .vb
4690: PetscErrorCode update(SNES snes, PetscInt iteration)
4691: {
4692: PetscFunctionBeginUser;
4693: if (iteration > 0) {
4694: // update the model parameters here
4695: Vec x,f;
4696: PetscCall(SNESGetSolution(snes,&x));
4697: PetcCall(SNESGetFunction(snes,&f,NULL,NULL));
4698: PetscCall(SNESComputeFunction(snes,x,f));
4699: }
4700: PetscFunctionReturn(PETSC_SUCCESS);
4701: }
4702: .ve
4704: There are a variety of function hooks one many set that are called at different stages of the nonlinear solution process, see the functions listed below.
4706: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESSetJacobian()`, `SNESLineSearchSetPreCheck()`, `SNESLineSearchSetPostCheck()`, `SNESNewtonTRSetPreCheck()`, `SNESNewtonTRSetPostCheck()`,
4707: `SNESMonitorSet()`
4708: @*/
4709: PetscErrorCode SNESSetUpdate(SNES snes, SNESUpdateFn *func)
4710: {
4711: PetscFunctionBegin;
4713: snes->ops->update = func;
4714: PetscFunctionReturn(PETSC_SUCCESS);
4715: }
4717: /*@
4718: SNESConvergedReasonView - Displays the reason a `SNES` solve converged or diverged to a viewer
4720: Collective
4722: Input Parameters:
4723: + snes - iterative context obtained from `SNESCreate()`
4724: - viewer - the viewer to display the reason
4726: Options Database Keys:
4727: + -snes_converged_reason - print reason for converged or diverged, also prints number of iterations
4728: - -snes_converged_reason ::failed - only print reason and number of iterations when diverged
4730: Level: beginner
4732: Note:
4733: To change the format of the output call `PetscViewerPushFormat`(viewer,format) before this call. Use `PETSC_VIEWER_DEFAULT` for the default,
4734: use `PETSC_VIEWER_FAILED` to only display a reason if it fails.
4736: .seealso: [](ch_snes), `SNESConvergedReason`, `PetscViewer`, `SNES`,
4737: `SNESCreate()`, `SNESSetUp()`, `SNESDestroy()`, `SNESSetTolerances()`, `SNESConvergedDefault()`, `SNESGetConvergedReason()`,
4738: `SNESConvergedReasonViewFromOptions()`,
4739: `PetscViewerPushFormat()`, `PetscViewerPopFormat()`
4740: @*/
4741: PetscErrorCode SNESConvergedReasonView(SNES snes, PetscViewer viewer)
4742: {
4743: PetscViewerFormat format;
4744: PetscBool isAscii;
4746: PetscFunctionBegin;
4747: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)snes));
4748: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isAscii));
4749: if (isAscii) {
4750: PetscCall(PetscViewerGetFormat(viewer, &format));
4751: PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)snes)->tablevel + 1));
4752: if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
4753: DM dm;
4754: Vec u;
4755: PetscDS prob;
4756: PetscInt Nf;
4757: PetscErrorCode (**exactSol)(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar[], void *);
4758: void **exactCtx;
4759: PetscReal error;
4761: PetscCall(SNESGetDM(snes, &dm));
4762: PetscCall(SNESGetSolution(snes, &u));
4763: PetscCall(DMGetDS(dm, &prob));
4764: PetscCall(PetscDSGetNumFields(prob, &Nf));
4765: PetscCall(PetscMalloc2(Nf, &exactSol, Nf, &exactCtx));
4766: for (PetscInt f = 0; f < Nf; ++f) PetscCall(PetscDSGetExactSolution(prob, f, &exactSol[f], &exactCtx[f]));
4767: PetscCall(DMComputeL2Diff(dm, 0.0, exactSol, exactCtx, u, &error));
4768: PetscCall(PetscFree2(exactSol, exactCtx));
4769: if (error < 1.0e-11) PetscCall(PetscViewerASCIIPrintf(viewer, "L_2 Error: < 1.0e-11\n"));
4770: else PetscCall(PetscViewerASCIIPrintf(viewer, "L_2 Error: %g\n", (double)error));
4771: }
4772: if (snes->reason > 0 && format != PETSC_VIEWER_FAILED) {
4773: if (((PetscObject)snes)->prefix) {
4774: PetscCall(PetscViewerASCIIPrintf(viewer, "Nonlinear %s solve converged due to %s iterations %" PetscInt_FMT "\n", ((PetscObject)snes)->prefix, SNESConvergedReasons[snes->reason], snes->iter));
4775: } else {
4776: PetscCall(PetscViewerASCIIPrintf(viewer, "Nonlinear solve converged due to %s iterations %" PetscInt_FMT "\n", SNESConvergedReasons[snes->reason], snes->iter));
4777: }
4778: } else if (snes->reason <= 0) {
4779: if (((PetscObject)snes)->prefix) {
4780: PetscCall(PetscViewerASCIIPrintf(viewer, "Nonlinear %s solve did not converge due to %s iterations %" PetscInt_FMT "\n", ((PetscObject)snes)->prefix, SNESConvergedReasons[snes->reason], snes->iter));
4781: } else {
4782: PetscCall(PetscViewerASCIIPrintf(viewer, "Nonlinear solve did not converge due to %s iterations %" PetscInt_FMT "\n", SNESConvergedReasons[snes->reason], snes->iter));
4783: }
4784: }
4785: PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)snes)->tablevel + 1));
4786: }
4787: PetscFunctionReturn(PETSC_SUCCESS);
4788: }
4790: /*@
4791: SNESConvergedReasonViewSet - Sets an ADDITIONAL function that is to be used at the
4792: end of the nonlinear solver to display the convergence reason of the nonlinear solver.
4794: Logically Collective
4796: Input Parameters:
4797: + snes - the `SNES` context
4798: . f - the `SNESConvergedReason` view function
4799: . vctx - [optional] user-defined context for private data for the `SNESConvergedReason` view function (use `NULL` if no context is desired)
4800: - reasonviewdestroy - [optional] routine that frees the context (may be `NULL`), see `PetscCtxDestroyFn` for the calling sequence
4802: Calling sequence of `f`:
4803: + snes - the `SNES` context
4804: - vctx - [optional] context for private data for the function
4806: Options Database Keys:
4807: + -snes_converged_reason - sets a default `SNESConvergedReasonView()`
4808: - -snes_converged_reason_view_cancel - cancels all converged reason viewers that have been hardwired into a code by
4809: calls to `SNESConvergedReasonViewSet()`, but does not cancel those set via the options database.
4811: Level: intermediate
4813: Note:
4814: Several different converged reason view routines may be set by calling
4815: `SNESConvergedReasonViewSet()` multiple times; all will be called in the
4816: order in which they were set.
4818: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESConvergedReason`, `SNESGetConvergedReason()`, `SNESConvergedReasonView()`, `SNESConvergedReasonViewCancel()`,
4819: `PetscCtxDestroyFn`
4820: @*/
4821: PetscErrorCode SNESConvergedReasonViewSet(SNES snes, PetscErrorCode (*f)(SNES snes, PetscCtx vctx), PetscCtx vctx, PetscCtxDestroyFn *reasonviewdestroy)
4822: {
4823: PetscFunctionBegin;
4825: for (PetscInt i = 0; i < snes->numberreasonviews; i++) {
4826: PetscBool identical;
4828: PetscCall(PetscMonitorCompare((PetscErrorCode (*)(void))(PetscVoidFn *)f, vctx, reasonviewdestroy, (PetscErrorCode (*)(void))(PetscVoidFn *)snes->reasonview[i], snes->reasonviewcontext[i], snes->reasonviewdestroy[i], &identical));
4829: if (identical) PetscFunctionReturn(PETSC_SUCCESS);
4830: }
4831: PetscCheck(snes->numberreasonviews < MAXSNESREASONVIEWS, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Too many SNES reasonview set");
4832: snes->reasonview[snes->numberreasonviews] = f;
4833: snes->reasonviewdestroy[snes->numberreasonviews] = reasonviewdestroy;
4834: snes->reasonviewcontext[snes->numberreasonviews++] = vctx;
4835: PetscFunctionReturn(PETSC_SUCCESS);
4836: }
4838: /*@
4839: SNESConvergedReasonViewFromOptions - Processes command line options to determine if/how a `SNESConvergedReason` is to be viewed at the end of `SNESSolve()`
4840: All the user-provided viewer routines set with `SNESConvergedReasonViewSet()` will be called, if they exist.
4842: Collective
4844: Input Parameter:
4845: . snes - the `SNES` object
4847: Level: advanced
4849: Note:
4850: This function has a different API and behavior than `PetscObjectViewFromOptions()`
4852: .seealso: [](ch_snes), `SNES`, `SNESConvergedReason`, `SNESConvergedReasonViewSet()`, `SNESCreate()`, `SNESSetUp()`, `SNESDestroy()`,
4853: `SNESSetTolerances()`, `SNESConvergedDefault()`, `SNESGetConvergedReason()`, `SNESConvergedReasonView()`
4854: @*/
4855: PetscErrorCode SNESConvergedReasonViewFromOptions(SNES snes)
4856: {
4857: static PetscBool incall = PETSC_FALSE;
4859: PetscFunctionBegin;
4860: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
4861: incall = PETSC_TRUE;
4863: /* All user-provided viewers are called first, if they exist. */
4864: for (PetscInt i = 0; i < snes->numberreasonviews; i++) PetscCall((*snes->reasonview[i])(snes, snes->reasonviewcontext[i]));
4866: /* Call PETSc default routine if users ask for it */
4867: if (snes->convergedreasonviewer) {
4868: PetscCall(PetscViewerPushFormat(snes->convergedreasonviewer, snes->convergedreasonformat));
4869: PetscCall(SNESConvergedReasonView(snes, snes->convergedreasonviewer));
4870: PetscCall(PetscViewerPopFormat(snes->convergedreasonviewer));
4871: }
4872: incall = PETSC_FALSE;
4873: PetscFunctionReturn(PETSC_SUCCESS);
4874: }
4876: /*@
4877: SNESSolve - Solves a nonlinear system $F(x) = b$ associated with a `SNES` object
4879: Collective
4881: Input Parameters:
4882: + snes - the `SNES` context
4883: . b - the constant part of the equation $F(x) = b$, or `NULL` to use zero.
4884: - x - the solution vector.
4886: Level: beginner
4888: Note:
4889: The user should initialize the vector, `x`, with the initial guess
4890: for the nonlinear solve prior to calling `SNESSolve()` .
4892: .seealso: [](ch_snes), `SNES`, `SNESCreate()`, `SNESDestroy()`, `SNESSetFunction()`, `SNESSetJacobian()`, `SNESSetGridSequence()`, `SNESGetSolution()`,
4893: `SNESNewtonTRSetPreCheck()`, `SNESNewtonTRGetPreCheck()`, `SNESNewtonTRSetPostCheck()`, `SNESNewtonTRGetPostCheck()`,
4894: `SNESLineSearchSetPostCheck()`, `SNESLineSearchGetPostCheck()`, `SNESLineSearchSetPreCheck()`, `SNESLineSearchGetPreCheck()`
4895: @*/
4896: PetscErrorCode SNESSolve(SNES snes, Vec b, Vec x)
4897: {
4898: PetscBool flg;
4899: Vec xcreated = NULL;
4900: DM dm;
4902: PetscFunctionBegin;
4905: if (x) PetscCheckSameComm(snes, 1, x, 3);
4907: if (b) PetscCheckSameComm(snes, 1, b, 2);
4909: /* High level operations using the nonlinear solver */
4910: {
4911: PetscViewer viewer;
4912: PetscViewerFormat format;
4913: PetscInt num;
4914: PetscBool flg;
4915: static PetscBool incall = PETSC_FALSE;
4917: if (!incall) {
4918: /* Estimate the convergence rate of the discretization */
4919: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_convergence_estimate", &viewer, &format, &flg));
4920: if (flg) {
4921: PetscConvEst conv;
4922: DM dm;
4923: PetscReal *alpha; /* Convergence rate of the solution error for each field in the L_2 norm */
4924: PetscInt Nf;
4926: incall = PETSC_TRUE;
4927: PetscCall(SNESGetDM(snes, &dm));
4928: PetscCall(DMGetNumFields(dm, &Nf));
4929: PetscCall(PetscCalloc1(Nf, &alpha));
4930: PetscCall(PetscConvEstCreate(PetscObjectComm((PetscObject)snes), &conv));
4931: PetscCall(PetscConvEstSetSolver(conv, (PetscObject)snes));
4932: PetscCall(PetscConvEstSetFromOptions(conv));
4933: PetscCall(PetscConvEstSetUp(conv));
4934: PetscCall(PetscConvEstGetConvRate(conv, alpha));
4935: PetscCall(PetscViewerPushFormat(viewer, format));
4936: PetscCall(PetscConvEstRateView(conv, alpha, viewer));
4937: PetscCall(PetscViewerPopFormat(viewer));
4938: PetscCall(PetscViewerDestroy(&viewer));
4939: PetscCall(PetscConvEstDestroy(&conv));
4940: PetscCall(PetscFree(alpha));
4941: incall = PETSC_FALSE;
4942: }
4943: /* Adaptively refine the initial grid */
4944: num = 1;
4945: PetscCall(PetscOptionsGetInt(NULL, ((PetscObject)snes)->prefix, "-snes_adapt_initial", &num, &flg));
4946: if (flg) {
4947: DMAdaptor adaptor;
4949: incall = PETSC_TRUE;
4950: PetscCall(DMAdaptorCreate(PetscObjectComm((PetscObject)snes), &adaptor));
4951: PetscCall(DMAdaptorSetSolver(adaptor, snes));
4952: PetscCall(DMAdaptorSetSequenceLength(adaptor, num));
4953: PetscCall(DMAdaptorSetFromOptions(adaptor));
4954: PetscCall(DMAdaptorSetUp(adaptor));
4955: PetscCall(DMAdaptorAdapt(adaptor, x, DM_ADAPTATION_INITIAL, &dm, &x));
4956: PetscCall(DMAdaptorDestroy(&adaptor));
4957: incall = PETSC_FALSE;
4958: }
4959: /* Use grid sequencing to adapt */
4960: num = 0;
4961: PetscCall(PetscOptionsGetInt(NULL, ((PetscObject)snes)->prefix, "-snes_adapt_sequence", &num, NULL));
4962: if (num) {
4963: DMAdaptor adaptor;
4964: const char *prefix;
4966: incall = PETSC_TRUE;
4967: PetscCall(DMAdaptorCreate(PetscObjectComm((PetscObject)snes), &adaptor));
4968: PetscCall(SNESGetOptionsPrefix(snes, &prefix));
4969: PetscCall(DMAdaptorSetOptionsPrefix(adaptor, prefix));
4970: PetscCall(DMAdaptorSetSolver(adaptor, snes));
4971: PetscCall(DMAdaptorSetSequenceLength(adaptor, num));
4972: PetscCall(DMAdaptorSetFromOptions(adaptor));
4973: PetscCall(DMAdaptorSetUp(adaptor));
4974: PetscCall(PetscObjectViewFromOptions((PetscObject)adaptor, NULL, "-snes_adapt_view"));
4975: PetscCall(DMAdaptorAdapt(adaptor, x, DM_ADAPTATION_SEQUENTIAL, &dm, &x));
4976: PetscCall(DMAdaptorDestroy(&adaptor));
4977: incall = PETSC_FALSE;
4978: }
4979: }
4980: }
4981: if (!x) x = snes->vec_sol;
4982: if (!x) {
4983: PetscCall(SNESGetDM(snes, &dm));
4984: PetscCall(DMCreateGlobalVector(dm, &xcreated));
4985: x = xcreated;
4986: }
4987: PetscCall(SNESViewFromOptions(snes, NULL, "-snes_view_pre"));
4989: for (PetscInt grid = 0; grid < snes->gridsequence; grid++) PetscCall(PetscViewerASCIIPushTab(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)snes))));
4990: for (PetscInt grid = 0; grid < snes->gridsequence + 1; grid++) {
4991: /* set solution vector */
4992: if (!grid) PetscCall(PetscObjectReference((PetscObject)x));
4993: PetscCall(VecDestroy(&snes->vec_sol));
4994: snes->vec_sol = x;
4995: PetscCall(SNESGetDM(snes, &dm));
4997: /* set affine vector if provided */
4998: PetscCall(PetscObjectReference((PetscObject)b));
4999: PetscCall(VecDestroy(&snes->vec_rhs));
5000: snes->vec_rhs = b;
5002: if (snes->vec_rhs) PetscCheck(snes->vec_func != snes->vec_rhs, PETSC_COMM_SELF, PETSC_ERR_ARG_IDN, "Right hand side vector cannot be function vector");
5003: PetscCheck(snes->vec_func != snes->vec_sol, PETSC_COMM_SELF, PETSC_ERR_ARG_IDN, "Solution vector cannot be function vector");
5004: PetscCheck(snes->vec_rhs != snes->vec_sol, PETSC_COMM_SELF, PETSC_ERR_ARG_IDN, "Solution vector cannot be right-hand side vector");
5005: if (!snes->vec_sol_update /* && snes->vec_sol */) PetscCall(VecDuplicate(snes->vec_sol, &snes->vec_sol_update));
5006: PetscCall(DMShellSetGlobalVector(dm, snes->vec_sol));
5007: PetscCall(SNESSetUp(snes));
5009: if (!grid) {
5010: if (snes->ops->computeinitialguess) PetscCallBack("SNES callback compute initial guess", (*snes->ops->computeinitialguess)(snes, snes->vec_sol, snes->initialguessP));
5011: }
5013: if (snes->conv_hist_reset) snes->conv_hist_len = 0;
5014: PetscCall(SNESResetCounters(snes));
5015: snes->reason = SNES_CONVERGED_ITERATING;
5016: PetscCall(PetscLogEventBegin(SNES_Solve, snes, 0, 0, 0));
5017: PetscUseTypeMethod(snes, solve);
5018: PetscCall(PetscLogEventEnd(SNES_Solve, snes, 0, 0, 0));
5019: PetscCheck(snes->reason, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Internal error, solver %s returned without setting converged reason", ((PetscObject)snes)->type_name);
5020: snes->functiondomainerror = PETSC_FALSE; /* clear the flag if it has been set */
5021: snes->objectivedomainerror = PETSC_FALSE; /* clear the flag if it has been set */
5022: snes->jacobiandomainerror = PETSC_FALSE; /* clear the flag if it has been set */
5024: if (snes->lagjac_persist) snes->jac_iter += snes->iter;
5025: if (snes->lagpre_persist) snes->pre_iter += snes->iter;
5027: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_test_local_min", NULL, NULL, &flg));
5028: if (flg && !PetscPreLoadingOn) PetscCall(SNESTestLocalMin(snes));
5029: /* Call converged reason views. This may involve user-provided viewers as well */
5030: PetscCall(SNESConvergedReasonViewFromOptions(snes));
5032: if (snes->errorifnotconverged) {
5033: if (snes->reason < 0) PetscCall(SNESMonitorCancel(snes));
5034: PetscCheck(snes->reason >= 0, PetscObjectComm((PetscObject)snes), PETSC_ERR_NOT_CONVERGED, "SNESSolve has not converged");
5035: }
5036: if (snes->reason < 0) break;
5037: if (grid < snes->gridsequence) {
5038: DM fine;
5039: Vec xnew;
5040: Mat interp;
5042: PetscCall(DMRefine(snes->dm, PetscObjectComm((PetscObject)snes), &fine));
5043: PetscCheck(fine, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_INCOMP, "DMRefine() did not perform any refinement, cannot continue grid sequencing");
5044: PetscCall(DMGetCoordinatesLocalSetUp(fine));
5045: PetscCall(DMCreateInterpolation(snes->dm, fine, &interp, NULL));
5046: PetscCall(DMCreateGlobalVector(fine, &xnew));
5047: PetscCall(MatInterpolate(interp, x, xnew));
5048: PetscCall(DMInterpolate(snes->dm, interp, fine));
5049: PetscCall(MatDestroy(&interp));
5050: x = xnew;
5052: PetscCall(SNESReset(snes));
5053: PetscCall(SNESSetDM(snes, fine));
5054: PetscCall(SNESResetFromOptions(snes));
5055: PetscCall(DMDestroy(&fine));
5056: PetscCall(PetscViewerASCIIPopTab(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)snes))));
5057: }
5058: }
5059: PetscCall(SNESViewFromOptions(snes, NULL, "-snes_view"));
5060: PetscCall(VecViewFromOptions(snes->vec_sol, (PetscObject)snes, "-snes_view_solution"));
5061: PetscCall(DMMonitor(snes->dm));
5062: PetscCall(SNESMonitorPauseFinal_Internal(snes));
5064: PetscCall(VecDestroy(&xcreated));
5065: PetscCall(PetscObjectSAWsBlock((PetscObject)snes));
5066: PetscFunctionReturn(PETSC_SUCCESS);
5067: }
5069: /* --------- Internal routines for SNES Package --------- */
5071: /*@
5072: SNESSetType - Sets the algorithm/method to be used to solve the nonlinear system with the given `SNES`
5074: Collective
5076: Input Parameters:
5077: + snes - the `SNES` context
5078: - type - a known method
5080: Options Database Key:
5081: . -snes_type type - Sets the method; see `SNESType`
5083: Level: intermediate
5085: Notes:
5086: See `SNESType` for available methods (for instance)
5087: + `SNESNEWTONLS` - Newton's method with line search
5088: (systems of nonlinear equations)
5089: - `SNESNEWTONTR` - Newton's method with trust region
5090: (systems of nonlinear equations)
5092: Normally, it is best to use the `SNESSetFromOptions()` command and then
5093: set the `SNES` solver type from the options database rather than by using
5094: this routine. Using the options database provides the user with
5095: maximum flexibility in evaluating the many nonlinear solvers.
5096: The `SNESSetType()` routine is provided for those situations where it
5097: is necessary to set the nonlinear solver independently of the command
5098: line or options database. This might be the case, for example, when
5099: the choice of solver changes during the execution of the program,
5100: and the user's application is taking responsibility for choosing the
5101: appropriate method.
5103: Developer Note:
5104: `SNESRegister()` adds a constructor for a new `SNESType` to `SNESList`, `SNESSetType()` locates
5105: the constructor in that list and calls it to create the specific object.
5107: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESType`, `SNESCreate()`, `SNESDestroy()`, `SNESGetType()`, `SNESSetFromOptions()`
5108: @*/
5109: PetscErrorCode SNESSetType(SNES snes, SNESType type)
5110: {
5111: PetscBool match;
5112: PetscErrorCode (*r)(SNES);
5114: PetscFunctionBegin;
5116: PetscAssertPointer(type, 2);
5118: PetscCall(PetscObjectTypeCompare((PetscObject)snes, type, &match));
5119: if (match) PetscFunctionReturn(PETSC_SUCCESS);
5121: PetscCall(PetscFunctionListFind(SNESList, type, &r));
5122: PetscCheck(r, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unable to find requested SNES type %s", type);
5123: /* Destroy the previous private SNES context */
5124: PetscTryTypeMethod(snes, destroy);
5125: /* Reinitialize type-specific function pointers in SNESOps structure */
5126: snes->ops->reset = NULL;
5127: snes->ops->setup = NULL;
5128: snes->ops->solve = NULL;
5129: snes->ops->view = NULL;
5130: snes->ops->setfromoptions = NULL;
5131: snes->ops->destroy = NULL;
5133: /* It may happen the user has customized the line search before calling SNESSetType */
5134: if (((PetscObject)snes)->type_name) PetscCall(SNESLineSearchDestroy(&snes->linesearch));
5136: /* Reinitialize default parameters */
5137: PetscCall(SNESParametersInitialize(snes));
5139: /* Call the SNESCreate_XXX routine for this particular Nonlinear solver */
5140: snes->setupcalled = PETSC_FALSE;
5141: PetscCall(PetscObjectChangeTypeName((PetscObject)snes, type));
5142: PetscCall((*r)(snes));
5143: PetscFunctionReturn(PETSC_SUCCESS);
5144: }
5146: /*@
5147: SNESGetType - Gets the `SNES` method type and name (as a string).
5149: Not Collective
5151: Input Parameter:
5152: . snes - nonlinear solver context
5154: Output Parameter:
5155: . type - `SNES` method (a character string)
5157: Level: intermediate
5159: Note:
5160: `type` should not be retained for later use as it will be an invalid pointer if the `SNESType` of `snes` is changed.
5162: .seealso: [](ch_snes), `SNESSetType()`, `SNESType`, `SNESSetFromOptions()`, `SNES`, `PetscObjectTypeCompare()`, `PetscObjectTypeCompareAny()`
5163: @*/
5164: PetscErrorCode SNESGetType(SNES snes, SNESType *type)
5165: {
5166: PetscFunctionBegin;
5168: PetscAssertPointer(type, 2);
5169: *type = ((PetscObject)snes)->type_name;
5170: PetscFunctionReturn(PETSC_SUCCESS);
5171: }
5173: /*@
5174: SNESSetSolution - Sets the solution vector for use by the `SNES` routines.
5176: Logically Collective
5178: Input Parameters:
5179: + snes - the `SNES` context obtained from `SNESCreate()`
5180: - u - the solution vector
5182: Level: beginner
5184: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetSolution()`, `Vec`
5185: @*/
5186: PetscErrorCode SNESSetSolution(SNES snes, Vec u)
5187: {
5188: DM dm;
5190: PetscFunctionBegin;
5193: PetscCall(PetscObjectReference((PetscObject)u));
5194: PetscCall(VecDestroy(&snes->vec_sol));
5196: snes->vec_sol = u;
5198: PetscCall(SNESGetDM(snes, &dm));
5199: PetscCall(DMShellSetGlobalVector(dm, u));
5200: PetscFunctionReturn(PETSC_SUCCESS);
5201: }
5203: /*@
5204: SNESGetSolution - Returns the vector where the approximate solution is
5205: stored. This is the fine grid solution when using `SNESSetGridSequence()`.
5207: Not Collective, but `x` is parallel if `snes` is parallel
5209: Input Parameter:
5210: . snes - the `SNES` context
5212: Output Parameter:
5213: . x - the solution
5215: Level: intermediate
5217: .seealso: [](ch_snes), `SNESSetSolution()`, `SNESSolve()`, `SNES`, `SNESGetSolutionUpdate()`, `SNESGetFunction()`
5218: @*/
5219: PetscErrorCode SNESGetSolution(SNES snes, Vec *x)
5220: {
5221: PetscFunctionBegin;
5223: PetscAssertPointer(x, 2);
5224: *x = snes->vec_sol;
5225: PetscFunctionReturn(PETSC_SUCCESS);
5226: }
5228: /*@
5229: SNESGetSolutionUpdate - Returns the vector where the solution update is
5230: stored.
5232: Not Collective, but `x` is parallel if `snes` is parallel
5234: Input Parameter:
5235: . snes - the `SNES` context
5237: Output Parameter:
5238: . x - the solution update
5240: Level: advanced
5242: .seealso: [](ch_snes), `SNES`, `SNESGetSolution()`, `SNESGetFunction()`
5243: @*/
5244: PetscErrorCode SNESGetSolutionUpdate(SNES snes, Vec *x)
5245: {
5246: PetscFunctionBegin;
5248: PetscAssertPointer(x, 2);
5249: *x = snes->vec_sol_update;
5250: PetscFunctionReturn(PETSC_SUCCESS);
5251: }
5253: /*@
5254: SNESGetFunction - Returns the function that defines the nonlinear system set with `SNESSetFunction()`
5256: Not Collective, but `r` is parallel if `snes` is parallel. Collective if `r` is requested, but has not been created yet.
5258: Input Parameter:
5259: . snes - the `SNES` context
5261: Output Parameters:
5262: + r - the vector that is used to store residuals (or `NULL` if you don't want it)
5263: . f - the function (or `NULL` if you don't want it); for calling sequence see `SNESFunctionFn`
5264: - ctx - the function context (or `NULL` if you don't want it)
5266: Level: advanced
5268: Note:
5269: The vector `r` DOES NOT, in general, contain the current value of the `SNES` nonlinear function
5271: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESSetFunction()`, `SNESGetSolution()`, `SNESFunctionFn`
5272: @*/
5273: PetscErrorCode SNESGetFunction(SNES snes, Vec *r, SNESFunctionFn **f, PetscCtxRt ctx)
5274: {
5275: DM dm;
5277: PetscFunctionBegin;
5279: if (r) {
5280: if (!snes->vec_func) {
5281: if (snes->vec_rhs) {
5282: PetscCall(VecDuplicate(snes->vec_rhs, &snes->vec_func));
5283: } else if (snes->vec_sol) {
5284: PetscCall(VecDuplicate(snes->vec_sol, &snes->vec_func));
5285: } else if (snes->dm) {
5286: PetscCall(DMCreateGlobalVector(snes->dm, &snes->vec_func));
5287: }
5288: }
5289: *r = snes->vec_func;
5290: }
5291: PetscCall(SNESGetDM(snes, &dm));
5292: PetscCall(DMSNESGetFunction(dm, f, ctx));
5293: PetscFunctionReturn(PETSC_SUCCESS);
5294: }
5296: /*@
5297: SNESGetNGS - Returns the function and context set with `SNESSetNGS()`
5299: Input Parameter:
5300: . snes - the `SNES` context
5302: Output Parameters:
5303: + f - the function (or `NULL`) see `SNESNGSFn` for calling sequence
5304: - ctx - the function context (or `NULL`)
5306: Level: advanced
5308: .seealso: [](ch_snes), `SNESSetNGS()`, `SNESGetFunction()`, `SNESNGSFn`
5309: @*/
5310: PetscErrorCode SNESGetNGS(SNES snes, SNESNGSFn **f, PetscCtxRt ctx)
5311: {
5312: DM dm;
5314: PetscFunctionBegin;
5316: PetscCall(SNESGetDM(snes, &dm));
5317: PetscCall(DMSNESGetNGS(dm, f, ctx));
5318: PetscFunctionReturn(PETSC_SUCCESS);
5319: }
5321: /*@
5322: SNESSetOptionsPrefix - Sets the prefix used for searching for all
5323: `SNES` options in the database.
5325: Logically Collective
5327: Input Parameters:
5328: + snes - the `SNES` context
5329: - prefix - the prefix to prepend to all option names
5331: Level: advanced
5333: Note:
5334: A hyphen (-) must NOT be given at the beginning of the prefix name.
5335: The first character of all runtime options is AUTOMATICALLY the hyphen.
5337: .seealso: [](ch_snes), `SNES`, `SNESSetFromOptions()`, `SNESAppendOptionsPrefix()`
5338: @*/
5339: PetscErrorCode SNESSetOptionsPrefix(SNES snes, const char prefix[])
5340: {
5341: PetscFunctionBegin;
5343: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)snes, prefix));
5344: if (!snes->ksp) PetscCall(SNESGetKSP(snes, &snes->ksp));
5345: if (snes->linesearch) {
5346: PetscCall(SNESGetLineSearch(snes, &snes->linesearch));
5347: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)snes->linesearch, prefix));
5348: }
5349: PetscCall(KSPSetOptionsPrefix(snes->ksp, prefix));
5350: PetscFunctionReturn(PETSC_SUCCESS);
5351: }
5353: /*@
5354: SNESAppendOptionsPrefix - Appends to the prefix used for searching for all
5355: `SNES` options in the database.
5357: Logically Collective
5359: Input Parameters:
5360: + snes - the `SNES` context
5361: - prefix - the prefix to prepend to all option names
5363: Level: advanced
5365: Note:
5366: A hyphen (-) must NOT be given at the beginning of the prefix name.
5367: The first character of all runtime options is AUTOMATICALLY the hyphen.
5369: .seealso: [](ch_snes), `SNESGetOptionsPrefix()`, `SNESSetOptionsPrefix()`
5370: @*/
5371: PetscErrorCode SNESAppendOptionsPrefix(SNES snes, const char prefix[])
5372: {
5373: PetscFunctionBegin;
5375: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)snes, prefix));
5376: if (!snes->ksp) PetscCall(SNESGetKSP(snes, &snes->ksp));
5377: if (snes->linesearch) {
5378: PetscCall(SNESGetLineSearch(snes, &snes->linesearch));
5379: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)snes->linesearch, prefix));
5380: }
5381: PetscCall(KSPAppendOptionsPrefix(snes->ksp, prefix));
5382: PetscFunctionReturn(PETSC_SUCCESS);
5383: }
5385: /*@
5386: SNESGetOptionsPrefix - Gets the prefix used for searching for all
5387: `SNES` options in the database.
5389: Not Collective
5391: Input Parameter:
5392: . snes - the `SNES` context
5394: Output Parameter:
5395: . prefix - pointer to the prefix string used
5397: Level: advanced
5399: .seealso: [](ch_snes), `SNES`, `SNESSetOptionsPrefix()`, `SNESAppendOptionsPrefix()`
5400: @*/
5401: PetscErrorCode SNESGetOptionsPrefix(SNES snes, const char *prefix[])
5402: {
5403: PetscFunctionBegin;
5405: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)snes, prefix));
5406: PetscFunctionReturn(PETSC_SUCCESS);
5407: }
5409: /*@
5410: SNESRegister - Adds a method to the nonlinear solver package.
5412: Not Collective
5414: Input Parameters:
5415: + sname - name of a new user-defined solver
5416: - function - routine to create method context
5418: Level: advanced
5420: Note:
5421: `SNESRegister()` may be called multiple times to add several user-defined solvers.
5423: Example Usage:
5424: .vb
5425: SNESRegister("my_solver", MySolverCreate);
5426: .ve
5428: Then, your solver can be chosen with the procedural interface via
5429: .vb
5430: SNESSetType(snes, "my_solver")
5431: .ve
5432: or at runtime via the option
5433: .vb
5434: -snes_type my_solver
5435: .ve
5437: .seealso: [](ch_snes), `SNESRegisterAll()`, `SNESRegisterDestroy()`
5438: @*/
5439: PetscErrorCode SNESRegister(const char sname[], PetscErrorCode (*function)(SNES))
5440: {
5441: PetscFunctionBegin;
5442: PetscCall(SNESInitializePackage());
5443: PetscCall(PetscFunctionListAdd(&SNESList, sname, function));
5444: PetscFunctionReturn(PETSC_SUCCESS);
5445: }
5447: /*@
5448: SNESTestLocalMin - Diagnostic that probes each entry of the current `SNES` solution to check whether the residual norm has a local minimum along the coordinate directions
5450: Collective
5452: Input Parameter:
5453: . snes - the `SNES` context
5455: Level: developer
5457: Note:
5458: Currently intended for serial runs. For each degree of freedom it perturbs the solution by increasing amounts and prints the resulting `SNESComputeFunction()` residual norms so the user can inspect local behavior.
5460: .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESComputeFunction()`
5461: @*/
5462: PetscErrorCode SNESTestLocalMin(SNES snes)
5463: {
5464: PetscInt N, i, j;
5465: Vec u, uh, fh;
5466: PetscScalar value;
5467: PetscReal norm;
5469: PetscFunctionBegin;
5470: PetscCall(SNESGetSolution(snes, &u));
5471: PetscCall(VecDuplicate(u, &uh));
5472: PetscCall(VecDuplicate(u, &fh));
5474: /* currently only works for sequential */
5475: PetscCall(PetscPrintf(PetscObjectComm((PetscObject)snes), "Testing FormFunction() for local min\n"));
5476: PetscCall(VecGetSize(u, &N));
5477: for (i = 0; i < N; i++) {
5478: PetscCall(VecCopy(u, uh));
5479: PetscCall(PetscPrintf(PetscObjectComm((PetscObject)snes), "i = %" PetscInt_FMT "\n", i));
5480: for (j = -10; j < 11; j++) {
5481: value = PetscSign(j) * PetscExpReal(PetscAbs(j) - 10.0);
5482: PetscCall(VecSetValue(uh, i, value, ADD_VALUES));
5483: PetscCall(SNESComputeFunction(snes, uh, fh));
5484: PetscCall(VecNorm(fh, NORM_2, &norm)); /* does not handle use of SNESSetFunctionDomainError() correctly */
5485: PetscCall(PetscPrintf(PetscObjectComm((PetscObject)snes), " j norm %" PetscInt_FMT " %18.16e\n", j, (double)norm));
5486: value = -value;
5487: PetscCall(VecSetValue(uh, i, value, ADD_VALUES));
5488: }
5489: }
5490: PetscCall(VecDestroy(&uh));
5491: PetscCall(VecDestroy(&fh));
5492: PetscFunctionReturn(PETSC_SUCCESS);
5493: }
5495: /*@
5496: SNESGetLineSearch - Returns the line search associated with the `SNES`.
5498: Not Collective
5500: Input Parameter:
5501: . snes - iterative context obtained from `SNESCreate()`
5503: Output Parameter:
5504: . linesearch - linesearch context
5506: Level: beginner
5508: Notes:
5509: It creates a default line search instance which can be configured as needed in case it has not been already set with `SNESSetLineSearch()`.
5511: You can also use the options database keys `-snes_linesearch_*` to configure the line search. See `SNESLineSearchSetFromOptions()` for the possible options.
5513: .seealso: [](ch_snes), `SNESLineSearch`, `SNESSetLineSearch()`, `SNESLineSearchCreate()`, `SNESLineSearchSetFromOptions()`
5514: @*/
5515: PetscErrorCode SNESGetLineSearch(SNES snes, SNESLineSearch *linesearch)
5516: {
5517: const char *optionsprefix;
5519: PetscFunctionBegin;
5521: PetscAssertPointer(linesearch, 2);
5522: if (!snes->linesearch) {
5523: PetscCall(SNESGetOptionsPrefix(snes, &optionsprefix));
5524: PetscCall(SNESLineSearchCreate(PetscObjectComm((PetscObject)snes), &snes->linesearch));
5525: PetscCall(SNESLineSearchSetSNES(snes->linesearch, snes));
5526: PetscCall(SNESLineSearchAppendOptionsPrefix(snes->linesearch, optionsprefix));
5527: PetscCall(PetscObjectIncrementTabLevel((PetscObject)snes->linesearch, (PetscObject)snes, 1));
5528: }
5529: *linesearch = snes->linesearch;
5530: PetscFunctionReturn(PETSC_SUCCESS);
5531: }
5533: /*@
5534: SNESKSPSetUseEW - Sets `SNES` to the use Eisenstat-Walker method for
5535: computing relative tolerance for linear solvers within an inexact
5536: Newton method.
5538: Logically Collective
5540: Input Parameters:
5541: + snes - `SNES` context
5542: - flag - `PETSC_TRUE` or `PETSC_FALSE`
5544: Options Database Keys:
5545: + -snes_ksp_ew - use Eisenstat-Walker method for determining linear system convergence
5546: . -snes_ksp_ew_version ver - version of Eisenstat-Walker method
5547: . -snes_ksp_ew_rtol0 rtol0 - Sets rtol0
5548: . -snes_ksp_ew_rtolmax rtolmax - Sets rtolmax
5549: . -snes_ksp_ew_gamma gamma - Sets gamma
5550: . -snes_ksp_ew_alpha alpha - Sets alpha
5551: . -snes_ksp_ew_alpha2 alpha2 - Sets alpha2
5552: - -snes_ksp_ew_threshold threshold - Sets threshold
5554: Level: advanced
5556: Note:
5557: The default is to use a constant relative tolerance for
5558: the inner linear solvers. Alternatively, one can use the
5559: Eisenstat-Walker method {cite}`ew96`, where the relative convergence tolerance
5560: is reset at each Newton iteration according progress of the nonlinear
5561: solver.
5563: .seealso: [](ch_snes), `KSP`, `SNES`, `SNESKSPGetUseEW()`, `SNESKSPGetParametersEW()`, `SNESKSPSetParametersEW()`
5564: @*/
5565: PetscErrorCode SNESKSPSetUseEW(SNES snes, PetscBool flag)
5566: {
5567: PetscFunctionBegin;
5570: snes->ksp_ewconv = flag;
5571: PetscFunctionReturn(PETSC_SUCCESS);
5572: }
5574: /*@
5575: SNESKSPGetUseEW - Gets if `SNES` is using Eisenstat-Walker method
5576: for computing relative tolerance for linear solvers within an
5577: inexact Newton method.
5579: Not Collective
5581: Input Parameter:
5582: . snes - `SNES` context
5584: Output Parameter:
5585: . flag - `PETSC_TRUE` or `PETSC_FALSE`
5587: Level: advanced
5589: .seealso: [](ch_snes), `SNESKSPSetUseEW()`, `SNESKSPGetParametersEW()`, `SNESKSPSetParametersEW()`
5590: @*/
5591: PetscErrorCode SNESKSPGetUseEW(SNES snes, PetscBool *flag)
5592: {
5593: PetscFunctionBegin;
5595: PetscAssertPointer(flag, 2);
5596: *flag = snes->ksp_ewconv;
5597: PetscFunctionReturn(PETSC_SUCCESS);
5598: }
5600: /*@
5601: SNESKSPSetParametersEW - Sets parameters for Eisenstat-Walker
5602: convergence criteria for the linear solvers within an inexact
5603: Newton method.
5605: Logically Collective
5607: Input Parameters:
5608: + snes - `SNES` context
5609: . version - version 1, 2 (default is 2), 3 or 4
5610: . rtol_0 - initial relative tolerance (0 <= rtol_0 < 1)
5611: . rtol_max - maximum relative tolerance (0 <= rtol_max < 1)
5612: . gamma - multiplicative factor for version 2 rtol computation
5613: (0 <= gamma2 <= 1)
5614: . alpha - power for version 2 rtol computation (1 < alpha <= 2)
5615: . alpha2 - power for safeguard
5616: - threshold - threshold for imposing safeguard (0 < threshold < 1)
5618: Level: advanced
5620: Notes:
5621: Version 3 was contributed by Luis Chacon, June 2006.
5623: Use `PETSC_CURRENT` to retain the default for any of the parameters.
5625: .seealso: [](ch_snes), `SNES`, `SNESKSPSetUseEW()`, `SNESKSPGetUseEW()`, `SNESKSPGetParametersEW()`
5626: @*/
5627: PetscErrorCode SNESKSPSetParametersEW(SNES snes, PetscInt version, PetscReal rtol_0, PetscReal rtol_max, PetscReal gamma, PetscReal alpha, PetscReal alpha2, PetscReal threshold)
5628: {
5629: SNESKSPEW *kctx;
5631: PetscFunctionBegin;
5633: kctx = (SNESKSPEW *)snes->kspconvctx;
5634: PetscCheck(kctx, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "No Eisenstat-Walker context existing");
5643: if (version != PETSC_CURRENT) kctx->version = version;
5644: if (rtol_0 != (PetscReal)PETSC_CURRENT) kctx->rtol_0 = rtol_0;
5645: if (rtol_max != (PetscReal)PETSC_CURRENT) kctx->rtol_max = rtol_max;
5646: if (gamma != (PetscReal)PETSC_CURRENT) kctx->gamma = gamma;
5647: if (alpha != (PetscReal)PETSC_CURRENT) kctx->alpha = alpha;
5648: if (alpha2 != (PetscReal)PETSC_CURRENT) kctx->alpha2 = alpha2;
5649: if (threshold != (PetscReal)PETSC_CURRENT) kctx->threshold = threshold;
5651: PetscCheck(kctx->version >= 1 && kctx->version <= 4, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Only versions 1 to 4 are supported: %" PetscInt_FMT, kctx->version);
5652: PetscCheck(kctx->rtol_0 >= 0.0 && kctx->rtol_0 < 1.0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "0.0 <= rtol_0 < 1.0: %g", (double)kctx->rtol_0);
5653: PetscCheck(kctx->rtol_max >= 0.0 && kctx->rtol_max < 1.0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "0.0 <= rtol_max (%g) < 1.0", (double)kctx->rtol_max);
5654: PetscCheck(kctx->gamma >= 0.0 && kctx->gamma <= 1.0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "0.0 <= gamma (%g) <= 1.0", (double)kctx->gamma);
5655: PetscCheck(kctx->alpha > 1.0 && kctx->alpha <= 2.0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "1.0 < alpha (%g) <= 2.0", (double)kctx->alpha);
5656: PetscCheck(kctx->threshold > 0.0 && kctx->threshold < 1.0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "0.0 < threshold (%g) < 1.0", (double)kctx->threshold);
5657: PetscFunctionReturn(PETSC_SUCCESS);
5658: }
5660: /*@
5661: SNESKSPGetParametersEW - Gets parameters for Eisenstat-Walker
5662: convergence criteria for the linear solvers within an inexact
5663: Newton method.
5665: Not Collective
5667: Input Parameter:
5668: . snes - `SNES` context
5670: Output Parameters:
5671: + version - version 1, 2 (default is 2), 3 or 4
5672: . rtol_0 - initial relative tolerance (0 <= rtol_0 < 1)
5673: . rtol_max - maximum relative tolerance (0 <= rtol_max < 1)
5674: . gamma - multiplicative factor for version 2 rtol computation (0 <= gamma2 <= 1)
5675: . alpha - power for version 2 rtol computation (1 < alpha <= 2)
5676: . alpha2 - power for safeguard
5677: - threshold - threshold for imposing safeguard (0 < threshold < 1)
5679: Level: advanced
5681: .seealso: [](ch_snes), `SNES`, `SNESKSPSetUseEW()`, `SNESKSPGetUseEW()`, `SNESKSPSetParametersEW()`
5682: @*/
5683: PetscErrorCode SNESKSPGetParametersEW(SNES snes, PetscInt *version, PetscReal *rtol_0, PetscReal *rtol_max, PetscReal *gamma, PetscReal *alpha, PetscReal *alpha2, PetscReal *threshold)
5684: {
5685: SNESKSPEW *kctx;
5687: PetscFunctionBegin;
5689: kctx = (SNESKSPEW *)snes->kspconvctx;
5690: PetscCheck(kctx, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "No Eisenstat-Walker context existing");
5691: if (version) *version = kctx->version;
5692: if (rtol_0) *rtol_0 = kctx->rtol_0;
5693: if (rtol_max) *rtol_max = kctx->rtol_max;
5694: if (gamma) *gamma = kctx->gamma;
5695: if (alpha) *alpha = kctx->alpha;
5696: if (alpha2) *alpha2 = kctx->alpha2;
5697: if (threshold) *threshold = kctx->threshold;
5698: PetscFunctionReturn(PETSC_SUCCESS);
5699: }
5701: PetscErrorCode KSPPreSolve_SNESEW(KSP ksp, Vec b, Vec x, PetscCtx ctx)
5702: {
5703: SNES snes = (SNES)ctx;
5704: SNESKSPEW *kctx = (SNESKSPEW *)snes->kspconvctx;
5705: PetscReal rtol = PETSC_CURRENT, stol;
5707: PetscFunctionBegin;
5708: if (!snes->ksp_ewconv) PetscFunctionReturn(PETSC_SUCCESS);
5709: if (!snes->iter) {
5710: rtol = kctx->rtol_0; /* first time in, so use the original user rtol */
5711: PetscCall(VecNorm(snes->vec_func, NORM_2, &kctx->norm_first));
5712: } else {
5713: PetscCheck(kctx->version >= 1 && kctx->version <= 4, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Only versions 1-4 are supported: %" PetscInt_FMT, kctx->version);
5714: if (kctx->version == 1) {
5715: rtol = PetscAbsReal(snes->norm - kctx->lresid_last) / kctx->norm_last;
5716: stol = PetscPowReal(kctx->rtol_last, kctx->alpha2);
5717: if (stol > kctx->threshold) rtol = PetscMax(rtol, stol);
5718: } else if (kctx->version == 2) {
5719: rtol = kctx->gamma * PetscPowReal(snes->norm / kctx->norm_last, kctx->alpha);
5720: stol = kctx->gamma * PetscPowReal(kctx->rtol_last, kctx->alpha);
5721: if (stol > kctx->threshold) rtol = PetscMax(rtol, stol);
5722: } else if (kctx->version == 3) { /* contributed by Luis Chacon, June 2006. */
5723: rtol = kctx->gamma * PetscPowReal(snes->norm / kctx->norm_last, kctx->alpha);
5724: /* safeguard: avoid sharp decrease of rtol */
5725: stol = kctx->gamma * PetscPowReal(kctx->rtol_last, kctx->alpha);
5726: stol = PetscMax(rtol, stol);
5727: rtol = PetscMin(kctx->rtol_0, stol);
5728: /* safeguard: avoid oversolving */
5729: stol = kctx->gamma * (kctx->norm_first * snes->rtol) / snes->norm;
5730: stol = PetscMax(rtol, stol);
5731: rtol = PetscMin(kctx->rtol_0, stol);
5732: } else /* if (kctx->version == 4) */ {
5733: /* H.-B. An et al. Journal of Computational and Applied Mathematics 200 (2007) 47-60 */
5734: PetscReal ared = PetscAbsReal(kctx->norm_last - snes->norm);
5735: PetscReal pred = PetscAbsReal(kctx->norm_last - kctx->lresid_last);
5736: PetscReal rk = ared / pred;
5737: if (rk < kctx->v4_p1) rtol = 1. - 2. * kctx->v4_p1;
5738: else if (rk < kctx->v4_p2) rtol = kctx->rtol_last;
5739: else if (rk < kctx->v4_p3) rtol = kctx->v4_m1 * kctx->rtol_last;
5740: else rtol = kctx->v4_m2 * kctx->rtol_last;
5742: if (kctx->rtol_last_2 > kctx->v4_m3 && kctx->rtol_last > kctx->v4_m3 && kctx->rk_last_2 < kctx->v4_p1 && kctx->rk_last < kctx->v4_p1) rtol = kctx->v4_m4 * kctx->rtol_last;
5743: kctx->rtol_last_2 = kctx->rtol_last;
5744: kctx->rk_last_2 = kctx->rk_last;
5745: kctx->rk_last = rk;
5746: }
5747: }
5748: /* safeguard: avoid rtol greater than rtol_max */
5749: rtol = PetscMin(rtol, kctx->rtol_max);
5750: PetscCall(KSPSetTolerances(ksp, rtol, PETSC_CURRENT, PETSC_CURRENT, PETSC_CURRENT));
5751: PetscCall(PetscInfo(snes, "iter %" PetscInt_FMT ", Eisenstat-Walker (version %" PetscInt_FMT ") KSP rtol=%g\n", snes->iter, kctx->version, (double)rtol));
5752: PetscFunctionReturn(PETSC_SUCCESS);
5753: }
5755: PetscErrorCode KSPPostSolve_SNESEW(KSP ksp, Vec b, Vec x, PetscCtx ctx)
5756: {
5757: SNES snes = (SNES)ctx;
5758: SNESKSPEW *kctx = (SNESKSPEW *)snes->kspconvctx;
5759: PCSide pcside;
5760: Vec lres;
5762: PetscFunctionBegin;
5763: if (!snes->ksp_ewconv) PetscFunctionReturn(PETSC_SUCCESS);
5764: PetscCall(KSPGetTolerances(ksp, &kctx->rtol_last, NULL, NULL, NULL));
5765: kctx->norm_last = snes->norm;
5766: if (kctx->version == 1 || kctx->version == 4) {
5767: PC pc;
5768: PetscBool getRes;
5770: PetscCall(KSPGetPC(ksp, &pc));
5771: PetscCall(PetscObjectTypeCompare((PetscObject)pc, PCNONE, &getRes));
5772: if (!getRes) {
5773: KSPNormType normtype;
5775: PetscCall(KSPGetNormType(ksp, &normtype));
5776: getRes = (PetscBool)(normtype == KSP_NORM_UNPRECONDITIONED);
5777: }
5778: PetscCall(KSPGetPCSide(ksp, &pcside));
5779: if (pcside == PC_RIGHT || getRes) { /* KSP residual is true linear residual */
5780: PetscCall(KSPGetResidualNorm(ksp, &kctx->lresid_last));
5781: } else {
5782: /* KSP residual is preconditioned residual */
5783: /* compute true linear residual norm */
5784: Mat J;
5785: PetscCall(KSPGetOperators(ksp, &J, NULL));
5786: PetscCall(VecDuplicate(b, &lres));
5787: PetscCall(MatMult(J, x, lres));
5788: PetscCall(VecAYPX(lres, -1.0, b));
5789: PetscCall(VecNorm(lres, NORM_2, &kctx->lresid_last));
5790: PetscCall(VecDestroy(&lres));
5791: }
5792: }
5793: PetscFunctionReturn(PETSC_SUCCESS);
5794: }
5796: #include <petsc/private/kspimpl.h>
5797: /*@
5798: SNESGetKSP - Returns the `KSP` context for a `SNES` solver.
5800: Not Collective, but if `snes` is parallel, then `ksp` is parallel
5802: Input Parameter:
5803: . snes - the `SNES` context
5805: Output Parameter:
5806: . ksp - the `KSP` context
5808: Level: beginner
5810: Notes:
5811: The user can then directly manipulate the `KSP` context to set various
5812: options, etc. Likewise, the user can then extract and manipulate the
5813: `PC` contexts as well.
5815: Some `SNESType`s do not use a `KSP` but a `KSP` is still returned by this function, changes to that `KSP` will have no effect.
5817: .seealso: [](ch_snes), `SNES`, `KSP`, `PC`, `KSPGetPC()`, `SNESCreate()`, `KSPCreate()`, `SNESSetKSP()`
5818: @*/
5819: PetscErrorCode SNESGetKSP(SNES snes, KSP *ksp)
5820: {
5821: PetscFunctionBegin;
5823: PetscAssertPointer(ksp, 2);
5825: if (!snes->ksp) {
5826: PetscCall(KSPCreate(PetscObjectComm((PetscObject)snes), &snes->ksp));
5827: PetscCall(PetscObjectIncrementTabLevel((PetscObject)snes->ksp, (PetscObject)snes, 1));
5829: snes->ksp->presolve_ew = KSPPreSolve_SNESEW;
5830: snes->ksp->prectx_ew = snes;
5831: snes->ksp->postsolve_ew = KSPPostSolve_SNESEW;
5832: snes->ksp->postctx_ew = snes;
5834: PetscCall(KSPMonitorSetFromOptions(snes->ksp, "-snes_monitor_ksp", "snes_preconditioned_residual", snes));
5835: PetscCall(PetscObjectSetOptions((PetscObject)snes->ksp, ((PetscObject)snes)->options));
5836: }
5837: *ksp = snes->ksp;
5838: PetscFunctionReturn(PETSC_SUCCESS);
5839: }
5841: #include <petsc/private/dmimpl.h>
5842: /*@
5843: SNESSetDM - Sets the `DM` that may be used by some `SNES` nonlinear solvers or their underlying preconditioners
5845: Logically Collective
5847: Input Parameters:
5848: + snes - the nonlinear solver context
5849: - dm - the `DM`, cannot be `NULL`
5851: Level: intermediate
5853: Note:
5854: A `DM` can only be used for solving one problem at a time because information about the problem is stored on the `DM`,
5855: even when not using interfaces like `DMSNESSetFunction()`. Use `DMClone()` to get a distinct `DM` when solving different
5856: problems using the same function space.
5858: .seealso: [](ch_snes), `DM`, `SNES`, `SNESGetDM()`, `KSPSetDM()`, `KSPGetDM()`
5859: @*/
5860: PetscErrorCode SNESSetDM(SNES snes, DM dm)
5861: {
5862: KSP ksp;
5863: DMSNES sdm;
5864: DM odm;
5866: PetscFunctionBegin;
5869: PetscCall(PetscObjectReference((PetscObject)dm));
5870: odm = snes->dm;
5871: if (snes->dm) { /* Move the DMSNES context over to the new DM unless the new DM already has one */
5872: if (snes->dm->dmsnes && !dm->dmsnes) {
5873: PetscCall(DMCopyDMSNES(snes->dm, dm));
5874: PetscCall(DMGetDMSNES(snes->dm, &sdm));
5875: if (sdm->originaldm == snes->dm) sdm->originaldm = dm; /* Grant write privileges to the replacement DM */
5876: }
5877: PetscCall(DMCoarsenHookRemove(snes->dm, DMCoarsenHook_SNESVecSol, DMRestrictHook_SNESVecSol, snes));
5878: PetscCall(DMDestroy(&snes->dm));
5879: }
5880: snes->dm = dm;
5881: snes->dmAuto = PETSC_FALSE;
5883: PetscCall(SNESGetKSP(snes, &ksp));
5884: PetscCall(KSPSetDM(ksp, dm));
5885: PetscCall(KSPSetDMActive(ksp, KSP_DMACTIVE_ALL, PETSC_FALSE));
5886: /* Propagate DM to NPC if npc does not have one yet or
5887: if it has the same DM SNES had before (like for gridsequencing) */
5888: if (snes->npc && (!snes->npc->dm || snes->npc->dm == odm)) PetscCall(SNESSetDM(snes->npc, snes->dm));
5889: PetscFunctionReturn(PETSC_SUCCESS);
5890: }
5892: /*@
5893: SNESGetDM - Gets the `DM` that may be used by some `SNES` nonlinear solvers/preconditioners
5895: Not Collective but `dm` obtained is parallel on `snes`
5897: Input Parameter:
5898: . snes - the `SNES` context
5900: Output Parameter:
5901: . dm - the `DM`
5903: Level: intermediate
5905: .seealso: [](ch_snes), `DM`, `SNES`, `SNESSetDM()`, `KSPSetDM()`, `KSPGetDM()`
5906: @*/
5907: PetscErrorCode SNESGetDM(SNES snes, DM *dm)
5908: {
5909: PetscFunctionBegin;
5911: if (!snes->dm) {
5912: PetscCall(DMShellCreate(PetscObjectComm((PetscObject)snes), &snes->dm));
5913: snes->dmAuto = PETSC_TRUE;
5914: }
5915: *dm = snes->dm;
5916: PetscFunctionReturn(PETSC_SUCCESS);
5917: }
5919: /*@
5920: SNESSetNPC - Sets the nonlinear preconditioner to be used.
5922: Collective
5924: Input Parameters:
5925: + snes - iterative context obtained from `SNESCreate()`
5926: - npc - the `SNES` nonlinear preconditioner object
5928: Level: developer
5930: Notes:
5931: This is rarely used, rather use `SNESGetNPC()` to retrieve the preconditioner and configure it using the API.
5933: Only some `SNESType` can use a nonlinear preconditioner
5935: .seealso: [](ch_snes), `SNES`, `SNESNGS`, `SNESFAS`, `SNESGetNPC()`, `SNESHasNPC()`
5936: @*/
5937: PetscErrorCode SNESSetNPC(SNES snes, SNES npc)
5938: {
5939: PetscFunctionBegin;
5942: PetscCheckSameComm(snes, 1, npc, 2);
5943: PetscCall(PetscObjectReference((PetscObject)npc));
5944: PetscCall(SNESDestroy(&snes->npc));
5945: snes->npc = npc;
5946: PetscFunctionReturn(PETSC_SUCCESS);
5947: }
5949: /*@
5950: SNESGetNPC - Gets a nonlinear preconditioning solver SNES` to be used to precondition the original nonlinear solver.
5952: Collective the first time it is called if the `SNES` has no NPC set.
5954: Input Parameter:
5955: . snes - iterative context obtained from `SNESCreate()`
5957: Output Parameter:
5958: . npc - the `SNES` preconditioner context
5960: Options Database Key:
5961: . -npc_snes_type type - set the type of the `SNES` to use as the nonlinear preconditioner
5963: Level: advanced
5965: Notes:
5966: If a `SNES` was previously set with `SNESSetNPC()` then that object is returned, otherwise a new `SNES` object is created that will
5967: be used as the nonlinear preconditioner for the current `SNES` if no nonlinear preconditioner is present.
5969: The (preconditioner) `SNES` returned automatically inherits the same nonlinear function and Jacobian supplied to the original
5970: `SNES`. These may be overwritten if needed by calling `SNESSetDM()` on the nonlinear preconditioner followed by `SNESSetFunction()`
5971: and `SNESSetJacobian()`.
5973: The default preconditioner uses the options database prefixes `-npc_snes`, `-npc_ksp`, etc., to control the configuration.
5975: .seealso: [](ch_snes), `SNESSetNPC()`, `SNESHasNPC()`, `SNES`, `SNESCreate()`
5976: @*/
5977: PetscErrorCode SNESGetNPC(SNES snes, SNES *npc)
5978: {
5979: const char *optionsprefix;
5981: PetscFunctionBegin;
5983: PetscAssertPointer(npc, 2);
5984: if (!snes->npc) {
5985: PetscCall(SNESCreate(PetscObjectComm((PetscObject)snes), &snes->npc));
5986: PetscCall(PetscObjectIncrementTabLevel((PetscObject)snes->npc, (PetscObject)snes, 1));
5987: PetscCall(SNESGetOptionsPrefix(snes, &optionsprefix));
5988: PetscCall(SNESSetOptionsPrefix(snes->npc, optionsprefix));
5989: PetscCall(SNESAppendOptionsPrefix(snes->npc, "npc_"));
5990: PetscCall(SNESSetCountersReset(snes->npc, PETSC_FALSE));
5991: PetscCall(SNESSetNormSchedule(snes->npc, SNES_NORM_DEFAULT));
5993: /* default to 1 iteration */
5994: PetscCall(SNESSetTolerances(snes->npc, 0.0, 0.0, 0.0, 1, snes->npc->max_funcs));
5995: }
5996: *npc = snes->npc;
5997: PetscFunctionReturn(PETSC_SUCCESS);
5998: }
6000: /*@
6001: SNESHasNPC - Returns whether a nonlinear preconditioner is associated with the given `SNES`
6003: Not Collective
6005: Input Parameter:
6006: . snes - iterative context obtained from `SNESCreate()`
6008: Output Parameter:
6009: . has_npc - whether the `SNES` has a nonlinear preconditioner or not
6011: Level: developer
6013: .seealso: [](ch_snes), `SNESSetNPC()`, `SNESGetNPC()`
6014: @*/
6015: PetscErrorCode SNESHasNPC(SNES snes, PetscBool *has_npc)
6016: {
6017: PetscFunctionBegin;
6019: PetscAssertPointer(has_npc, 2);
6020: *has_npc = snes->npc ? PETSC_TRUE : PETSC_FALSE;
6021: PetscFunctionReturn(PETSC_SUCCESS);
6022: }
6024: /*@
6025: SNESSetNPCSide - Sets the nonlinear preconditioning side used by the nonlinear preconditioner inside `SNES`.
6027: Logically Collective
6029: Input Parameter:
6030: . snes - iterative context obtained from `SNESCreate()`
6032: Output Parameter:
6033: . side - the preconditioning side, where side is one of
6034: .vb
6035: PC_LEFT - left preconditioning
6036: PC_RIGHT - right preconditioning (default for most nonlinear solvers)
6037: .ve
6039: Options Database Key:
6040: . -snes_npc_side (right|left) - nonlinear preconditioner side
6042: Level: intermediate
6044: Note:
6045: `SNESNRICHARDSON` and `SNESNCG` only support left preconditioning.
6047: .seealso: [](ch_snes), `SNES`, `SNESGetNPC()`, `SNESNRICHARDSON`, `SNESNCG`, `SNESType`, `SNESGetNPCSide()`, `KSPSetPCSide()`, `PC_LEFT`, `PC_RIGHT`, `PCSide`
6048: @*/
6049: PetscErrorCode SNESSetNPCSide(SNES snes, PCSide side)
6050: {
6051: PetscFunctionBegin;
6054: if (side == PC_SIDE_DEFAULT) side = PC_RIGHT;
6055: PetscCheck((side == PC_LEFT) || (side == PC_RIGHT), PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONG, "Only PC_LEFT and PC_RIGHT are supported");
6056: snes->npcside = side;
6057: PetscFunctionReturn(PETSC_SUCCESS);
6058: }
6060: /*@
6061: SNESGetNPCSide - Gets the preconditioning side used by the nonlinear preconditioner inside `SNES`.
6063: Not Collective
6065: Input Parameter:
6066: . snes - iterative context obtained from `SNESCreate()`
6068: Output Parameter:
6069: . side - the preconditioning side, where side is one of
6070: .vb
6071: `PC_LEFT` - left preconditioning
6072: `PC_RIGHT` - right preconditioning (default for most nonlinear solvers)
6073: .ve
6075: Level: intermediate
6077: .seealso: [](ch_snes), `SNES`, `SNESGetNPC()`, `SNESSetNPCSide()`, `KSPGetPCSide()`, `PC_LEFT`, `PC_RIGHT`, `PCSide`
6078: @*/
6079: PetscErrorCode SNESGetNPCSide(SNES snes, PCSide *side)
6080: {
6081: PetscFunctionBegin;
6083: PetscAssertPointer(side, 2);
6084: *side = snes->npcside;
6085: PetscFunctionReturn(PETSC_SUCCESS);
6086: }
6088: /*@
6089: SNESSetLineSearch - Sets the `SNESLineSearch` to be used for a given `SNES`
6091: Collective
6093: Input Parameters:
6094: + snes - iterative context obtained from `SNESCreate()`
6095: - linesearch - the linesearch object
6097: Level: developer
6099: Note:
6100: This is almost never used, rather one uses `SNESGetLineSearch()` to retrieve the line search and set options on it
6101: to configure it using the API).
6103: .seealso: [](ch_snes), `SNES`, `SNESLineSearch`, `SNESGetLineSearch()`
6104: @*/
6105: PetscErrorCode SNESSetLineSearch(SNES snes, SNESLineSearch linesearch)
6106: {
6107: PetscFunctionBegin;
6110: PetscCheckSameComm(snes, 1, linesearch, 2);
6111: PetscCall(PetscObjectReference((PetscObject)linesearch));
6112: PetscCall(SNESLineSearchDestroy(&snes->linesearch));
6114: snes->linesearch = linesearch;
6115: PetscFunctionReturn(PETSC_SUCCESS);
6116: }