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