Actual source code: petscsnes.h

  1: /*
  2:     User interface for the nonlinear solvers package.
  3: */
 6:  #include petscksp.h
  7: PETSC_EXTERN_CXX_BEGIN

  9: /*S
 10:      SNES - Abstract PETSc object that manages all nonlinear solves

 12:    Level: beginner

 14:   Concepts: nonlinear solvers

 16: .seealso:  SNESCreate(), SNESSetType(), SNESType, TS, KSP, KSP, PC
 17: S*/
 18: typedef struct _p_SNES* SNES;

 20: /*J
 21:     SNESType - String with the name of a PETSc SNES method or the creation function
 22:        with an optional dynamic library name, for example
 23:        http://www.mcs.anl.gov/petsc/lib.a:mysnescreate()

 25:    Level: beginner

 27: .seealso: SNESSetType(), SNES
 28: J*/
 29: #define SNESType char*
 30: #define SNESLS      "ls"
 31: #define SNESTR      "tr"
 32: #define SNESPYTHON  "python"
 33: #define SNESTEST    "test"
 34: #define SNESPICARD  "picard"
 35: #define SNESKSPONLY "ksponly"
 36: #define SNESVI      "vi"
 37: #define SNESNGMRES  "ngmres"
 38: #define SNESSORQN   "sorqn"

 40: /* Logging support */
 41: extern PetscClassId  SNES_CLASSID;

 43: extern PetscErrorCode  SNESInitializePackage(const char[]);

 45: extern PetscErrorCode  SNESCreate(MPI_Comm,SNES*);
 46: extern PetscErrorCode  SNESReset(SNES);
 47: extern PetscErrorCode  SNESDestroy(SNES*);
 48: extern PetscErrorCode  SNESSetType(SNES,const SNESType);
 49: extern PetscErrorCode  SNESMonitor(SNES,PetscInt,PetscReal);
 50: extern PetscErrorCode  SNESMonitorSet(SNES,PetscErrorCode(*)(SNES,PetscInt,PetscReal,void*),void *,PetscErrorCode (*)(void**));
 51: extern PetscErrorCode  SNESMonitorCancel(SNES);
 52: extern PetscErrorCode  SNESSetConvergenceHistory(SNES,PetscReal[],PetscInt[],PetscInt,PetscBool );
 53: extern PetscErrorCode  SNESGetConvergenceHistory(SNES,PetscReal*[],PetscInt *[],PetscInt *);
 54: extern PetscErrorCode  SNESSetUp(SNES);
 55: extern PetscErrorCode  SNESSolve(SNES,Vec,Vec);
 56: extern PetscErrorCode  SNESSetErrorIfNotConverged(SNES,PetscBool );
 57: extern PetscErrorCode  SNESGetErrorIfNotConverged(SNES,PetscBool  *);


 60: extern PetscErrorCode  SNESAddOptionsChecker(PetscErrorCode (*)(SNES));

 62: extern PetscErrorCode  SNESSetUpdate(SNES, PetscErrorCode (*)(SNES, PetscInt));
 63: extern PetscErrorCode  SNESDefaultUpdate(SNES, PetscInt);

 65: extern PetscFList SNESList;
 66: extern PetscErrorCode  SNESRegisterDestroy(void);
 67: extern PetscErrorCode  SNESRegisterAll(const char[]);

 69: extern PetscErrorCode  SNESRegister(const char[],const char[],const char[],PetscErrorCode (*)(SNES));

 71: /*MC
 72:    SNESRegisterDynamic - Adds a method to the nonlinear solver package.

 74:    Synopsis:
 75:    PetscErrorCode SNESRegisterDynamic(const char *name_solver,const char *path,const char *name_create,PetscErrorCode (*routine_create)(SNES))

 77:    Not collective

 79:    Input Parameters:
 80: +  name_solver - name of a new user-defined solver
 81: .  path - path (either absolute or relative) the library containing this solver
 82: .  name_create - name of routine to create method context
 83: -  routine_create - routine to create method context

 85:    Notes:
 86:    SNESRegisterDynamic() may be called multiple times to add several user-defined solvers.

 88:    If dynamic libraries are used, then the fourth input argument (routine_create)
 89:    is ignored.

 91:    Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR},
 92:    and others of the form ${any_environmental_variable} occuring in pathname will be 
 93:    replaced with appropriate values.

 95:    Sample usage:
 96: .vb
 97:    SNESRegisterDynamic("my_solver",/home/username/my_lib/lib/libg/solaris/mylib.a,
 98:                 "MySolverCreate",MySolverCreate);
 99: .ve

101:    Then, your solver can be chosen with the procedural interface via
102: $     SNESSetType(snes,"my_solver")
103:    or at runtime via the option
104: $     -snes_type my_solver

106:    Level: advanced

108:     Note: If your function is not being put into a shared library then use SNESRegister() instead

110: .keywords: SNES, nonlinear, register

112: .seealso: SNESRegisterAll(), SNESRegisterDestroy()
113: M*/
114: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
115: #define SNESRegisterDynamic(a,b,c,d) SNESRegister(a,b,c,0)
116: #else
117: #define SNESRegisterDynamic(a,b,c,d) SNESRegister(a,b,c,d)
118: #endif

