Actual source code: snes.c
1: /*$Id: snes.c,v 1.235 2001/08/21 21:03:49 bsmith Exp $*/
3: #include src/snes/snesimpl.h
5: PetscTruth SNESRegisterAllCalled = PETSC_FALSE;
6: PetscFList SNESList = PETSC_NULL;
8: /* Logging support */
9: int SNES_COOKIE;
10: int MATSNESMFCTX_COOKIE;
11: int SNES_Solve, SNES_LineSearch, SNES_FunctionEval, SNES_JacobianEval, SNES_MinimizationFunctionEval, SNES_GradientEval;
12: int SNES_HessianEval;
14: /*@C
15: SNESGetProblemType -Indicates if SNES is solving a nonlinear system or a minimization
17: Not Collective
19: Input Parameter:
20: . SNES - the SNES context
22: Output Parameter:
23: . type - SNES_NONLINEAR_EQUATIONS (for systems of nonlinear equations)
24: or SNES_UNCONSTRAINED_MINIMIZATION (for unconstrained minimization)
26: Level: intermediate
28: .keywords: SNES, problem type
30: .seealso: SNESCreate()
31: @*/
32: int SNESGetProblemType(SNES snes,SNESProblemType *type)
33: {
36: *type = snes->method_class;
37: return(0);
38: }
40: /*@C
41: SNESView - Prints the SNES data structure.
43: Collective on SNES
45: Input Parameters:
46: + SNES - the SNES context
47: - viewer - visualization context
49: Options Database Key:
50: . -snes_view - Calls SNESView() at end of SNESSolve()
52: Notes:
53: The available visualization contexts include
54: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
55: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
56: output where only the first processor opens
57: the file. All other processors send their
58: data to the first processor to print.
60: The user can open an alternative visualization context with
61: PetscViewerASCIIOpen() - output to a specified file.
63: Level: beginner
65: .keywords: SNES, view
67: .seealso: PetscViewerASCIIOpen()
68: @*/
69: int SNESView(SNES snes,PetscViewer viewer)
70: {
71: SNES_KSP_EW_ConvCtx *kctx;
72: int ierr;
73: SLES sles;
74: char *type;
75: PetscTruth isascii,isstring;
79: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(snes->comm);
83: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
84: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
85: if (isascii) {
86: if (snes->prefix) {
87: PetscViewerASCIIPrintf(viewer,"SNES Object:(%s)n",snes->prefix);
88: } else {
89: PetscViewerASCIIPrintf(viewer,"SNES Object:n");
90: }
91: SNESGetType(snes,&type);
92: if (type) {
93: PetscViewerASCIIPrintf(viewer," type: %sn",type);
94: } else {
95: PetscViewerASCIIPrintf(viewer," type: not set yetn");
96: }
97: if (snes->view) {
98: PetscViewerASCIIPushTab(viewer);
99: (*snes->view)(snes,viewer);
100: PetscViewerASCIIPopTab(viewer);
101: }
102: PetscViewerASCIIPrintf(viewer," maximum iterations=%d, maximum function evaluations=%dn",snes->max_its,snes->max_funcs);
103: PetscViewerASCIIPrintf(viewer," tolerances: relative=%g, absolute=%g, solution=%gn",
104: snes->rtol,snes->atol,snes->xtol);
105: PetscViewerASCIIPrintf(viewer," total number of linear solver iterations=%dn",snes->linear_its);
106: PetscViewerASCIIPrintf(viewer," total number of function evaluations=%dn",snes->nfuncs);
107: if (snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) {
108: PetscViewerASCIIPrintf(viewer," min function tolerance=%gn",snes->fmin);
109: }
110: if (snes->ksp_ewconv) {
111: kctx = (SNES_KSP_EW_ConvCtx *)snes->kspconvctx;
112: if (kctx) {
113: PetscViewerASCIIPrintf(viewer," Eisenstat-Walker computation of KSP relative tolerance (version %d)n",kctx->version);
114: PetscViewerASCIIPrintf(viewer," rtol_0=%g, rtol_max=%g, threshold=%gn",kctx->rtol_0,kctx->rtol_max,kctx->threshold);
115: PetscViewerASCIIPrintf(viewer," gamma=%g, alpha=%g, alpha2=%gn",kctx->gamma,kctx->alpha,kctx->alpha2);
116: }
117: }
118: } else if (isstring) {
119: SNESGetType(snes,&type);
120: PetscViewerStringSPrintf(viewer," %-3.3s",type);
121: }
122: SNESGetSLES(snes,&sles);
123: PetscViewerASCIIPushTab(viewer);
124: SLESView(sles,viewer);
125: PetscViewerASCIIPopTab(viewer);
126: return(0);
127: }
129: /*
130: We retain a list of functions that also take SNES command
131: line options. These are called at the end SNESSetFromOptions()
132: */
133: #define MAXSETFROMOPTIONS 5
134: static int numberofsetfromoptions;
135: static int (*othersetfromoptions[MAXSETFROMOPTIONS])(SNES);
137: /*@
138: SNESAddOptionsChecker - Adds an additional function to check for SNES options.
140: Not Collective
142: Input Parameter:
143: . snescheck - function that checks for options
145: Level: developer
147: .seealso: SNESSetFromOptions()
148: @*/
149: int SNESAddOptionsChecker(int (*snescheck)(SNES))
150: {
152: if (numberofsetfromoptions >= MAXSETFROMOPTIONS) {
153: SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE, "Too many options checkers, only %d allowed", MAXSETFROMOPTIONS);
154: }
156: othersetfromoptions[numberofsetfromoptions++] = snescheck;
157: return(0);
158: }
160: /*@
161: SNESSetFromOptions - Sets various SNES and SLES parameters from user options.
163: Collective on SNES
165: Input Parameter:
166: . snes - the SNES context
168: Options Database Keys:
169: + -snes_type <type> - ls, tr, umls, umtr, test
170: . -snes_stol - convergence tolerance in terms of the norm
171: of the change in the solution between steps
172: . -snes_atol <atol> - absolute tolerance of residual norm
173: . -snes_rtol <rtol> - relative decrease in tolerance norm from initial
174: . -snes_max_it <max_it> - maximum number of iterations
175: . -snes_max_funcs <max_funcs> - maximum number of function evaluations
176: . -snes_max_fail <max_fail> - maximum number of failures
177: . -snes_trtol <trtol> - trust region tolerance
178: . -snes_no_convergence_test - skip convergence test in nonlinear or minimization
179: solver; hence iterations will continue until max_it
180: or some other criterion is reached. Saves expense
181: of convergence test
182: . -snes_monitor - prints residual norm at each iteration
183: . -snes_vecmonitor - plots solution at each iteration
184: . -snes_vecmonitor_update - plots update to solution at each iteration
185: . -snes_xmonitor - plots residual norm at each iteration
186: . -snes_fd - use finite differences to compute Jacobian; very slow, only for testing
187: - -snes_mf_ksp_monitor - if using matrix-free multiply then print h at each KSP iteration
189: Options Database for Eisenstat-Walker method:
190: + -snes_ksp_eq_conv - use Eisenstat-Walker method for determining linear system convergence
191: . -snes_ksp_eq_version ver - version of Eisenstat-Walker method
192: . -snes_ksp_ew_rtol0 <rtol0> - Sets rtol0
193: . -snes_ksp_ew_rtolmax <rtolmax> - Sets rtolmax
194: . -snes_ksp_ew_gamma <gamma> - Sets gamma
195: . -snes_ksp_ew_alpha <alpha> - Sets alpha
196: . -snes_ksp_ew_alpha2 <alpha2> - Sets alpha2
197: - -snes_ksp_ew_threshold <threshold> - Sets threshold
199: Notes:
200: To see all options, run your program with the -help option or consult
201: the users manual.
203: Level: beginner
205: .keywords: SNES, nonlinear, set, options, database
207: .seealso: SNESSetOptionsPrefix()
208: @*/
209: int SNESSetFromOptions(SNES snes)
210: {
211: SLES sles;
212: SNES_KSP_EW_ConvCtx *kctx = (SNES_KSP_EW_ConvCtx *)snes->kspconvctx;
213: PetscTruth flg;
214: int ierr, i;
215: char *deft,type[256];
220: PetscOptionsBegin(snes->comm,snes->prefix,"Nonlinear solver (SNES) options","SNES");
221: if (snes->type_name) {
222: deft = snes->type_name;
223: } else {
224: if (snes->method_class == SNES_NONLINEAR_EQUATIONS) {
225: deft = SNESEQLS;
226: } else {
227: deft = SNESUMTR;
228: }
229: }
231: if (!SNESRegisterAllCalled) {SNESRegisterAll(PETSC_NULL);}
232: PetscOptionsList("-snes_type","Nonlinear solver method","SNESSetType",SNESList,deft,type,256,&flg);
233: if (flg) {
234: SNESSetType(snes,type);
235: } else if (!snes->type_name) {
236: SNESSetType(snes,deft);
237: }
239: PetscOptionsReal("-snes_stol","Stop if step length less then","SNESSetTolerances",snes->xtol,&snes->xtol,0);
240: PetscOptionsReal("-snes_atol","Stop if function norm less then","SNESSetTolerances",snes->atol,&snes->atol,0);
242: PetscOptionsReal("-snes_rtol","Stop if decrease in function norm less then","SNESSetTolerances",snes->rtol,&snes->rtol,0);
243: PetscOptionsInt("-snes_max_it","Maximum iterations","SNESSetTolerances",snes->max_its,&snes->max_its,PETSC_NULL);
244: PetscOptionsInt("-snes_max_funcs","Maximum function evaluations","SNESSetTolerances",snes->max_funcs,&snes->max_funcs,PETSC_NULL);
245: PetscOptionsInt("-snes_max_fail","Maximum failures","SNESSetTolerances",snes->maxFailures,&snes->maxFailures,PETSC_NULL);
246: PetscOptionsReal("-snes_fmin","Minimization function tolerance","SNESSetMinimizationFunctionTolerance",snes->fmin,&snes->fmin,0);
248: PetscOptionsName("-snes_ksp_ew_conv","Use Eisentat-Walker linear system convergence test","SNES_KSP_SetParametersEW",&snes->ksp_ewconv);
250: PetscOptionsInt("-snes_ksp_ew_version","Version 1 or 2","SNES_KSP_SetParametersEW",kctx->version,&kctx->version,0);
251: PetscOptionsReal("-snes_ksp_ew_rtol0","0 <= rtol0 < 1","SNES_KSP_SetParametersEW",kctx->rtol_0,&kctx->rtol_0,0);
252: PetscOptionsReal("-snes_ksp_ew_rtolmax","0 <= rtolmax < 1","SNES_KSP_SetParametersEW",kctx->rtol_max,&kctx->rtol_max,0);
253: PetscOptionsReal("-snes_ksp_ew_gamma","0 <= gamma <= 1","SNES_KSP_SetParametersEW",kctx->gamma,&kctx->gamma,0);
254: PetscOptionsReal("-snes_ksp_ew_alpha","1 < alpha <= 2","SNES_KSP_SetParametersEW",kctx->alpha,&kctx->alpha,0);
255: PetscOptionsReal("-snes_ksp_ew_alpha2","alpha2","SNES_KSP_SetParametersEW",kctx->alpha2,&kctx->alpha2,0);
256: PetscOptionsReal("-snes_ksp_ew_threshold","0 < threshold < 1","SNES_KSP_SetParametersEW",kctx->threshold,&kctx->threshold,0);
258: PetscOptionsName("-snes_no_convergence_test","Don't test for convergence","None",&flg);
259: if (flg) {snes->converged = 0;}
260: PetscOptionsName("-snes_cancelmonitors","Remove all monitors","SNESClearMonitor",&flg);
261: if (flg) {SNESClearMonitor(snes);}
262: PetscOptionsName("-snes_monitor","Monitor norm of function","SNESDefaultMonitor",&flg);
263: if (flg) {SNESSetMonitor(snes,SNESDefaultMonitor,0,0);}
264: PetscOptionsName("-snes_ratiomonitor","Monitor norm of function","SNESSetRatioMonitor",&flg);
265: if (flg) {SNESSetRatioMonitor(snes);}
266: PetscOptionsName("-snes_smonitor","Monitor norm of function (fewer digits)","SNESDefaultSMonitor",&flg);
267: if (flg) {SNESSetMonitor(snes,SNESDefaultSMonitor,0,0);}
268: PetscOptionsName("-snes_vecmonitor","Plot solution at each iteration","SNESVecViewMonitor",&flg);
269: if (flg) {SNESSetMonitor(snes,SNESVecViewMonitor,0,0);}
270: PetscOptionsName("-snes_vecmonitor_update","Plot correction at each iteration","SNESVecViewUpdateMonitor",&flg);
271: if (flg) {SNESSetMonitor(snes,SNESVecViewUpdateMonitor,0,0);}
272: PetscOptionsName("-snes_xmonitor","Plot function norm at each iteration","SNESLGMonitor",&flg);
273: if (flg) {SNESSetMonitor(snes,SNESLGMonitor,PETSC_NULL,PETSC_NULL);}
275: PetscOptionsName("-snes_fd","Use finite differences (slow) to compute Jacobian","SNESDefaultComputeJacobian",&flg);
276: if (flg && snes->method_class == SNES_NONLINEAR_EQUATIONS) {
277: SNESSetJacobian(snes,snes->jacobian,snes->jacobian_pre,SNESDefaultComputeJacobian,snes->funP);
278: PetscLogInfo(snes,"SNESSetFromOptions: Setting default finite difference Jacobian matrixn");
279: } else if (flg && snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) {
280: SNESSetHessian(snes,snes->jacobian,snes->jacobian_pre,SNESDefaultComputeHessian,snes->funP);
281: PetscLogInfo(snes,"SNESSetFromOptions: Setting default finite difference Hessian matrixn");
282: }
284: for(i = 0; i < numberofsetfromoptions; i++) {
285: (*othersetfromoptions[i])(snes);
286: }
288: if (snes->setfromoptions) {
289: (*snes->setfromoptions)(snes);
290: }
292: PetscOptionsEnd();
294: SNESGetSLES(snes,&sles);
295: SLESSetFromOptions(sles);
297: return(0);
298: }
301: /*@
302: SNESSetApplicationContext - Sets the optional user-defined context for
303: the nonlinear solvers.
305: Collective on SNES
307: Input Parameters:
308: + snes - the SNES context
309: - usrP - optional user context
311: Level: intermediate
313: .keywords: SNES, nonlinear, set, application, context
315: .seealso: SNESGetApplicationContext()
316: @*/
317: int SNESSetApplicationContext(SNES snes,void *usrP)
318: {
321: snes->user = usrP;
322: return(0);
323: }
325: /*@C
326: SNESGetApplicationContext - Gets the user-defined context for the
327: nonlinear solvers.
329: Not Collective
331: Input Parameter:
332: . snes - SNES context
334: Output Parameter:
335: . usrP - user context
337: Level: intermediate
339: .keywords: SNES, nonlinear, get, application, context
341: .seealso: SNESSetApplicationContext()
342: @*/
343: int SNESGetApplicationContext(SNES snes,void **usrP)
344: {
347: *usrP = snes->user;
348: return(0);
349: }
351: /*@
352: SNESGetIterationNumber - Gets the number of nonlinear iterations completed
353: at this time.
355: Not Collective
357: Input Parameter:
358: . snes - SNES context
360: Output Parameter:
361: . iter - iteration number
363: Notes:
364: For example, during the computation of iteration 2 this would return 1.
366: This is useful for using lagged Jacobians (where one does not recompute the
367: Jacobian at each SNES iteration). For example, the code
368: .vb
369: SNESGetIterationNumber(snes,&it);
370: if (!(it % 2)) {
371: [compute Jacobian here]
372: }
373: .ve
374: can be used in your ComputeJacobian() function to cause the Jacobian to be
375: recomputed every second SNES iteration.
377: Level: intermediate
379: .keywords: SNES, nonlinear, get, iteration, number
380: @*/
381: int SNESGetIterationNumber(SNES snes,int* iter)
382: {
386: *iter = snes->iter;
387: return(0);
388: }
390: /*@
391: SNESGetFunctionNorm - Gets the norm of the current function that was set
392: with SNESSSetFunction().
394: Collective on SNES
396: Input Parameter:
397: . snes - SNES context
399: Output Parameter:
400: . fnorm - 2-norm of function
402: Note:
403: SNESGetFunctionNorm() is valid for SNES_NONLINEAR_EQUATIONS methods only.
404: A related routine for SNES_UNCONSTRAINED_MINIMIZATION methods is
405: SNESGetGradientNorm().
407: Level: intermediate
409: .keywords: SNES, nonlinear, get, function, norm
411: .seealso: SNESGetFunction()
412: @*/
413: int SNESGetFunctionNorm(SNES snes,PetscScalar *fnorm)
414: {
418: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
419: SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"For SNES_NONLINEAR_EQUATIONS only");
420: }
421: *fnorm = snes->norm;
422: return(0);
423: }
425: /*@
426: SNESGetGradientNorm - Gets the norm of the current gradient that was set
427: with SNESSSetGradient().
429: Collective on SNES
431: Input Parameter:
432: . snes - SNES context
434: Output Parameter:
435: . fnorm - 2-norm of gradient
437: Note:
438: SNESGetGradientNorm() is valid for SNES_UNCONSTRAINED_MINIMIZATION
439: methods only. A related routine for SNES_NONLINEAR_EQUATIONS methods
440: is SNESGetFunctionNorm().
442: Level: intermediate
444: .keywords: SNES, nonlinear, get, gradient, norm
446: .seelso: SNESSetGradient()
447: @*/
448: int SNESGetGradientNorm(SNES snes,PetscScalar *gnorm)
449: {
453: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
454: SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"For SNES_UNCONSTRAINED_MINIMIZATION only");
455: }
456: *gnorm = snes->norm;
457: return(0);
458: }
460: /*@
461: SNESGetNumberUnsuccessfulSteps - Gets the number of unsuccessful steps
462: attempted by the nonlinear solver.
464: Not Collective
466: Input Parameter:
467: . snes - SNES context
469: Output Parameter:
470: . nfails - number of unsuccessful steps attempted
472: Notes:
473: This counter is reset to zero for each successive call to SNESSolve().
475: Level: intermediate
477: .keywords: SNES, nonlinear, get, number, unsuccessful, steps
478: @*/
479: int SNESGetNumberUnsuccessfulSteps(SNES snes,int* nfails)
480: {
484: *nfails = snes->numFailures;
485: return(0);
486: }
488: /*@
489: SNESSetMaximumUnsuccessfulSteps - Sets the maximum number of unsuccessful steps
490: attempted by the nonlinear solver before it gives up.
492: Not Collective
494: Input Parameters:
495: + snes - SNES context
496: - maxFails - maximum of unsuccessful steps
498: Level: intermediate
500: .keywords: SNES, nonlinear, set, maximum, unsuccessful, steps
501: @*/
502: int SNESSetMaximumUnsuccessfulSteps(SNES snes, int maxFails)
503: {
506: snes->maxFailures = maxFails;
507: return(0);
508: }
510: /*@
511: SNESGetMaximumUnsuccessfulSteps - Gets the maximum number of unsuccessful steps
512: attempted by the nonlinear solver before it gives up.
514: Not Collective
516: Input Parameter:
517: . snes - SNES context
519: Output Parameter:
520: . maxFails - maximum of unsuccessful steps
522: Level: intermediate
524: .keywords: SNES, nonlinear, get, maximum, unsuccessful, steps
525: @*/
526: int SNESGetMaximumUnsuccessfulSteps(SNES snes, int *maxFails)
527: {
531: *maxFails = snes->maxFailures;
532: return(0);
533: }
535: /*@
536: SNESGetNumberLinearIterations - Gets the total number of linear iterations
537: used by the nonlinear solver.
539: Not Collective
541: Input Parameter:
542: . snes - SNES context
544: Output Parameter:
545: . lits - number of linear iterations
547: Notes:
548: This counter is reset to zero for each successive call to SNESSolve().
550: Level: intermediate
552: .keywords: SNES, nonlinear, get, number, linear, iterations
553: @*/
554: int SNESGetNumberLinearIterations(SNES snes,int* lits)
555: {
559: *lits = snes->linear_its;
560: return(0);
561: }
563: /*@C
564: SNESGetSLES - Returns the SLES context for a SNES solver.
566: Not Collective, but if SNES object is parallel, then SLES object is parallel
568: Input Parameter:
569: . snes - the SNES context
571: Output Parameter:
572: . sles - the SLES context
574: Notes:
575: The user can then directly manipulate the SLES context to set various
576: options, etc. Likewise, the user can then extract and manipulate the
577: KSP and PC contexts as well.
579: Level: beginner
581: .keywords: SNES, nonlinear, get, SLES, context
583: .seealso: SLESGetPC(), SLESGetKSP()
584: @*/
585: int SNESGetSLES(SNES snes,SLES *sles)
586: {
589: *sles = snes->sles;
590: return(0);
591: }
593: static int SNESPublish_Petsc(PetscObject obj)
594: {
595: #if defined(PETSC_HAVE_AMS)
596: SNES v = (SNES) obj;
597: int ierr;
598: #endif
602: #if defined(PETSC_HAVE_AMS)
603: /* if it is already published then return */
604: if (v->amem >=0) return(0);
606: PetscObjectPublishBaseBegin(obj);
607: AMS_Memory_add_field((AMS_Memory)v->amem,"Iteration",&v->iter,1,AMS_INT,AMS_READ,
608: AMS_COMMON,AMS_REDUCT_UNDEF);
609: AMS_Memory_add_field((AMS_Memory)v->amem,"Residual",&v->norm,1,AMS_DOUBLE,AMS_READ,
610: AMS_COMMON,AMS_REDUCT_UNDEF);
611: PetscObjectPublishBaseEnd(obj);
612: #endif
613: return(0);
614: }
616: /* -----------------------------------------------------------*/
617: /*@C
618: SNESCreate - Creates a nonlinear solver context.
620: Collective on MPI_Comm
622: Input Parameters:
623: + comm - MPI communicator
624: - type - type of method, either
625: SNES_NONLINEAR_EQUATIONS (for systems of nonlinear equations)
626: or SNES_UNCONSTRAINED_MINIMIZATION (for unconstrained minimization)
628: Output Parameter:
629: . outsnes - the new SNES context
631: Options Database Keys:
632: + -snes_mf - Activates default matrix-free Jacobian-vector products,
633: and no preconditioning matrix
634: . -snes_mf_operator - Activates default matrix-free Jacobian-vector
635: products, and a user-provided preconditioning matrix
636: as set by SNESSetJacobian()
637: - -snes_fd - Uses (slow!) finite differences to compute Jacobian
639: Level: beginner
641: .keywords: SNES, nonlinear, create, context
643: .seealso: SNESSolve(), SNESDestroy(), SNESProblemType, SNES
644: @*/
645: int SNESCreate(MPI_Comm comm,SNESProblemType type,SNES *outsnes)
646: {
647: int ierr;
648: SNES snes;
649: SNES_KSP_EW_ConvCtx *kctx;
653: *outsnes = PETSC_NULL;
654: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
655: SNESInitializePackage(PETSC_NULL);
656: #endif
658: if (type != SNES_UNCONSTRAINED_MINIMIZATION && type != SNES_NONLINEAR_EQUATIONS){
659: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"incorrect method type");
660: }
661: PetscHeaderCreate(snes,_p_SNES,int,SNES_COOKIE,0,"SNES",comm,SNESDestroy,SNESView);
662: PetscLogObjectCreate(snes);
663: snes->bops->publish = SNESPublish_Petsc;
664: snes->max_its = 50;
665: snes->max_funcs = 10000;
666: snes->norm = 0.0;
667: if (type == SNES_UNCONSTRAINED_MINIMIZATION) {
668: snes->rtol = 1.e-8;
669: snes->ttol = 0.0;
670: snes->atol = 1.e-10;
671: } else {
672: snes->rtol = 1.e-8;
673: snes->ttol = 0.0;
674: snes->atol = 1.e-50;
675: }
676: snes->xtol = 1.e-8;
677: snes->nfuncs = 0;
678: snes->numFailures = 0;
679: snes->maxFailures = 1;
680: snes->linear_its = 0;
681: snes->numbermonitors = 0;
682: snes->data = 0;
683: snes->view = 0;
684: snes->computeumfunction = 0;
685: snes->umfunP = 0;
686: snes->fc = 0;
687: snes->deltatol = 1.e-12;
688: snes->fmin = -1.e30;
689: snes->method_class = type;
690: snes->set_method_called = 0;
691: snes->setupcalled = 0;
692: snes->ksp_ewconv = PETSC_FALSE;
693: snes->vwork = 0;
694: snes->nwork = 0;
695: snes->conv_hist_len = 0;
696: snes->conv_hist_max = 0;
697: snes->conv_hist = PETSC_NULL;
698: snes->conv_hist_its = PETSC_NULL;
699: snes->conv_hist_reset = PETSC_TRUE;
700: snes->reason = SNES_CONVERGED_ITERATING;
702: /* Create context to compute Eisenstat-Walker relative tolerance for KSP */
703: PetscNew(SNES_KSP_EW_ConvCtx,&kctx);
704: PetscLogObjectMemory(snes,sizeof(SNES_KSP_EW_ConvCtx));
705: snes->kspconvctx = (void*)kctx;
706: kctx->version = 2;
707: kctx->rtol_0 = .3; /* Eisenstat and Walker suggest rtol_0=.5, but
708: this was too large for some test cases */
709: kctx->rtol_last = 0;
710: kctx->rtol_max = .9;
711: kctx->gamma = 1.0;
712: kctx->alpha2 = .5*(1.0 + sqrt(5.0));
713: kctx->alpha = kctx->alpha2;
714: kctx->threshold = .1;
715: kctx->lresid_last = 0;
716: kctx->norm_last = 0;
718: SLESCreate(comm,&snes->sles);
719: PetscLogObjectParent(snes,snes->sles)
721: *outsnes = snes;
722: PetscPublishAll(snes);
723: return(0);
724: }
726: /* --------------------------------------------------------------- */
727: /*@C
728: SNESSetFunction - Sets the function evaluation routine and function
729: vector for use by the SNES routines in solving systems of nonlinear
730: equations.
732: Collective on SNES
734: Input Parameters:
735: + snes - the SNES context
736: . func - function evaluation routine
737: . r - vector to store function value
738: - ctx - [optional] user-defined context for private data for the
739: function evaluation routine (may be PETSC_NULL)
741: Calling sequence of func:
742: $ func (SNES snes,Vec x,Vec f,void *ctx);
744: . f - function vector
745: - ctx - optional user-defined function context
747: Notes:
748: The Newton-like methods typically solve linear systems of the form
749: $ f'(x) x = -f(x),
750: where f'(x) denotes the Jacobian matrix and f(x) is the function.
752: SNESSetFunction() is valid for SNES_NONLINEAR_EQUATIONS methods only.
753: Analogous routines for SNES_UNCONSTRAINED_MINIMIZATION methods are
754: SNESSetMinimizationFunction() and SNESSetGradient();
756: Level: beginner
758: .keywords: SNES, nonlinear, set, function
760: .seealso: SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian()
761: @*/
762: int SNESSetFunction(SNES snes,Vec r,int (*func)(SNES,Vec,Vec,void*),void *ctx)
763: {
768: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
769: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
770: }
772: snes->computefunction = func;
773: snes->vec_func = snes->vec_func_always = r;
774: snes->funP = ctx;
775: return(0);
776: }
778: /*@
779: SNESComputeFunction - Calls the function that has been set with
780: SNESSetFunction().
782: Collective on SNES
784: Input Parameters:
785: + snes - the SNES context
786: - x - input vector
788: Output Parameter:
789: . y - function vector, as set by SNESSetFunction()
791: Notes:
792: SNESComputeFunction() is valid for SNES_NONLINEAR_EQUATIONS methods only.
793: Analogous routines for SNES_UNCONSTRAINED_MINIMIZATION methods are
794: SNESComputeMinimizationFunction() and SNESComputeGradient();
796: SNESComputeFunction() is typically used within nonlinear solvers
797: implementations, so most users would not generally call this routine
798: themselves.
800: Level: developer
802: .keywords: SNES, nonlinear, compute, function
804: .seealso: SNESSetFunction(), SNESGetFunction()
805: @*/
806: int SNESComputeFunction(SNES snes,Vec x,Vec y)
807: {
808: int ierr;
816: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
817: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
818: }
820: PetscLogEventBegin(SNES_FunctionEval,snes,x,y,0);
821: PetscStackPush("SNES user function");
822: (*snes->computefunction)(snes,x,y,snes->funP);
823: PetscStackPop;
824: snes->nfuncs++;
825: PetscLogEventEnd(SNES_FunctionEval,snes,x,y,0);
826: return(0);
827: }
829: /*@C
830: SNESSetMinimizationFunction - Sets the function evaluation routine for
831: unconstrained minimization.
833: Collective on SNES
835: Input Parameters:
836: + snes - the SNES context
837: . func - function evaluation routine
838: - ctx - [optional] user-defined context for private data for the
839: function evaluation routine (may be PETSC_NULL)
841: Calling sequence of func:
842: $ func (SNES snes,Vec x,PetscReal *f,void *ctx);
844: + x - input vector
845: . f - function
846: - ctx - [optional] user-defined function context
848: Level: beginner
850: Notes:
851: SNESSetMinimizationFunction() is valid for SNES_UNCONSTRAINED_MINIMIZATION
852: methods only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is
853: SNESSetFunction().
855: .keywords: SNES, nonlinear, set, minimization, function
857: .seealso: SNESGetMinimizationFunction(), SNESComputeMinimizationFunction(),
858: SNESSetHessian(), SNESSetGradient()
859: @*/
860: int SNESSetMinimizationFunction(SNES snes,int (*func)(SNES,Vec,PetscReal*,void*),void *ctx)
861: {
864: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
865: SETERRQ(PETSC_ERR_ARG_WRONG,"Only for SNES_UNCONSTRAINED_MINIMIZATION");
866: }
867: snes->computeumfunction = func;
868: snes->umfunP = ctx;
869: return(0);
870: }
872: /*@
873: SNESComputeMinimizationFunction - Computes the function that has been
874: set with SNESSetMinimizationFunction().
876: Collective on SNES
878: Input Parameters:
879: + snes - the SNES context
880: - x - input vector
882: Output Parameter:
883: . y - function value
885: Notes:
886: SNESComputeMinimizationFunction() is valid only for
887: SNES_UNCONSTRAINED_MINIMIZATION methods. An analogous routine for
888: SNES_NONLINEAR_EQUATIONS methods is SNESComputeFunction().
890: SNESComputeMinimizationFunction() is typically used within minimization
891: implementations, so most users would not generally call this routine
892: themselves.
894: Level: developer
896: .keywords: SNES, nonlinear, compute, minimization, function
898: .seealso: SNESSetMinimizationFunction(), SNESGetMinimizationFunction(),
899: SNESComputeGradient(), SNESComputeHessian()
900: @*/
901: int SNESComputeMinimizationFunction(SNES snes,Vec x,PetscReal *y)
902: {
903: int ierr;
909: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
910: SETERRQ(PETSC_ERR_ARG_WRONG,"Only for SNES_UNCONSTRAINED_MINIMIZATION");
911: }
913: PetscLogEventBegin(SNES_MinimizationFunctionEval,snes,x,y,0);
914: PetscStackPush("SNES user minimzation function");
915: (*snes->computeumfunction)(snes,x,y,snes->umfunP);
916: PetscStackPop;
917: snes->nfuncs++;
918: PetscLogEventEnd(SNES_MinimizationFunctionEval,snes,x,y,0);
919: return(0);
920: }
922: /*@C
923: SNESSetGradient - Sets the gradient evaluation routine and gradient
924: vector for use by the SNES routines.
926: Collective on SNES
928: Input Parameters:
929: + snes - the SNES context
930: . func - function evaluation routine
931: . ctx - optional user-defined context for private data for the
932: gradient evaluation routine (may be PETSC_NULL)
933: - r - vector to store gradient value
935: Calling sequence of func:
936: $ func (SNES, Vec x, Vec g, void *ctx);
938: + x - input vector
939: . g - gradient vector
940: - ctx - optional user-defined gradient context
942: Notes:
943: SNESSetMinimizationFunction() is valid for SNES_UNCONSTRAINED_MINIMIZATION
944: methods only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is
945: SNESSetFunction().
947: Level: beginner
949: .keywords: SNES, nonlinear, set, function
951: .seealso: SNESGetGradient(), SNESComputeGradient(), SNESSetHessian(),
952: SNESSetMinimizationFunction(),
953: @*/
954: int SNESSetGradient(SNES snes,Vec r,int (*func)(SNES,Vec,Vec,void*),void *ctx)
955: {
960: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
961: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
962: }
963: snes->computefunction = func;
964: snes->vec_func = snes->vec_func_always = r;
965: snes->funP = ctx;
966: return(0);
967: }
969: /*@
970: SNESComputeGradient - Computes the gradient that has been set with
971: SNESSetGradient().
973: Collective on SNES
975: Input Parameters:
976: + snes - the SNES context
977: - x - input vector
979: Output Parameter:
980: . y - gradient vector
982: Notes:
983: SNESComputeGradient() is valid only for
984: SNES_UNCONSTRAINED_MINIMIZATION methods. An analogous routine for
985: SNES_NONLINEAR_EQUATIONS methods is SNESComputeFunction().
987: SNESComputeGradient() is typically used within minimization
988: implementations, so most users would not generally call this routine
989: themselves.
991: Level: developer
993: .keywords: SNES, nonlinear, compute, gradient
995: .seealso: SNESSetGradient(), SNESGetGradient(),
996: SNESComputeMinimizationFunction(), SNESComputeHessian()
997: @*/
998: int SNESComputeGradient(SNES snes,Vec x,Vec y)
999: {
1000: int ierr;
1008: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
1009: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
1010: }
1011: PetscLogEventBegin(SNES_GradientEval,snes,x,y,0);
1012: PetscStackPush("SNES user gradient function");
1013: (*snes->computefunction)(snes,x,y,snes->funP);
1014: PetscStackPop;
1015: PetscLogEventEnd(SNES_GradientEval,snes,x,y,0);
1016: return(0);
1017: }
1019: /*@
1020: SNESComputeJacobian - Computes the Jacobian matrix that has been
1021: set with SNESSetJacobian().
1023: Collective on SNES and Mat
1025: Input Parameters:
1026: + snes - the SNES context
1027: - x - input vector
1029: Output Parameters:
1030: + A - Jacobian matrix
1031: . B - optional preconditioning matrix
1032: - flag - flag indicating matrix structure
1034: Notes:
1035: Most users should not need to explicitly call this routine, as it
1036: is used internally within the nonlinear solvers.
1038: See SLESSetOperators() for important information about setting the
1039: flag parameter.
1041: SNESComputeJacobian() is valid only for SNES_NONLINEAR_EQUATIONS
1042: methods. An analogous routine for SNES_UNCONSTRAINED_MINIMIZATION
1043: methods is SNESComputeHessian().
1045: Level: developer
1047: .keywords: SNES, compute, Jacobian, matrix
1049: .seealso: SNESSetJacobian(), SLESSetOperators()
1050: @*/
1051: int SNESComputeJacobian(SNES snes,Vec X,Mat *A,Mat *B,MatStructure *flg)
1052: {
1053: int ierr;
1059: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
1060: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
1061: }
1062: if (!snes->computejacobian) return(0);
1063: PetscLogEventBegin(SNES_JacobianEval,snes,X,*A,*B);
1064: *flg = DIFFERENT_NONZERO_PATTERN;
1065: PetscStackPush("SNES user Jacobian function");
1066: (*snes->computejacobian)(snes,X,A,B,flg,snes->jacP);
1067: PetscStackPop;
1068: PetscLogEventEnd(SNES_JacobianEval,snes,X,*A,*B);
1069: /* make sure user returned a correct Jacobian and preconditioner */
1072: return(0);
1073: }
1075: /*@
1076: SNESComputeHessian - Computes the Hessian matrix that has been
1077: set with SNESSetHessian().
1079: Collective on SNES and Mat
1081: Input Parameters:
1082: + snes - the SNES context
1083: - x - input vector
1085: Output Parameters:
1086: + A - Hessian matrix
1087: . B - optional preconditioning matrix
1088: - flag - flag indicating matrix structure
1090: Notes:
1091: Most users should not need to explicitly call this routine, as it
1092: is used internally within the nonlinear solvers.
1094: See SLESSetOperators() for important information about setting the
1095: flag parameter.
1097: SNESComputeHessian() is valid only for
1098: SNES_UNCONSTRAINED_MINIMIZATION methods. An analogous routine for
1099: SNES_NONLINEAR_EQUATIONS methods is SNESComputeJacobian().
1101: SNESComputeHessian() is typically used within minimization
1102: implementations, so most users would not generally call this routine
1103: themselves.
1105: Level: developer
1107: .keywords: SNES, compute, Hessian, matrix
1109: .seealso: SNESSetHessian(), SLESSetOperators(), SNESComputeGradient(),
1110: SNESComputeMinimizationFunction()
1111: @*/
1112: int SNESComputeHessian(SNES snes,Vec x,Mat *A,Mat *B,MatStructure *flag)
1113: {
1114: int ierr;
1120: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
1121: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
1122: }
1123: if (!snes->computejacobian) return(0);
1124: PetscLogEventBegin(SNES_HessianEval,snes,x,*A,*B);
1125: *flag = DIFFERENT_NONZERO_PATTERN;
1126: PetscStackPush("SNES user Hessian function");
1127: (*snes->computejacobian)(snes,x,A,B,flag,snes->jacP);
1128: PetscStackPop;
1129: PetscLogEventEnd(SNES_HessianEval,snes,x,*A,*B);
1130: /* make sure user returned a correct Jacobian and preconditioner */
1133: return(0);
1134: }
1136: /*@C
1137: SNESSetJacobian - Sets the function to compute Jacobian as well as the
1138: location to store the matrix.
1140: Collective on SNES and Mat
1142: Input Parameters:
1143: + snes - the SNES context
1144: . A - Jacobian matrix
1145: . B - preconditioner matrix (usually same as the Jacobian)
1146: . func - Jacobian evaluation routine
1147: - ctx - [optional] user-defined context for private data for the
1148: Jacobian evaluation routine (may be PETSC_NULL)
1150: Calling sequence of func:
1151: $ func (SNES snes,Vec x,Mat *A,Mat *B,int *flag,void *ctx);
1153: + x - input vector
1154: . A - Jacobian matrix
1155: . B - preconditioner matrix, usually the same as A
1156: . flag - flag indicating information about the preconditioner matrix
1157: structure (same as flag in SLESSetOperators())
1158: - ctx - [optional] user-defined Jacobian context
1160: Notes:
1161: See SLESSetOperators() for important information about setting the flag
1162: output parameter in the routine func(). Be sure to read this information!
1164: The routine func() takes Mat * as the matrix arguments rather than Mat.
1165: This allows the Jacobian evaluation routine to replace A and/or B with a
1166: completely new new matrix structure (not just different matrix elements)
1167: when appropriate, for instance, if the nonzero structure is changing
1168: throughout the global iterations.
1170: Level: beginner
1172: .keywords: SNES, nonlinear, set, Jacobian, matrix
1174: .seealso: SLESSetOperators(), SNESSetFunction()
1175: @*/
1176: int SNESSetJacobian(SNES snes,Mat A,Mat B,int (*func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
1177: {
1186: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
1187: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
1188: }
1190: if (func) snes->computejacobian = func;
1191: if (ctx) snes->jacP = ctx;
1192: if (A) {
1193: if (snes->jacobian) {MatDestroy(snes->jacobian);}
1194: snes->jacobian = A;
1195: ierr = PetscObjectReference((PetscObject)A);
1196: }
1197: if (B) {
1198: if (snes->jacobian_pre) {MatDestroy(snes->jacobian_pre);}
1199: snes->jacobian_pre = B;
1200: ierr = PetscObjectReference((PetscObject)B);
1201: }
1202: return(0);
1203: }
1205: /*@C
1206: SNESGetJacobian - Returns the Jacobian matrix and optionally the user
1207: provided context for evaluating the Jacobian.
1209: Not Collective, but Mat object will be parallel if SNES object is
1211: Input Parameter:
1212: . snes - the nonlinear solver context
1214: Output Parameters:
1215: + A - location to stash Jacobian matrix (or PETSC_NULL)
1216: . B - location to stash preconditioner matrix (or PETSC_NULL)
1217: . ctx - location to stash Jacobian ctx (or PETSC_NULL)
1218: - func - location to put Jacobian function (or PETSC_NULL)
1220: Level: advanced
1222: .seealso: SNESSetJacobian(), SNESComputeJacobian()
1223: @*/
1224: int SNESGetJacobian(SNES snes,Mat *A,Mat *B,void **ctx,int (**func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*))
1225: {
1228: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
1229: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
1230: }
1231: if (A) *A = snes->jacobian;
1232: if (B) *B = snes->jacobian_pre;
1233: if (ctx) *ctx = snes->jacP;
1234: if (func) *func = snes->computejacobian;
1235: return(0);
1236: }
1238: /*@C
1239: SNESSetHessian - Sets the function to compute Hessian as well as the
1240: location to store the matrix.
1242: Collective on SNES and Mat
1244: Input Parameters:
1245: + snes - the SNES context
1246: . A - Hessian matrix
1247: . B - preconditioner matrix (usually same as the Hessian)
1248: . func - Jacobian evaluation routine
1249: - ctx - [optional] user-defined context for private data for the
1250: Hessian evaluation routine (may be PETSC_NULL)
1252: Calling sequence of func:
1253: $ func (SNES snes,Vec x,Mat *A,Mat *B,int *flag,void *ctx);
1255: + x - input vector
1256: . A - Hessian matrix
1257: . B - preconditioner matrix, usually the same as A
1258: . flag - flag indicating information about the preconditioner matrix
1259: structure (same as flag in SLESSetOperators())
1260: - ctx - [optional] user-defined Hessian context
1262: Notes:
1263: See SLESSetOperators() for important information about setting the flag
1264: output parameter in the routine func(). Be sure to read this information!
1266: The function func() takes Mat * as the matrix arguments rather than Mat.
1267: This allows the Hessian evaluation routine to replace A and/or B with a
1268: completely new new matrix structure (not just different matrix elements)
1269: when appropriate, for instance, if the nonzero structure is changing
1270: throughout the global iterations.
1272: Level: beginner
1274: .keywords: SNES, nonlinear, set, Hessian, matrix
1276: .seealso: SNESSetMinimizationFunction(), SNESSetGradient(), SLESSetOperators()
1277: @*/
1278: int SNESSetHessian(SNES snes,Mat A,Mat B,int (*func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
1279: {
1288: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
1289: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
1290: }
1291: if (func) snes->computejacobian = func;
1292: if (ctx) snes->jacP = ctx;
1293: if (A) {
1294: if (snes->jacobian) {MatDestroy(snes->jacobian);}
1295: snes->jacobian = A;
1296: ierr = PetscObjectReference((PetscObject)A);
1297: }
1298: if (B) {
1299: if (snes->jacobian_pre) {MatDestroy(snes->jacobian_pre);}
1300: snes->jacobian_pre = B;
1301: ierr = PetscObjectReference((PetscObject)B);
1302: }
1303: return(0);
1304: }
1306: /*@
1307: SNESGetHessian - Returns the Hessian matrix and optionally the user
1308: provided context for evaluating the Hessian.
1310: Not Collective, but Mat object is parallel if SNES object is parallel
1312: Input Parameter:
1313: . snes - the nonlinear solver context
1315: Output Parameters:
1316: + A - location to stash Hessian matrix (or PETSC_NULL)
1317: . B - location to stash preconditioner matrix (or PETSC_NULL)
1318: - ctx - location to stash Hessian ctx (or PETSC_NULL)
1320: Level: advanced
1322: .seealso: SNESSetHessian(), SNESComputeHessian()
1324: .keywords: SNES, get, Hessian
1325: @*/
1326: int SNESGetHessian(SNES snes,Mat *A,Mat *B,void **ctx)
1327: {
1330: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION){
1331: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
1332: }
1333: if (A) *A = snes->jacobian;
1334: if (B) *B = snes->jacobian_pre;
1335: if (ctx) *ctx = snes->jacP;
1336: return(0);
1337: }
1339: /* ----- Routines to initialize and destroy a nonlinear solver ---- */
1341: /*@
1342: SNESSetUp - Sets up the internal data structures for the later use
1343: of a nonlinear solver.
1345: Collective on SNES
1347: Input Parameters:
1348: + snes - the SNES context
1349: - x - the solution vector
1351: Notes:
1352: For basic use of the SNES solvers the user need not explicitly call
1353: SNESSetUp(), since these actions will automatically occur during
1354: the call to SNESSolve(). However, if one wishes to control this
1355: phase separately, SNESSetUp() should be called after SNESCreate()
1356: and optional routines of the form SNESSetXXX(), but before SNESSolve().
1358: Level: advanced
1360: .keywords: SNES, nonlinear, setup
1362: .seealso: SNESCreate(), SNESSolve(), SNESDestroy()
1363: @*/
1364: int SNESSetUp(SNES snes,Vec x)
1365: {
1366: int ierr;
1367: PetscTruth flg;
1373: snes->vec_sol = snes->vec_sol_always = x;
1375: PetscOptionsHasName(snes->prefix,"-snes_mf_operator",&flg);
1376: /*
1377: This version replaces the user provided Jacobian matrix with a
1378: matrix-free version but still employs the user-provided preconditioner matrix
1379: */
1380: if (flg) {
1381: Mat J;
1382: MatCreateSNESMF(snes,snes->vec_sol,&J);
1383: MatSNESMFSetFromOptions(J);
1384: PetscLogInfo(snes,"SNESSetUp: Setting default matrix-free operator routinesn");
1385: SNESSetJacobian(snes,J,0,0,0);
1386: MatDestroy(J);
1387: }
1388: PetscOptionsHasName(snes->prefix,"-snes_mf",&flg);
1389: /*
1390: This version replaces both the user-provided Jacobian and the user-
1391: provided preconditioner matrix with the default matrix free version.
1392: */
1393: if (flg) {
1394: Mat J;
1395: SLES sles;
1396: PC pc;
1398: MatCreateSNESMF(snes,snes->vec_sol,&J);
1399: MatSNESMFSetFromOptions(J);
1400: PetscLogInfo(snes,"SNESSetUp: Setting default matrix-free operator and preconditioner routinesn");
1401: if (snes->method_class == SNES_NONLINEAR_EQUATIONS) {
1402: SNESSetJacobian(snes,J,J,MatSNESMFComputeJacobian,snes->funP);
1403: } else if (snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) {
1404: SNESSetHessian(snes,J,J,MatSNESMFComputeJacobian,snes->funP);
1405: } else {
1406: SETERRQ(PETSC_ERR_SUP,"Method class doesn't support matrix-free option");
1407: }
1408: MatDestroy(J);
1410: /* force no preconditioner */
1411: SNESGetSLES(snes,&sles);
1412: SLESGetPC(sles,&pc);
1413: PCSetType(pc,PCNONE);
1414: }
1416: if ((snes->method_class == SNES_NONLINEAR_EQUATIONS)) {
1417: PetscTruth iseqtr;
1419: if (!snes->vec_func) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetFunction() first");
1420: if (!snes->computefunction) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetFunction() first");
1421: if (!snes->jacobian) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetJacobian() first n or use -snes_mf option");
1422: if (snes->vec_func == snes->vec_sol) {
1423: SETERRQ(PETSC_ERR_ARG_IDN,"Solution vector cannot be function vector");
1424: }
1426: /* Set the KSP stopping criterion to use the Eisenstat-Walker method */
1427: PetscTypeCompare((PetscObject)snes,SNESEQTR,&iseqtr);
1428: if (snes->ksp_ewconv && !iseqtr) {
1429: SLES sles; KSP ksp;
1430: SNESGetSLES(snes,&sles);
1431: SLESGetKSP(sles,&ksp);
1432: KSPSetConvergenceTest(ksp,SNES_KSP_EW_Converged_Private,snes);
1433: }
1434: } else if ((snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION)) {
1435: if (!snes->vec_func) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetGradient() first");
1436: if (!snes->computefunction) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetGradient() first");
1437: if (!snes->computeumfunction) {
1438: SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetMinimizationFunction() first");
1439: }
1440: if (!snes->jacobian) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetHessian()");
1441: } else {
1442: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Unknown method class");
1443: }
1444: if (snes->setup) {(*snes->setup)(snes);}
1445: snes->setupcalled = 1;
1446: return(0);
1447: }
1449: /*@C
1450: SNESDestroy - Destroys the nonlinear solver context that was created
1451: with SNESCreate().
1453: Collective on SNES
1455: Input Parameter:
1456: . snes - the SNES context
1458: Level: beginner
1460: .keywords: SNES, nonlinear, destroy
1462: .seealso: SNESCreate(), SNESSolve()
1463: @*/
1464: int SNESDestroy(SNES snes)
1465: {
1466: int i,ierr;
1470: if (--snes->refct > 0) return(0);
1472: /* if memory was published with AMS then destroy it */
1473: PetscObjectDepublish(snes);
1475: if (snes->destroy) {(*(snes)->destroy)(snes);}
1476: if (snes->kspconvctx) {PetscFree(snes->kspconvctx);}
1477: if (snes->jacobian) {MatDestroy(snes->jacobian);}
1478: if (snes->jacobian_pre) {MatDestroy(snes->jacobian_pre);}
1479: SLESDestroy(snes->sles);
1480: if (snes->vwork) {VecDestroyVecs(snes->vwork,snes->nvwork);}
1481: for (i=0; i<snes->numbermonitors; i++) {
1482: if (snes->monitordestroy[i]) {
1483: (*snes->monitordestroy[i])(snes->monitorcontext[i]);
1484: }
1485: }
1486: PetscLogObjectDestroy((PetscObject)snes);
1487: PetscHeaderDestroy((PetscObject)snes);
1488: return(0);
1489: }
1491: /* ----------- Routines to set solver parameters ---------- */
1493: /*@
1494: SNESSetTolerances - Sets various parameters used in convergence tests.
1496: Collective on SNES
1498: Input Parameters:
1499: + snes - the SNES context
1500: . atol - absolute convergence tolerance
1501: . rtol - relative convergence tolerance
1502: . stol - convergence tolerance in terms of the norm
1503: of the change in the solution between steps
1504: . maxit - maximum number of iterations
1505: - maxf - maximum number of function evaluations
1507: Options Database Keys:
1508: + -snes_atol <atol> - Sets atol
1509: . -snes_rtol <rtol> - Sets rtol
1510: . -snes_stol <stol> - Sets stol
1511: . -snes_max_it <maxit> - Sets maxit
1512: - -snes_max_funcs <maxf> - Sets maxf
1514: Notes:
1515: The default maximum number of iterations is 50.
1516: The default maximum number of function evaluations is 1000.
1518: Level: intermediate
1520: .keywords: SNES, nonlinear, set, convergence, tolerances
1522: .seealso: SNESSetTrustRegionTolerance(), SNESSetMinimizationFunctionTolerance()
1523: @*/
1524: int SNESSetTolerances(SNES snes,PetscReal atol,PetscReal rtol,PetscReal stol,int maxit,int maxf)
1525: {
1528: if (atol != PETSC_DEFAULT) snes->atol = atol;
1529: if (rtol != PETSC_DEFAULT) snes->rtol = rtol;
1530: if (stol != PETSC_DEFAULT) snes->xtol = stol;
1531: if (maxit != PETSC_DEFAULT) snes->max_its = maxit;
1532: if (maxf != PETSC_DEFAULT) snes->max_funcs = maxf;
1533: return(0);
1534: }
1536: /*@
1537: SNESGetTolerances - Gets various parameters used in convergence tests.
1539: Not Collective
1541: Input Parameters:
1542: + snes - the SNES context
1543: . atol - absolute convergence tolerance
1544: . rtol - relative convergence tolerance
1545: . stol - convergence tolerance in terms of the norm
1546: of the change in the solution between steps
1547: . maxit - maximum number of iterations
1548: - maxf - maximum number of function evaluations
1550: Notes:
1551: The user can specify PETSC_NULL for any parameter that is not needed.
1553: Level: intermediate
1555: .keywords: SNES, nonlinear, get, convergence, tolerances
1557: .seealso: SNESSetTolerances()
1558: @*/
1559: int SNESGetTolerances(SNES snes,PetscReal *atol,PetscReal *rtol,PetscReal *stol,int *maxit,int *maxf)
1560: {
1563: if (atol) *atol = snes->atol;
1564: if (rtol) *rtol = snes->rtol;
1565: if (stol) *stol = snes->xtol;
1566: if (maxit) *maxit = snes->max_its;
1567: if (maxf) *maxf = snes->max_funcs;
1568: return(0);
1569: }
1571: /*@
1572: SNESSetTrustRegionTolerance - Sets the trust region parameter tolerance.
1574: Collective on SNES
1576: Input Parameters:
1577: + snes - the SNES context
1578: - tol - tolerance
1579:
1580: Options Database Key:
1581: . -snes_trtol <tol> - Sets tol
1583: Level: intermediate
1585: .keywords: SNES, nonlinear, set, trust region, tolerance
1587: .seealso: SNESSetTolerances(), SNESSetMinimizationFunctionTolerance()
1588: @*/
1589: int SNESSetTrustRegionTolerance(SNES snes,PetscReal tol)
1590: {
1593: snes->deltatol = tol;
1594: return(0);
1595: }
1597: /*@
1598: SNESSetMinimizationFunctionTolerance - Sets the minimum allowable function tolerance
1599: for unconstrained minimization solvers.
1600:
1601: Collective on SNES
1603: Input Parameters:
1604: + snes - the SNES context
1605: - ftol - minimum function tolerance
1607: Options Database Key:
1608: . -snes_fmin <ftol> - Sets ftol
1610: Note:
1611: SNESSetMinimizationFunctionTolerance() is valid for SNES_UNCONSTRAINED_MINIMIZATION
1612: methods only.
1614: Level: intermediate
1616: .keywords: SNES, nonlinear, set, minimum, convergence, function, tolerance
1618: .seealso: SNESSetTolerances(), SNESSetTrustRegionTolerance()
1619: @*/
1620: int SNESSetMinimizationFunctionTolerance(SNES snes,PetscReal ftol)
1621: {
1624: snes->fmin = ftol;
1625: return(0);
1626: }
1627: /*
1628: Duplicate the lg monitors for SNES from KSP; for some reason with
1629: dynamic libraries things don't work under Sun4 if we just use
1630: macros instead of functions
1631: */
1632: int SNESLGMonitor(SNES snes,int it,PetscReal norm,void *ctx)
1633: {
1638: KSPLGMonitor((KSP)snes,it,norm,ctx);
1639: return(0);
1640: }
1642: int SNESLGMonitorCreate(char *host,char *label,int x,int y,int m,int n,PetscDrawLG *draw)
1643: {
1647: KSPLGMonitorCreate(host,label,x,y,m,n,draw);
1648: return(0);
1649: }
1651: int SNESLGMonitorDestroy(PetscDrawLG draw)
1652: {
1656: KSPLGMonitorDestroy(draw);
1657: return(0);
1658: }
1660: /* ------------ Routines to set performance monitoring options ----------- */
1662: /*@C
1663: SNESSetMonitor - Sets an ADDITIONAL function that is to be used at every
1664: iteration of the nonlinear solver to display the iteration's
1665: progress.
1667: Collective on SNES
1669: Input Parameters:
1670: + snes - the SNES context
1671: . func - monitoring routine
1672: . mctx - [optional] user-defined context for private data for the
1673: monitor routine (use PETSC_NULL if no context is desitre)
1674: - monitordestroy - [optional] routine that frees monitor context
1675: (may be PETSC_NULL)
1677: Calling sequence of func:
1678: $ int func(SNES snes,int its, PetscReal norm,void *mctx)
1680: + snes - the SNES context
1681: . its - iteration number
1682: . norm - 2-norm function value (may be estimated)
1683: for SNES_NONLINEAR_EQUATIONS methods
1684: . norm - 2-norm gradient value (may be estimated)
1685: for SNES_UNCONSTRAINED_MINIMIZATION methods
1686: - mctx - [optional] monitoring context
1688: Options Database Keys:
1689: + -snes_monitor - sets SNESDefaultMonitor()
1690: . -snes_xmonitor - sets line graph monitor,
1691: uses SNESLGMonitorCreate()
1692: _ -snes_cancelmonitors - cancels all monitors that have
1693: been hardwired into a code by
1694: calls to SNESSetMonitor(), but
1695: does not cancel those set via
1696: the options database.
1698: Notes:
1699: Several different monitoring routines may be set by calling
1700: SNESSetMonitor() multiple times; all will be called in the
1701: order in which they were set.
1703: Level: intermediate
1705: .keywords: SNES, nonlinear, set, monitor
1707: .seealso: SNESDefaultMonitor(), SNESClearMonitor()
1708: @*/
1709: int SNESSetMonitor(SNES snes,int (*func)(SNES,int,PetscReal,void*),void *mctx,int (*monitordestroy)(void *))
1710: {
1713: if (snes->numbermonitors >= MAXSNESMONITORS) {
1714: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1715: }
1717: snes->monitor[snes->numbermonitors] = func;
1718: snes->monitordestroy[snes->numbermonitors] = monitordestroy;
1719: snes->monitorcontext[snes->numbermonitors++] = (void*)mctx;
1720: return(0);
1721: }
1723: /*@C
1724: SNESClearMonitor - Clears all the monitor functions for a SNES object.
1726: Collective on SNES
1728: Input Parameters:
1729: . snes - the SNES context
1731: Options Database:
1732: . -snes_cancelmonitors - cancels all monitors that have been hardwired
1733: into a code by calls to SNESSetMonitor(), but does not cancel those
1734: set via the options database
1736: Notes:
1737: There is no way to clear one specific monitor from a SNES object.
1739: Level: intermediate
1741: .keywords: SNES, nonlinear, set, monitor
1743: .seealso: SNESDefaultMonitor(), SNESSetMonitor()
1744: @*/
1745: int SNESClearMonitor(SNES snes)
1746: {
1749: snes->numbermonitors = 0;
1750: return(0);
1751: }
1753: /*@C
1754: SNESSetConvergenceTest - Sets the function that is to be used
1755: to test for convergence of the nonlinear iterative solution.
1757: Collective on SNES
1759: Input Parameters:
1760: + snes - the SNES context
1761: . func - routine to test for convergence
1762: - cctx - [optional] context for private data for the convergence routine
1763: (may be PETSC_NULL)
1765: Calling sequence of func:
1766: $ int func (SNES snes,PetscReal xnorm,PetscReal gnorm,PetscReal f,SNESConvergedReason *reason,void *cctx)
1768: + snes - the SNES context
1769: . cctx - [optional] convergence context
1770: . reason - reason for convergence/divergence
1771: . xnorm - 2-norm of current iterate
1772: . gnorm - 2-norm of current step (SNES_NONLINEAR_EQUATIONS methods)
1773: . f - 2-norm of function (SNES_NONLINEAR_EQUATIONS methods)
1774: . gnorm - 2-norm of current gradient (SNES_UNCONSTRAINED_MINIMIZATION methods)
1775: - f - function value (SNES_UNCONSTRAINED_MINIMIZATION methods)
1777: Level: advanced
1779: .keywords: SNES, nonlinear, set, convergence, test
1781: .seealso: SNESConverged_EQ_LS(), SNESConverged_EQ_TR(),
1782: SNESConverged_UM_LS(), SNESConverged_UM_TR()
1783: @*/
1784: int SNESSetConvergenceTest(SNES snes,int (*func)(SNES,PetscReal,PetscReal,PetscReal,SNESConvergedReason*,void*),void *cctx)
1785: {
1788: (snes)->converged = func;
1789: (snes)->cnvP = cctx;
1790: return(0);
1791: }
1793: /*@C
1794: SNESGetConvergedReason - Gets the reason the SNES iteration was stopped.
1796: Not Collective
1798: Input Parameter:
1799: . snes - the SNES context
1801: Output Parameter:
1802: . reason - negative value indicates diverged, positive value converged, see petscsnes.h or the
1803: manual pages for the individual convergence tests for complete lists
1805: Level: intermediate
1807: Notes: Can only be called after the call the SNESSolve() is complete.
1809: .keywords: SNES, nonlinear, set, convergence, test
1811: .seealso: SNESSetConvergenceTest(), SNESConverged_EQ_LS(), SNESConverged_EQ_TR(),
1812: SNESConverged_UM_LS(), SNESConverged_UM_TR(), SNESConvergedReason
1813: @*/
1814: int SNESGetConvergedReason(SNES snes,SNESConvergedReason *reason)
1815: {
1818: *reason = snes->reason;
1819: return(0);
1820: }
1822: /*@
1823: SNESSetConvergenceHistory - Sets the array used to hold the convergence history.
1825: Collective on SNES
1827: Input Parameters:
1828: + snes - iterative context obtained from SNESCreate()
1829: . a - array to hold history
1830: . its - integer array holds the number of linear iterations for each solve.
1831: . na - size of a and its
1832: - reset - PETSC_TRUE indicates each new nonlinear solve resets the history counter to zero,
1833: else it continues storing new values for new nonlinear solves after the old ones
1835: Notes:
1836: If set, this array will contain the function norms (for
1837: SNES_NONLINEAR_EQUATIONS methods) or gradient norms
1838: (for SNES_UNCONSTRAINED_MINIMIZATION methods) computed
1839: at each step.
1841: This routine is useful, e.g., when running a code for purposes
1842: of accurate performance monitoring, when no I/O should be done
1843: during the section of code that is being timed.
1845: Level: intermediate
1847: .keywords: SNES, set, convergence, history
1849: .seealso: SNESGetConvergenceHistory()
1851: @*/
1852: int SNESSetConvergenceHistory(SNES snes,PetscReal *a,int *its,int na,PetscTruth reset)
1853: {
1857: snes->conv_hist = a;
1858: snes->conv_hist_its = its;
1859: snes->conv_hist_max = na;
1860: snes->conv_hist_reset = reset;
1861: return(0);
1862: }
1864: /*@C
1865: SNESGetConvergenceHistory - Gets the array used to hold the convergence history.
1867: Collective on SNES
1869: Input Parameter:
1870: . snes - iterative context obtained from SNESCreate()
1872: Output Parameters:
1873: . a - array to hold history
1874: . its - integer array holds the number of linear iterations (or
1875: negative if not converged) for each solve.
1876: - na - size of a and its
1878: Notes:
1879: The calling sequence for this routine in Fortran is
1880: $ call SNESGetConvergenceHistory(SNES snes, integer na, integer ierr)
1882: This routine is useful, e.g., when running a code for purposes
1883: of accurate performance monitoring, when no I/O should be done
1884: during the section of code that is being timed.
1886: Level: intermediate
1888: .keywords: SNES, get, convergence, history
1890: .seealso: SNESSetConvergencHistory()
1892: @*/
1893: int SNESGetConvergenceHistory(SNES snes,PetscReal **a,int **its,int *na)
1894: {
1897: if (a) *a = snes->conv_hist;
1898: if (its) *its = snes->conv_hist_its;
1899: if (na) *na = snes->conv_hist_len;
1900: return(0);
1901: }
1903: /*@
1904: SNESSetRhsBC - Sets the function which applies boundary conditions
1905: to the Rhs of each system.
1907: Collective on SNES
1909: Input Parameters:
1910: . snes - The nonlinear solver context
1911: . func - The function
1913: Calling sequence of func:
1914: . func (SNES snes, Vec rhs, void *ctx);
1916: . rhs - The current rhs vector
1917: . ctx - The user-context
1919: Level: intermediate
1921: .keywords: SNES, Rhs, boundary conditions
1922: .seealso SNESDefaultRhsBC(), SNESSetSolutionBC(), SNESSetUpdate()
1923: @*/
1924: int SNESSetRhsBC(SNES snes, int (*func)(SNES, Vec, void *))
1925: {
1928: snes->applyrhsbc = func;
1929: return(0);
1930: }
1932: /*@
1933: SNESDefaultRhsBC - The default boundary condition function which does nothing.
1935: Not collective
1937: Input Parameters:
1938: . snes - The nonlinear solver context
1939: . rhs - The Rhs
1940: . ctx - The user-context
1942: Level: beginner
1944: .keywords: SNES, Rhs, boundary conditions
1945: .seealso SNESSetRhsBC(), SNESDefaultSolutionBC(), SNESDefaultUpdate()
1946: @*/
1947: int SNESDefaultRhsBC(SNES snes, Vec rhs, void *ctx)
1948: {
1950: return(0);
1951: }
1953: /*@
1954: SNESSetSolutionBC - Sets the function which applies boundary conditions
1955: to the solution of each system.
1957: Collective on SNES
1959: Input Parameters:
1960: . snes - The nonlinear solver context
1961: . func - The function
1963: Calling sequence of func:
1964: . func (TS ts, Vec rsol, void *ctx);
1966: . sol - The current solution vector
1967: . ctx - The user-context
1969: Level: intermediate
1971: .keywords: SNES, solution, boundary conditions
1972: .seealso SNESDefaultSolutionBC(), SNESSetRhsBC(), SNESSetUpdate()
1973: @*/
1974: int SNESSetSolutionBC(SNES snes, int (*func)(SNES, Vec, void *))
1975: {
1978: snes->applysolbc = func;
1979: return(0);
1980: }
1982: /*@
1983: SNESDefaultSolutionBC - The default boundary condition function which does nothing.
1985: Not collective
1987: Input Parameters:
1988: . snes - The nonlinear solver context
1989: . sol - The solution
1990: . ctx - The user-context
1992: Level: beginner
1994: .keywords: SNES, solution, boundary conditions
1995: .seealso SNESSetSolutionBC(), SNESDefaultRhsBC(), SNESDefaultUpdate()
1996: @*/
1997: int SNESDefaultSolutionBC(SNES snes, Vec sol, void *ctx)
1998: {
2000: return(0);
2001: }
2003: /*@
2004: SNESSetUpdate - Sets the general-purpose update function called
2005: at the beginning of every step of the iteration.
2007: Collective on SNES
2009: Input Parameters:
2010: . snes - The nonlinear solver context
2011: . func - The function
2013: Calling sequence of func:
2014: . func (TS ts, int step);
2016: . step - The current step of the iteration
2018: Level: intermediate
2020: .keywords: SNES, update
2021: .seealso SNESDefaultUpdate(), SNESSetRhsBC(), SNESSetSolutionBC()
2022: @*/
2023: int SNESSetUpdate(SNES snes, int (*func)(SNES, int))
2024: {
2027: snes->update = func;
2028: return(0);
2029: }
2031: /*@
2032: SNESDefaultUpdate - The default update function which does nothing.
2034: Not collective
2036: Input Parameters:
2037: . snes - The nonlinear solver context
2038: . step - The current step of the iteration
2040: Level: intermediate
2042: .keywords: SNES, update
2043: .seealso SNESSetUpdate(), SNESDefaultRhsBC(), SNESDefaultSolutionBC()
2044: @*/
2045: int SNESDefaultUpdate(SNES snes, int step)
2046: {
2048: return(0);
2049: }
2051: /*
2052: SNESScaleStep_Private - Scales a step so that its length is less than the
2053: positive parameter delta.
2055: Input Parameters:
2056: + snes - the SNES context
2057: . y - approximate solution of linear system
2058: . fnorm - 2-norm of current function
2059: - delta - trust region size
2061: Output Parameters:
2062: + gpnorm - predicted function norm at the new point, assuming local
2063: linearization. The value is zero if the step lies within the trust
2064: region, and exceeds zero otherwise.
2065: - ynorm - 2-norm of the step
2067: Note:
2068: For non-trust region methods such as SNESEQLS, the parameter delta
2069: is set to be the maximum allowable step size.
2071: .keywords: SNES, nonlinear, scale, step
2072: */
2073: int SNESScaleStep_Private(SNES snes,Vec y,PetscReal *fnorm,PetscReal *delta,PetscReal *gpnorm,PetscReal *ynorm)
2074: {
2075: PetscReal nrm;
2076: PetscScalar cnorm;
2077: int ierr;
2084: VecNorm(y,NORM_2,&nrm);
2085: if (nrm > *delta) {
2086: nrm = *delta/nrm;
2087: *gpnorm = (1.0 - nrm)*(*fnorm);
2088: cnorm = nrm;
2089: VecScale(&cnorm,y);
2090: *ynorm = *delta;
2091: } else {
2092: *gpnorm = 0.0;
2093: *ynorm = nrm;
2094: }
2095: return(0);
2096: }
2098: /*@
2099: SNESSolve - Solves a nonlinear system. Call SNESSolve after calling
2100: SNESCreate() and optional routines of the form SNESSetXXX().
2102: Collective on SNES
2104: Input Parameters:
2105: + snes - the SNES context
2106: - x - the solution vector
2108: Output Parameter:
2109: . its - number of iterations until termination
2111: Notes:
2112: The user should initialize the vector,x, with the initial guess
2113: for the nonlinear solve prior to calling SNESSolve. In particular,
2114: to employ an initial guess of zero, the user should explicitly set
2115: this vector to zero by calling VecSet().
2117: Level: beginner
2119: .keywords: SNES, nonlinear, solve
2121: .seealso: SNESCreate(), SNESDestroy()
2122: @*/
2123: int SNESSolve(SNES snes,Vec x,int *its)
2124: {
2125: int ierr;
2126: PetscTruth flg;
2133: if (!snes->solve) SETERRQ(1,"SNESSetType() or SNESSetFromOptions() must be called before SNESSolve()");
2135: if (!snes->setupcalled) {SNESSetUp(snes,x);}
2136: else {snes->vec_sol = snes->vec_sol_always = x;}
2137: if (snes->conv_hist_reset == PETSC_TRUE) snes->conv_hist_len = 0;
2138: PetscLogEventBegin(SNES_Solve,snes,0,0,0);
2139: snes->nfuncs = 0; snes->linear_its = 0; snes->numFailures = 0;
2140: (*(snes)->solve)(snes,its);
2141: PetscLogEventEnd(SNES_Solve,snes,0,0,0);
2142: PetscOptionsHasName(snes->prefix,"-snes_view",&flg);
2143: if (flg && !PetscPreLoadingOn) { SNESView(snes,PETSC_VIEWER_STDOUT_WORLD); }
2144: return(0);
2145: }
2147: /* --------- Internal routines for SNES Package --------- */
2149: /*@C
2150: SNESSetType - Sets the method for the nonlinear solver.
2152: Collective on SNES
2154: Input Parameters:
2155: + snes - the SNES context
2156: - type - a known method
2158: Options Database Key:
2159: . -snes_type <type> - Sets the method; use -help for a list
2160: of available methods (for instance, ls or tr)
2162: Notes:
2163: See "petsc/include/petscsnes.h" for available methods (for instance)
2164: + SNESEQLS - Newton's method with line search
2165: (systems of nonlinear equations)
2166: . SNESEQTR - Newton's method with trust region
2167: (systems of nonlinear equations)
2168: . SNESUMTR - Newton's method with trust region
2169: (unconstrained minimization)
2170: - SNESUMLS - Newton's method with line search
2171: (unconstrained minimization)
2173: Normally, it is best to use the SNESSetFromOptions() command and then
2174: set the SNES solver type from the options database rather than by using
2175: this routine. Using the options database provides the user with
2176: maximum flexibility in evaluating the many nonlinear solvers.
2177: The SNESSetType() routine is provided for those situations where it
2178: is necessary to set the nonlinear solver independently of the command
2179: line or options database. This might be the case, for example, when
2180: the choice of solver changes during the execution of the program,
2181: and the user's application is taking responsibility for choosing the
2182: appropriate method.
2184: Level: intermediate
2186: .keywords: SNES, set, type
2188: .seealso: SNESType, SNESCreate()
2190: @*/
2191: int SNESSetType(SNES snes,SNESType type)
2192: {
2193: int ierr,(*r)(SNES);
2194: PetscTruth match;
2200: PetscTypeCompare((PetscObject)snes,type,&match);
2201: if (match) return(0);
2203: if (snes->setupcalled) {
2204: ierr = (*(snes)->destroy)(snes);
2205: snes->data = 0;
2206: }
2208: /* Get the function pointers for the iterative method requested */
2209: if (!SNESRegisterAllCalled) {SNESRegisterAll(PETSC_NULL);}
2211: PetscFListFind(snes->comm,SNESList,type,(void (**)(void)) &r);
2213: if (!r) SETERRQ1(1,"Unable to find requested SNES type %s",type);
2215: if (snes->data) {PetscFree(snes->data);}
2216: snes->data = 0;
2217: (*r)(snes);
2219: PetscObjectChangeTypeName((PetscObject)snes,type);
2220: snes->set_method_called = 1;
2222: return(0);
2223: }
2226: /* --------------------------------------------------------------------- */
2227: /*@C
2228: SNESRegisterDestroy - Frees the list of nonlinear solvers that were
2229: registered by SNESRegisterDynamic().
2231: Not Collective
2233: Level: advanced
2235: .keywords: SNES, nonlinear, register, destroy
2237: .seealso: SNESRegisterAll(), SNESRegisterAll()
2238: @*/
2239: int SNESRegisterDestroy(void)
2240: {
2244: if (SNESList) {
2245: PetscFListDestroy(&SNESList);
2246: SNESList = 0;
2247: }
2248: SNESRegisterAllCalled = PETSC_FALSE;
2249: return(0);
2250: }
2252: /*@C
2253: SNESGetType - Gets the SNES method type and name (as a string).
2255: Not Collective
2257: Input Parameter:
2258: . snes - nonlinear solver context
2260: Output Parameter:
2261: . type - SNES method (a character string)
2263: Level: intermediate
2265: .keywords: SNES, nonlinear, get, type, name
2266: @*/
2267: int SNESGetType(SNES snes,SNESType *type)
2268: {
2271: *type = snes->type_name;
2272: return(0);
2273: }
2275: /*@C
2276: SNESGetSolution - Returns the vector where the approximate solution is
2277: stored.
2279: Not Collective, but Vec is parallel if SNES is parallel
2281: Input Parameter:
2282: . snes - the SNES context
2284: Output Parameter:
2285: . x - the solution
2287: Level: advanced
2289: .keywords: SNES, nonlinear, get, solution
2291: .seealso: SNESGetFunction(), SNESGetGradient(), SNESGetSolutionUpdate()
2292: @*/
2293: int SNESGetSolution(SNES snes,Vec *x)
2294: {
2297: *x = snes->vec_sol_always;
2298: return(0);
2299: }
2301: /*@C
2302: SNESGetSolutionUpdate - Returns the vector where the solution update is
2303: stored.
2305: Not Collective, but Vec is parallel if SNES is parallel
2307: Input Parameter:
2308: . snes - the SNES context
2310: Output Parameter:
2311: . x - the solution update
2313: Level: advanced
2315: .keywords: SNES, nonlinear, get, solution, update
2317: .seealso: SNESGetSolution(), SNESGetFunction
2318: @*/
2319: int SNESGetSolutionUpdate(SNES snes,Vec *x)
2320: {
2323: *x = snes->vec_sol_update_always;
2324: return(0);
2325: }
2327: /*@C
2328: SNESGetFunction - Returns the vector where the function is stored.
2330: Not Collective, but Vec is parallel if SNES is parallel
2332: Input Parameter:
2333: . snes - the SNES context
2335: Output Parameter:
2336: + r - the function (or PETSC_NULL)
2337: . ctx - the function context (or PETSC_NULL)
2338: - func - the function (or PETSC_NULL)
2340: Notes:
2341: SNESGetFunction() is valid for SNES_NONLINEAR_EQUATIONS methods only
2342: Analogous routines for SNES_UNCONSTRAINED_MINIMIZATION methods are
2343: SNESGetMinimizationFunction() and SNESGetGradient();
2345: Level: advanced
2347: .keywords: SNES, nonlinear, get, function
2349: .seealso: SNESSetFunction(), SNESGetSolution(), SNESGetMinimizationFunction(),
2350: SNESGetGradient()
2352: @*/
2353: int SNESGetFunction(SNES snes,Vec *r,void **ctx,int (**func)(SNES,Vec,Vec,void*))
2354: {
2357: if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
2358: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
2359: }
2360: if (r) *r = snes->vec_func_always;
2361: if (ctx) *ctx = snes->funP;
2362: if (func) *func = snes->computefunction;
2363: return(0);
2364: }
2366: /*@C
2367: SNESGetGradient - Returns the vector where the gradient is stored.
2369: Not Collective, but Vec is parallel if SNES is parallel
2371: Input Parameter:
2372: . snes - the SNES context
2374: Output Parameter:
2375: + r - the gradient (or PETSC_NULL)
2376: - ctx - the gradient context (or PETSC_NULL)
2378: Notes:
2379: SNESGetGradient() is valid for SNES_UNCONSTRAINED_MINIMIZATION methods
2380: only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is
2381: SNESGetFunction().
2383: Level: advanced
2385: .keywords: SNES, nonlinear, get, gradient
2387: .seealso: SNESGetMinimizationFunction(), SNESGetSolution(), SNESGetFunction(),
2388: SNESSetGradient(), SNESSetFunction()
2390: @*/
2391: int SNESGetGradient(SNES snes,Vec *r,void **ctx)
2392: {
2395: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
2396: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
2397: }
2398: if (r) *r = snes->vec_func_always;
2399: if (ctx) *ctx = snes->funP;
2400: return(0);
2401: }
2403: /*@C
2404: SNESGetMinimizationFunction - Returns the scalar function value for
2405: unconstrained minimization problems.
2407: Not Collective
2409: Input Parameter:
2410: . snes - the SNES context
2412: Output Parameter:
2413: + r - the function (or PETSC_NULL)
2414: - ctx - the function context (or PETSC_NULL)
2416: Notes:
2417: SNESGetMinimizationFunction() is valid for SNES_UNCONSTRAINED_MINIMIZATION
2418: methods only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is
2419: SNESGetFunction().
2421: Level: advanced
2423: .keywords: SNES, nonlinear, get, function
2425: .seealso: SNESGetGradient(), SNESGetSolution(), SNESGetFunction(), SNESSetFunction()
2427: @*/
2428: int SNESGetMinimizationFunction(SNES snes,PetscReal *r,void **ctx)
2429: {
2433: if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
2434: SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
2435: }
2436: if (r) *r = snes->fc;
2437: if (ctx) *ctx = snes->umfunP;
2438: return(0);
2439: }
2441: /*@C
2442: SNESSetOptionsPrefix - Sets the prefix used for searching for all
2443: SNES options in the database.
2445: Collective on SNES
2447: Input Parameter:
2448: + snes - the SNES context
2449: - prefix - the prefix to prepend to all option names
2451: Notes:
2452: A hyphen (-) must NOT be given at the beginning of the prefix name.
2453: The first character of all runtime options is AUTOMATICALLY the hyphen.
2455: Level: advanced
2457: .keywords: SNES, set, options, prefix, database
2459: .seealso: SNESSetFromOptions()
2460: @*/
2461: int SNESSetOptionsPrefix(SNES snes,char *prefix)
2462: {
2467: PetscObjectSetOptionsPrefix((PetscObject)snes,prefix);
2468: SLESSetOptionsPrefix(snes->sles,prefix);
2469: return(0);
2470: }
2472: /*@C
2473: SNESAppendOptionsPrefix - Appends to the prefix used for searching for all
2474: SNES options in the database.
2476: Collective on SNES
2478: Input Parameters:
2479: + snes - the SNES context
2480: - prefix - the prefix to prepend to all option names
2482: Notes:
2483: A hyphen (-) must NOT be given at the beginning of the prefix name.
2484: The first character of all runtime options is AUTOMATICALLY the hyphen.
2486: Level: advanced
2488: .keywords: SNES, append, options, prefix, database
2490: .seealso: SNESGetOptionsPrefix()
2491: @*/
2492: int SNESAppendOptionsPrefix(SNES snes,char *prefix)
2493: {
2495:
2498: PetscObjectAppendOptionsPrefix((PetscObject)snes,prefix);
2499: SLESAppendOptionsPrefix(snes->sles,prefix);
2500: return(0);
2501: }
2503: /*@C
2504: SNESGetOptionsPrefix - Sets the prefix used for searching for all
2505: SNES options in the database.
2507: Not Collective
2509: Input Parameter:
2510: . snes - the SNES context
2512: Output Parameter:
2513: . prefix - pointer to the prefix string used
2515: Notes: On the fortran side, the user should pass in a string 'prifix' of
2516: sufficient length to hold the prefix.
2518: Level: advanced
2520: .keywords: SNES, get, options, prefix, database
2522: .seealso: SNESAppendOptionsPrefix()
2523: @*/
2524: int SNESGetOptionsPrefix(SNES snes,char **prefix)
2525: {
2530: PetscObjectGetOptionsPrefix((PetscObject)snes,prefix);
2531: return(0);
2532: }
2534: /*MC
2535: SNESRegisterDynamic - Adds a method to the nonlinear solver package.
2537: Synopsis:
2538: int SNESRegisterDynamic(char *name_solver,char *path,char *name_create,int (*routine_create)(SNES))
2540: Not collective
2542: Input Parameters:
2543: + name_solver - name of a new user-defined solver
2544: . path - path (either absolute or relative) the library containing this solver
2545: . name_create - name of routine to create method context
2546: - routine_create - routine to create method context
2548: Notes:
2549: SNESRegisterDynamic() may be called multiple times to add several user-defined solvers.
2551: If dynamic libraries are used, then the fourth input argument (routine_create)
2552: is ignored.
2554: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR}, ${BOPT},
2555: and others of the form ${any_environmental_variable} occuring in pathname will be
2556: replaced with appropriate values.
2558: Sample usage:
2559: .vb
2560: SNESRegisterDynamic("my_solver",/home/username/my_lib/lib/libg/solaris/mylib.a,
2561: "MySolverCreate",MySolverCreate);
2562: .ve
2564: Then, your solver can be chosen with the procedural interface via
2565: $ SNESSetType(snes,"my_solver")
2566: or at runtime via the option
2567: $ -snes_type my_solver
2569: Level: advanced
2571: .keywords: SNES, nonlinear, register
2573: .seealso: SNESRegisterAll(), SNESRegisterDestroy()
2574: M*/
2576: int SNESRegister(char *sname,char *path,char *name,int (*function)(SNES))
2577: {
2578: char fullname[256];
2579: int ierr;
2582: PetscFListConcat(path,name,fullname);
2583: PetscFListAdd(&SNESList,sname,fullname,(void (*)(void))function);
2584: return(0);
2585: }