1 /*
2 User interface for the nonlinear solvers package.
3 */
4 #pragma once
5
6 #include <petscksp.h>
7 #include <petscdmtypes.h>
8 #include <petscfvtypes.h>
9 #include <petscdmdatypes.h>
10 #include <petscsnestypes.h>
11
12 /* SUBMANSEC = SNES */
13
14 /*J
15 SNESType - String with the name of a PETSc `SNES` method. These are all the nonlinear solvers that PETSc provides.
16
17 Level: beginner
18
19 Note:
20 Use `SNESSetType()` or the options database key `-snes_type` to set the specific nonlinear solver algorithm to use with a given `SNES` object
21
22 .seealso: [](doc_nonlinsolve), [](ch_snes), `SNESSetType()`, `SNES`, `SNESCreate()`, `SNESDestroy()`, `SNESSetFromOptions()`
23 J*/
24 typedef const char *SNESType;
25 #define SNESNEWTONLS "newtonls"
26 #define SNESNEWTONTR "newtontr"
27 #define SNESNEWTONTRDC "newtontrdc"
28 #define SNESPYTHON "python"
29 #define SNESNRICHARDSON "nrichardson"
30 #define SNESKSPONLY "ksponly"
31 #define SNESKSPTRANSPOSEONLY "ksptransposeonly"
32 #define SNESVINEWTONRSLS "vinewtonrsls"
33 #define SNESVINEWTONSSLS "vinewtonssls"
34 #define SNESNGMRES "ngmres"
35 #define SNESQN "qn"
36 #define SNESSHELL "shell"
37 #define SNESNGS "ngs"
38 #define SNESNCG "ncg"
39 #define SNESFAS "fas"
40 #define SNESMS "ms"
41 #define SNESNASM "nasm"
42 #define SNESANDERSON "anderson"
43 #define SNESASPIN "aspin"
44 #define SNESCOMPOSITE "composite"
45 #define SNESPATCH "patch"
46 #define SNESNEWTONAL "newtonal"
47
48 /* Logging support */
49 PETSC_EXTERN PetscClassId SNES_CLASSID;
50 PETSC_EXTERN PetscClassId DMSNES_CLASSID;
51
52 PETSC_EXTERN PetscErrorCode SNESInitializePackage(void);
53 PETSC_EXTERN PetscErrorCode SNESFinalizePackage(void);
54
55 PETSC_EXTERN PetscErrorCode SNESCreate(MPI_Comm, SNES *);
56 PETSC_EXTERN PetscErrorCode SNESParametersInitialize(SNES);
57 PETSC_EXTERN PetscErrorCode SNESReset(SNES);
58 PETSC_EXTERN PetscErrorCode SNESDestroy(SNES *);
59 PETSC_EXTERN PetscErrorCode SNESSetType(SNES, SNESType);
60 PETSC_EXTERN PetscErrorCode SNESMonitor(SNES, PetscInt, PetscReal);
61 PETSC_EXTERN PetscErrorCode SNESMonitorSet(SNES, PetscErrorCode (*)(SNES, PetscInt, PetscReal, void *), void *, PetscCtxDestroyFn *);
62 PETSC_EXTERN PetscErrorCode SNESMonitorSetFromOptions(SNES, const char[], const char[], const char[], PetscErrorCode (*)(SNES, PetscInt, PetscReal, PetscViewerAndFormat *), PetscErrorCode (*)(SNES, PetscViewerAndFormat *));
63 PETSC_EXTERN PetscErrorCode SNESMonitorCancel(SNES);
64 PETSC_EXTERN PetscErrorCode SNESMonitorSAWs(SNES, PetscInt, PetscReal, PetscCtx);
65 PETSC_EXTERN PetscErrorCode SNESMonitorSAWsCreate(SNES, void **);
66 PETSC_EXTERN PetscErrorCode SNESMonitorSAWsDestroy(PetscCtxRt);
67 PETSC_EXTERN PetscErrorCode SNESSetConvergenceHistory(SNES, PetscReal[], PetscInt[], PetscInt, PetscBool);
68 PETSC_EXTERN PetscErrorCode SNESGetConvergenceHistory(SNES, PetscReal *[], PetscInt *[], PetscInt *);
69 PETSC_EXTERN PetscErrorCode SNESSetUp(SNES);
70 PETSC_EXTERN PetscErrorCode SNESSolve(SNES, Vec, Vec);
71 PETSC_EXTERN PetscErrorCode SNESSetErrorIfNotConverged(SNES, PetscBool);
72 PETSC_EXTERN PetscErrorCode SNESGetErrorIfNotConverged(SNES, PetscBool *);
73 PETSC_EXTERN PetscErrorCode SNESConverged(SNES, PetscInt, PetscReal, PetscReal, PetscReal);
74
75 PETSC_EXTERN PetscErrorCode SNESSetWorkVecs(SNES, PetscInt);
76
77 PETSC_EXTERN PetscErrorCode SNESAddOptionsChecker(PetscErrorCode (*)(SNES));
78
79 PETSC_EXTERN PetscErrorCode SNESRegister(const char[], PetscErrorCode (*)(SNES));
80
81 PETSC_EXTERN PetscErrorCode SNESGetKSP(SNES, KSP *);
82 PETSC_EXTERN PetscErrorCode SNESSetKSP(SNES, KSP);
83 PETSC_EXTERN PetscErrorCode SNESSetSolution(SNES, Vec);
84 PETSC_EXTERN PetscErrorCode SNESGetSolution(SNES, Vec *);
85 PETSC_EXTERN PetscErrorCode SNESGetSolutionUpdate(SNES, Vec *);
86 PETSC_EXTERN PetscErrorCode SNESGetRhs(SNES, Vec *);
87 PETSC_EXTERN PetscErrorCode SNESView(SNES, PetscViewer);
88 PETSC_EXTERN PetscErrorCode SNESLoad(SNES, PetscViewer);
89 PETSC_EXTERN PetscErrorCode SNESConvergedReasonViewSet(SNES, PetscErrorCode (*)(SNES, void *), void *, PetscCtxDestroyFn *);
90 PETSC_EXTERN PetscErrorCode SNESViewFromOptions(SNES, PetscObject, const char[]);
91 PETSC_EXTERN PetscErrorCode SNESConvergedReasonView(SNES, PetscViewer);
92 PETSC_EXTERN PetscErrorCode SNESConvergedReasonViewFromOptions(SNES);
93 PETSC_EXTERN PetscErrorCode SNESConvergedReasonViewCancel(SNES);
94
SNESReasonView(SNES snes,PetscViewer v)95 PETSC_DEPRECATED_FUNCTION(3, 14, 0, "SNESConvergedReasonView()", ) static inline PetscErrorCode SNESReasonView(SNES snes, PetscViewer v)
96 {
97 return SNESConvergedReasonView(snes, v);
98 }
SNESReasonViewFromOptions(SNES snes)99 PETSC_DEPRECATED_FUNCTION(3, 14, 0, "SNESConvergedReasonViewFromOptions()", ) static inline PetscErrorCode SNESReasonViewFromOptions(SNES snes)
100 {
101 return SNESConvergedReasonViewFromOptions(snes);
102 }
103
104 #define SNES_FILE_CLASSID 1211224
105
106 PETSC_EXTERN PetscErrorCode SNESSetOptionsPrefix(SNES, const char[]);
107 PETSC_EXTERN PetscErrorCode SNESAppendOptionsPrefix(SNES, const char[]);
108 PETSC_EXTERN PetscErrorCode SNESGetOptionsPrefix(SNES, const char *[]);
109 PETSC_EXTERN PetscErrorCode SNESSetFromOptions(SNES);
110 PETSC_EXTERN PetscErrorCode SNESResetFromOptions(SNES);
111
112 PETSC_EXTERN PetscErrorCode SNESSetUseMatrixFree(SNES, PetscBool, PetscBool);
113 PETSC_EXTERN PetscErrorCode SNESGetUseMatrixFree(SNES, PetscBool *, PetscBool *);
114 PETSC_EXTERN PetscErrorCode MatCreateSNESMF(SNES, Mat *);
115 PETSC_EXTERN PetscErrorCode MatSNESMFGetSNES(Mat, SNES *);
116 PETSC_EXTERN PetscErrorCode MatSNESMFSetReuseBase(Mat, PetscBool);
117 PETSC_EXTERN PetscErrorCode MatSNESMFGetReuseBase(Mat, PetscBool *);
118 PETSC_EXTERN PetscErrorCode MatMFFDComputeJacobian(SNES, Vec, Mat, Mat, void *);
119 PETSC_EXTERN PetscErrorCode MatCreateSNESMFMore(SNES, Vec, Mat *);
120 PETSC_EXTERN PetscErrorCode MatSNESMFMoreSetParameters(Mat, PetscReal, PetscReal, PetscReal);
121
122 PETSC_EXTERN PetscErrorCode SNESGetType(SNES, SNESType *);
123 PETSC_EXTERN PetscErrorCode SNESMonitorDefaultSetUp(SNES, PetscViewerAndFormat *);
124 PETSC_EXTERN PetscErrorCode SNESMonitorDefault(SNES, PetscInt, PetscReal, PetscViewerAndFormat *);
125 PETSC_EXTERN PetscErrorCode SNESMonitorScaling(SNES, PetscInt, PetscReal, PetscViewerAndFormat *);
126 PETSC_EXTERN PetscErrorCode SNESMonitorRange(SNES, PetscInt, PetscReal, PetscViewerAndFormat *);
127 PETSC_EXTERN PetscErrorCode SNESMonitorRatio(SNES, PetscInt, PetscReal, PetscViewerAndFormat *);
128 PETSC_EXTERN PetscErrorCode SNESMonitorRatioSetUp(SNES, PetscViewerAndFormat *);
129 PETSC_EXTERN PetscErrorCode SNESMonitorSolution(SNES, PetscInt, PetscReal, PetscViewerAndFormat *);
130 PETSC_EXTERN PetscErrorCode SNESMonitorResidual(SNES, PetscInt, PetscReal, PetscViewerAndFormat *);
131 PETSC_EXTERN PetscErrorCode SNESMonitorSolutionUpdate(SNES, PetscInt, PetscReal, PetscViewerAndFormat *);
132 PETSC_EXTERN PetscErrorCode SNESMonitorDefaultShort(SNES, PetscInt, PetscReal, PetscViewerAndFormat *);
133 PETSC_EXTERN PetscErrorCode SNESMonitorDefaultField(SNES, PetscInt, PetscReal, PetscViewerAndFormat *);
134 PETSC_EXTERN PetscErrorCode SNESMonitorJacUpdateSpectrum(SNES, PetscInt, PetscReal, PetscViewerAndFormat *);
135 PETSC_EXTERN PetscErrorCode SNESMonitorFields(SNES, PetscInt, PetscReal, PetscViewerAndFormat *);
136 PETSC_EXTERN PetscErrorCode KSPMonitorSNESResidual(KSP, PetscInt, PetscReal, PetscViewerAndFormat *);
137 PETSC_EXTERN PetscErrorCode KSPMonitorSNESResidualDrawLG(KSP, PetscInt, PetscReal, PetscViewerAndFormat *);
138 PETSC_EXTERN PetscErrorCode KSPMonitorSNESResidualDrawLGCreate(PetscViewer, PetscViewerFormat, void *, PetscViewerAndFormat **);
139
140 PETSC_EXTERN PetscErrorCode SNESSetTolerances(SNES, PetscReal, PetscReal, PetscReal, PetscInt, PetscInt);
141 PETSC_EXTERN PetscErrorCode SNESSetDivergenceTolerance(SNES, PetscReal);
142 PETSC_EXTERN PetscErrorCode SNESGetTolerances(SNES, PetscReal *, PetscReal *, PetscReal *, PetscInt *, PetscInt *);
143 PETSC_EXTERN PetscErrorCode SNESGetDivergenceTolerance(SNES, PetscReal *);
144 PETSC_EXTERN PetscErrorCode SNESGetForceIteration(SNES, PetscBool *);
145 PETSC_EXTERN PetscErrorCode SNESSetForceIteration(SNES, PetscBool);
146 PETSC_EXTERN PetscErrorCode SNESGetIterationNumber(SNES, PetscInt *);
147 PETSC_EXTERN PetscErrorCode SNESSetIterationNumber(SNES, PetscInt);
148
149 /*E
150 SNESNewtonTRFallbackType - type of fallback in case the solution of the trust-region subproblem is outside of the radius
151
152 Values:
153 + `SNES_TR_FALLBACK_NEWTON` - use scaled Newton step
154 . `SNES_TR_FALLBACK_CAUCHY` - use Cauchy direction
155 - `SNES_TR_FALLBACK_DOGLEG` - use dogleg method
156
157 Level: intermediate
158
159 .seealso: [](ch_snes), `SNES`, `SNESNEWTONTR`, `SNESNEWTONTRDC`
160 E*/
161 typedef enum {
162 SNES_TR_FALLBACK_NEWTON,
163 SNES_TR_FALLBACK_CAUCHY,
164 SNES_TR_FALLBACK_DOGLEG,
165 } SNESNewtonTRFallbackType;
166
167 PETSC_EXTERN const char *const SNESNewtonTRFallbackTypes[];
168
169 PETSC_EXTERN PetscErrorCode SNESNewtonTRSetPreCheck(SNES, PetscErrorCode (*)(SNES, Vec, Vec, PetscBool *, void *), PetscCtx);
170 PETSC_EXTERN PetscErrorCode SNESNewtonTRGetPreCheck(SNES, PetscErrorCode (**)(SNES, Vec, Vec, PetscBool *, void *), PetscCtxRt);
171 PETSC_EXTERN PetscErrorCode SNESNewtonTRSetPostCheck(SNES, PetscErrorCode (*)(SNES, Vec, Vec, Vec, PetscBool *, PetscBool *, void *), PetscCtx);
172 PETSC_EXTERN PetscErrorCode SNESNewtonTRGetPostCheck(SNES, PetscErrorCode (**)(SNES, Vec, Vec, Vec, PetscBool *, PetscBool *, void *), PetscCtxRt);
173 PETSC_EXTERN PetscErrorCode SNESNewtonTRSetFallbackType(SNES, SNESNewtonTRFallbackType);
174 PETSC_EXTERN PetscErrorCode SNESNewtonTRPreCheck(SNES, Vec, Vec, PetscBool *);
175 PETSC_EXTERN PetscErrorCode SNESNewtonTRPostCheck(SNES, Vec, Vec, Vec, PetscBool *, PetscBool *);
176 PETSC_EXTERN PetscErrorCode SNESNewtonTRSetNormType(SNES, NormType);
177
178 /*E
179 SNESNewtonTRQNType - type of quasi-Newton model to use
180
181 Values:
182 + `SNES_TR_QN_NONE` - do not use a quasi-Newton model
183 . `SNES_TR_QN_SAME` - use the same quasi-Newton model for matrix and preconditioner
184 - `SNES_TR_QN_DIFFERENT` - use different quasi-Newton models for matrix and preconditioner
185
186 Level: intermediate
187
188 .seealso: [](ch_snes), `SNES`, `SNESNEWTONTR`
189 E*/
190 typedef enum {
191 SNES_TR_QN_NONE,
192 SNES_TR_QN_SAME,
193 SNES_TR_QN_DIFFERENT,
194 } SNESNewtonTRQNType;
195
196 PETSC_EXTERN const char *const SNESNewtonTRQNTypes[];
197
198 PETSC_EXTERN PetscErrorCode SNESNewtonTRSetQNType(SNES, SNESNewtonTRQNType);
199
200 PETSC_EXTERN PETSC_DEPRECATED_FUNCTION(3, 22, 0, "SNESNewtonTRSetTolerances()", ) PetscErrorCode SNESSetTrustRegionTolerance(SNES, PetscReal);
201 PETSC_EXTERN PetscErrorCode SNESNewtonTRSetTolerances(SNES, PetscReal, PetscReal, PetscReal);
202 PETSC_EXTERN PetscErrorCode SNESNewtonTRGetTolerances(SNES, PetscReal *, PetscReal *, PetscReal *);
203 PETSC_EXTERN PetscErrorCode SNESNewtonTRSetUpdateParameters(SNES, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal);
204 PETSC_EXTERN PetscErrorCode SNESNewtonTRGetUpdateParameters(SNES, PetscReal *, PetscReal *, PetscReal *, PetscReal *, PetscReal *);
205
206 PETSC_EXTERN PetscErrorCode SNESNewtonTRDCGetRhoFlag(SNES, PetscBool *);
207 PETSC_EXTERN PetscErrorCode SNESNewtonTRDCSetPreCheck(SNES, PetscErrorCode (*)(SNES, Vec, Vec, PetscBool *, PetscCtx), PetscCtx);
208 PETSC_EXTERN PetscErrorCode SNESNewtonTRDCGetPreCheck(SNES, PetscErrorCode (**)(SNES, Vec, Vec, PetscBool *, PetscCtx), PetscCtxRt);
209 PETSC_EXTERN PetscErrorCode SNESNewtonTRDCSetPostCheck(SNES, PetscErrorCode (*)(SNES, Vec, Vec, Vec, PetscBool *, PetscBool *, PetscCtx), PetscCtx);
210 PETSC_EXTERN PetscErrorCode SNESNewtonTRDCGetPostCheck(SNES, PetscErrorCode (**)(SNES, Vec, Vec, Vec, PetscBool *, PetscBool *, PetscCtx), PetscCtxRt);
211
212 PETSC_EXTERN PetscErrorCode SNESGetNonlinearStepFailures(SNES, PetscInt *);
213 PETSC_EXTERN PetscErrorCode SNESSetMaxNonlinearStepFailures(SNES, PetscInt);
214 PETSC_EXTERN PetscErrorCode SNESGetMaxNonlinearStepFailures(SNES, PetscInt *);
215 PETSC_EXTERN PetscErrorCode SNESGetNumberFunctionEvals(SNES, PetscInt *);
216
217 PETSC_EXTERN PetscErrorCode SNESSetLagPreconditioner(SNES, PetscInt);
218 PETSC_EXTERN PetscErrorCode SNESGetLagPreconditioner(SNES, PetscInt *);
219 PETSC_EXTERN PetscErrorCode SNESSetLagJacobian(SNES, PetscInt);
220 PETSC_EXTERN PetscErrorCode SNESGetLagJacobian(SNES, PetscInt *);
221 PETSC_EXTERN PetscErrorCode SNESSetLagPreconditionerPersists(SNES, PetscBool);
222 PETSC_EXTERN PetscErrorCode SNESSetLagJacobianPersists(SNES, PetscBool);
223 PETSC_EXTERN PetscErrorCode SNESSetGridSequence(SNES, PetscInt);
224 PETSC_EXTERN PetscErrorCode SNESGetGridSequence(SNES, PetscInt *);
225
226 PETSC_EXTERN PetscErrorCode SNESGetLinearSolveIterations(SNES, PetscInt *);
227 PETSC_EXTERN PetscErrorCode SNESGetLinearSolveFailures(SNES, PetscInt *);
228 PETSC_EXTERN PetscErrorCode SNESSetMaxLinearSolveFailures(SNES, PetscInt);
229 PETSC_EXTERN PetscErrorCode SNESGetMaxLinearSolveFailures(SNES, PetscInt *);
230 PETSC_EXTERN PetscErrorCode SNESSetCountersReset(SNES, PetscBool);
231 PETSC_EXTERN PetscErrorCode SNESResetCounters(SNES);
232
233 PETSC_EXTERN PetscErrorCode SNESKSPSetUseEW(SNES, PetscBool);
234 PETSC_EXTERN PetscErrorCode SNESKSPGetUseEW(SNES, PetscBool *);
235 PETSC_EXTERN PetscErrorCode SNESKSPSetParametersEW(SNES, PetscInt, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal);
236 PETSC_EXTERN PetscErrorCode SNESKSPGetParametersEW(SNES, PetscInt *, PetscReal *, PetscReal *, PetscReal *, PetscReal *, PetscReal *, PetscReal *);
237
238 PETSC_EXTERN PetscErrorCode SNESMonitorLGRange(SNES, PetscInt, PetscReal, void *);
239
240 PETSC_EXTERN PetscErrorCode SNESSetApplicationContext(SNES, PetscCtx);
241 PETSC_EXTERN PetscErrorCode SNESGetApplicationContext(SNES, PetscCtxRt);
242 PETSC_EXTERN PetscErrorCode SNESSetComputeApplicationContext(SNES, PetscErrorCode (*)(SNES, PetscCtxRt), PetscCtxDestroyFn *);
243
244 PETSC_EXTERN PetscErrorCode SNESPythonSetType(SNES, const char[]);
245 PETSC_EXTERN PetscErrorCode SNESPythonGetType(SNES, const char *[]);
246
247 PETSC_EXTERN PetscErrorCode SNESSetFunctionDomainError(SNES);
248 PETSC_EXTERN PetscErrorCode SNESSetObjectiveDomainError(SNES);
249 PETSC_EXTERN PetscErrorCode SNESSetJacobianDomainError(SNES);
250 PETSC_EXTERN PetscErrorCode SNESSetCheckJacobianDomainError(SNES, PetscBool);
251 PETSC_EXTERN PetscErrorCode SNESGetCheckJacobianDomainError(SNES, PetscBool *);
252
253 #define SNES_CONVERGED_TR_DELTA_DEPRECATED SNES_CONVERGED_TR_DELTA PETSC_DEPRECATED_ENUM(3, 12, 0, "SNES_DIVERGED_TR_DELTA", )
254 #define SNES_DIVERGED_FNORM_NAN_DEPRECATED SNES_DIVERGED_FNORM_NAN PETSC_DEPRECATED_ENUM(3, 25, 0, "SNES_DIVERGED_FUNCTION_NANORINF", )
255 /*E
256 SNESConvergedReason - reason a `SNESSolve()` was determined to have converged or diverged
257
258 Values:
259 + `SNES_CONVERGED_FNORM_ABS` - $ ||F|| \le abstol $
260 . `SNES_CONVERGED_FNORM_RELATIVE` - $ ||F|| <= rtol*||F(x_0))|| $ where $x_0 $ is the initial guess
261 . `SNES_CONVERGED_SNORM_RELATIVE` - The 2-norm of the last step $ \le stol * ||x|| $ where $ x $ is the current solution
262 . `SNES_CONVERGED_USER` - The user has indicated convergence for an arbitrary reason
263 . `SNES_DIVERGED_FUNCTION_COUNT` - The user provided function has been called more times than the maximum set in `SNESSetTolerances()`
264 . `SNES_DIVERGED_DTOL` - The norm of the function has increased by a factor of divtol set with `SNESSetDivergenceTolerance()`
265 . `SNES_DIVERGED_FUNCTION_NANORINF` - the 2-norm of the current function evaluation is not-a-number (NaN) or infinity, (this
266 is usually caused by a division of 0 by 0) and the solver could not recover from this (by, for example, cutting the step size)
267 . `SNES_DIVERGED_OBJECTIVE_NANORINF` - the object function evaluation is not-a-number (NaN) or infinity, (this
268 is usually caused by a division of 0 by 0) and the solver could not recover from this (by, for example, cutting the step size)
269 . `SNES_DIVERGED_FUNCTION_DOMAIN` - the function evaluation occurred outside the function's domain (function callback provided by
270 `SNESSetFunction()` called `SNESSetObjectiveDomainError()`) and the solver could not recover from this (by, for example, cutting the step size)
271 . `SNES_DIVERGED_OBJECTIVE_DOMAIN` - the object function evaluation occurred outside the function's domain (function callback provided by
272 `SNESSetObjective()` called `SNESSetObjectiveDomainError()`) and the solver could not recover from this (by, for example, cutting the step size)
273 . `SNES_DIVERGED_JACOBIAN_DOMAIN` - the Jacobian evaluation occurred outside the function's domain (function callback provided by
274 `SNESSetJacobian()` called `SNESSetJacobianDomainError()`)
275 . `SNES_DIVERGED_MAX_IT` - `SNESSolve()` has reached the maximum number of iterations requested
276 . `SNES_DIVERGED_LINE_SEARCH` - The line search has failed. This only occurs for `SNES` solvers that use a line search
277 . `SNES_DIVERGED_LOCAL_MIN` - the algorithm seems to have stagnated at a local minimum that is not zero.
278 - `SNES_CONVERGED_ITERATING - this only occurs if `SNESGetConvergedReason()` is called during the `SNESSolve()`
279
280 Level: beginner
281
282 Notes:
283 The two most common reasons for divergence are an incorrectly coded or computed Jacobian or failure or lack of convergence in the linear system
284 (in this case we recommend
285 testing with `-pc_type lu` to eliminate the linear solver as the cause of the problem).
286
287 `SNES_DIVERGED_LOCAL_MIN` can only occur when using a `SNES` solver that uses a line search (`SNESLineSearch`).
288 The line search wants to $ \min Q(\alpha) = 1/2 || F(x + \alpha s) ||^2_2 $ this occurs
289 at $ Q'(\alpha) = s^T F'(x+\alpha s)^T F(x+\alpha s) = 0$. If $s$ is the Newton direction $ - F'(x)^(-1)F(x)$ then
290 $ Q'(\alpha) = -F(x)^T F'(x)^(-1)^T F'(x+\alpha s)F(x+\alpha s)$; when $\alpha = 0$
291 $Q'(0) = - ||F(x)||^2_2 $ which is always NEGATIVE if $F'(x)$ is invertible. This means the Newton
292 direction is a descent direction and the line search should succeed if $\alpha $ is small enough.
293
294 If $F'(x)$ is NOT invertible AND $F'(x)^T F(x) = 0 $ then $Q'(0) = 0 $ and the Newton direction
295 is NOT a descent direction so the line search will fail. All one can do at this point
296 is change the initial guess and try again.
297
298 An alternative explanation: Newton's method can be regarded as replacing the function with
299 its linear approximation and minimizing the 2-norm of that. That is $F(x+s) \approx F(x) + F'(x)s$
300 so we minimize $ || F(x) + F'(x) s ||^2_2$ using Least Squares. If $F'(x)$ is invertible then
301 $s = - F'(x)^(-1)F(x)$ otherwise $F'(x)^T F'(x) s = -F'(x)^T F(x)$. If $F'(x)^T F(x)$ is NOT zero then there
302 exists a nontrivial (that is $F'(x)s \ne 0$) solution to the equation and this direction is
303 $s = - [F'(x)^T F'(x)]^(-1) F'(x)^T F(x)$ so $Q'(0) = - F(x)^T F'(x) [F'(x)^T F'(x)]^(-T) F'(x)^T F(x)
304 = - (F'(x)^T F(x)) [F'(x)^T F'(x)]^(-T) (F'(x)^T F(x))$. Since we are assuming $(F'(x)^T F(x)) \ne 0$
305 and $F'(x)^T F'(x)$ has no negative eigenvalues $Q'(0) < 0$ so $s$ is a descent direction and the line
306 search should succeed for small enough $\alpha$.
307
308 Note that this RARELY happens in practice. Far more likely the linear system is not being solved
309 (well enough?) or the Jacobian is wrong.
310
311 `SNES_DIVERGED_MAX_IT` means that the solver reached the maximum number of iterations without satisfying any
312 convergence criteria. `SNES_CONVERGED_ITS` means that `SNESConvergedSkip()` was chosen as the convergence test;
313 thus the usual convergence criteria have not been checked and may or may not be satisfied.
314
315 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `KSPConvergedReason`, `SNESSetConvergenceTest()`, `SNESSetTolerances()`
316 E*/
317 typedef enum { /* converged */
318 SNES_CONVERGED_FNORM_ABS = 2, /* ||F|| < atol */
319 SNES_CONVERGED_FNORM_RELATIVE = 3, /* ||F|| < rtol*||F_initial|| */
320 SNES_CONVERGED_SNORM_RELATIVE = 4, /* Newton computed step size small; || delta x || < stol || x || */
321 SNES_CONVERGED_ITS = 5, /* maximum iterations reached */
322 SNES_BREAKOUT_INNER_ITER = 6, /* Flag to break out of inner loop after checking custom convergence, used in multi-phase flow when state changes */
323 SNES_CONVERGED_USER = 7, /* The user has indicated convergence for an arbitrary reason */
324 /* diverged */
325 SNES_DIVERGED_FUNCTION_DOMAIN = -1, /* the new x location passed the function is not in the domain of F */
326 SNES_DIVERGED_FUNCTION_COUNT = -2,
327 SNES_DIVERGED_LINEAR_SOLVE = -3, /* the linear solve failed */
328 SNES_DIVERGED_FUNCTION_NANORINF = -4,
329 SNES_DIVERGED_FNORM_NAN_DEPRECATED = -4,
330 SNES_DIVERGED_MAX_IT = -5,
331 SNES_DIVERGED_LINE_SEARCH = -6, /* the line search failed */
332 SNES_DIVERGED_INNER = -7, /* inner solve failed */
333 SNES_DIVERGED_LOCAL_MIN = -8, /* || J^T b || is small, implies converged to local minimum of F() */
334 SNES_DIVERGED_DTOL = -9, /* || F || > divtol*||F_initial|| */
335 SNES_DIVERGED_JACOBIAN_DOMAIN = -10, /* Jacobian calculation does not make sense */
336 SNES_DIVERGED_TR_DELTA = -11,
337 SNES_CONVERGED_TR_DELTA_DEPRECATED = -11,
338 SNES_DIVERGED_USER = -12, /* The user has indicated divergence for an arbitrary reason */
339 SNES_DIVERGED_OBJECTIVE_DOMAIN = -13,
340 SNES_DIVERGED_OBJECTIVE_NANORINF = -14,
341
342 SNES_CONVERGED_ITERATING = 0
343 } SNESConvergedReason;
344 PETSC_EXTERN const char *const *SNESConvergedReasons;
345
346 /*MC
347 SNES_CONVERGED_FNORM_ABS - $||F|| \le abstol$
348
349 Level: beginner
350
351 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `SNESConvergedReason`, `SNESSetTolerances()`
352 M*/
353
354 /*MC
355 SNES_CONVERGED_FNORM_RELATIVE - $||F|| \le rtol*||F(x_0)||$ where $x_0$ is the initial guess
356
357 Level: beginner
358
359 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `SNESConvergedReason`, `SNESSetTolerances()`
360 M*/
361
362 /*MC
363 SNES_CONVERGED_SNORM_RELATIVE - The 2-norm of the last step $\le stol * ||x||$ where $x$ is the current
364 solution and $stol$ is the 4th argument to `SNESSetTolerances()`
365
366 Options Database Key:
367 -snes_stol <stol> - the step tolerance
368
369 Level: beginner
370
371 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `SNESConvergedReason`, `SNESSetTolerances()`
372 M*/
373
374 /*MC
375 SNES_DIVERGED_FUNCTION_COUNT - The user provided function has been called more times then the final
376 argument to `SNESSetTolerances()`
377
378 Level: beginner
379
380 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `SNESConvergedReason`, `SNESSetTolerances()`
381 M*/
382
383 /*MC
384 SNES_DIVERGED_DTOL - The norm of the function has increased by a factor of divtol set with `SNESSetDivergenceTolerance()`
385
386 Level: beginner
387
388 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `SNESConvergedReason`, `SNESSetTolerances()`, `SNESSetDivergenceTolerance()`
389 M*/
390
391 /*MC
392 SNES_DIVERGED_FUNCTION_NANORINF - the 2-norm of the current function evaluation is not-a-number (NaN) or infinity, this
393 is usually caused by a division of 0 by 0, or infinity. See `SNESSetFunctionDomainError()`
394
395 Level: beginner
396
397 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `SNESConvergedReason`, `SNESSetTolerances()`
398 M*/
399
400 /*MC
401 SNES_DIVERGED_FUNCTION_DOMAIN - the function provided with `SNESSetFunction()` called `SNESSetFunctionDomainError()` and
402 the solver could not recoverer.
403
404 Level: beginner
405
406 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `SNESConvergedReason`, `SNESSetTolerances()`
407 M*/
408
409 /*MC
410 SNES_DIVERGED_OBJECTIVE_DOMAIN - the function provided with `SNESSetObjective()` called `SNESSetObjectiveDomainError()` and
411 the solver could not recoverer.
412
413 Level: beginner
414
415 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `SNESConvergedReason`, `SNESSetTolerances()`
416 M*/
417
418 /*MC
419 SNES_DIVERGED_JACOBIAN_DOMAIN - the function provided with `SNESSetJacobian()` called `SNESSetJacobianDomainError()`
420
421 Level: beginner
422
423 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `SNESConvergedReason`, `SNESSetTolerances()`
424 M*/
425
426 /*MC
427 SNES_DIVERGED_MAX_IT - SNESSolve() has reached the maximum number of iterations requested
428
429 Level: beginner
430
431 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `SNESConvergedReason`, `SNESSetTolerances()`
432 M*/
433
434 /*MC
435 SNES_DIVERGED_LINE_SEARCH - The line search has failed. This only occurs for a `SNES` solvers that use a line search
436
437 Level: beginner
438
439 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `SNESConvergedReason`, `SNESSetTolerances()`, `SNESLineSearch`
440 M*/
441
442 /*MC
443 SNES_DIVERGED_LOCAL_MIN - the algorithm seems to have stagnated at a local minimum that is not zero.
444 See the manual page for `SNESConvergedReason` for more details
445
446 Level: beginner
447
448 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `SNESConvergedReason`, `SNESSetTolerances()`
449 M*/
450
451 /*MC
452 SNES_CONERGED_ITERATING - this only occurs if `SNESGetConvergedReason()` is called during the `SNESSolve()`
453
454 Level: beginner
455
456 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `SNESConvergedReason`, `SNESSetTolerances()`
457 M*/
458
459 PETSC_EXTERN PetscErrorCode SNESSetConvergenceTest(SNES, PetscErrorCode (*)(SNES, PetscInt, PetscReal, PetscReal, PetscReal, SNESConvergedReason *, void *), void *, PetscCtxDestroyFn *);
460 PETSC_EXTERN PetscErrorCode SNESConvergedDefault(SNES, PetscInt, PetscReal, PetscReal, PetscReal, SNESConvergedReason *, void *);
461 PETSC_EXTERN PetscErrorCode SNESConvergedSkip(SNES, PetscInt, PetscReal, PetscReal, PetscReal, SNESConvergedReason *, void *);
462 PETSC_EXTERN PetscErrorCode SNESConvergedCorrectPressure(SNES, PetscInt, PetscReal, PetscReal, PetscReal, SNESConvergedReason *, void *);
463 PETSC_EXTERN PetscErrorCode SNESGetConvergedReason(SNES, SNESConvergedReason *);
464 PETSC_EXTERN PetscErrorCode SNESGetConvergedReasonString(SNES, const char **);
465 PETSC_EXTERN PetscErrorCode SNESSetConvergedReason(SNES, SNESConvergedReason);
466
SNESSkipConverged(void)467 PETSC_DEPRECATED_FUNCTION(3, 5, 0, "SNESConvergedSkip()", ) static inline void SNESSkipConverged(void)
468 { /* never called */
469 }
470 #define SNESSkipConverged (SNESSkipConverged, SNESConvergedSkip)
471
472 /*S
473 SNESInitialGuessFn - A prototype of a `SNES` compute initial guess function that would be passed to `SNESSetComputeInitialGuess()`
474
475 Calling Sequence:
476 + snes - `SNES` context
477 . u - output vector to contain initial guess
478 - ctx - [optional] user-defined function context
479
480 Level: beginner
481
482 .seealso: [](ch_snes), `SNES`, `SNESSetComputeInitialGuess()`, `SNESSetFunction()`, `SNESGetFunction()`, `SNESJacobianFn`, `SNESFunctionFn`
483 S*/
484 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode SNESInitialGuessFn(SNES snes, Vec u, PetscCtx ctx);
485
486 /*S
487 SNESFunctionFn - A prototype of a `SNES` evaluation function that would be passed to `SNESSetFunction()`
488
489 Calling Sequence:
490 + snes - `SNES` context
491 . u - input vector
492 . F - function vector
493 - ctx - [optional] user-defined function context
494
495 Level: beginner
496
497 .seealso: [](ch_snes), `SNES`, `SNESSetFunction()`, `SNESGetFunction()`, `SNESJacobianFn`, `SNESNGSFn`
498 S*/
499 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode SNESFunctionFn(SNES snes, Vec u, Vec F, PetscCtx ctx);
500
501 /*S
502 SNESObjectiveFn - A prototype of a `SNES` objective evaluation function that would be passed to `SNESSetObjective()`
503
504 Calling Sequence:
505 + snes - `SNES` context
506 . u - input vector
507 . o - output value
508 - ctx - [optional] user-defined function context
509
510 Level: beginner
511
512 .seealso: [](ch_snes), `SNES`, `SNESSetFunction()`, `SNESGetFunction()`, `SNESJacobianFn`, `SNESNGSFn`
513 S*/
514 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode SNESObjectiveFn(SNES snes, Vec u, PetscReal *o, PetscCtx ctx);
515
516 /*S
517 SNESJacobianFn - A prototype of a `SNES` Jacobian evaluation function that would be passed to `SNESSetJacobian()`
518
519 Calling Sequence:
520 + snes - the `SNES` context obtained from `SNESCreate()`
521 . u - input vector
522 . Amat - (approximate) Jacobian matrix
523 . Pmat - matrix used to construct the preconditioner, often the same as `Amat`
524 - ctx - [optional] user-defined context for matrix evaluation routine
525
526 Level: beginner
527
528 .seealso: [](ch_snes), `SNES`, `SNESSetJacobian()`, `SNESGetJacobian()`, `SNESFunctionFn`, `SNESNGSFn`
529 S*/
530 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode SNESJacobianFn(SNES snes, Vec u, Mat Amat, Mat Pmat, PetscCtx ctx);
531
532 /*S
533 SNESNGSFn - A prototype of a `SNES` nonlinear Gauss-Seidel function that would be passed to `SNESSetNGS()`
534
535 Calling Sequence:
536 + snes - the `SNES` context obtained from `SNESCreate()`
537 . u - the current solution, updated in place
538 . b - the right-hand side vector (which may be `NULL`)
539 - ctx - [optional] user-defined context for matrix evaluation routine
540
541 Level: beginner
542
543 .seealso: [](ch_snes), `SNES`, `SNESSetJacobian()`, `SNESGetJacobian()`, `SNESFunctionFn`, `SNESSetFunction()`, `SNESGetFunction()`, `SNESJacobianFn`
544 S*/
545 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode SNESNGSFn(SNES snes, Vec u, Vec b, PetscCtx ctx);
546
547 /*S
548 SNESUpdateFn - A prototype of a `SNES` update function that would be passed to `SNESSetUpdate()`
549
550 Calling Sequence:
551 + snes - `SNES` context
552 - step - the current iteration index
553
554 Level: advanced
555
556 .seealso: [](ch_snes), `SNES`, `SNESSetUpdate()`
557 S*/
558 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode SNESUpdateFn(SNES snes, PetscInt step);
559
560 /* --------- Solving systems of nonlinear equations --------------- */
561 PETSC_EXTERN PetscErrorCode SNESSetFunction(SNES, Vec, SNESFunctionFn *, PetscCtx);
562 PETSC_EXTERN PetscErrorCode SNESGetFunction(SNES, Vec *, SNESFunctionFn **, PetscCtxRt);
563 PETSC_EXTERN PetscErrorCode SNESComputeFunction(SNES, Vec, Vec);
564 PETSC_EXTERN PetscErrorCode SNESComputeMFFunction(SNES, Vec, Vec);
565 PETSC_EXTERN PetscErrorCode SNESSetInitialFunction(SNES, Vec);
566
567 PETSC_EXTERN PetscErrorCode SNESSetJacobian(SNES, Mat, Mat, SNESJacobianFn *, PetscCtx);
568 PETSC_EXTERN PetscErrorCode SNESGetJacobian(SNES, Mat *, Mat *, SNESJacobianFn **, PetscCtxRt);
569 PETSC_EXTERN SNESFunctionFn SNESObjectiveComputeFunctionDefaultFD;
570 PETSC_EXTERN SNESJacobianFn SNESComputeJacobianDefault;
571 PETSC_EXTERN SNESJacobianFn SNESComputeJacobianDefaultColor;
572 PETSC_EXTERN PetscErrorCode SNESPruneJacobianColor(SNES, Mat, Mat);
573 PETSC_EXTERN PetscErrorCode SNESSetComputeInitialGuess(SNES, SNESInitialGuessFn *, PetscCtx);
574 PETSC_EXTERN PetscErrorCode SNESSetPicard(SNES, Vec, SNESFunctionFn *, Mat, Mat, SNESJacobianFn *, PetscCtx);
575 PETSC_EXTERN PetscErrorCode SNESGetPicard(SNES, Vec *, SNESFunctionFn **, Mat *, Mat *, SNESJacobianFn **, PetscCtxRt);
576 PETSC_EXTERN SNESFunctionFn SNESPicardComputeFunction;
577 PETSC_EXTERN SNESFunctionFn SNESPicardComputeMFFunction;
578 PETSC_EXTERN SNESJacobianFn SNESPicardComputeJacobian;
579
580 PETSC_EXTERN PetscErrorCode SNESSetObjective(SNES, SNESObjectiveFn *, PetscCtx);
581 PETSC_EXTERN PetscErrorCode SNESGetObjective(SNES, SNESObjectiveFn **, PetscCtxRt);
582 PETSC_EXTERN PetscErrorCode SNESComputeObjective(SNES, Vec, PetscReal *);
583
584 PETSC_EXTERN PetscErrorCode SNESSetUpdate(SNES, SNESUpdateFn *);
585
586 /*E
587 SNESNormSchedule - Frequency with which the norm is computed during a nonliner solve
588
589 Values:
590 + `SNES_NORM_DEFAULT` - use the default behavior for the current `SNESType`
591 . `SNES_NORM_NONE` - avoid all norm computations
592 . `SNES_NORM_ALWAYS` - compute the norms whenever possible
593 . `SNES_NORM_INITIAL_ONLY` - compute the norm only when the algorithm starts
594 . `SNES_NORM_FINAL_ONLY` - compute the norm only when the algorithm finishes
595 - `SNES_NORM_INITIAL_FINAL_ONLY` - compute the norm at the start and end of the algorithm
596
597 Level: advanced
598
599 Notes:
600 Support for these is highly dependent on the solver.
601
602 Some options limit the convergence tests that can be used.
603
604 The `SNES_NORM_NONE` option is most commonly used when the nonlinear solver is being used as a smoother, for example for `SNESFAS`
605
606 This is primarily used to turn off extra norm and function computation
607 when the solvers are composed.
608
609 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `KSPSetNormType()`,
610 `KSPSetConvergenceTest()`, `KSPSetPCSide()`
611 E*/
612 typedef enum {
613 SNES_NORM_DEFAULT = -1,
614 SNES_NORM_NONE = 0,
615 SNES_NORM_ALWAYS = 1,
616 SNES_NORM_INITIAL_ONLY = 2,
617 SNES_NORM_FINAL_ONLY = 3,
618 SNES_NORM_INITIAL_FINAL_ONLY = 4
619 } SNESNormSchedule;
620 PETSC_EXTERN const char *const *const SNESNormSchedules;
621
622 /*MC
623 SNES_NORM_NONE - Don't compute function and its L2 norm when possible
624
625 Level: advanced
626
627 Note:
628 This is most useful for stationary solvers with a fixed number of iterations used as smoothers.
629
630 .seealso: [](ch_snes), `SNESNormSchedule`, `SNES`, `SNESSetNormSchedule()`, `SNES_NORM_DEFAULT`
631 M*/
632
633 /*MC
634 SNES_NORM_ALWAYS - Compute the function and its L2 norm at each iteration.
635
636 Level: advanced
637
638 Note:
639 Most solvers will use this no matter what norm type is passed to them.
640
641 .seealso: [](ch_snes), `SNESNormSchedule`, `SNES`, `SNESSetNormSchedule()`, `SNES_NORM_NONE`
642 M*/
643
644 /*MC
645 SNES_NORM_INITIAL_ONLY - Compute the function and its L2 at iteration 0, but do not update it.
646
647 Level: advanced
648
649 Notes:
650 This method is useful in composed methods, when a true solution might actually be found before `SNESSolve()` is called.
651 This option enables the solve to abort on the zeroth iteration if this is the case.
652
653 For solvers that require the computation of the L2 norm of the function as part of the method, this merely cancels
654 the norm computation at the last iteration (if possible).
655
656 .seealso: [](ch_snes), `SNESNormSchedule`, `SNES`, `SNESSetNormSchedule()`, `SNES_NORM_FINAL_ONLY`, `SNES_NORM_INITIAL_FINAL_ONLY`
657 M*/
658
659 /*MC
660 SNES_NORM_FINAL_ONLY - Compute the function and its L2 norm on only the final iteration.
661
662 Level: advanced
663
664 Note:
665 For solvers that require the computation of the L2 norm of the function as part of the method, behaves
666 exactly as `SNES_NORM_DEFAULT`. This method is useful when the function is gotten after `SNESSolve()` and
667 used in subsequent computation for methods that do not need the norm computed during the rest of the
668 solution procedure.
669
670 .seealso: [](ch_snes), `SNESNormSchedule`, `SNES`, `SNESSetNormSchedule()`, `SNES_NORM_INITIAL_ONLY`, `SNES_NORM_INITIAL_FINAL_ONLY`
671 M*/
672
673 /*MC
674 SNES_NORM_INITIAL_FINAL_ONLY - Compute the function and its L2 norm on only the initial and final iterations.
675
676 Level: advanced
677
678 Note:
679 This method combines the benefits of `SNES_NORM_INITIAL_ONLY` and `SNES_NORM_FINAL_ONLY`.
680
681 .seealso: [](ch_snes), `SNESNormSchedule`, `SNES`, `SNESSetNormSchedule()`, `SNES_NORM_SNES_NORM_INITIAL_ONLY`, `SNES_NORM_FINAL_ONLY`
682 M*/
683
684 PETSC_EXTERN PetscErrorCode SNESSetNormSchedule(SNES, SNESNormSchedule);
685 PETSC_EXTERN PetscErrorCode SNESGetNormSchedule(SNES, SNESNormSchedule *);
686 PETSC_EXTERN PetscErrorCode SNESSetFunctionNorm(SNES, PetscReal);
687 PETSC_EXTERN PetscErrorCode SNESGetFunctionNorm(SNES, PetscReal *);
688 PETSC_EXTERN PetscErrorCode SNESGetUpdateNorm(SNES, PetscReal *);
689 PETSC_EXTERN PetscErrorCode SNESGetSolutionNorm(SNES, PetscReal *);
690
691 /*E
692 SNESFunctionType - Type of function computed
693
694 Values:
695 + `SNES_FUNCTION_DEFAULT` - the default behavior for the current `SNESType`
696 . `SNES_FUNCTION_UNPRECONDITIONED` - the original function provided
697 - `SNES_FUNCTION_PRECONDITIONED` - the modification of the function by the preconditioner
698
699 Level: advanced
700
701 Note:
702 Support for these is dependent on the solver.
703
704 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `KSPSetNormType()`,
705 `KSPSetConvergenceTest()`, `KSPSetPCSide()`
706 E*/
707 typedef enum {
708 SNES_FUNCTION_DEFAULT = -1,
709 SNES_FUNCTION_UNPRECONDITIONED = 0,
710 SNES_FUNCTION_PRECONDITIONED = 1
711 } SNESFunctionType;
712 PETSC_EXTERN const char *const *const SNESFunctionTypes;
713
714 PETSC_EXTERN PetscErrorCode SNESSetFunctionType(SNES, SNESFunctionType);
715 PETSC_EXTERN PetscErrorCode SNESGetFunctionType(SNES, SNESFunctionType *);
716
717 PETSC_EXTERN PetscErrorCode SNESSetNGS(SNES, SNESNGSFn *, PetscCtx);
718 PETSC_EXTERN PetscErrorCode SNESGetNGS(SNES, SNESNGSFn **, PetscCtxRt);
719 PETSC_EXTERN PetscErrorCode SNESComputeNGS(SNES, Vec, Vec);
720
721 PETSC_EXTERN PetscErrorCode SNESNGSSetSweeps(SNES, PetscInt);
722 PETSC_EXTERN PetscErrorCode SNESNGSGetSweeps(SNES, PetscInt *);
723 PETSC_EXTERN PetscErrorCode SNESNGSSetTolerances(SNES, PetscReal, PetscReal, PetscReal, PetscInt);
724 PETSC_EXTERN PetscErrorCode SNESNGSGetTolerances(SNES, PetscReal *, PetscReal *, PetscReal *, PetscInt *);
725
726 PETSC_EXTERN PetscErrorCode SNESSetAlwaysComputesFinalResidual(SNES, PetscBool);
727 PETSC_EXTERN PetscErrorCode SNESGetAlwaysComputesFinalResidual(SNES, PetscBool *);
728
729 PETSC_EXTERN PetscErrorCode SNESShellGetContext(SNES, PetscCtxRt);
730 PETSC_EXTERN PetscErrorCode SNESShellSetContext(SNES, PetscCtx);
731 PETSC_EXTERN PetscErrorCode SNESShellSetSolve(SNES, PetscErrorCode (*)(SNES, Vec));
732
733 /* --------- Routines specifically for line search methods --------------- */
734
735 /*S
736 SNESLineSearch - Abstract PETSc object that manages line-search operations for nonlinear solvers
737
738 Level: beginner
739
740 .seealso: [](ch_snes), `SNESLineSearchType`, `SNESLineSearchCreate()`, `SNESLineSearchSetType()`, `SNES`
741 S*/
742 typedef struct _p_LineSearch *SNESLineSearch;
743
744 /*J
745 SNESLineSearchType - String with the name of a PETSc line search method `SNESLineSearch`. Provides all the linesearches for the nonlinear solvers, `SNES`,
746 in PETSc.
747
748 Values:
749 + `SNESLINESEARCHBASIC` - (or equivalently `SNESLINESEARCHNONE`) Simple damping line search, defaults to using the full Newton step
750 . `SNESLINESEARCHBT` - Backtracking line search over the L2 norm of the function or an objective function
751 . `SNESLINESEARCHSECANT` - Secant line search over the L2 norm of the function or an objective function
752 . `SNESLINESEARCHCP` - Critical point secant line search assuming $F(x) = \nabla G(x)$ for some unknown $G(x)$
753 . `SNESLINESEARCHNLEQERR` - Affine-covariant error-oriented linesearch
754 - `SNESLINESEARCHBISECTION` - bisection line search for a root in the directional derivative
755 - `SNESLINESEARCHSHELL` - User provided `SNESLineSearch` implementation
756
757 Level: beginner
758
759 Note:
760 Use `SNESLineSearchSetType()` or the options database key `-snes_linesearch_type` to set
761 the specific line search algorithm to use with a given `SNES` object. Not all `SNESType` can utilize a line search.
762
763 .seealso: [](ch_snes), `SNESLineSearch`, `SNESLineSearchSetType()`, `SNES`
764 J*/
765 typedef const char *SNESLineSearchType;
766 #define SNESLINESEARCHBT "bt"
767 #define SNESLINESEARCHNLEQERR "nleqerr"
768 #define SNESLINESEARCHBASIC "basic"
769 #define SNESLINESEARCHNONE "none"
770 #define SNESLINESEARCHSECANT "secant"
771 #define SNESLINESEARCHL2 PETSC_DEPRECATED_MACRO(3, 24, 0, "SNESLINESEARCHSECANT", ) "secant"
772 #define SNESLINESEARCHCP "cp"
773 #define SNESLINESEARCHSHELL "shell"
774 #define SNESLINESEARCHNCGLINEAR "ncglinear"
775 #define SNESLINESEARCHBISECTION "bisection"
776
777 PETSC_EXTERN PetscFunctionList SNESList;
778 PETSC_EXTERN PetscClassId SNESLINESEARCH_CLASSID;
779 PETSC_EXTERN PetscFunctionList SNESLineSearchList;
780
781 #define SNES_LINESEARCH_ORDER_LINEAR 1
782 #define SNES_LINESEARCH_ORDER_QUADRATIC 2
783 #define SNES_LINESEARCH_ORDER_CUBIC 3
784
785 /*S
786 SNESLineSearchVIProjectFn - A prototype of a `SNES` function that projects a vector onto the VI bounds, passed to `SNESLineSearchSetVIFunctions()`
787
788 Calling Sequence:
789 + snes - `SNES` context
790 - u - the vector to project to the bounds
791
792 Level: advanced
793
794 Note:
795 The deprecated `SNESLineSearchVIProjectFunc` still works as a replacement for `SNESLineSearchVIProjectFn` *.
796
797 .seealso: [](ch_snes), `SNES`, `SNESLineSearch`
798 S*/
799 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode SNESLineSearchVIProjectFn(SNES snes, Vec u);
800 PETSC_EXTERN_TYPEDEF typedef SNESLineSearchVIProjectFn *SNESLineSearchVIProjectFunc PETSC_DEPRECATED_TYPEDEF(3, 21, 0, "SNESLineSearchVIProjectFn*", );
801
802 /*S
803 SNESLineSearchVIProjectFn - A prototype of a `SNES` function that computes the norm of the active set variables in a vector in a VI solve,
804 passed to `SNESLineSearchSetVIFunctions()`
805
806 Calling Sequence:
807 + snes - `SNES` context
808 . f - the vector to compute the norm of
809 . u - the current solution, entries that are on the VI bounds are ignored
810 - fnorm - the resulting norm
811
812 Level: advanced
813
814 Note:
815 The deprecated `SNESLineSearchVINormFunc` still works as a replacement for `SNESLineSearchVINormFn` *.
816
817 .seealso: [](ch_snes), `SNES`, `SNESLineSearch`
818 S*/
819 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode SNESLineSearchVINormFn(SNES snes, Vec f, Vec u, PetscReal *fnorm);
820 PETSC_EXTERN_TYPEDEF typedef SNESLineSearchVINormFn *SNESLineSearchVINormFunc PETSC_DEPRECATED_TYPEDEF(3, 21, 0, "SNESLineSearchVINormFnn*", );
821
822 /*S
823 SNESLineSearchVIDirDerivFn - A prototype of a `SNES` function that computes the directional derivative considering the VI bounds, passed to `SNESLineSearchSetVIFunctions()`
824
825 Calling Sequence:
826 + snes - `SNES` context
827 . f - the function vector to compute the directional derivative with
828 . u - the current solution, entries that are on the VI bounds are ignored
829 . y - the direction to compute the directional derivative
830 - fty - the resulting directional derivative
831
832 Level: advanced
833
834 .seealso: [](ch_snes), `SNES`, `SNESLineSearch`, `SNESLineSearchVIProjectFn`, `SNESLineSearchVIProjectFn`, `SNESLineSearchSetVIFunctions()`, `SNESLineSearchGetVIFunctions()`
835 S*/
836 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode SNESLineSearchVIDirDerivFn(SNES snes, Vec f, Vec u, Vec y, PetscScalar *fty);
837
838 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode SNESLineSearchApplyFn(SNESLineSearch);
839 PETSC_EXTERN_TYPEDEF typedef SNESLineSearchApplyFn *SNESLineSearchApplyFunc PETSC_DEPRECATED_TYPEDEF(3, 21, 0, "SNESLineSearchApplyFn*", );
840 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode SNESLineSearchShellApplyFn(SNESLineSearch, PetscCtx);
841 PETSC_EXTERN_TYPEDEF typedef SNESLineSearchShellApplyFn *SNESLineSearchUserFunc PETSC_DEPRECATED_TYPEDEF(3, 21, 0, "SNESLineSearchApplyFn*", );
842
843 PETSC_EXTERN PetscErrorCode SNESLineSearchCreate(MPI_Comm, SNESLineSearch *);
844 PETSC_EXTERN PetscErrorCode SNESLineSearchReset(SNESLineSearch);
845 PETSC_EXTERN PetscErrorCode SNESLineSearchView(SNESLineSearch, PetscViewer);
846 PETSC_EXTERN PetscErrorCode SNESLineSearchDestroy(SNESLineSearch *);
847 PETSC_EXTERN PetscErrorCode SNESLineSearchGetType(SNESLineSearch, SNESLineSearchType *);
848 PETSC_EXTERN PetscErrorCode SNESLineSearchSetType(SNESLineSearch, SNESLineSearchType);
849 PETSC_EXTERN PetscErrorCode SNESLineSearchSetFromOptions(SNESLineSearch);
850 PETSC_EXTERN PetscErrorCode SNESLineSearchSetFunction(SNESLineSearch, PetscErrorCode (*)(SNES, Vec, Vec));
851 PETSC_EXTERN PetscErrorCode SNESLineSearchSetUp(SNESLineSearch);
852 PETSC_EXTERN PetscErrorCode SNESLineSearchApply(SNESLineSearch, Vec, Vec, PetscReal *, Vec);
853 PETSC_EXTERN PetscErrorCode SNESLineSearchPreCheck(SNESLineSearch, Vec, Vec, PetscBool *);
854 PETSC_EXTERN PetscErrorCode SNESLineSearchPostCheck(SNESLineSearch, Vec, Vec, Vec, PetscBool *, PetscBool *);
855 PETSC_EXTERN PetscErrorCode SNESLineSearchSetWorkVecs(SNESLineSearch, PetscInt);
856
857 /* set the functions for precheck and postcheck */
858
859 PETSC_EXTERN PetscErrorCode SNESLineSearchSetPreCheck(SNESLineSearch, PetscErrorCode (*)(SNESLineSearch, Vec, Vec, PetscBool *, PetscCtx), PetscCtx);
860 PETSC_EXTERN PetscErrorCode SNESLineSearchSetPostCheck(SNESLineSearch, PetscErrorCode (*)(SNESLineSearch, Vec, Vec, Vec, PetscBool *, PetscBool *, PetscCtx), PetscCtx);
861
862 PETSC_EXTERN PetscErrorCode SNESLineSearchGetPreCheck(SNESLineSearch, PetscErrorCode (**)(SNESLineSearch, Vec, Vec, PetscBool *, PetscCtx), PetscCtxRt);
863 PETSC_EXTERN PetscErrorCode SNESLineSearchGetPostCheck(SNESLineSearch, PetscErrorCode (**)(SNESLineSearch, Vec, Vec, Vec, PetscBool *, PetscBool *, PetscCtx), PetscCtxRt);
864
865 /* set the functions for VI-specific line search operations */
866
867 PETSC_EXTERN PetscErrorCode SNESLineSearchSetVIFunctions(SNESLineSearch, SNESLineSearchVIProjectFn *, SNESLineSearchVINormFn *, SNESLineSearchVIDirDerivFn *);
868 PETSC_EXTERN PetscErrorCode SNESLineSearchGetVIFunctions(SNESLineSearch, SNESLineSearchVIProjectFn **, SNESLineSearchVINormFn **, SNESLineSearchVIDirDerivFn **);
869
870 /* pointers to the associated SNES in order to be able to get the function evaluation out */
871 PETSC_EXTERN PetscErrorCode SNESLineSearchSetSNES(SNESLineSearch, SNES);
872 PETSC_EXTERN PetscErrorCode SNESLineSearchGetSNES(SNESLineSearch, SNES *);
873
874 /* set and get the parameters and vectors */
875 PETSC_EXTERN PetscErrorCode SNESLineSearchGetTolerances(SNESLineSearch, PetscReal *, PetscReal *, PetscReal *, PetscReal *, PetscReal *, PetscInt *);
876 PETSC_EXTERN PetscErrorCode SNESLineSearchSetTolerances(SNESLineSearch, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal, PetscInt);
877
878 PETSC_EXTERN PetscErrorCode SNESLineSearchPreCheckPicard(SNESLineSearch, Vec, Vec, PetscBool *, PetscCtx);
879
880 PETSC_EXTERN PetscErrorCode SNESLineSearchGetLambda(SNESLineSearch, PetscReal *);
881 PETSC_EXTERN PetscErrorCode SNESLineSearchSetLambda(SNESLineSearch, PetscReal);
882
883 PETSC_EXTERN PetscErrorCode SNESLineSearchGetDamping(SNESLineSearch, PetscReal *);
884 PETSC_EXTERN PetscErrorCode SNESLineSearchSetDamping(SNESLineSearch, PetscReal);
885
886 PETSC_EXTERN PetscErrorCode SNESLineSearchGetOrder(SNESLineSearch, PetscInt *);
887 PETSC_EXTERN PetscErrorCode SNESLineSearchSetOrder(SNESLineSearch, PetscInt);
888
889 /*E
890 SNESLineSearchReason - indication if the line search has succeeded or failed and why
891
892 Values:
893 + `SNES_LINESEARCH_SUCCEEDED` - the line search succeeded
894 . `SNES_LINESEARCH_FAILED_NANORINF` - a not a number of infinity appeared in the computions
895 . `SNES_LINESEARCH_FAILED_FUNCTION_DOMAIN` - the function was evaluated outside of its domain, see `SNESSetFunctionDomainError()`
896 . `SNES_LINESEARCH_FAILED_OBJECTIVE_DOMAIN`- the objective function was evaluated outside of its domain, see `SNESSetObjectiveDomainError()`
897 . `SNES_LINESEARCH_FAILED_JACOBIAN_DOMAIN` - the Jacobian was evaluated outside of its domain, see `SNESSetJacobianDomainError()`
898 . `SNES_LINESEARCH_FAILED_REDUCT` - the linear search failed to get the requested decrease in its norm or objective
899 . `SNES_LINESEARCH_FAILED_USER` - used by `SNESLINESEARCHNLEQERR` to indicate the user changed the search direction inappropriately
900 - `SNES_LINESEARCH_FAILED_FUNCTION` - indicates the maximum number of function evaluations allowed has been surpassed, `SNESConvergedReason` is also
901 set to `SNES_DIVERGED_FUNCTION_COUNT`
902
903 Level: intermediate
904
905 Developer Note:
906 Some of these reasons overlap with values of `SNESConvergedReason`. It is possibly a better design to have `SNESConvergedReaon` alone used also for indicating line
907 search failures.
908
909 .seealso: [](ch_snes), `SNES`, `SNESSolve()`, `SNESGetConvergedReason()`, `KSPConvergedReason`, `SNESSetConvergenceTest()`,
910 `SNESSetFunctionDomainError()` and `SNESSetJacobianDomainError()`
911 E*/
912 typedef enum {
913 SNES_LINESEARCH_SUCCEEDED,
914 SNES_LINESEARCH_FAILED_NANORINF,
915 SNES_LINESEARCH_FAILED_FUNCTION_DOMAIN,
916 SNES_LINESEARCH_FAILED_OBJECTIVE_DOMAIN,
917 SNES_LINESEARCH_FAILED_JACOBIAN_DOMAIN,
918 SNES_LINESEARCH_FAILED_REDUCT, /* INSUFFICIENT REDUCTION */
919 SNES_LINESEARCH_FAILED_USER,
920 SNES_LINESEARCH_FAILED_FUNCTION
921 } SNESLineSearchReason;
922
923 PETSC_EXTERN PetscErrorCode SNESLineSearchGetReason(SNESLineSearch, SNESLineSearchReason *);
924 PETSC_EXTERN PetscErrorCode SNESLineSearchSetReason(SNESLineSearch, SNESLineSearchReason);
925
926 PETSC_EXTERN PetscErrorCode SNESLineSearchGetVecs(SNESLineSearch, Vec *, Vec *, Vec *, Vec *, Vec *);
927 PETSC_EXTERN PetscErrorCode SNESLineSearchSetVecs(SNESLineSearch, Vec, Vec, Vec, Vec, Vec);
928
929 PETSC_EXTERN PetscErrorCode SNESLineSearchGetNorms(SNESLineSearch, PetscReal *, PetscReal *, PetscReal *);
930 PETSC_EXTERN PetscErrorCode SNESLineSearchSetNorms(SNESLineSearch, PetscReal, PetscReal, PetscReal);
931 PETSC_EXTERN PetscErrorCode SNESLineSearchComputeNorms(SNESLineSearch);
932 PETSC_EXTERN PetscErrorCode SNESLineSearchSetComputeNorms(SNESLineSearch, PetscBool);
933
934 PETSC_EXTERN PetscErrorCode SNESLineSearchMonitor(SNESLineSearch);
935 PETSC_EXTERN PetscErrorCode SNESLineSearchMonitorSet(SNESLineSearch, PetscErrorCode (*)(SNESLineSearch, PetscCtx), PetscCtx, PetscCtxDestroyFn *);
936 PETSC_EXTERN PetscErrorCode SNESLineSearchMonitorSetFromOptions(SNESLineSearch, const char[], const char[], const char[], PetscErrorCode (*)(SNESLineSearch, PetscViewerAndFormat *), PetscErrorCode (*)(SNESLineSearch, PetscViewerAndFormat *));
937 PETSC_EXTERN PetscErrorCode SNESLineSearchMonitorCancel(SNESLineSearch);
938 PETSC_EXTERN PetscErrorCode SNESLineSearchSetDefaultMonitor(SNESLineSearch, PetscViewer);
939 PETSC_EXTERN PetscErrorCode SNESLineSearchGetDefaultMonitor(SNESLineSearch, PetscViewer *);
940 PETSC_EXTERN PetscErrorCode SNESLineSearchMonitorSolutionUpdate(SNESLineSearch, PetscViewerAndFormat *);
941
942 PETSC_EXTERN PetscErrorCode SNESLineSearchAppendOptionsPrefix(SNESLineSearch, const char[]);
943 PETSC_EXTERN PetscErrorCode SNESLineSearchGetOptionsPrefix(SNESLineSearch, const char *[]);
944
945 /* Shell interface functions */
946 PETSC_EXTERN PetscErrorCode SNESLineSearchShellSetApply(SNESLineSearch, SNESLineSearchShellApplyFn *, PetscCtx);
947 PETSC_EXTERN PetscErrorCode SNESLineSearchShellGetApply(SNESLineSearch, SNESLineSearchShellApplyFn **, PetscCtxRt);
948
SNESLineSearchShellSetUserFunc(SNESLineSearch ls,SNESLineSearchShellApplyFn * f,PetscCtx ctx)949 PETSC_DEPRECATED_FUNCTION(3, 21, 0, "SNESLinesearchShellSetApply()", ) static inline PetscErrorCode SNESLineSearchShellSetUserFunc(SNESLineSearch ls, SNESLineSearchShellApplyFn *f, PetscCtx ctx)
950 {
951 return SNESLineSearchShellSetApply(ls, f, ctx);
952 }
953
SNESLineSearchShellGetUserFunc(SNESLineSearch ls,SNESLineSearchShellApplyFn ** f,PetscCtxRt ctx)954 PETSC_DEPRECATED_FUNCTION(3, 21, 0, "SNESLinesearchShellGetApply()", ) static inline PetscErrorCode SNESLineSearchShellGetUserFunc(SNESLineSearch ls, SNESLineSearchShellApplyFn **f, PetscCtxRt ctx)
955 {
956 return SNESLineSearchShellGetApply(ls, f, ctx);
957 }
958
959 /* BT interface functions */
960 PETSC_EXTERN PetscErrorCode SNESLineSearchBTSetAlpha(SNESLineSearch, PetscReal);
961 PETSC_EXTERN PetscErrorCode SNESLineSearchBTGetAlpha(SNESLineSearch, PetscReal *);
962
963 /*register line search types */
964 PETSC_EXTERN PetscErrorCode SNESLineSearchRegister(const char[], PetscErrorCode (*)(SNESLineSearch));
965
966 /* Routines for VI solver */
967 PETSC_EXTERN PetscErrorCode SNESVISetVariableBounds(SNES, Vec, Vec);
968 PETSC_EXTERN PetscErrorCode SNESVIGetVariableBounds(SNES, Vec *, Vec *);
969 PETSC_EXTERN PetscErrorCode SNESVISetComputeVariableBounds(SNES, PetscErrorCode (*)(SNES, Vec, Vec));
970 PETSC_EXTERN PetscErrorCode SNESVIGetInactiveSet(SNES, IS *);
971 PETSC_EXTERN PetscErrorCode SNESVIGetActiveSetIS(SNES, Vec, Vec, IS *);
972 PETSC_EXTERN PetscErrorCode SNESVIComputeInactiveSetFnorm(SNES, Vec, Vec, PetscReal *);
973 PETSC_EXTERN PetscErrorCode SNESVIComputeInactiveSetFtY(SNES, Vec, Vec, Vec, PetscScalar *);
974 PETSC_EXTERN PetscErrorCode SNESVISetRedundancyCheck(SNES, PetscErrorCode (*)(SNES, IS, IS *, PetscCtx), PetscCtx);
975 PETSC_EXTERN PetscErrorCode SNESVIComputeMeritFunction(Vec, PetscReal *, PetscReal *);
976 PETSC_EXTERN PetscErrorCode SNESVIComputeFunction(SNES, Vec, Vec, void *);
977 PETSC_EXTERN PetscErrorCode DMSetVI(DM, IS);
978 PETSC_EXTERN PetscErrorCode DMDestroyVI(DM);
979
980 PETSC_EXTERN PetscErrorCode SNESTestLocalMin(SNES);
981
982 /* Should this routine be private? */
983 PETSC_EXTERN PetscErrorCode SNESComputeJacobian(SNES, Vec, Mat, Mat);
984 PETSC_EXTERN PetscErrorCode SNESTestJacobian(SNES, PetscReal *, PetscReal *);
985 PETSC_EXTERN PetscErrorCode SNESTestFunction(SNES);
986
987 PETSC_EXTERN PetscErrorCode SNESSetDM(SNES, DM);
988 PETSC_EXTERN PetscErrorCode SNESGetDM(SNES, DM *);
989 PETSC_EXTERN PetscErrorCode SNESSetNPC(SNES, SNES);
990 PETSC_EXTERN PetscErrorCode SNESGetNPC(SNES, SNES *);
991 PETSC_EXTERN PetscErrorCode SNESHasNPC(SNES, PetscBool *);
992 PETSC_EXTERN PetscErrorCode SNESApplyNPC(SNES, Vec, Vec, Vec);
993 PETSC_EXTERN PetscErrorCode SNESGetNPCFunction(SNES, Vec, PetscReal *);
994 PETSC_EXTERN PetscErrorCode SNESComputeFunctionDefaultNPC(SNES, Vec, Vec);
995 PETSC_EXTERN PetscErrorCode SNESSetNPCSide(SNES, PCSide);
996 PETSC_EXTERN PetscErrorCode SNESGetNPCSide(SNES, PCSide *);
997 PETSC_EXTERN PetscErrorCode SNESSetLineSearch(SNES, SNESLineSearch);
998 PETSC_EXTERN PetscErrorCode SNESGetLineSearch(SNES, SNESLineSearch *);
999
SNESGetSNESLineSearch(SNES snes,SNESLineSearch * ls)1000 PETSC_DEPRECATED_FUNCTION(3, 4, 0, "SNESGetLineSearch()", ) static inline PetscErrorCode SNESGetSNESLineSearch(SNES snes, SNESLineSearch *ls)
1001 {
1002 return SNESGetLineSearch(snes, ls);
1003 }
SNESSetSNESLineSearch(SNES snes,SNESLineSearch ls)1004 PETSC_DEPRECATED_FUNCTION(3, 4, 0, "SNESSetLineSearch()", ) static inline PetscErrorCode SNESSetSNESLineSearch(SNES snes, SNESLineSearch ls)
1005 {
1006 return SNESSetLineSearch(snes, ls);
1007 }
1008
1009 PETSC_EXTERN PetscErrorCode SNESSetUpMatrices(SNES);
1010 PETSC_EXTERN PetscErrorCode DMSNESSetFunction(DM, SNESFunctionFn *, PetscCtx);
1011 PETSC_EXTERN PetscErrorCode DMSNESGetFunction(DM, SNESFunctionFn **, PetscCtxRt);
1012 PETSC_EXTERN PetscErrorCode DMSNESSetFunctionContextDestroy(DM, PetscCtxDestroyFn *);
1013 PETSC_EXTERN PetscErrorCode DMSNESSetMFFunction(DM, SNESFunctionFn *, PetscCtx);
1014 PETSC_EXTERN PetscErrorCode DMSNESSetNGS(DM, SNESNGSFn *, PetscCtx);
1015 PETSC_EXTERN PetscErrorCode DMSNESGetNGS(DM, SNESNGSFn **, PetscCtxRt);
1016 PETSC_EXTERN PetscErrorCode DMSNESSetJacobian(DM, SNESJacobianFn *, PetscCtx);
1017 PETSC_EXTERN PetscErrorCode DMSNESGetJacobian(DM, SNESJacobianFn **, PetscCtxRt);
1018 PETSC_EXTERN PetscErrorCode DMSNESSetJacobianContextDestroy(DM, PetscCtxDestroyFn *);
1019 PETSC_EXTERN PetscErrorCode DMSNESSetPicard(DM, SNESFunctionFn *, SNESJacobianFn *, PetscCtx);
1020 PETSC_EXTERN PetscErrorCode DMSNESGetPicard(DM, SNESFunctionFn **, SNESJacobianFn **, PetscCtxRt);
1021 PETSC_EXTERN PetscErrorCode DMSNESSetObjective(DM, SNESObjectiveFn *, PetscCtx);
1022 PETSC_EXTERN PetscErrorCode DMSNESGetObjective(DM, SNESObjectiveFn **, PetscCtxRt);
1023 PETSC_EXTERN PetscErrorCode DMCopyDMSNES(DM, DM);
1024
1025 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode DMDASNESFunctionFn(DMDALocalInfo *, void *, void *, void *);
1026 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode DMDASNESJacobianFn(DMDALocalInfo *, void *, Mat, Mat, void *);
1027 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode DMDASNESObjectiveFn(DMDALocalInfo *, void *, PetscReal *, void *);
1028
1029 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode DMDASNESFunctionVecFn(DMDALocalInfo *, Vec, Vec, void *);
1030 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode DMDASNESJacobianVecFn(DMDALocalInfo *, Vec, Mat, Mat, void *);
1031 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode DMDASNESObjectiveVecFn(DMDALocalInfo *, Vec, PetscReal *, void *);
1032
1033 PETSC_EXTERN PetscErrorCode DMDASNESSetFunctionLocal(DM, InsertMode, DMDASNESFunctionFn *, void *);
1034 PETSC_EXTERN PetscErrorCode DMDASNESSetJacobianLocal(DM, DMDASNESJacobianFn *, void *);
1035 PETSC_EXTERN PetscErrorCode DMDASNESSetObjectiveLocal(DM, DMDASNESObjectiveFn *, void *);
1036 PETSC_EXTERN PetscErrorCode DMDASNESSetPicardLocal(DM, InsertMode, DMDASNESFunctionFn *, DMDASNESJacobianFn, void *);
1037
1038 PETSC_EXTERN PetscErrorCode DMDASNESSetFunctionLocalVec(DM, InsertMode, DMDASNESFunctionVecFn *, void *);
1039 PETSC_EXTERN PetscErrorCode DMDASNESSetJacobianLocalVec(DM, DMDASNESJacobianVecFn *, void *);
1040 PETSC_EXTERN PetscErrorCode DMDASNESSetObjectiveLocalVec(DM, DMDASNESObjectiveVecFn *, void *);
1041
1042 PETSC_EXTERN PetscErrorCode DMSNESSetBoundaryLocal(DM, PetscErrorCode (*)(DM, Vec, PetscCtx), PetscCtx);
1043 PETSC_EXTERN PetscErrorCode DMSNESSetObjectiveLocal(DM, PetscErrorCode (*)(DM, Vec, PetscReal *, PetscCtx), PetscCtx);
1044 PETSC_EXTERN PetscErrorCode DMSNESSetFunctionLocal(DM, PetscErrorCode (*)(DM, Vec, Vec, PetscCtx), PetscCtx);
1045 PETSC_EXTERN PetscErrorCode DMSNESSetJacobianLocal(DM, PetscErrorCode (*)(DM, Vec, Mat, Mat, PetscCtx), PetscCtx);
1046 PETSC_EXTERN PetscErrorCode DMSNESGetBoundaryLocal(DM, PetscErrorCode (**)(DM, Vec, PetscCtx), PetscCtxRt);
1047 PETSC_EXTERN PetscErrorCode DMSNESGetObjectiveLocal(DM, PetscErrorCode (**)(DM, Vec, PetscReal *, PetscCtx), PetscCtxRt);
1048 PETSC_EXTERN PetscErrorCode DMSNESGetFunctionLocal(DM, PetscErrorCode (**)(DM, Vec, Vec, PetscCtx), PetscCtxRt);
1049 PETSC_EXTERN PetscErrorCode DMSNESGetJacobianLocal(DM, PetscErrorCode (**)(DM, Vec, Mat, Mat, PetscCtx), PetscCtxRt);
1050
1051 /* Routines for Multiblock solver */
1052 PETSC_EXTERN PetscErrorCode SNESMultiblockSetFields(SNES, const char[], PetscInt, const PetscInt *);
1053 PETSC_EXTERN PetscErrorCode SNESMultiblockSetIS(SNES, const char[], IS);
1054 PETSC_EXTERN PetscErrorCode SNESMultiblockSetBlockSize(SNES, PetscInt);
1055 PETSC_EXTERN PetscErrorCode SNESMultiblockSetType(SNES, PCCompositeType);
1056 PETSC_EXTERN PetscErrorCode SNESMultiblockGetSubSNES(SNES, PetscInt *, SNES *[]);
1057
1058 /*J
1059 SNESMSType - String with the name of a PETSc `SNESMS` method.
1060
1061 Level: intermediate
1062
1063 .seealso: [](ch_snes), `SNESMS`, `SNESMSGetType()`, `SNESMSSetType()`, `SNES`
1064 J*/
1065 typedef const char *SNESMSType;
1066 #define SNESMSM62 "m62"
1067 #define SNESMSEULER "euler"
1068 #define SNESMSJAMESON83 "jameson83"
1069 #define SNESMSVLTP11 "vltp11"
1070 #define SNESMSVLTP21 "vltp21"
1071 #define SNESMSVLTP31 "vltp31"
1072 #define SNESMSVLTP41 "vltp41"
1073 #define SNESMSVLTP51 "vltp51"
1074 #define SNESMSVLTP61 "vltp61"
1075
1076 PETSC_EXTERN PetscErrorCode SNESMSRegister(SNESMSType, PetscInt, PetscInt, PetscReal, const PetscReal[], const PetscReal[], const PetscReal[]);
1077 PETSC_EXTERN PetscErrorCode SNESMSRegisterAll(void);
1078 PETSC_EXTERN PetscErrorCode SNESMSGetType(SNES, SNESMSType *);
1079 PETSC_EXTERN PetscErrorCode SNESMSSetType(SNES, SNESMSType);
1080 PETSC_EXTERN PetscErrorCode SNESMSGetDamping(SNES, PetscReal *);
1081 PETSC_EXTERN PetscErrorCode SNESMSSetDamping(SNES, PetscReal);
1082 PETSC_EXTERN PetscErrorCode SNESMSFinalizePackage(void);
1083 PETSC_EXTERN PetscErrorCode SNESMSInitializePackage(void);
1084 PETSC_EXTERN PetscErrorCode SNESMSRegisterDestroy(void);
1085
1086 /*MC
1087 SNESNGMRESRestartType - the restart approach used by `SNESNGMRES`
1088
1089 Values:
1090 + `SNES_NGMRES_RESTART_NONE` - never restart
1091 . `SNES_NGMRES_RESTART_DIFFERENCE` - restart based upon difference criteria
1092 - `SNES_NGMRES_RESTART_PERIODIC` - restart after a fixed number of iterations
1093
1094 Options Database Keys:
1095 + -snes_ngmres_restart_type <difference,periodic,none> - set the restart type
1096 - -snes_ngmres_restart <30> - sets the number of iterations before restart for periodic
1097
1098 Level: intermediate
1099
1100 .seealso: `SNES, `SNESNGMRES`, `SNESNGMRESSetSelectType()`, `SNESNGMRESGetSelectType()`, `SNESNGMRESSetRestartType()`,
1101 `SNESNGMRESGetRestartType()`, `SNESNGMRESSelectType`
1102 M*/
1103 typedef enum {
1104 SNES_NGMRES_RESTART_NONE = 0,
1105 SNES_NGMRES_RESTART_PERIODIC = 1,
1106 SNES_NGMRES_RESTART_DIFFERENCE = 2
1107 } SNESNGMRESRestartType;
1108 PETSC_EXTERN const char *const SNESNGMRESRestartTypes[];
1109
1110 /*MC
1111 SNESNGMRESSelectType - the approach used by `SNESNGMRES` to determine how the candidate solution and
1112 combined solution are used to create the next iterate.
1113
1114 Values:
1115 + `SNES_NGMRES_SELECT_NONE` - choose the combined solution all the time
1116 . `SNES_NGMRES_SELECT_DIFFERENCE` - choose based upon the selection criteria
1117 - `SNES_NGMRES_SELECT_LINESEARCH` - choose based upon line search combination
1118
1119 Options Database Key:
1120 . -snes_ngmres_select_type<difference,none,linesearch> - select type
1121
1122 Level: intermediate
1123
1124 .seealso: `SNES, `SNESNGMRES`, `SNESNGMRESSetSelectType()`, `SNESNGMRESGetSelectType()`, `SNESNGMRESSetRestartType()`,
1125 `SNESNGMRESGetRestartType()`, `SNESNGMRESRestartType`
1126 M*/
1127 typedef enum {
1128 SNES_NGMRES_SELECT_NONE = 0,
1129 SNES_NGMRES_SELECT_DIFFERENCE = 1,
1130 SNES_NGMRES_SELECT_LINESEARCH = 2
1131 } SNESNGMRESSelectType;
1132 PETSC_EXTERN const char *const SNESNGMRESSelectTypes[];
1133
1134 PETSC_EXTERN PetscErrorCode SNESNGMRESSetRestartType(SNES, SNESNGMRESRestartType);
1135 PETSC_EXTERN PetscErrorCode SNESNGMRESSetSelectType(SNES, SNESNGMRESSelectType);
1136 PETSC_EXTERN PetscErrorCode SNESNGMRESSetRestartFmRise(SNES, PetscBool);
1137 PETSC_EXTERN PetscErrorCode SNESNGMRESGetRestartFmRise(SNES, PetscBool *);
1138
1139 /*MC
1140 SNESNCGType - the conjugate update approach for `SNESNCG`
1141
1142 Values:
1143 + `SNES_NCG_FR` - Fletcher-Reeves update
1144 . `SNES_NCG_PRP` - Polak-Ribiere-Polyak update, the default and the only one that tolerates generalized search directions
1145 . `SNES_NCG_HS` - Hestenes-Steifel update
1146 . `SNES_NCG_DY` - Dai-Yuan update
1147 - `SNES_NCG_CD` - Conjugate Descent update
1148
1149 Options Database Key:
1150 . -snes_ncg_type<fr,prp,hs,dy,cd> - select type
1151
1152 Level: intermediate
1153
1154 .seealso: `SNES, `SNESNCG`, `SNESNCGSetType()`
1155 M*/
1156 typedef enum {
1157 SNES_NCG_FR = 0,
1158 SNES_NCG_PRP = 1,
1159 SNES_NCG_HS = 2,
1160 SNES_NCG_DY = 3,
1161 SNES_NCG_CD = 4
1162 } SNESNCGType;
1163 PETSC_EXTERN const char *const SNESNCGTypes[];
1164
1165 PETSC_EXTERN PetscErrorCode SNESNCGSetType(SNES, SNESNCGType);
1166
1167 /*MC
1168 SNESQNScaleType - the scaling type used by `SNESQN`
1169
1170 Values:
1171 + `SNES_QN_SCALE_NONE` - don't scale the problem
1172 . `SNES_QN_SCALE_SCALAR` - use Shanno scaling
1173 . `SNES_QN_SCALE_DIAGONAL` - scale with a diagonalized BFGS formula (see Gilbert and Lemarechal 1989), available
1174 - `SNES_QN_SCALE_JACOBIAN` - scale by solving a linear system coming from the Jacobian you provided with `SNESSetJacobian()`
1175 computed at the first iteration of `SNESQN` and at ever restart.
1176
1177 Options Database Key:
1178 . -snes_qn_scale_type <diagonal,none,scalar,jacobian> - Scaling type
1179
1180 Level: intermediate
1181
1182 .seealso: `SNES, `SNESQN`, `SNESQNSetScaleType()`, `SNESQNType`, `SNESQNSetType()`, `SNESQNSetRestartType()`, `SNESQNRestartType`
1183 M*/
1184 typedef enum {
1185 SNES_QN_SCALE_DEFAULT = 0,
1186 SNES_QN_SCALE_NONE = 1,
1187 SNES_QN_SCALE_SCALAR = 2,
1188 SNES_QN_SCALE_DIAGONAL = 3,
1189 SNES_QN_SCALE_JACOBIAN = 4
1190 } SNESQNScaleType;
1191 PETSC_EXTERN const char *const SNESQNScaleTypes[];
1192
1193 /*MC
1194 SNESQNRestartType - the restart approached used by `SNESQN`
1195
1196 Values:
1197 + `SNES_QN_RESTART_NONE` - never restart
1198 . `SNES_QN_RESTART_POWELL` - restart based upon descent criteria
1199 - `SNES_QN_RESTART_PERIODIC` - restart after a fixed number of iterations
1200
1201 Options Database Keys:
1202 + -snes_qn_restart_type <powell,periodic,none> - set the restart type
1203 - -snes_qn_m <m> - sets the number of stored updates and the restart period for periodic
1204
1205 Level: intermediate
1206
1207 .seealso: `SNES, `SNESQN`, `SNESQNSetScaleType()`, `SNESQNType`, `SNESQNSetType()`, `SNESQNSetRestartType()`, `SNESQNScaleType`
1208 M*/
1209 typedef enum {
1210 SNES_QN_RESTART_DEFAULT = 0,
1211 SNES_QN_RESTART_NONE = 1,
1212 SNES_QN_RESTART_POWELL = 2,
1213 SNES_QN_RESTART_PERIODIC = 3
1214 } SNESQNRestartType;
1215 PETSC_EXTERN const char *const SNESQNRestartTypes[];
1216
1217 /*MC
1218 SNESQNType - the type used by `SNESQN`
1219
1220 Values:
1221 + `SNES_QN_LBFGS` - LBFGS variant
1222 . `SNES_QN_BROYDEN` - Broyden variant
1223 - `SNES_QN_BADBROYDEN` - Bad Broyden variant
1224
1225 Options Database Key:
1226 . -snes_qn_type <lbfgs,broyden,badbroyden> - quasi-Newton type
1227
1228 Level: intermediate
1229
1230 .seealso: `SNES, `SNESQN`, `SNESQNSetScaleType()`, `SNESQNSetType()`, `SNESQNScaleType`, `SNESQNRestartType`, `SNESQNSetRestartType()`
1231 M*/
1232 typedef enum {
1233 SNES_QN_LBFGS = 0,
1234 SNES_QN_BROYDEN = 1,
1235 SNES_QN_BADBROYDEN = 2
1236 } SNESQNType;
1237 PETSC_EXTERN const char *const SNESQNTypes[];
1238
1239 PETSC_EXTERN PetscErrorCode SNESQNSetType(SNES, SNESQNType);
1240 PETSC_EXTERN PetscErrorCode SNESQNSetScaleType(SNES, SNESQNScaleType);
1241 PETSC_EXTERN PetscErrorCode SNESQNSetRestartType(SNES, SNESQNRestartType);
1242
1243 PETSC_EXTERN PetscErrorCode SNESNASMGetType(SNES, PCASMType *);
1244 PETSC_EXTERN PetscErrorCode SNESNASMSetType(SNES, PCASMType);
1245 PETSC_EXTERN PetscErrorCode SNESNASMGetSubdomains(SNES, PetscInt *, SNES *[], VecScatter *[], VecScatter *[], VecScatter *[]);
1246 PETSC_EXTERN PetscErrorCode SNESNASMSetSubdomains(SNES, PetscInt, SNES[], VecScatter[], VecScatter[], VecScatter[]);
1247 PETSC_EXTERN PetscErrorCode SNESNASMSetDamping(SNES, PetscReal);
1248 PETSC_EXTERN PetscErrorCode SNESNASMGetDamping(SNES, PetscReal *);
1249 PETSC_EXTERN PetscErrorCode SNESNASMGetSubdomainVecs(SNES, PetscInt *, Vec *[], Vec *[], Vec *[], Vec *[]);
1250 PETSC_EXTERN PetscErrorCode SNESNASMSetComputeFinalJacobian(SNES, PetscBool);
1251 PETSC_EXTERN PetscErrorCode SNESNASMGetSNES(SNES, PetscInt, SNES *);
1252 PETSC_EXTERN PetscErrorCode SNESNASMGetNumber(SNES, PetscInt *);
1253 PETSC_EXTERN PetscErrorCode SNESNASMSetWeight(SNES, Vec);
1254
1255 /*E
1256 SNESCompositeType - Determines how two or more preconditioners are composed with the `SNESType` of `SNESCOMPOSITE`
1257
1258 Values:
1259 + `SNES_COMPOSITE_ADDITIVE` - results from application of all preconditioners are added together
1260 . `SNES_COMPOSITE_MULTIPLICATIVE` - preconditioners are applied sequentially to the residual freshly
1261 computed after the previous preconditioner application
1262 - `SNES_COMPOSITE_ADDITIVEOPTIMAL` - uses a linear combination of the solutions obtained with each preconditioner that approximately minimize the function
1263 value at the new iteration.
1264
1265 Level: beginner
1266
1267 .seealso: [](sec_pc), `PCCOMPOSITE`, `PCFIELDSPLIT`, `PC`, `PCCompositeSetType()`, `PCCompositeType`
1268 E*/
1269 typedef enum {
1270 SNES_COMPOSITE_ADDITIVE,
1271 SNES_COMPOSITE_MULTIPLICATIVE,
1272 SNES_COMPOSITE_ADDITIVEOPTIMAL
1273 } SNESCompositeType;
1274 PETSC_EXTERN const char *const SNESCompositeTypes[];
1275
1276 PETSC_EXTERN PetscErrorCode SNESCompositeSetType(SNES, SNESCompositeType);
1277 PETSC_EXTERN PetscErrorCode SNESCompositeAddSNES(SNES, SNESType);
1278 PETSC_EXTERN PetscErrorCode SNESCompositeGetSNES(SNES, PetscInt, SNES *);
1279 PETSC_EXTERN PetscErrorCode SNESCompositeGetNumber(SNES, PetscInt *);
1280 PETSC_EXTERN PetscErrorCode SNESCompositeSetDamping(SNES, PetscInt, PetscReal);
1281
1282 PETSC_EXTERN PetscErrorCode SNESPatchSetDiscretisationInfo(SNES, PetscInt, DM *, PetscInt *, PetscInt *, const PetscInt **, const PetscInt *, PetscInt, const PetscInt *, PetscInt, const PetscInt *);
1283 PETSC_EXTERN PetscErrorCode SNESPatchSetComputeOperator(SNES, PetscErrorCode (*func)(PC, PetscInt, Vec, Mat, IS, PetscInt, const PetscInt *, const PetscInt *, PetscCtx), PetscCtx);
1284 PETSC_EXTERN PetscErrorCode SNESPatchSetComputeFunction(SNES, PetscErrorCode (*func)(PC, PetscInt, Vec, Vec, IS, PetscInt, const PetscInt *, const PetscInt *, PetscCtx), PetscCtx);
1285 PETSC_EXTERN PetscErrorCode SNESPatchSetConstructType(SNES, PCPatchConstructType, PetscErrorCode (*func)(PC, PetscInt *, IS **, IS *, PetscCtx), PetscCtx);
1286 PETSC_EXTERN PetscErrorCode SNESPatchSetCellNumbering(SNES, PetscSection);
1287
1288 /*E
1289 SNESFASType - Determines the type of nonlinear multigrid method that is run.
1290
1291 Values:
1292 + `SNES_FAS_MULTIPLICATIVE` (default) - traditional V or W cycle as determined by `SNESFASSetCycles()`
1293 . `SNES_FAS_ADDITIVE` - additive FAS cycle
1294 . `SNES_FAS_FULL` - full FAS cycle
1295 - `SNES_FAS_KASKADE` - Kaskade FAS cycle
1296
1297 Level: beginner
1298
1299 .seealso: [](ch_snes), `SNESFAS`, `PCMGSetType()`, `PCMGType`
1300 E*/
1301 typedef enum {
1302 SNES_FAS_MULTIPLICATIVE,
1303 SNES_FAS_ADDITIVE,
1304 SNES_FAS_FULL,
1305 SNES_FAS_KASKADE
1306 } SNESFASType;
1307 PETSC_EXTERN const char *const SNESFASTypes[];
1308
1309 /* called on the finest level FAS instance*/
1310 PETSC_EXTERN PetscErrorCode SNESFASSetType(SNES, SNESFASType);
1311 PETSC_EXTERN PetscErrorCode SNESFASGetType(SNES, SNESFASType *);
1312 PETSC_EXTERN PetscErrorCode SNESFASSetLevels(SNES, PetscInt, MPI_Comm *);
1313 PETSC_EXTERN PetscErrorCode SNESFASGetLevels(SNES, PetscInt *);
1314 PETSC_EXTERN PetscErrorCode SNESFASGetCycleSNES(SNES, PetscInt, SNES *);
1315 PETSC_EXTERN PetscErrorCode SNESFASSetNumberSmoothUp(SNES, PetscInt);
1316 PETSC_EXTERN PetscErrorCode SNESFASSetNumberSmoothDown(SNES, PetscInt);
1317 PETSC_EXTERN PetscErrorCode SNESFASSetCycles(SNES, PetscInt);
1318 PETSC_EXTERN PetscErrorCode SNESFASSetMonitor(SNES, PetscViewerAndFormat *, PetscBool);
1319 PETSC_EXTERN PetscErrorCode SNESFASSetLog(SNES, PetscBool);
1320
1321 PETSC_EXTERN PetscErrorCode SNESFASSetGalerkin(SNES, PetscBool);
1322 PETSC_EXTERN PetscErrorCode SNESFASGetGalerkin(SNES, PetscBool *);
1323 PETSC_EXTERN PetscErrorCode SNESFASGalerkinFunctionDefault(SNES, Vec, Vec, PetscCtx);
1324
1325 /* called on any level -- "Cycle" FAS instance */
1326 PETSC_EXTERN PetscErrorCode SNESFASCycleGetSmoother(SNES, SNES *);
1327 PETSC_EXTERN PetscErrorCode SNESFASCycleGetSmootherUp(SNES, SNES *);
1328 PETSC_EXTERN PetscErrorCode SNESFASCycleGetSmootherDown(SNES, SNES *);
1329 PETSC_EXTERN PetscErrorCode SNESFASCycleGetCorrection(SNES, SNES *);
1330 PETSC_EXTERN PetscErrorCode SNESFASCycleGetInterpolation(SNES, Mat *);
1331 PETSC_EXTERN PetscErrorCode SNESFASCycleGetRestriction(SNES, Mat *);
1332 PETSC_EXTERN PetscErrorCode SNESFASCycleGetInjection(SNES, Mat *);
1333 PETSC_EXTERN PetscErrorCode SNESFASCycleGetRScale(SNES, Vec *);
1334 PETSC_EXTERN PetscErrorCode SNESFASCycleSetCycles(SNES, PetscInt);
1335 PETSC_EXTERN PetscErrorCode SNESFASCycleIsFine(SNES, PetscBool *);
1336
1337 /* called on the (outer) finest level FAS to set/get parameters on any level instance */
1338 PETSC_EXTERN PetscErrorCode SNESFASSetInterpolation(SNES, PetscInt, Mat);
1339 PETSC_EXTERN PetscErrorCode SNESFASGetInterpolation(SNES, PetscInt, Mat *);
1340 PETSC_EXTERN PetscErrorCode SNESFASSetRestriction(SNES, PetscInt, Mat);
1341 PETSC_EXTERN PetscErrorCode SNESFASGetRestriction(SNES, PetscInt, Mat *);
1342 PETSC_EXTERN PetscErrorCode SNESFASSetInjection(SNES, PetscInt, Mat);
1343 PETSC_EXTERN PetscErrorCode SNESFASGetInjection(SNES, PetscInt, Mat *);
1344 PETSC_EXTERN PetscErrorCode SNESFASSetRScale(SNES, PetscInt, Vec);
1345 PETSC_EXTERN PetscErrorCode SNESFASGetRScale(SNES, PetscInt, Vec *);
1346 PETSC_EXTERN PetscErrorCode SNESFASSetContinuation(SNES, PetscBool);
1347
1348 PETSC_EXTERN PetscErrorCode SNESFASGetSmoother(SNES, PetscInt, SNES *);
1349 PETSC_EXTERN PetscErrorCode SNESFASGetSmootherUp(SNES, PetscInt, SNES *);
1350 PETSC_EXTERN PetscErrorCode SNESFASGetSmootherDown(SNES, PetscInt, SNES *);
1351 PETSC_EXTERN PetscErrorCode SNESFASGetCoarseSolve(SNES, SNES *);
1352
1353 /* parameters for full FAS */
1354 PETSC_EXTERN PetscErrorCode SNESFASFullSetDownSweep(SNES, PetscBool);
1355 PETSC_EXTERN PetscErrorCode SNESFASCreateCoarseVec(SNES, Vec *);
1356 PETSC_EXTERN PetscErrorCode SNESFASRestrict(SNES, Vec, Vec);
1357 PETSC_EXTERN PetscErrorCode SNESFASFullSetTotal(SNES, PetscBool);
1358 PETSC_EXTERN PetscErrorCode SNESFASFullGetTotal(SNES, PetscBool *);
1359
1360 PETSC_EXTERN PetscErrorCode DMPlexSetSNESVariableBounds(DM, SNES);
1361 PETSC_EXTERN PetscErrorCode DMSNESCheckDiscretization(SNES, DM, PetscReal, Vec, PetscReal, PetscReal[]);
1362 PETSC_EXTERN PetscErrorCode DMSNESCheckResidual(SNES, DM, Vec, PetscReal, PetscReal *);
1363 PETSC_EXTERN PetscErrorCode DMSNESCheckJacobian(SNES, DM, Vec, PetscReal, PetscBool *, PetscReal *);
1364 PETSC_EXTERN PetscErrorCode DMSNESCheckFromOptions(SNES, Vec);
1365 PETSC_EXTERN PetscErrorCode DMSNESComputeJacobianAction(DM, Vec, Vec, Vec, PetscCtx);
1366 PETSC_EXTERN PetscErrorCode DMSNESCreateJacobianMF(DM, Vec, PetscCtx, Mat *);
1367
1368 PETSC_EXTERN PetscErrorCode SNESNewtonALSetFunction(SNES, SNESFunctionFn *, PetscCtx);
1369 PETSC_EXTERN PetscErrorCode SNESNewtonALGetFunction(SNES, SNESFunctionFn **, PetscCtxRt);
1370 PETSC_EXTERN PetscErrorCode SNESNewtonALComputeFunction(SNES, Vec, Vec);
1371 PETSC_EXTERN PetscErrorCode SNESNewtonALGetLoadParameter(SNES, PetscReal *);
1372
1373 /*MC
1374 SNESNewtonALCorrectionType - the approach used by `SNESNEWTONAL` to determine
1375 the correction to the current increment. While the exact correction satisfies
1376 the constraint surface at every iteration, it also requires solving a quadratic
1377 equation which may not have real roots. Conversely, the normal correction is more
1378 efficient and always yields a real correction and is the default.
1379
1380 Values:
1381 + `SNES_NEWTONAL_CORRECTION_EXACT` - choose the correction which exactly satisfies the constraint
1382 - `SNES_NEWTONAL_CORRECTION_NORMAL` - choose the correction in the updated normal hyper-surface to the constraint surface
1383
1384 Options Database Key:
1385 . -snes_newtonal_correction_type <exact> - select type from <exact,normal>
1386
1387 Level: intermediate
1388
1389 .seealso: `SNES`, `SNESNEWTONAL`, `SNESNewtonALSetCorrectionType()`
1390 M*/
1391 typedef enum {
1392 SNES_NEWTONAL_CORRECTION_EXACT = 0,
1393 SNES_NEWTONAL_CORRECTION_NORMAL = 1,
1394 } SNESNewtonALCorrectionType;
1395 PETSC_EXTERN const char *const SNESNewtonALCorrectionTypes[];
1396
1397 PETSC_EXTERN PetscErrorCode SNESNewtonALSetCorrectionType(SNES, SNESNewtonALCorrectionType);
1398