120: extern PetscErrorCode  SNESGetKSP(SNES,KSP*);
121: extern PetscErrorCode  SNESSetKSP(SNES,KSP);
122: extern PetscErrorCode  SNESGetSolution(SNES,Vec*);
123: extern PetscErrorCode  SNESGetSolutionUpdate(SNES,Vec*);
124: extern PetscErrorCode  SNESGetRhs(SNES,Vec*);
125: extern PetscErrorCode  SNESView(SNES,PetscViewer);

127: extern PetscErrorCode  SNESSetOptionsPrefix(SNES,const char[]);
128: extern PetscErrorCode  SNESAppendOptionsPrefix(SNES,const char[]);
129: extern PetscErrorCode  SNESGetOptionsPrefix(SNES,const char*[]);
130: extern PetscErrorCode  SNESSetFromOptions(SNES);
131: extern PetscErrorCode  SNESDefaultGetWork(SNES,PetscInt);

133: extern PetscErrorCode  MatCreateSNESMF(SNES,Mat*);
134: extern PetscErrorCode  MatMFFDComputeJacobian(SNES,Vec,Mat*,Mat*,MatStructure*,void*);

136: extern PetscErrorCode  MatDAADSetSNES(Mat,SNES);

138: extern PetscErrorCode  SNESGetType(SNES,const SNESType*);
139: extern PetscErrorCode  SNESMonitorDefault(SNES,PetscInt,PetscReal,void *);
140: extern PetscErrorCode  SNESMonitorRange(SNES,PetscInt,PetscReal,void *);
141: extern PetscErrorCode  SNESMonitorRatio(SNES,PetscInt,PetscReal,void *);
142: extern PetscErrorCode  SNESMonitorSetRatio(SNES,PetscViewer);
143: extern PetscErrorCode  SNESMonitorSolution(SNES,PetscInt,PetscReal,void *);
144: extern PetscErrorCode  SNESMonitorResidual(SNES,PetscInt,PetscReal,void *);
145: extern PetscErrorCode  SNESMonitorSolutionUpdate(SNES,PetscInt,PetscReal,void *);
146: extern PetscErrorCode  SNESMonitorDefaultShort(SNES,PetscInt,PetscReal,void *);
147: extern PetscErrorCode  SNESSetTolerances(SNES,PetscReal,PetscReal,PetscReal,PetscInt,PetscInt);
148: extern PetscErrorCode  SNESGetTolerances(SNES,PetscReal*,PetscReal*,PetscReal*,PetscInt*,PetscInt*);
149: extern PetscErrorCode  SNESSetTrustRegionTolerance(SNES,PetscReal);
150: extern PetscErrorCode  SNESGetFunctionNorm(SNES,PetscReal*);
151: extern PetscErrorCode  SNESGetIterationNumber(SNES,PetscInt*);

153: extern PetscErrorCode  SNESGetNonlinearStepFailures(SNES,PetscInt*);
154: extern PetscErrorCode  SNESSetMaxNonlinearStepFailures(SNES,PetscInt);
155: extern PetscErrorCode  SNESGetMaxNonlinearStepFailures(SNES,PetscInt*);
156: extern PetscErrorCode  SNESGetNumberFunctionEvals(SNES,PetscInt*);

158: extern PetscErrorCode  SNESSetLagPreconditioner(SNES,PetscInt);
159: extern PetscErrorCode  SNESGetLagPreconditioner(SNES,PetscInt*);
160: extern PetscErrorCode  SNESSetLagJacobian(SNES,PetscInt);
161: extern PetscErrorCode  SNESGetLagJacobian(SNES,PetscInt*);
162: extern PetscErrorCode  SNESSetGridSequence(SNES,PetscInt);

164: extern PetscErrorCode  SNESGetLinearSolveIterations(SNES,PetscInt*);
165: extern PetscErrorCode  SNESGetLinearSolveFailures(SNES,PetscInt*);
166: extern PetscErrorCode  SNESSetMaxLinearSolveFailures(SNES,PetscInt);
167: extern PetscErrorCode  SNESGetMaxLinearSolveFailures(SNES,PetscInt*);

