Actual source code: snes.c
1: /*$Id: snes.c,v 1.229 2001/04/10 19:36:48 bsmith Exp bsmith $*/
3: #include src/snes/snesimpl.h
5: PetscTruth SNESRegisterAllCalled = PETSC_FALSE;
6: PetscFList SNESList = 0;
8: /*@C
9: SNESView - Prints the SNES data structure.
11: Collective on SNES
13: Input Parameters:
14: + SNES - the SNES context
15: - viewer - visualization context
17: Options Database Key:
18: . -snes_view - Calls SNESView() at end of SNESSolve()
20: Notes:
21: The available visualization contexts include
22: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
23: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
24: output where only the first processor opens
25: the file. All other processors send their
26: data to the first processor to print.
28: The user can open an alternative visualization context with
29: PetscViewerASCIIOpen() - output to a specified file.
31: Level: beginner
33: .keywords: SNES, view
35: .seealso: PetscViewerASCIIOpen()
36: @*/
37: int SNESView(SNES snes,PetscViewer viewer)
38: {
39: SNES_KSP_EW_ConvCtx *kctx;
40: int ierr;
41: SLES sles;
42: char *type;
43: PetscTruth isascii,isstring;
47: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(snes->comm);
51: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
52: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
53: if (isascii) {
54: PetscViewerASCIIPrintf(viewer,"SNES Object:n");
55: SNESGetType(snes,&type);
56: if (type) {
57: PetscViewerASCIIPrintf(viewer," type: %sn",type);
58: } else {
59: PetscViewerASCIIPrintf(viewer," type: not set yetn");
60: }
61: if (snes->view) {
62: PetscViewerASCIIPushTab(viewer);
63: (*snes->view)(snes,viewer);
64: PetscViewerASCIIPopTab(viewer);
65: }
66: PetscViewerASCIIPrintf(viewer," maximum iterations=%d, maximum function evaluations=%dn",snes->max_its,snes->max_funcs);
67: PetscViewerASCIIPrintf(viewer," tolerances: relative=%g, absolute=%g, truncation=%g, solution=%gn",
68: snes->rtol,snes->atol,snes->trunctol,snes->xtol);
69: PetscViewerASCIIPrintf(viewer," total number of linear solver iterations=%dn",snes->linear_its);
70: PetscViewerASCIIPrintf(viewer," total number of function evaluations=%dn",snes->nfuncs);
71: if (snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) {
72: PetscViewerASCIIPrintf(viewer," min function tolerance=%gn",snes->fmin);
73: }
74: if (snes->ksp_ewconv) {
75: kctx = (SNES_KSP_EW_ConvCtx *)snes->kspconvctx;
76: if (kctx) {
77: PetscViewerASCIIPrintf(viewer," Eisenstat-Walker computation of KSP relative tolerance (version %d)n",kctx->version);
78: PetscViewerASCIIPrintf(viewer," rtol_0=%g, rtol_max=%g, threshold=%gn",kctx->rtol_0,kctx->rtol_max,kctx->threshold);
79: PetscViewerASCIIPrintf(viewer," gamma=%g, alpha=%g, alpha2=%gn",kctx->gamma,kctx->alpha,kctx->alpha2);
80: }
81: }
82: } else if (isstring) {
83: SNESGetType(snes,&type);
84: PetscViewerStringSPrintf(viewer," %-3.3s",type);
85: }
86: SNESGetSLES(snes,&sles);
87: PetscViewerASCIIPushTab(viewer);
88: SLESView(sles,viewer);
89: PetscViewerASCIIPopTab(viewer);
90: return(0);
91: }
93: /*@
94: SNESSetFromOptions - Sets various SNES and SLES parameters from user options.
96: Collective on SNES
98: Input Parameter:
99: . snes - the SNES context
101: Options Database Keys:
102: + -snes_type <type> - ls, tr, umls, umtr, test
103: . -snes_stol - convergence tolerance in terms of the norm
104: of the change in the solution between steps
105: . -snes_atol <atol> - absolute tolerance of residual norm
106: . -snes_rtol <rtol> - relative decrease in tolerance norm from initial
107: . -snes_max_it <max_it> - maximum number of iterations
108: . -snes_max_funcs <max_funcs> - maximum number of function evaluations
109: . -snes_trtol <trtol> - trust region tolerance
110: . -snes_no_convergence_test - skip convergence test in nonlinear or minimization
111: solver; hence iterations will continue until max_it
112: or some other criterion is reached. Saves expense
113: of convergence test
114: . -snes_monitor - prints residual norm at each iteration
115: . -snes_vecmonitor - plots solution at each iteration
116: . -snes_vecmonitor_update - plots update to solution at each iteration
117: . -snes_xmonitor - plots residual norm at each iteration
118: . -snes_fd - use finite differences to compute Jacobian; very slow, only for testing
119: - -snes_mf_ksp_monitor - if using matrix-free multiply then print h at each KSP iteration
121: Options Database for Eisenstat-Walker method:
122: + -snes_ksp_eq_conv - use Eisenstat-Walker method for determining linear system convergence
123: . -snes_ksp_eq_version ver - version of Eisenstat-Walker method
124: . -snes_ksp_ew_rtol0 <rtol0> - Sets rtol0
125: . -snes_ksp_ew_rtolmax <rtolmax> - Sets rtolmax
126: . -snes_ksp_ew_gamma <gamma> - Sets gamma
127: . -snes_ksp_ew_alpha <alpha> - Sets alpha
128: . -snes_ksp_ew_alpha2 <alpha2> - Sets alpha2
129: - -snes_ksp_ew_threshold <threshold> - Sets threshold
131: Notes:
132: To see all options, run your program with the -help option or consult
133: the users manual.
135: Level: beginner
137: .keywords: SNES, nonlinear, set, options, database
139: .seealso: SNESSetOptionsPrefix(), SNESSetTypeFromOptions()
140: @*/
141: int SNESSetFromOptions(SNES snes)
142: {
143: SLES sles;
144: SNES_KSP_EW_ConvCtx *kctx = (SNES_KSP_EW_ConvCtx *)snes->kspconvctx;
145: PetscTruth flg;
146: int ierr;
147: char *deft,type[256];
152: PetscOptionsBegin(snes->comm,snes->prefix,"Nonlinear solver (SNES) options","SNES");
153: if (snes->type_name) {
154: deft = snes->type_name;
155: } else {
156: if (snes->method_class == SNES_NONLINEAR_EQUATIONS) {
157: deft = SNESEQLS;
158: } else {
159: deft = SNESUMTR;
160: }
161: }
163: if (!SNESRegisterAllCalled) {SNESRegisterAll(PETSC_NULL);}
164: PetscOptionsList("-snes_type","Nonlinear solver method","SNESSetType",SNESList,deft,type,256,&flg);
165: if (flg) {
166: SNESSetType(snes,type);
167: } else if (!snes->type_name) {
168: SNESSetType(snes,deft);
169: }
171: PetscOptionsDouble("-snes_stol","Stop if step length less then","SNESSetTolerances",snes->xtol,&snes->xtol,0);
172: PetscOptionsDouble("-snes_atol","Stop if function norm less then","SNESSetTolerances",snes->atol,&snes->atol,0);
174: PetscOptionsDouble("-snes_rtol","Stop if decrease in function norm less then","SNESSetTolerances",snes->rtol,&snes->rtol,0);
175: PetscOptionsInt("-snes_max_it","Maximum iterations","SNESSetTolerances",snes->max_its,&snes->max_its,PETSC_NULL);
176: PetscOptionsInt("-snes_max_funcs","Maximum function evaluations","SNESSetTolerances",snes->max_funcs,&snes->max_funcs,PETSC_NULL);
177: PetscOptionsDouble("-snes_fmin","Minimization function tolerance","SNESSetMinimizationFunctionTolerance",snes->fmin,&snes->fmin,0);
179: PetscOptionsName("-snes_ksp_ew_conv","Use Eisentat-Walker linear system convergence test","SNES_KSP_SetParametersEW",&snes->ksp_ewconv);
181: PetscOptionsInt("-snes_ksp_ew_version","Version 1 or 2","SNES_KSP_SetParametersEW",kctx->version,&kctx->version,0);
182: PetscOptionsDouble("-snes_ksp_ew_rtol0","0 <= rtol0 < 1","SNES_KSP_SetParametersEW",kctx->rtol_0,&kctx->rtol_0,0);
183: PetscOptionsDouble("-snes_ksp_ew_rtolmax","0 <= rtolmax < 1","SNES_KSP_SetParametersEW",kctx->rtol_max,&kctx->rtol_max,0);
184: PetscOptionsDouble("-snes_ksp_ew_gamma","0 <= gamma <= 1","SNES_KSP_SetParametersEW",kctx->gamma,&kctx->gamma,0);
185: PetscOptionsDouble("-snes_ksp_ew_alpha","1 < alpha <= 2","SNES_KSP_SetParametersEW",kctx->alpha,&kctx->alpha,0);
186: PetscOptionsDouble("-snes_ksp_ew_alpha2","alpha2","SNES_KSP_SetParametersEW",kctx->alpha2,&kctx->alpha2,0);
187: PetscOptionsDouble("-snes_ksp_ew_threshold","0 < threshold < 1","SNES_KSP_SetParametersEW",kctx->threshold,&kctx->threshold,0);
189: PetscOptionsName("-snes_no_convergence_test","Don't test for convergence","None",&flg);
190: if (flg) {snes->converged = 0;}
191: PetscOptionsName("-snes_cancelmonitors","Remove all monitors","SNESClearMonitor",&flg);
192: if (flg) {SNESClearMonitor(snes);}
193: PetscOptionsName("-snes_monitor","Monitor norm of function","SNESDefaultMonitor",&flg);
194: if (flg) {SNESSetMonitor(snes,SNESDefaultMonitor,0,0);}
195: PetscOptionsName("-snes_smonitor","Monitor norm of function (fewer digits)","SNESDefaultSMonitor",&flg);
196: if (flg) {SNESSetMonitor(snes,SNESDefaultSMonitor,0,0);}
197: PetscOptionsName("-snes_vecmonitor","Plot solution at each iteration","SNESVecViewMonitor",&flg);
198: if (flg) {SNESSetMonitor(snes,SNESVecViewMonitor,0,0);}
199: PetscOptionsName("-snes_vecmonitor_update","Plot correction at each iteration","SNESVecViewUpdateMonitor",&flg);
200: if (flg) {SNESSetMonitor(snes,SNESVecViewUpdateMonitor,0,0);}
201: PetscOptionsName("-snes_xmonitor","Plot function norm at each iteration","SNESLGMonitor",&flg);
202: if (flg) {SNESSetMonitor(snes,SNESLGMonitor,PETSC_NULL,PETSC_NULL);}
204: PetscOptionsName("-snes_fd","Use finite differences (slow) to compute Jacobian","SNESDefaultComputeJacobian",&flg);
205: if (flg && snes->method_class == SNES_NONLINEAR_EQUATIONS) {
206: SNESSetJacobian(snes,snes->jacobian,snes->jacobian_pre,SNESDefaultComputeJacobian,snes->funP);
207: PetscLogInfo(snes,"SNESSetFromOptions: Setting default finite difference Jacobian matrixn");
208: } else if (flg && snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) {
209: SNESSetHessian(snes,snes->jacobian,snes->jacobian_pre,SNESDefaultComputeHessian,snes->funP);
210: PetscLogInfo(snes,"SNESSetFromOptions: Setting default finite difference Hessian matrixn");
211: }
213: if (snes->setfromoptions) {
214: (*snes->setfromoptions)(snes);
215: }
217: PetscOptionsEnd();
219: SNESGetSLES(snes,&sles);
220: SLESSetFromOptions(sles);
222: return(0);
223: }
226: /*@
227: SNESSetApplicationContext - Sets the optional user-defined context for
228: the nonlinear solvers.
230: Collective on SNES
232: Input Parameters:
233: + snes - the SNES context
234: - usrP - optional user context
236: Level: intermediate
238: .keywords: SNES, nonlinear, set, application, context
240: .seealso: SNESGetApplicationContext()
241: @*/
242: int SNESSetApplicationContext(SNES snes,void *usrP)
243: {
246: snes->user = usrP;
247: return(0);
248: }
250: /*@C
251: SNESGetApplicationContext - Gets the user-defined context for the
252: nonlinear solvers.
254: Not Collective
256: Input Parameter:
257: . snes - SNES context
259: Output Parameter:
260: . usrP - user context
262: Level: intermediate
264: .keywords: SNES, nonlinear, get, application, context
266: .seealso: SNESSetApplicationContext()
267: @*/
268: int SNESGetApplicationContext(SNES snes,void **usrP)
269: {
272: *usrP = snes->user;
273: return(0);
274: }
276: /*@
277: SNESGetIterationNumber - Gets the number of nonlinear iterations completed
278: at this time.
280: Not Collective
282: Input Parameter:
283: . snes - SNES context
285: Output Parameter:
286: . iter - iteration number
288: Notes:
289: For example, during the computation of iteration 2 this would return 1.
291: This is useful for using lagged Jacobians (where one does not recompute the
292: Jacobian at each SNES iteration). For example, the code
293: .vb
294: SNESGetIterationNumber(snes,&it);
295: if (!(it % 2)) {
296: [compute Jacobian here]
297: }
298: .ve
299: can be used in your ComputeJacobian() function to cause the Jacobian to be
300: recomputed every second SNES iteration.
302: Level: intermediate
304: .keywords: SNES, nonlinear, get, iteration, number
305: @*/
306: int SNESGetIterationNumber(SNES snes,int* iter)
307: {
311: *iter = snes->iter;
312: return(0);
313: }
315: /*@
316: SNESGetFunctionNorm - Gets the norm of the current function that was set
317: with SNESSSetFunction().
319: Collective on SNES
321: Input Parameter:
322: . snes - SNES context
324: Output Parameter:
325: . fnorm - 2-norm of function
327: Note:
328: SNESGetFunctionNorm() is valid for SNES_NONLINEAR_EQUATIONS methods only.
329: A related routine for SNES_UNCONSTRAINED_MINIMIZATION methods is
330: SNESGetGradientNorm().
332: Level: intermediate
334: .keywords: SNES, nonlinear, get, function, norm
336: .seealso: SNESGetFunction()
337: @*/
338: int SNESGetFunctionNorm(SNES snes,Scalar *fnorm)
339: {
343: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
344: SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"For SNES_NONLINEAR_EQUATIONS only");
345: }
346: *fnorm = snes->norm;
347: return(0);
348: }
350: /*@
351: SNESGetGradientNorm - Gets the norm of the current gradient that was set
352: with SNESSSetGradient().
354: Collective on SNES
356: Input Parameter:
357: . snes - SNES context
359: Output Parameter:
360: . fnorm - 2-norm of gradient
362: Note:
363: SNESGetGradientNorm() is valid for SNES_UNCONSTRAINED_MINIMIZATION
364: methods only. A related routine for SNES_NONLINEAR_EQUATIONS methods
365: is SNESGetFunctionNorm().
367: Level: intermediate
369: .keywords: SNES, nonlinear, get, gradient, norm
371: .seelso: SNESSetGradient()
372: @*/
373: int SNESGetGradientNorm(SNES snes,Scalar *gnorm)
374: {
378: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
379: SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"For SNES_UNCONSTRAINED_MINIMIZATION only");
380: }
381: *gnorm = snes->norm;
382: return(0);
383: }
385: /*@
386: SNESGetNumberUnsuccessfulSteps - Gets the number of unsuccessful steps
387: attempted by the nonlinear solver.
389: Not Collective
391: Input Parameter:
392: . snes - SNES context
394: Output Parameter:
395: . nfails - number of unsuccessful steps attempted
397: Notes:
398: This counter is reset to zero for each successive call to SNESSolve().
400: Level: intermediate
402: .keywords: SNES, nonlinear, get, number, unsuccessful, steps
403: @*/
404: int SNESGetNumberUnsuccessfulSteps(SNES snes,int* nfails)
405: {
409: *nfails = snes->nfailures;
410: return(0);
411: }
413: /*@
414: SNESGetNumberLinearIterations - Gets the total number of linear iterations
415: used by the nonlinear solver.
417: Not Collective
419: Input Parameter:
420: . snes - SNES context
422: Output Parameter:
423: . lits - number of linear iterations
425: Notes:
426: This counter is reset to zero for each successive call to SNESSolve().
428: Level: intermediate
430: .keywords: SNES, nonlinear, get, number, linear, iterations
431: @*/
432: int SNESGetNumberLinearIterations(SNES snes,int* lits)
433: {
437: *lits = snes->linear_its;
438: return(0);
439: }
441: /*@C
442: SNESGetSLES - Returns the SLES context for a SNES solver.
444: Not Collective, but if SNES object is parallel, then SLES object is parallel
446: Input Parameter:
447: . snes - the SNES context
449: Output Parameter:
450: . sles - the SLES context
452: Notes:
453: The user can then directly manipulate the SLES context to set various
454: options, etc. Likewise, the user can then extract and manipulate the
455: KSP and PC contexts as well.
457: Level: beginner
459: .keywords: SNES, nonlinear, get, SLES, context
461: .seealso: SLESGetPC(), SLESGetKSP()
462: @*/
463: int SNESGetSLES(SNES snes,SLES *sles)
464: {
467: *sles = snes->sles;
468: return(0);
469: }
471: static int SNESPublish_Petsc(PetscObject obj)
472: {
473: #if defined(PETSC_HAVE_AMS)
474: SNES v = (SNES) obj;
475: int ierr;
476: #endif
480: #if defined(PETSC_HAVE_AMS)
481: /* if it is already published then return */
482: if (v->amem >=0) return(0);
484: PetscObjectPublishBaseBegin(obj);
485: AMS_Memory_add_field((AMS_Memory)v->amem,"Iteration",&v->iter,1,AMS_INT,AMS_READ,
486: AMS_COMMON,AMS_REDUCT_UNDEF);
487: AMS_Memory_add_field((AMS_Memory)v->amem,"Residual",&v->norm,1,AMS_DOUBLE,AMS_READ,
488: AMS_COMMON,AMS_REDUCT_UNDEF);
489: PetscObjectPublishBaseEnd(obj);
490: #endif
491: return(0);
492: }
494: /* -----------------------------------------------------------*/
495: /*@C
496: SNESCreate - Creates a nonlinear solver context.
498: Collective on MPI_Comm
500: Input Parameters:
501: + comm - MPI communicator
502: - type - type of method, either
503: SNES_NONLINEAR_EQUATIONS (for systems of nonlinear equations)
504: or SNES_UNCONSTRAINED_MINIMIZATION (for unconstrained minimization)
506: Output Parameter:
507: . outsnes - the new SNES context
509: Options Database Keys:
510: + -snes_mf - Activates default matrix-free Jacobian-vector products,
511: and no preconditioning matrix
512: . -snes_mf_operator - Activates default matrix-free Jacobian-vector
513: products, and a user-provided preconditioning matrix
514: as set by SNESSetJacobian()
515: - -snes_fd - Uses (slow!) finite differences to compute Jacobian
517: Level: beginner
519: .keywords: SNES, nonlinear, create, context
521: .seealso: SNESSolve(), SNESDestroy(), SNESProblemType, SNES
522: @*/
523: int SNESCreate(MPI_Comm comm,SNESProblemType type,SNES *outsnes)
524: {
525: int ierr;
526: SNES snes;
527: SNES_KSP_EW_ConvCtx *kctx;
530: *outsnes = 0;
531: if (type != SNES_UNCONSTRAINED_MINIMIZATION && type != SNES_NONLINEAR_EQUATIONS){
532: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"incorrect method type");
533: }
534: PetscHeaderCreate(snes,_p_SNES,int,SNES_COOKIE,0,"SNES",comm,SNESDestroy,SNESView);
535: PetscLogObjectCreate(snes);
536: snes->bops->publish = SNESPublish_Petsc;
537: snes->max_its = 50;
538: snes->max_funcs = 10000;
539: snes->norm = 0.0;
540: if (type == SNES_UNCONSTRAINED_MINIMIZATION) {
541: snes->rtol = 1.e-8;
542: snes->ttol = 0.0;
543: snes->atol = 1.e-10;
544: } else {
545: snes->rtol = 1.e-8;
546: snes->ttol = 0.0;
547: snes->atol = 1.e-50;
548: }
549: snes->xtol = 1.e-8;
550: snes->trunctol = 1.e-12; /* no longer used */
551: snes->nfuncs = 0;
552: snes->nfailures = 0;
553: snes->linear_its = 0;
554: snes->numbermonitors = 0;
555: snes->data = 0;
556: snes->view = 0;
557: snes->computeumfunction = 0;
558: snes->umfunP = 0;
559: snes->fc = 0;
560: snes->deltatol = 1.e-12;
561: snes->fmin = -1.e30;
562: snes->method_class = type;
563: snes->set_method_called = 0;
564: snes->setupcalled = 0;
565: snes->ksp_ewconv = PETSC_FALSE;
566: snes->vwork = 0;
567: snes->nwork = 0;
568: snes->conv_hist_len = 0;
569: snes->conv_hist_max = 0;
570: snes->conv_hist = PETSC_NULL;
571: snes->conv_hist_its = PETSC_NULL;
572: snes->conv_hist_reset = PETSC_TRUE;
573: snes->reason = SNES_CONVERGED_ITERATING;
575: /* Create context to compute Eisenstat-Walker relative tolerance for KSP */
576: PetscNew(SNES_KSP_EW_ConvCtx,&kctx);
577: PetscLogObjectMemory(snes,sizeof(SNES_KSP_EW_ConvCtx));
578: snes->kspconvctx = (void*)kctx;
579: kctx->version = 2;
580: kctx->rtol_0 = .3; /* Eisenstat and Walker suggest rtol_0=.5, but
581: this was too large for some test cases */
582: kctx->rtol_last = 0;
583: kctx->rtol_max = .9;
584: kctx->gamma = 1.0;
585: kctx->alpha2 = .5*(1.0 + sqrt(5.0));
586: kctx->alpha = kctx->alpha2;
587: kctx->threshold = .1;
588: kctx->lresid_last = 0;
589: kctx->norm_last = 0;
591: SLESCreate(comm,&snes->sles);
592: PetscLogObjectParent(snes,snes->sles)
594: *outsnes = snes;
595: PetscPublishAll(snes);
596: return(0);
597: }
599: /* --------------------------------------------------------------- */
600: /*@C
601: SNESSetFunction - Sets the function evaluation routine and function
602: vector for use by the SNES routines in solving systems of nonlinear
603: equations.
605: Collective on SNES
607: Input Parameters:
608: + snes - the SNES context
609: . func - function evaluation routine
610: . r - vector to store function value
611: - ctx - [optional] user-defined context for private data for the
612: function evaluation routine (may be PETSC_NULL)
614: Calling sequence of func:
615: $ func (SNES snes,Vec x,Vec f,void *ctx);
617: . f - function vector
618: - ctx - optional user-defined function context
620: Notes:
621: The Newton-like methods typically solve linear systems of the form
622: $ f'(x) x = -f(x),
623: where f'(x) denotes the Jacobian matrix and f(x) is the function.
625: SNESSetFunction() is valid for SNES_NONLINEAR_EQUATIONS methods only.
626: Analogous routines for SNES_UNCONSTRAINED_MINIMIZATION methods are
627: SNESSetMinimizationFunction() and SNESSetGradient();
629: Level: beginner
631: .keywords: SNES, nonlinear, set, function
633: .seealso: SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian()
634: @*/
635: int SNESSetFunction(SNES snes,Vec r,int (*func)(SNES,Vec,Vec,void*),void *ctx)
636: {
641: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
642: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
643: }
645: snes->computefunction = func;
646: snes->vec_func = snes->vec_func_always = r;
647: snes->funP = ctx;
648: return(0);
649: }
651: /*@
652: SNESComputeFunction - Calls the function that has been set with
653: SNESSetFunction().
655: Collective on SNES
657: Input Parameters:
658: + snes - the SNES context
659: - x - input vector
661: Output Parameter:
662: . y - function vector, as set by SNESSetFunction()
664: Notes:
665: SNESComputeFunction() is valid for SNES_NONLINEAR_EQUATIONS methods only.
666: Analogous routines for SNES_UNCONSTRAINED_MINIMIZATION methods are
667: SNESComputeMinimizationFunction() and SNESComputeGradient();
669: SNESComputeFunction() is typically used within nonlinear solvers
670: implementations, so most users would not generally call this routine
671: themselves.
673: Level: developer
675: .keywords: SNES, nonlinear, compute, function
677: .seealso: SNESSetFunction(), SNESGetFunction()
678: @*/
679: int SNESComputeFunction(SNES snes,Vec x,Vec y)
680: {
681: int ierr;
689: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
690: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
691: }
693: PetscLogEventBegin(SNES_FunctionEval,snes,x,y,0);
694: PetscStackPush("SNES user function");
695: (*snes->computefunction)(snes,x,y,snes->funP);
696: PetscStackPop;
697: snes->nfuncs++;
698: PetscLogEventEnd(SNES_FunctionEval,snes,x,y,0);
699: return(0);
700: }
702: /*@C
703: SNESSetMinimizationFunction - Sets the function evaluation routine for
704: unconstrained minimization.
706: Collective on SNES
708: Input Parameters:
709: + snes - the SNES context
710: . func - function evaluation routine
711: - ctx - [optional] user-defined context for private data for the
712: function evaluation routine (may be PETSC_NULL)
714: Calling sequence of func:
715: $ func (SNES snes,Vec x,PetscReal *f,void *ctx);
717: + x - input vector
718: . f - function
719: - ctx - [optional] user-defined function context
721: Level: beginner
723: Notes:
724: SNESSetMinimizationFunction() is valid for SNES_UNCONSTRAINED_MINIMIZATION
725: methods only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is
726: SNESSetFunction().
728: .keywords: SNES, nonlinear, set, minimization, function
730: .seealso: SNESGetMinimizationFunction(), SNESComputeMinimizationFunction(),
731: SNESSetHessian(), SNESSetGradient()
732: @*/
733: int SNESSetMinimizationFunction(SNES snes,int (*func)(SNES,Vec,PetscReal*,void*),void *ctx)
734: {
737: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
738: SETERRQ(PETSC_ERR_ARG_WRONG,"Only for SNES_UNCONSTRAINED_MINIMIZATION");
739: }
740: snes->computeumfunction = func;
741: snes->umfunP = ctx;
742: return(0);
743: }
745: /*@
746: SNESComputeMinimizationFunction - Computes the function that has been
747: set with SNESSetMinimizationFunction().
749: Collective on SNES
751: Input Parameters:
752: + snes - the SNES context
753: - x - input vector
755: Output Parameter:
756: . y - function value
758: Notes:
759: SNESComputeMinimizationFunction() is valid only for
760: SNES_UNCONSTRAINED_MINIMIZATION methods. An analogous routine for
761: SNES_NONLINEAR_EQUATIONS methods is SNESComputeFunction().
763: SNESComputeMinimizationFunction() is typically used within minimization
764: implementations, so most users would not generally call this routine
765: themselves.
767: Level: developer
769: .keywords: SNES, nonlinear, compute, minimization, function
771: .seealso: SNESSetMinimizationFunction(), SNESGetMinimizationFunction(),
772: SNESComputeGradient(), SNESComputeHessian()
773: @*/
774: int SNESComputeMinimizationFunction(SNES snes,Vec x,PetscReal *y)
775: {
776: int ierr;
782: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
783: SETERRQ(PETSC_ERR_ARG_WRONG,"Only for SNES_UNCONSTRAINED_MINIMIZATION");
784: }
786: PetscLogEventBegin(SNES_MinimizationFunctionEval,snes,x,y,0);
787: PetscStackPush("SNES user minimzation function");
788: (*snes->computeumfunction)(snes,x,y,snes->umfunP);
789: PetscStackPop;
790: snes->nfuncs++;
791: PetscLogEventEnd(SNES_MinimizationFunctionEval,snes,x,y,0);
792: return(0);
793: }
795: /*@C
796: SNESSetGradient - Sets the gradient evaluation routine and gradient
797: vector for use by the SNES routines.
799: Collective on SNES
801: Input Parameters:
802: + snes - the SNES context
803: . func - function evaluation routine
804: . ctx - optional user-defined context for private data for the
805: gradient evaluation routine (may be PETSC_NULL)
806: - r - vector to store gradient value
808: Calling sequence of func:
809: $ func (SNES, Vec x, Vec g, void *ctx);
811: + x - input vector
812: . g - gradient vector
813: - ctx - optional user-defined gradient context
815: Notes:
816: SNESSetMinimizationFunction() is valid for SNES_UNCONSTRAINED_MINIMIZATION
817: methods only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is
818: SNESSetFunction().
820: Level: beginner
822: .keywords: SNES, nonlinear, set, function
824: .seealso: SNESGetGradient(), SNESComputeGradient(), SNESSetHessian(),
825: SNESSetMinimizationFunction(),
826: @*/
827: int SNESSetGradient(SNES snes,Vec r,int (*func)(SNES,Vec,Vec,void*),void *ctx)
828: {
833: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
834: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
835: }
836: snes->computefunction = func;
837: snes->vec_func = snes->vec_func_always = r;
838: snes->funP = ctx;
839: return(0);
840: }
842: /*@
843: SNESComputeGradient - Computes the gradient that has been set with
844: SNESSetGradient().
846: Collective on SNES
848: Input Parameters:
849: + snes - the SNES context
850: - x - input vector
852: Output Parameter:
853: . y - gradient vector
855: Notes:
856: SNESComputeGradient() is valid only for
857: SNES_UNCONSTRAINED_MINIMIZATION methods. An analogous routine for
858: SNES_NONLINEAR_EQUATIONS methods is SNESComputeFunction().
860: SNESComputeGradient() is typically used within minimization
861: implementations, so most users would not generally call this routine
862: themselves.
864: Level: developer
866: .keywords: SNES, nonlinear, compute, gradient
868: .seealso: SNESSetGradient(), SNESGetGradient(),
869: SNESComputeMinimizationFunction(), SNESComputeHessian()
870: @*/
871: int SNESComputeGradient(SNES snes,Vec x,Vec y)
872: {
873: int ierr;
881: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
882: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
883: }
884: PetscLogEventBegin(SNES_GradientEval,snes,x,y,0);
885: PetscStackPush("SNES user gradient function");
886: (*snes->computefunction)(snes,x,y,snes->funP);
887: PetscStackPop;
888: PetscLogEventEnd(SNES_GradientEval,snes,x,y,0);
889: return(0);
890: }
892: /*@
893: SNESComputeJacobian - Computes the Jacobian matrix that has been
894: set with SNESSetJacobian().
896: Collective on SNES and Mat
898: Input Parameters:
899: + snes - the SNES context
900: - x - input vector
902: Output Parameters:
903: + A - Jacobian matrix
904: . B - optional preconditioning matrix
905: - flag - flag indicating matrix structure
907: Notes:
908: Most users should not need to explicitly call this routine, as it
909: is used internally within the nonlinear solvers.
911: See SLESSetOperators() for important information about setting the
912: flag parameter.
914: SNESComputeJacobian() is valid only for SNES_NONLINEAR_EQUATIONS
915: methods. An analogous routine for SNES_UNCONSTRAINED_MINIMIZATION
916: methods is SNESComputeHessian().
918: Level: developer
920: .keywords: SNES, compute, Jacobian, matrix
922: .seealso: SNESSetJacobian(), SLESSetOperators()
923: @*/
924: int SNESComputeJacobian(SNES snes,Vec X,Mat *A,Mat *B,MatStructure *flg)
925: {
926: int ierr;
932: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
933: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
934: }
935: if (!snes->computejacobian) return(0);
936: PetscLogEventBegin(SNES_JacobianEval,snes,X,*A,*B);
937: *flg = DIFFERENT_NONZERO_PATTERN;
938: PetscStackPush("SNES user Jacobian function");
939: (*snes->computejacobian)(snes,X,A,B,flg,snes->jacP);
940: PetscStackPop;
941: PetscLogEventEnd(SNES_JacobianEval,snes,X,*A,*B);
942: /* make sure user returned a correct Jacobian and preconditioner */
945: return(0);
946: }
948: /*@
949: SNESComputeHessian - Computes the Hessian matrix that has been
950: set with SNESSetHessian().
952: Collective on SNES and Mat
954: Input Parameters:
955: + snes - the SNES context
956: - x - input vector
958: Output Parameters:
959: + A - Hessian matrix
960: . B - optional preconditioning matrix
961: - flag - flag indicating matrix structure
963: Notes:
964: Most users should not need to explicitly call this routine, as it
965: is used internally within the nonlinear solvers.
967: See SLESSetOperators() for important information about setting the
968: flag parameter.
970: SNESComputeHessian() is valid only for
971: SNES_UNCONSTRAINED_MINIMIZATION methods. An analogous routine for
972: SNES_NONLINEAR_EQUATIONS methods is SNESComputeJacobian().
974: SNESComputeHessian() is typically used within minimization
975: implementations, so most users would not generally call this routine
976: themselves.
978: Level: developer
980: .keywords: SNES, compute, Hessian, matrix
982: .seealso: SNESSetHessian(), SLESSetOperators(), SNESComputeGradient(),
983: SNESComputeMinimizationFunction()
984: @*/
985: int SNESComputeHessian(SNES snes,Vec x,Mat *A,Mat *B,MatStructure *flag)
986: {
987: int ierr;
993: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
994: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
995: }
996: if (!snes->computejacobian) return(0);
997: PetscLogEventBegin(SNES_HessianEval,snes,x,*A,*B);
998: *flag = DIFFERENT_NONZERO_PATTERN;
999: PetscStackPush("SNES user Hessian function");
1000: (*snes->computejacobian)(snes,x,A,B,flag,snes->jacP);
1001: PetscStackPop;
1002: PetscLogEventEnd(SNES_HessianEval,snes,x,*A,*B);
1003: /* make sure user returned a correct Jacobian and preconditioner */
1006: return(0);
1007: }
1009: /*@C
1010: SNESSetJacobian - Sets the function to compute Jacobian as well as the
1011: location to store the matrix.
1013: Collective on SNES and Mat
1015: Input Parameters:
1016: + snes - the SNES context
1017: . A - Jacobian matrix
1018: . B - preconditioner matrix (usually same as the Jacobian)
1019: . func - Jacobian evaluation routine
1020: - ctx - [optional] user-defined context for private data for the
1021: Jacobian evaluation routine (may be PETSC_NULL)
1023: Calling sequence of func:
1024: $ func (SNES snes,Vec x,Mat *A,Mat *B,int *flag,void *ctx);
1026: + x - input vector
1027: . A - Jacobian matrix
1028: . B - preconditioner matrix, usually the same as A
1029: . flag - flag indicating information about the preconditioner matrix
1030: structure (same as flag in SLESSetOperators())
1031: - ctx - [optional] user-defined Jacobian context
1033: Notes:
1034: See SLESSetOperators() for important information about setting the flag
1035: output parameter in the routine func(). Be sure to read this information!
1037: The routine func() takes Mat * as the matrix arguments rather than Mat.
1038: This allows the Jacobian evaluation routine to replace A and/or B with a
1039: completely new new matrix structure (not just different matrix elements)
1040: when appropriate, for instance, if the nonzero structure is changing
1041: throughout the global iterations.
1043: Level: beginner
1045: .keywords: SNES, nonlinear, set, Jacobian, matrix
1047: .seealso: SLESSetOperators(), SNESSetFunction()
1048: @*/
1049: int SNESSetJacobian(SNES snes,Mat A,Mat B,int (*func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
1050: {
1057: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
1058: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
1059: }
1061: snes->computejacobian = func;
1062: snes->jacP = ctx;
1063: snes->jacobian = A;
1064: snes->jacobian_pre = B;
1065: return(0);
1066: }
1068: /*@C
1069: SNESGetJacobian - Returns the Jacobian matrix and optionally the user
1070: provided context for evaluating the Jacobian.
1072: Not Collective, but Mat object will be parallel if SNES object is
1074: Input Parameter:
1075: . snes - the nonlinear solver context
1077: Output Parameters:
1078: + A - location to stash Jacobian matrix (or PETSC_NULL)
1079: . B - location to stash preconditioner matrix (or PETSC_NULL)
1080: . ctx - location to stash Jacobian ctx (or PETSC_NULL)
1081: - func - location to put Jacobian function (or PETSC_NULL)
1083: Level: advanced
1085: .seealso: SNESSetJacobian(), SNESComputeJacobian()
1086: @*/
1087: int SNESGetJacobian(SNES snes,Mat *A,Mat *B,void **ctx,int (**func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*))
1088: {
1091: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
1092: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
1093: }
1094: if (A) *A = snes->jacobian;
1095: if (B) *B = snes->jacobian_pre;
1096: if (ctx) *ctx = snes->jacP;
1097: if (func) *func = snes->computejacobian;
1098: return(0);
1099: }
1101: /*@C
1102: SNESSetHessian - Sets the function to compute Hessian as well as the
1103: location to store the matrix.
1105: Collective on SNES and Mat
1107: Input Parameters:
1108: + snes - the SNES context
1109: . A - Hessian matrix
1110: . B - preconditioner matrix (usually same as the Hessian)
1111: . func - Jacobian evaluation routine
1112: - ctx - [optional] user-defined context for private data for the
1113: Hessian evaluation routine (may be PETSC_NULL)
1115: Calling sequence of func:
1116: $ func (SNES snes,Vec x,Mat *A,Mat *B,int *flag,void *ctx);
1118: + x - input vector
1119: . A - Hessian matrix
1120: . B - preconditioner matrix, usually the same as A
1121: . flag - flag indicating information about the preconditioner matrix
1122: structure (same as flag in SLESSetOperators())
1123: - ctx - [optional] user-defined Hessian context
1125: Notes:
1126: See SLESSetOperators() for important information about setting the flag
1127: output parameter in the routine func(). Be sure to read this information!
1129: The function func() takes Mat * as the matrix arguments rather than Mat.
1130: This allows the Hessian evaluation routine to replace A and/or B with a
1131: completely new new matrix structure (not just different matrix elements)
1132: when appropriate, for instance, if the nonzero structure is changing
1133: throughout the global iterations.
1135: Level: beginner
1137: .keywords: SNES, nonlinear, set, Hessian, matrix
1139: .seealso: SNESSetMinimizationFunction(), SNESSetGradient(), SLESSetOperators()
1140: @*/
1141: int SNESSetHessian(SNES snes,Mat A,Mat B,int (*func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
1142: {
1149: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
1150: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
1151: }
1152: snes->computejacobian = func;
1153: snes->jacP = ctx;
1154: snes->jacobian = A;
1155: snes->jacobian_pre = B;
1156: return(0);
1157: }
1159: /*@
1160: SNESGetHessian - Returns the Hessian matrix and optionally the user
1161: provided context for evaluating the Hessian.
1163: Not Collective, but Mat object is parallel if SNES object is parallel
1165: Input Parameter:
1166: . snes - the nonlinear solver context
1168: Output Parameters:
1169: + A - location to stash Hessian matrix (or PETSC_NULL)
1170: . B - location to stash preconditioner matrix (or PETSC_NULL)
1171: - ctx - location to stash Hessian ctx (or PETSC_NULL)
1173: Level: advanced
1175: .seealso: SNESSetHessian(), SNESComputeHessian()
1177: .keywords: SNES, get, Hessian
1178: @*/
1179: int SNESGetHessian(SNES snes,Mat *A,Mat *B,void **ctx)
1180: {
1183: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION){
1184: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
1185: }
1186: if (A) *A = snes->jacobian;
1187: if (B) *B = snes->jacobian_pre;
1188: if (ctx) *ctx = snes->jacP;
1189: return(0);
1190: }
1192: /* ----- Routines to initialize and destroy a nonlinear solver ---- */
1194: /*@
1195: SNESSetUp - Sets up the internal data structures for the later use
1196: of a nonlinear solver.
1198: Collective on SNES
1200: Input Parameters:
1201: + snes - the SNES context
1202: - x - the solution vector
1204: Notes:
1205: For basic use of the SNES solvers the user need not explicitly call
1206: SNESSetUp(), since these actions will automatically occur during
1207: the call to SNESSolve(). However, if one wishes to control this
1208: phase separately, SNESSetUp() should be called after SNESCreate()
1209: and optional routines of the form SNESSetXXX(), but before SNESSolve().
1211: Level: advanced
1213: .keywords: SNES, nonlinear, setup
1215: .seealso: SNESCreate(), SNESSolve(), SNESDestroy()
1216: @*/
1217: int SNESSetUp(SNES snes,Vec x)
1218: {
1219: int ierr;
1220: PetscTruth flg;
1226: snes->vec_sol = snes->vec_sol_always = x;
1228: PetscOptionsHasName(snes->prefix,"-snes_mf_operator",&flg);
1229: /*
1230: This version replaces the user provided Jacobian matrix with a
1231: matrix-free version but still employs the user-provided preconditioner matrix
1232: */
1233: if (flg) {
1234: Mat J;
1235: MatCreateSNESMF(snes,snes->vec_sol,&J);
1236: PetscLogObjectParent(snes,J);
1237: snes->mfshell = J;
1238: snes->jacobian = J;
1239: PetscLogInfo(snes,"SNESSetUp: Setting default matrix-free operator routinesn");
1240: MatSNESMFSetFromOptions(J);
1241: }
1242: PetscOptionsHasName(snes->prefix,"-snes_mf",&flg);
1243: /*
1244: This version replaces both the user-provided Jacobian and the user-
1245: provided preconditioner matrix with the default matrix free version.
1246: */
1247: if (flg) {
1248: Mat J;
1249: SLES sles;
1250: PC pc;
1252: MatCreateSNESMF(snes,snes->vec_sol,&J);
1253: PetscLogObjectParent(snes,J);
1254: snes->mfshell = J;
1255: PetscLogInfo(snes,"SNESSetUp: Setting default matrix-free operator and preconditioner routinesn");
1256: if (snes->method_class == SNES_NONLINEAR_EQUATIONS) {
1257: SNESSetJacobian(snes,J,J,MatSNESMFFormJacobian,snes->funP);
1258: } else if (snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) {
1259: SNESSetHessian(snes,J,J,MatSNESMFFormJacobian,snes->funP);
1260: } else {
1261: SETERRQ(PETSC_ERR_SUP,"Method class doesn't support matrix-free option");
1262: }
1263: MatSNESMFSetFromOptions(J);
1264: /* force no preconditioner */
1265: SNESGetSLES(snes,&sles);
1266: SLESGetPC(sles,&pc);
1267: PCSetType(pc,PCNONE);
1268: }
1270: if ((snes->method_class == SNES_NONLINEAR_EQUATIONS)) {
1271: PetscTruth iseqtr;
1273: if (!snes->vec_func) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetFunction() first");
1274: if (!snes->computefunction) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetFunction() first");
1275: if (!snes->jacobian) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetJacobian() first n or use -snes_mf option");
1276: if (snes->vec_func == snes->vec_sol) {
1277: SETERRQ(PETSC_ERR_ARG_IDN,"Solution vector cannot be function vector");
1278: }
1280: /* Set the KSP stopping criterion to use the Eisenstat-Walker method */
1281: PetscTypeCompare((PetscObject)snes,SNESEQTR,&iseqtr);
1282: if (snes->ksp_ewconv && !iseqtr) {
1283: SLES sles; KSP ksp;
1284: SNESGetSLES(snes,&sles);
1285: SLESGetKSP(sles,&ksp);
1286: KSPSetConvergenceTest(ksp,SNES_KSP_EW_Converged_Private,snes);
1287: }
1288: } else if ((snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION)) {
1289: if (!snes->vec_func) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetGradient() first");
1290: if (!snes->computefunction) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetGradient() first");
1291: if (!snes->computeumfunction) {
1292: SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetMinimizationFunction() first");
1293: }
1294: if (!snes->jacobian) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetHessian()");
1295: } else {
1296: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Unknown method class");
1297: }
1298: if (snes->setup) {(*snes->setup)(snes);}
1299: snes->setupcalled = 1;
1300: return(0);
1301: }
1303: /*@C
1304: SNESDestroy - Destroys the nonlinear solver context that was created
1305: with SNESCreate().
1307: Collective on SNES
1309: Input Parameter:
1310: . snes - the SNES context
1312: Level: beginner
1314: .keywords: SNES, nonlinear, destroy
1316: .seealso: SNESCreate(), SNESSolve()
1317: @*/
1318: int SNESDestroy(SNES snes)
1319: {
1320: int i,ierr;
1324: if (--snes->refct > 0) return(0);
1326: /* if memory was published with AMS then destroy it */
1327: PetscObjectDepublish(snes);
1329: if (snes->destroy) {(*(snes)->destroy)(snes);}
1330: if (snes->kspconvctx) {PetscFree(snes->kspconvctx);}
1331: if (snes->mfshell) {MatDestroy(snes->mfshell);}
1332: SLESDestroy(snes->sles);
1333: if (snes->vwork) {VecDestroyVecs(snes->vwork,snes->nvwork);}
1334: for (i=0; i<snes->numbermonitors; i++) {
1335: if (snes->monitordestroy[i]) {
1336: (*snes->monitordestroy[i])(snes->monitorcontext[i]);
1337: }
1338: }
1339: PetscLogObjectDestroy((PetscObject)snes);
1340: PetscHeaderDestroy((PetscObject)snes);
1341: return(0);
1342: }
1344: /* ----------- Routines to set solver parameters ---------- */
1346: /*@
1347: SNESSetTolerances - Sets various parameters used in convergence tests.
1349: Collective on SNES
1351: Input Parameters:
1352: + snes - the SNES context
1353: . atol - absolute convergence tolerance
1354: . rtol - relative convergence tolerance
1355: . stol - convergence tolerance in terms of the norm
1356: of the change in the solution between steps
1357: . maxit - maximum number of iterations
1358: - maxf - maximum number of function evaluations
1360: Options Database Keys:
1361: + -snes_atol <atol> - Sets atol
1362: . -snes_rtol <rtol> - Sets rtol
1363: . -snes_stol <stol> - Sets stol
1364: . -snes_max_it <maxit> - Sets maxit
1365: - -snes_max_funcs <maxf> - Sets maxf
1367: Notes:
1368: The default maximum number of iterations is 50.
1369: The default maximum number of function evaluations is 1000.
1371: Level: intermediate
1373: .keywords: SNES, nonlinear, set, convergence, tolerances
1375: .seealso: SNESSetTrustRegionTolerance(), SNESSetMinimizationFunctionTolerance()
1376: @*/
1377: int SNESSetTolerances(SNES snes,PetscReal atol,PetscReal rtol,PetscReal stol,int maxit,int maxf)
1378: {
1381: if (atol != PETSC_DEFAULT) snes->atol = atol;
1382: if (rtol != PETSC_DEFAULT) snes->rtol = rtol;
1383: if (stol != PETSC_DEFAULT) snes->xtol = stol;
1384: if (maxit != PETSC_DEFAULT) snes->max_its = maxit;
1385: if (maxf != PETSC_DEFAULT) snes->max_funcs = maxf;
1386: return(0);
1387: }
1389: /*@
1390: SNESGetTolerances - Gets various parameters used in convergence tests.
1392: Not Collective
1394: Input Parameters:
1395: + snes - the SNES context
1396: . atol - absolute convergence tolerance
1397: . rtol - relative convergence tolerance
1398: . stol - convergence tolerance in terms of the norm
1399: of the change in the solution between steps
1400: . maxit - maximum number of iterations
1401: - maxf - maximum number of function evaluations
1403: Notes:
1404: The user can specify PETSC_NULL for any parameter that is not needed.
1406: Level: intermediate
1408: .keywords: SNES, nonlinear, get, convergence, tolerances
1410: .seealso: SNESSetTolerances()
1411: @*/
1412: int SNESGetTolerances(SNES snes,PetscReal *atol,PetscReal *rtol,PetscReal *stol,int *maxit,int *maxf)
1413: {
1416: if (atol) *atol = snes->atol;
1417: if (rtol) *rtol = snes->rtol;
1418: if (stol) *stol = snes->xtol;
1419: if (maxit) *maxit = snes->max_its;
1420: if (maxf) *maxf = snes->max_funcs;
1421: return(0);
1422: }
1424: /*@
1425: SNESSetTrustRegionTolerance - Sets the trust region parameter tolerance.
1427: Collective on SNES
1429: Input Parameters:
1430: + snes - the SNES context
1431: - tol - tolerance
1432:
1433: Options Database Key:
1434: . -snes_trtol <tol> - Sets tol
1436: Level: intermediate
1438: .keywords: SNES, nonlinear, set, trust region, tolerance
1440: .seealso: SNESSetTolerances(), SNESSetMinimizationFunctionTolerance()
1441: @*/
1442: int SNESSetTrustRegionTolerance(SNES snes,PetscReal tol)
1443: {
1446: snes->deltatol = tol;
1447: return(0);
1448: }
1450: /*@
1451: SNESSetMinimizationFunctionTolerance - Sets the minimum allowable function tolerance
1452: for unconstrained minimization solvers.
1453:
1454: Collective on SNES
1456: Input Parameters:
1457: + snes - the SNES context
1458: - ftol - minimum function tolerance
1460: Options Database Key:
1461: . -snes_fmin <ftol> - Sets ftol
1463: Note:
1464: SNESSetMinimizationFunctionTolerance() is valid for SNES_UNCONSTRAINED_MINIMIZATION
1465: methods only.
1467: Level: intermediate
1469: .keywords: SNES, nonlinear, set, minimum, convergence, function, tolerance
1471: .seealso: SNESSetTolerances(), SNESSetTrustRegionTolerance()
1472: @*/
1473: int SNESSetMinimizationFunctionTolerance(SNES snes,PetscReal ftol)
1474: {
1477: snes->fmin = ftol;
1478: return(0);
1479: }
1480: /*
1481: Duplicate the lg monitors for SNES from KSP; for some reason with
1482: dynamic libraries things don't work under Sun4 if we just use
1483: macros instead of functions
1484: */
1485: int SNESLGMonitor(SNES snes,int it,PetscReal norm,void *ctx)
1486: {
1491: KSPLGMonitor((KSP)snes,it,norm,ctx);
1492: return(0);
1493: }
1495: int SNESLGMonitorCreate(char *host,char *label,int x,int y,int m,int n,PetscDrawLG *draw)
1496: {
1500: KSPLGMonitorCreate(host,label,x,y,m,n,draw);
1501: return(0);
1502: }
1504: int SNESLGMonitorDestroy(PetscDrawLG draw)
1505: {
1509: KSPLGMonitorDestroy(draw);
1510: return(0);
1511: }
1513: /* ------------ Routines to set performance monitoring options ----------- */
1515: /*@C
1516: SNESSetMonitor - Sets an ADDITIONAL function that is to be used at every
1517: iteration of the nonlinear solver to display the iteration's
1518: progress.
1520: Collective on SNES
1522: Input Parameters:
1523: + snes - the SNES context
1524: . func - monitoring routine
1525: . mctx - [optional] user-defined context for private data for the
1526: monitor routine (use PETSC_NULL if no context is desitre)
1527: - monitordestroy - [optional] routine that frees monitor context
1528: (may be PETSC_NULL)
1530: Calling sequence of func:
1531: $ int func(SNES snes,int its, PetscReal norm,void *mctx)
1533: + snes - the SNES context
1534: . its - iteration number
1535: . norm - 2-norm function value (may be estimated)
1536: for SNES_NONLINEAR_EQUATIONS methods
1537: . norm - 2-norm gradient value (may be estimated)
1538: for SNES_UNCONSTRAINED_MINIMIZATION methods
1539: - mctx - [optional] monitoring context
1541: Options Database Keys:
1542: + -snes_monitor - sets SNESDefaultMonitor()
1543: . -snes_xmonitor - sets line graph monitor,
1544: uses SNESLGMonitorCreate()
1545: _ -snes_cancelmonitors - cancels all monitors that have
1546: been hardwired into a code by
1547: calls to SNESSetMonitor(), but
1548: does not cancel those set via
1549: the options database.
1551: Notes:
1552: Several different monitoring routines may be set by calling
1553: SNESSetMonitor() multiple times; all will be called in the
1554: order in which they were set.
1556: Level: intermediate
1558: .keywords: SNES, nonlinear, set, monitor
1560: .seealso: SNESDefaultMonitor(), SNESClearMonitor()
1561: @*/
1562: int SNESSetMonitor(SNES snes,int (*func)(SNES,int,PetscReal,void*),void *mctx,int (*monitordestroy)(void *))
1563: {
1566: if (snes->numbermonitors >= MAXSNESMONITORS) {
1567: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1568: }
1570: snes->monitor[snes->numbermonitors] = func;
1571: snes->monitordestroy[snes->numbermonitors] = monitordestroy;
1572: snes->monitorcontext[snes->numbermonitors++] = (void*)mctx;
1573: return(0);
1574: }
1576: /*@C
1577: SNESClearMonitor - Clears all the monitor functions for a SNES object.
1579: Collective on SNES
1581: Input Parameters:
1582: . snes - the SNES context
1584: Options Database:
1585: . -snes_cancelmonitors - cancels all monitors that have been hardwired
1586: into a code by calls to SNESSetMonitor(), but does not cancel those
1587: set via the options database
1589: Notes:
1590: There is no way to clear one specific monitor from a SNES object.
1592: Level: intermediate
1594: .keywords: SNES, nonlinear, set, monitor
1596: .seealso: SNESDefaultMonitor(), SNESSetMonitor()
1597: @*/
1598: int SNESClearMonitor(SNES snes)
1599: {
1602: snes->numbermonitors = 0;
1603: return(0);
1604: }
1606: /*@C
1607: SNESSetConvergenceTest - Sets the function that is to be used
1608: to test for convergence of the nonlinear iterative solution.
1610: Collective on SNES
1612: Input Parameters:
1613: + snes - the SNES context
1614: . func - routine to test for convergence
1615: - cctx - [optional] context for private data for the convergence routine
1616: (may be PETSC_NULL)
1618: Calling sequence of func:
1619: $ int func (SNES snes,PetscReal xnorm,PetscReal gnorm,PetscReal f,SNESConvergedReason *reason,void *cctx)
1621: + snes - the SNES context
1622: . cctx - [optional] convergence context
1623: . reason - reason for convergence/divergence
1624: . xnorm - 2-norm of current iterate
1625: . gnorm - 2-norm of current step (SNES_NONLINEAR_EQUATIONS methods)
1626: . f - 2-norm of function (SNES_NONLINEAR_EQUATIONS methods)
1627: . gnorm - 2-norm of current gradient (SNES_UNCONSTRAINED_MINIMIZATION methods)
1628: - f - function value (SNES_UNCONSTRAINED_MINIMIZATION methods)
1630: Level: advanced
1632: .keywords: SNES, nonlinear, set, convergence, test
1634: .seealso: SNESConverged_EQ_LS(), SNESConverged_EQ_TR(),
1635: SNESConverged_UM_LS(), SNESConverged_UM_TR()
1636: @*/
1637: int SNESSetConvergenceTest(SNES snes,int (*func)(SNES,PetscReal,PetscReal,PetscReal,SNESConvergedReason*,void*),void *cctx)
1638: {
1641: (snes)->converged = func;
1642: (snes)->cnvP = cctx;
1643: return(0);
1644: }
1646: /*@C
1647: SNESGetConvergedReason - Gets the reason the SNES iteration was stopped.
1649: Not Collective
1651: Input Parameter:
1652: . snes - the SNES context
1654: Output Parameter:
1655: . reason - negative value indicates diverged, positive value converged, see petscsnes.h or the
1656: manual pages for the individual convergence tests for complete lists
1658: Level: intermediate
1660: Notes: Can only be called after the call the SNESSolve() is complete.
1662: .keywords: SNES, nonlinear, set, convergence, test
1664: .seealso: SNESSetConvergenceTest(), SNESConverged_EQ_LS(), SNESConverged_EQ_TR(),
1665: SNESConverged_UM_LS(), SNESConverged_UM_TR(), SNESConvergedReason
1666: @*/
1667: int SNESGetConvergedReason(SNES snes,SNESConvergedReason *reason)
1668: {
1671: *reason = snes->reason;
1672: return(0);
1673: }
1675: /*@
1676: SNESSetConvergenceHistory - Sets the array used to hold the convergence history.
1678: Collective on SNES
1680: Input Parameters:
1681: + snes - iterative context obtained from SNESCreate()
1682: . a - array to hold history
1683: . its - integer array holds the number of linear iterations for each solve.
1684: . na - size of a and its
1685: - reset - PETSC_TRUE indicates each new nonlinear solve resets the history counter to zero,
1686: else it continues storing new values for new nonlinear solves after the old ones
1688: Notes:
1689: If set, this array will contain the function norms (for
1690: SNES_NONLINEAR_EQUATIONS methods) or gradient norms
1691: (for SNES_UNCONSTRAINED_MINIMIZATION methods) computed
1692: at each step.
1694: This routine is useful, e.g., when running a code for purposes
1695: of accurate performance monitoring, when no I/O should be done
1696: during the section of code that is being timed.
1698: Level: intermediate
1700: .keywords: SNES, set, convergence, history
1702: .seealso: SNESGetConvergenceHistory()
1704: @*/
1705: int SNESSetConvergenceHistory(SNES snes,PetscReal *a,int *its,int na,PetscTruth reset)
1706: {
1710: snes->conv_hist = a;
1711: snes->conv_hist_its = its;
1712: snes->conv_hist_max = na;
1713: snes->conv_hist_reset = reset;
1714: return(0);
1715: }
1717: /*@C
1718: SNESGetConvergenceHistory - Gets the array used to hold the convergence history.
1720: Collective on SNES
1722: Input Parameter:
1723: . snes - iterative context obtained from SNESCreate()
1725: Output Parameters:
1726: . a - array to hold history
1727: . its - integer array holds the number of linear iterations (or
1728: negative if not converged) for each solve.
1729: - na - size of a and its
1731: Notes:
1732: The calling sequence for this routine in Fortran is
1733: $ call SNESGetConvergenceHistory(SNES snes, integer na, integer ierr)
1735: This routine is useful, e.g., when running a code for purposes
1736: of accurate performance monitoring, when no I/O should be done
1737: during the section of code that is being timed.
1739: Level: intermediate
1741: .keywords: SNES, get, convergence, history
1743: .seealso: SNESSetConvergencHistory()
1745: @*/
1746: int SNESGetConvergenceHistory(SNES snes,PetscReal **a,int **its,int *na)
1747: {
1750: if (a) *a = snes->conv_hist;
1751: if (its) *its = snes->conv_hist_its;
1752: if (na) *na = snes->conv_hist_len;
1753: return(0);
1754: }
1756: /*
1757: SNESScaleStep_Private - Scales a step so that its length is less than the
1758: positive parameter delta.
1760: Input Parameters:
1761: + snes - the SNES context
1762: . y - approximate solution of linear system
1763: . fnorm - 2-norm of current function
1764: - delta - trust region size
1766: Output Parameters:
1767: + gpnorm - predicted function norm at the new point, assuming local
1768: linearization. The value is zero if the step lies within the trust
1769: region, and exceeds zero otherwise.
1770: - ynorm - 2-norm of the step
1772: Note:
1773: For non-trust region methods such as SNESEQLS, the parameter delta
1774: is set to be the maximum allowable step size.
1776: .keywords: SNES, nonlinear, scale, step
1777: */
1778: int SNESScaleStep_Private(SNES snes,Vec y,PetscReal *fnorm,PetscReal *delta,
1779: PetscReal *gpnorm,PetscReal *ynorm)
1780: {
1781: PetscReal norm;
1782: Scalar cnorm;
1783: int ierr;
1790: VecNorm(y,NORM_2,&norm);
1791: if (norm > *delta) {
1792: norm = *delta/norm;
1793: *gpnorm = (1.0 - norm)*(*fnorm);
1794: cnorm = norm;
1795: VecScale(&cnorm,y);
1796: *ynorm = *delta;
1797: } else {
1798: *gpnorm = 0.0;
1799: *ynorm = norm;
1800: }
1801: return(0);
1802: }
1804: /*@
1805: SNESSolve - Solves a nonlinear system. Call SNESSolve after calling
1806: SNESCreate() and optional routines of the form SNESSetXXX().
1808: Collective on SNES
1810: Input Parameters:
1811: + snes - the SNES context
1812: - x - the solution vector
1814: Output Parameter:
1815: . its - number of iterations until termination
1817: Notes:
1818: The user should initialize the vector,x, with the initial guess
1819: for the nonlinear solve prior to calling SNESSolve. In particular,
1820: to employ an initial guess of zero, the user should explicitly set
1821: this vector to zero by calling VecSet().
1823: Level: beginner
1825: .keywords: SNES, nonlinear, solve
1827: .seealso: SNESCreate(), SNESDestroy()
1828: @*/
1829: int SNESSolve(SNES snes,Vec x,int *its)
1830: {
1831: int ierr;
1832: PetscTruth flg;
1839: if (!snes->solve) SETERRQ(1,"SNESSetType() or SNESSetFromOptions() must be called before SNESSolve()");
1841: if (!snes->setupcalled) {SNESSetUp(snes,x);}
1842: else {snes->vec_sol = snes->vec_sol_always = x;}
1843: if (snes->conv_hist_reset == PETSC_TRUE) snes->conv_hist_len = 0;
1844: PetscLogEventBegin(SNES_Solve,snes,0,0,0);
1845: snes->nfuncs = 0; snes->linear_its = 0; snes->nfailures = 0;
1846: (*(snes)->solve)(snes,its);
1847: PetscLogEventEnd(SNES_Solve,snes,0,0,0);
1848: PetscOptionsHasName(snes->prefix,"-snes_view",&flg);
1849: if (flg && !PetscPreLoadingOn) { SNESView(snes,PETSC_VIEWER_STDOUT_WORLD); }
1850: return(0);
1851: }
1853: /* --------- Internal routines for SNES Package --------- */
1855: /*@C
1856: SNESSetType - Sets the method for the nonlinear solver.
1858: Collective on SNES
1860: Input Parameters:
1861: + snes - the SNES context
1862: - type - a known method
1864: Options Database Key:
1865: . -snes_type <type> - Sets the method; use -help for a list
1866: of available methods (for instance, ls or tr)
1868: Notes:
1869: See "petsc/include/petscsnes.h" for available methods (for instance)
1870: + SNESEQLS - Newton's method with line search
1871: (systems of nonlinear equations)
1872: . SNESEQTR - Newton's method with trust region
1873: (systems of nonlinear equations)
1874: . SNESUMTR - Newton's method with trust region
1875: (unconstrained minimization)
1876: - SNESUMLS - Newton's method with line search
1877: (unconstrained minimization)
1879: Normally, it is best to use the SNESSetFromOptions() command and then
1880: set the SNES solver type from the options database rather than by using
1881: this routine. Using the options database provides the user with
1882: maximum flexibility in evaluating the many nonlinear solvers.
1883: The SNESSetType() routine is provided for those situations where it
1884: is necessary to set the nonlinear solver independently of the command
1885: line or options database. This might be the case, for example, when
1886: the choice of solver changes during the execution of the program,
1887: and the user's application is taking responsibility for choosing the
1888: appropriate method.
1890: Level: intermediate
1892: .keywords: SNES, set, type
1894: .seealso: SNESType, SNESCreate()
1896: @*/
1897: int SNESSetType(SNES snes,SNESType type)
1898: {
1899: int ierr,(*r)(SNES);
1900: PetscTruth match;
1906: PetscTypeCompare((PetscObject)snes,type,&match);
1907: if (match) return(0);
1909: if (snes->setupcalled) {
1910: ierr = (*(snes)->destroy)(snes);
1911: snes->data = 0;
1912: }
1914: /* Get the function pointers for the iterative method requested */
1915: if (!SNESRegisterAllCalled) {SNESRegisterAll(PETSC_NULL);}
1917: PetscFListFind(snes->comm,SNESList,type,(void (**)(void)) &r);
1919: if (!r) SETERRQ1(1,"Unable to find requested SNES type %s",type);
1921: if (snes->data) {PetscFree(snes->data);}
1922: snes->data = 0;
1923: (*r)(snes);
1925: PetscObjectChangeTypeName((PetscObject)snes,type);
1926: snes->set_method_called = 1;
1928: return(0);
1929: }
1932: /* --------------------------------------------------------------------- */
1933: /*@C
1934: SNESRegisterDestroy - Frees the list of nonlinear solvers that were
1935: registered by SNESRegisterDynamic().
1937: Not Collective
1939: Level: advanced
1941: .keywords: SNES, nonlinear, register, destroy
1943: .seealso: SNESRegisterAll(), SNESRegisterAll()
1944: @*/
1945: int SNESRegisterDestroy(void)
1946: {
1950: if (SNESList) {
1951: PetscFListDestroy(&SNESList);
1952: SNESList = 0;
1953: }
1954: SNESRegisterAllCalled = PETSC_FALSE;
1955: return(0);
1956: }
1958: /*@C
1959: SNESGetType - Gets the SNES method type and name (as a string).
1961: Not Collective
1963: Input Parameter:
1964: . snes - nonlinear solver context
1966: Output Parameter:
1967: . type - SNES method (a charactor string)
1969: Level: intermediate
1971: .keywords: SNES, nonlinear, get, type, name
1972: @*/
1973: int SNESGetType(SNES snes,SNESType *type)
1974: {
1977: *type = snes->type_name;
1978: return(0);
1979: }
1981: /*@C
1982: SNESGetSolution - Returns the vector where the approximate solution is
1983: stored.
1985: Not Collective, but Vec is parallel if SNES is parallel
1987: Input Parameter:
1988: . snes - the SNES context
1990: Output Parameter:
1991: . x - the solution
1993: Level: advanced
1995: .keywords: SNES, nonlinear, get, solution
1997: .seealso: SNESGetFunction(), SNESGetGradient(), SNESGetSolutionUpdate()
1998: @*/
1999: int SNESGetSolution(SNES snes,Vec *x)
2000: {
2003: *x = snes->vec_sol_always;
2004: return(0);
2005: }
2007: /*@C
2008: SNESGetSolutionUpdate - Returns the vector where the solution update is
2009: stored.
2011: Not Collective, but Vec is parallel if SNES is parallel
2013: Input Parameter:
2014: . snes - the SNES context
2016: Output Parameter:
2017: . x - the solution update
2019: Level: advanced
2021: .keywords: SNES, nonlinear, get, solution, update
2023: .seealso: SNESGetSolution(), SNESGetFunction
2024: @*/
2025: int SNESGetSolutionUpdate(SNES snes,Vec *x)
2026: {
2029: *x = snes->vec_sol_update_always;
2030: return(0);
2031: }
2033: /*@C
2034: SNESGetFunction - Returns the vector where the function is stored.
2036: Not Collective, but Vec is parallel if SNES is parallel
2038: Input Parameter:
2039: . snes - the SNES context
2041: Output Parameter:
2042: + r - the function (or PETSC_NULL)
2043: . ctx - the function context (or PETSC_NULL)
2044: - func - the function (or PETSC_NULL)
2046: Notes:
2047: SNESGetFunction() is valid for SNES_NONLINEAR_EQUATIONS methods only
2048: Analogous routines for SNES_UNCONSTRAINED_MINIMIZATION methods are
2049: SNESGetMinimizationFunction() and SNESGetGradient();
2051: Level: advanced
2053: .keywords: SNES, nonlinear, get, function
2055: .seealso: SNESSetFunction(), SNESGetSolution(), SNESGetMinimizationFunction(),
2056: SNESGetGradient()
2058: @*/
2059: int SNESGetFunction(SNES snes,Vec *r,void **ctx,int (**func)(SNES,Vec,Vec,void*))
2060: {
2063: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
2064: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
2065: }
2066: if (r) *r = snes->vec_func_always;
2067: if (ctx) *ctx = snes->funP;
2068: if (func) *func = snes->computefunction;
2069: return(0);
2070: }
2072: /*@C
2073: SNESGetGradient - Returns the vector where the gradient is stored.
2075: Not Collective, but Vec is parallel if SNES is parallel
2077: Input Parameter:
2078: . snes - the SNES context
2080: Output Parameter:
2081: + r - the gradient (or PETSC_NULL)
2082: - ctx - the gradient context (or PETSC_NULL)
2084: Notes:
2085: SNESGetGradient() is valid for SNES_UNCONSTRAINED_MINIMIZATION methods
2086: only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is
2087: SNESGetFunction().
2089: Level: advanced
2091: .keywords: SNES, nonlinear, get, gradient
2093: .seealso: SNESGetMinimizationFunction(), SNESGetSolution(), SNESGetFunction(),
2094: SNESSetGradient(), SNESSetFunction()
2096: @*/
2097: int SNESGetGradient(SNES snes,Vec *r,void **ctx)
2098: {
2101: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
2102: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
2103: }
2104: if (r) *r = snes->vec_func_always;
2105: if (ctx) *ctx = snes->funP;
2106: return(0);
2107: }
2109: /*@C
2110: SNESGetMinimizationFunction - Returns the scalar function value for
2111: unconstrained minimization problems.
2113: Not Collective
2115: Input Parameter:
2116: . snes - the SNES context
2118: Output Parameter:
2119: + r - the function (or PETSC_NULL)
2120: - ctx - the function context (or PETSC_NULL)
2122: Notes:
2123: SNESGetMinimizationFunction() is valid for SNES_UNCONSTRAINED_MINIMIZATION
2124: methods only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is
2125: SNESGetFunction().
2127: Level: advanced
2129: .keywords: SNES, nonlinear, get, function
2131: .seealso: SNESGetGradient(), SNESGetSolution(), SNESGetFunction(), SNESSetFunction()
2133: @*/
2134: int SNESGetMinimizationFunction(SNES snes,PetscReal *r,void **ctx)
2135: {
2139: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
2140: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
2141: }
2142: if (r) *r = snes->fc;
2143: if (ctx) *ctx = snes->umfunP;
2144: return(0);
2145: }
2147: /*@C
2148: SNESSetOptionsPrefix - Sets the prefix used for searching for all
2149: SNES options in the database.
2151: Collective on SNES
2153: Input Parameter:
2154: + snes - the SNES context
2155: - prefix - the prefix to prepend to all option names
2157: Notes:
2158: A hyphen (-) must NOT be given at the beginning of the prefix name.
2159: The first character of all runtime options is AUTOMATICALLY the hyphen.
2161: Level: advanced
2163: .keywords: SNES, set, options, prefix, database
2165: .seealso: SNESSetFromOptions()
2166: @*/
2167: int SNESSetOptionsPrefix(SNES snes,char *prefix)
2168: {
2173: PetscObjectSetOptionsPrefix((PetscObject)snes,prefix);
2174: SLESSetOptionsPrefix(snes->sles,prefix);
2175: return(0);
2176: }
2178: /*@C
2179: SNESAppendOptionsPrefix - Appends to the prefix used for searching for all
2180: SNES options in the database.
2182: Collective on SNES
2184: Input Parameters:
2185: + snes - the SNES context
2186: - prefix - the prefix to prepend to all option names
2188: Notes:
2189: A hyphen (-) must NOT be given at the beginning of the prefix name.
2190: The first character of all runtime options is AUTOMATICALLY the hyphen.
2192: Level: advanced
2194: .keywords: SNES, append, options, prefix, database
2196: .seealso: SNESGetOptionsPrefix()
2197: @*/
2198: int SNESAppendOptionsPrefix(SNES snes,char *prefix)
2199: {
2201:
2204: PetscObjectAppendOptionsPrefix((PetscObject)snes,prefix);
2205: SLESAppendOptionsPrefix(snes->sles,prefix);
2206: return(0);
2207: }
2209: /*@C
2210: SNESGetOptionsPrefix - Sets the prefix used for searching for all
2211: SNES options in the database.
2213: Not Collective
2215: Input Parameter:
2216: . snes - the SNES context
2218: Output Parameter:
2219: . prefix - pointer to the prefix string used
2221: Notes: On the fortran side, the user should pass in a string 'prifix' of
2222: sufficient length to hold the prefix.
2224: Level: advanced
2226: .keywords: SNES, get, options, prefix, database
2228: .seealso: SNESAppendOptionsPrefix()
2229: @*/
2230: int SNESGetOptionsPrefix(SNES snes,char **prefix)
2231: {
2236: PetscObjectGetOptionsPrefix((PetscObject)snes,prefix);
2237: return(0);
2238: }
2240: /*MC
2241: SNESRegisterDynamic - Adds a method to the nonlinear solver package.
2243: Synopsis:
2244: int SNESRegisterDynamic(char *name_solver,char *path,char *name_create,int (*routine_create)(SNES))
2246: Not collective
2248: Input Parameters:
2249: + name_solver - name of a new user-defined solver
2250: . path - path (either absolute or relative) the library containing this solver
2251: . name_create - name of routine to create method context
2252: - routine_create - routine to create method context
2254: Notes:
2255: SNESRegisterDynamic() may be called multiple times to add several user-defined solvers.
2257: If dynamic libraries are used, then the fourth input argument (routine_create)
2258: is ignored.
2260: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR}, ${BOPT},
2261: and others of the form ${any_environmental_variable} occuring in pathname will be
2262: replaced with appropriate values.
2264: Sample usage:
2265: .vb
2266: SNESRegisterDynamic("my_solver",/home/username/my_lib/lib/libg/solaris/mylib.a,
2267: "MySolverCreate",MySolverCreate);
2268: .ve
2270: Then, your solver can be chosen with the procedural interface via
2271: $ SNESSetType(snes,"my_solver")
2272: or at runtime via the option
2273: $ -snes_type my_solver
2275: Level: advanced
2277: .keywords: SNES, nonlinear, register
2279: .seealso: SNESRegisterAll(), SNESRegisterDestroy()
2280: M*/
2282: int SNESRegister(char *sname,char *path,char *name,int (*function)(SNES))
2283: {
2284: char fullname[256];
2285: int ierr;
2288: PetscFListConcat(path,name,fullname);
2289: PetscFListAdd(&SNESList,sname,fullname,(void (*)())function);
2290: return(0);
2291: }