169: extern PetscErrorCode  SNESKSPSetUseEW(SNES,PetscBool );
170: extern PetscErrorCode  SNESKSPGetUseEW(SNES,PetscBool *);
171: extern PetscErrorCode  SNESKSPSetParametersEW(SNES,PetscInt,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal);
172: extern PetscErrorCode  SNESKSPGetParametersEW(SNES,PetscInt*,PetscReal*,PetscReal*,PetscReal*,PetscReal*,PetscReal*,PetscReal*);

174: extern PetscErrorCode  SNESMonitorLGCreate(const char[],const char[],int,int,int,int,PetscDrawLG*);
175: extern PetscErrorCode  SNESMonitorLG(SNES,PetscInt,PetscReal,void*);
176: extern PetscErrorCode  SNESMonitorLGDestroy(PetscDrawLG*);
177: extern PetscErrorCode  SNESMonitorLGRangeCreate(const char[],const char[],int,int,int,int,PetscDrawLG*);
178: extern PetscErrorCode  SNESMonitorLGRange(SNES,PetscInt,PetscReal,void*);
179: extern PetscErrorCode  SNESMonitorLGRangeDestroy(PetscDrawLG*);

181: extern PetscErrorCode  SNESSetApplicationContext(SNES,void *);
182: extern PetscErrorCode  SNESGetApplicationContext(SNES,void *);
183: extern PetscErrorCode  SNESSetComputeApplicationContext(SNES,PetscErrorCode (*)(SNES,void**),PetscErrorCode (*)(void**));

185: extern PetscErrorCode  SNESPythonSetType(SNES,const char[]);

187: extern PetscErrorCode  SNESSetFunctionDomainError(SNES);
188: /*E
189:     SNESConvergedReason - reason a SNES method was said to 
190:          have converged or diverged

192:    Level: beginner

194:    The two most common reasons for divergence are 
195: $   1) an incorrectly coded or computed Jacobian or 
196: $   2) failure or lack of convergence in the linear system (in this case we recommend
197: $      testing with -pc_type lu to eliminate the linear solver as the cause of the problem).

199:    Diverged Reasons:
200: .    SNES_DIVERGED_LOCAL_MIN - this can only occur when using the line-search variant of SNES.
201:        The line search wants to minimize Q(alpha) = 1/2 || F(x + alpha s) ||^2_2  this occurs
202:        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
203:        you get Q'(alpha) = -F(x)^T F'(x)^(-1)^T F'(x+alpha s)F(x+alpha s); when alpha = 0
204:        Q'(0) = - ||F(x)||^2_2 which is always NEGATIVE if F'(x) is invertible. This means the Newton
205:        direction is a descent direction and the line search should succeed if alpha is small enough.

207:        If F'(x) is NOT invertible AND F'(x)^T F(x) = 0 then Q'(0) = 0 and the Newton direction 
208:        is NOT a descent direction so the line search will fail. All one can do at this point
209:        is change the initial guess and try again.

211:        An alternative explanation: Newton's method can be regarded as replacing the function with
212:        its linear approximation and minimizing the 2-norm of that. That is F(x+s) approx F(x) + F'(x)s
213:        so we minimize || F(x) + F'(x) s ||^2_2; do this using Least Squares. If F'(x) is invertible then
214:        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
215:        exists a nontrival (that is F'(x)s != 0) solution to the equation and this direction is 
216:        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)
217:        = - (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
218:        and F'(x)^T F'(x) has no negative eigenvalues Q'(0) < 0 so s is a descent direction and the line
219:        search should succeed for small enough alpha.

221:        Note that this RARELY happens in practice. Far more likely the linear system is not being solved
222:        (well enough?) or the Jacobian is wrong.
223:      
224:    SNES_DIVERGED_MAX_IT means that the solver reached the maximum number of iterations without satisfying any
225:    convergence criteria. SNES_CONVERGED_ITS means that SNESSkipConverged() was chosen as the convergence test;
226:    thus the usual convergence criteria have not been checked and may or may not be satisfied.

228:    Developer Notes: this must match finclude/petscsnes.h 

230:        The string versions of these are in SNESConvergedReason, if you change any value here you must
231:      also adjust that array.

233:    Each reason has its own manual page.

235: .seealso: SNESSolve(), SNESGetConvergedReason(), KSPConvergedReason, SNESSetConvergenceTest()
236: E*/
237: typedef enum {/* converged */
238:               SNES_CONVERGED_FNORM_ABS         =  2, /* ||F|| < atol */
239:               SNES_CONVERGED_FNORM_RELATIVE    =  3, /* ||F|| < rtol*||F_initial|| */
240:               SNES_CONVERGED_PNORM_RELATIVE    =  4, /* Newton computed step size small; || delta x || < stol */
241:               SNES_CONVERGED_ITS               =  5, /* maximum iterations reached */
242:               SNES_CONVERGED_TR_DELTA            =  7,
243:               /* diverged */
244:               SNES_DIVERGED_FUNCTION_DOMAIN     = -1, /* the new x location passed the function is not in the domain of F */
245:               SNES_DIVERGED_FUNCTION_COUNT      = -2,
246:               SNES_DIVERGED_LINEAR_SOLVE        = -3, /* the linear solve failed */
247:               SNES_DIVERGED_FNORM_NAN           = -4,
248:               SNES_DIVERGED_MAX_IT              = -5,
249:               SNES_DIVERGED_LINE_SEARCH         = -6, /* the line search failed */
250:               SNES_DIVERGED_LOCAL_MIN           = -8, /* || J^T b || is small, implies converged to local minimum of F() */
251:               SNES_CONVERGED_ITERATING          =  0} SNESConvergedReason;
252: extern const char *const*SNESConvergedReasons;

254: /*MC
255:      SNES_CONVERGED_FNORM_ABS - 2-norm(F) <= abstol

257:    Level: beginner

259: .seealso:  SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()

261: M*/

263: /*MC
264:      SNES_CONVERGED_FNORM_RELATIVE - 2-norm(F) <= rtol*2-norm(F(x_0)) where x_0 is the initial guess

266:    Level: beginner

268: .seealso:  SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()

270: M*/

272: /*MC
273:      SNES_CONVERGED_PNORM_RELATIVE - The 2-norm of the last step <= stol * 2-norm(x) where x is the current
274:           solution and stol is the 4th argument to SNESSetTolerances()

276:    Level: beginner

278: .seealso:  SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()

280: M*/

282: /*MC
283:      SNES_DIVERGED_FUNCTION_COUNT - The user provided function has been called more times then the final
284:          argument to SNESSetTolerances()

286:    Level: beginner

288: .seealso:  SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()

290: M*/

292: /*MC
293:      SNES_DIVERGED_FNORM_NAN - the 2-norm of the current function evaluation is not-a-number (NaN), this
294:       is usually caused by a division of 0 by 0.

296:    Level: beginner

298: .seealso:  SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()

300: M*/

302: /*MC
303:      SNES_DIVERGED_MAX_IT - SNESSolve() has reached the maximum number of iterations requested

305:    Level: beginner

307: .seealso:  SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()

309: M*/

311: /*MC
312:      SNES_DIVERGED_LINE_SEARCH - The line search has failed. This only occurs for a SNESType of SNESLS

314:    Level: beginner

316: .seealso:  SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()

318: M*/

320: /*MC
321:      SNES_DIVERGED_LOCAL_MIN - the algorithm seems to have stagnated at a local minimum that is not zero. 
322:         See the manual page for SNESConvergedReason for more details

324:    Level: beginner

326: .seealso:  SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()

328: M*/

330: /*MC
331:      SNES_CONERGED_ITERATING - this only occurs if SNESGetConvergedReason() is called during the SNESSolve()

333:    Level: beginner

335: .seealso:  SNESSolve(), SNESGetConvergedReason(), SNESConvergedReason, SNESSetTolerances()

337: M*/

339: extern PetscErrorCode  SNESSetConvergenceTest(SNES,PetscErrorCode (*)(SNES,PetscInt,PetscReal,PetscReal,PetscReal,SNESConvergedReason*,void*),void*,PetscErrorCode (*)(void*));
340: extern PetscErrorCode  SNESDefaultConverged(SNES,PetscInt,PetscReal,PetscReal,PetscReal,SNESConvergedReason*,void*);
341: extern PetscErrorCode  SNESSkipConverged(SNES,PetscInt,PetscReal,PetscReal,PetscReal,SNESConvergedReason*,void*);
342: extern PetscErrorCode  SNESGetConvergedReason(SNES,SNESConvergedReason*);

344: extern PetscErrorCode  SNESDAFormFunction(SNES,Vec,Vec,void*);
345: extern PetscErrorCode  SNESDAComputeJacobianWithAdic(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
346: extern PetscErrorCode  SNESDAComputeJacobianWithAdifor(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
347: extern PetscErrorCode  SNESDAComputeJacobian(SNES,Vec,Mat*,Mat*,MatStructure*,void*);

349: extern PetscErrorCode  SNESMeshFormFunction(SNES,Vec,Vec,void*);
350: extern PetscErrorCode SNESMeshFormJacobian(SNES,Vec,Mat*,Mat*,MatStructure*,void*);

352: /* --------- Solving systems of nonlinear equations --------------- */
353: typedef PetscErrorCode (*SNESFunction)(SNES,Vec,Vec,void*);
354: typedef PetscErrorCode (*SNESJacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
355: extern PetscErrorCode  SNESSetFunction(SNES,Vec,SNESFunction,void*);
356: extern PetscErrorCode  SNESGetFunction(SNES,Vec*,SNESFunction*,void**);
357: extern PetscErrorCode  SNESComputeFunction(SNES,Vec,Vec);
358: extern PetscErrorCode  SNESSetJacobian(SNES,Mat,Mat,SNESJacobian,void*);
359: extern PetscErrorCode  SNESGetJacobian(SNES,Mat*,Mat*,SNESJacobian*,void**);
360: extern PetscErrorCode  SNESDefaultComputeJacobian(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
361: extern PetscErrorCode  SNESDefaultComputeJacobianColor(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
362: extern PetscErrorCode  SNESSetComputeInitialGuess(SNES,PetscErrorCode (*)(SNES,Vec,void*),void*);

364: /* --------- Routines specifically for line search methods --------------- */
365: extern PetscErrorCode  SNESLineSearchSet(SNES,PetscErrorCode(*)(SNES,void*,Vec,Vec,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal*,PetscReal*,PetscBool *),void*);
366: extern PetscErrorCode  SNESLineSearchNo(SNES,void*,Vec,Vec,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal*,PetscReal*,PetscBool *);
367: extern PetscErrorCode  SNESLineSearchNoNorms(SNES,void*,Vec,Vec,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal*,PetscReal*,PetscBool *);
368: extern PetscErrorCode  SNESLineSearchCubic(SNES,void*,Vec,Vec,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal*,PetscReal*,PetscBool *);
369: extern PetscErrorCode  SNESLineSearchQuadratic(SNES,void*,Vec,Vec,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal*,PetscReal*,PetscBool *);

371: extern PetscErrorCode  SNESLineSearchSetPostCheck(SNES,PetscErrorCode(*)(SNES,Vec,Vec,Vec,void*,PetscBool *,PetscBool *),void*);
372: extern PetscErrorCode  SNESLineSearchSetPreCheck(SNES,PetscErrorCode(*)(SNES,Vec,Vec,void*,PetscBool *),void*);
373: extern PetscErrorCode  SNESLineSearchSetParams(SNES,PetscReal,PetscReal,PetscReal);
374: extern PetscErrorCode  SNESLineSearchGetParams(SNES,PetscReal*,PetscReal*,PetscReal*);
375: extern PetscErrorCode  SNESLineSearchSetMonitor(SNES,PetscBool );

377: /* Routines for VI solver */
378: extern PetscErrorCode  SNESVISetVariableBounds(SNES,Vec,Vec);
379: extern PetscErrorCode  SNESVISetComputeVariableBounds(SNES, PetscErrorCode (*)(SNES,Vec,Vec));
380: extern PetscErrorCode  SNESVIGetInactiveSet(SNES,IS*);
381: extern PetscErrorCode  SNESVISetRedundancyCheck(SNES,PetscErrorCode(*)(SNES,IS,IS*,void*),void*);
382: #define SNES_VI_INF   1.0e20
383: #define SNES_VI_NINF -1.0e20

385: extern PetscErrorCode  SNESTestLocalMin(SNES);

387: /* Should this routine be private? */
388: extern PetscErrorCode  SNESComputeJacobian(SNES,Vec,Mat*,Mat*,MatStructure*);

390: extern PetscErrorCode SNESSetDM(SNES,DM);
391: extern PetscErrorCode SNESGetDM(SNES,DM*);
392: extern PetscErrorCode SNESSetPC(SNES,SNES);
393: extern PetscErrorCode SNESGetPC(SNES,SNES*);

395: /* Routines for Multiblock solver */
396: PetscErrorCode SNESMultiblockSetFields(SNES, const char [], PetscInt, const PetscInt *);
397: PetscErrorCode SNESMultiblockSetIS(SNES, const char [], IS);
398: PetscErrorCode SNESMultiblockSetBlockSize(SNES, PetscInt);
399: PetscErrorCode SNESMultiblockSetType(SNES, PCCompositeType);

401: PETSC_EXTERN_CXX_END
402: #endif