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 SNES_Solve, SNES_LineSearch, SNES_FunctionEval, SNES_JacobianEval;
12: /*@C
13: SNESView - Prints the SNES data structure.
15: Collective on SNES
17: Input Parameters:
18: + SNES - the SNES context
19: - viewer - visualization context
21: Options Database Key:
22: . -snes_view - Calls SNESView() at end of SNESSolve()
24: Notes:
25: The available visualization contexts include
26: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
27: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
28: output where only the first processor opens
29: the file. All other processors send their
30: data to the first processor to print.
32: The user can open an alternative visualization context with
33: PetscViewerASCIIOpen() - output to a specified file.
35: Level: beginner
37: .keywords: SNES, view
39: .seealso: PetscViewerASCIIOpen()
40: @*/
41: int SNESView(SNES snes,PetscViewer viewer)
42: {
43: SNES_KSP_EW_ConvCtx *kctx;
44: int ierr;
45: SLES sles;
46: char *type;
47: PetscTruth isascii,isstring;
51: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(snes->comm);
55: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
56: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
57: if (isascii) {
58: if (snes->prefix) {
59: PetscViewerASCIIPrintf(viewer,"SNES Object:(%s)n",snes->prefix);
60: } else {
61: PetscViewerASCIIPrintf(viewer,"SNES Object:n");
62: }
63: SNESGetType(snes,&type);
64: if (type) {
65: PetscViewerASCIIPrintf(viewer," type: %sn",type);
66: } else {
67: PetscViewerASCIIPrintf(viewer," type: not set yetn");
68: }
69: if (snes->view) {
70: PetscViewerASCIIPushTab(viewer);
71: (*snes->view)(snes,viewer);
72: PetscViewerASCIIPopTab(viewer);
73: }
74: PetscViewerASCIIPrintf(viewer," maximum iterations=%d, maximum function evaluations=%dn",snes->max_its,snes->max_funcs);
75: PetscViewerASCIIPrintf(viewer," tolerances: relative=%g, absolute=%g, solution=%gn",
76: snes->rtol,snes->atol,snes->xtol);
77: PetscViewerASCIIPrintf(viewer," total number of linear solver iterations=%dn",snes->linear_its);
78: PetscViewerASCIIPrintf(viewer," total number of function evaluations=%dn",snes->nfuncs);
79: if (snes->ksp_ewconv) {
80: kctx = (SNES_KSP_EW_ConvCtx *)snes->kspconvctx;
81: if (kctx) {
82: PetscViewerASCIIPrintf(viewer," Eisenstat-Walker computation of KSP relative tolerance (version %d)n",kctx->version);
83: PetscViewerASCIIPrintf(viewer," rtol_0=%g, rtol_max=%g, threshold=%gn",kctx->rtol_0,kctx->rtol_max,kctx->threshold);
84: PetscViewerASCIIPrintf(viewer," gamma=%g, alpha=%g, alpha2=%gn",kctx->gamma,kctx->alpha,kctx->alpha2);
85: }
86: }
87: } else if (isstring) {
88: SNESGetType(snes,&type);
89: PetscViewerStringSPrintf(viewer," %-3.3s",type);
90: }
91: SNESGetSLES(snes,&sles);
92: PetscViewerASCIIPushTab(viewer);
93: SLESView(sles,viewer);
94: PetscViewerASCIIPopTab(viewer);
95: return(0);
96: }
98: /*
99: We retain a list of functions that also take SNES command
100: line options. These are called at the end SNESSetFromOptions()
101: */
102: #define MAXSETFROMOPTIONS 5
103: static int numberofsetfromoptions;
104: static int (*othersetfromoptions[MAXSETFROMOPTIONS])(SNES);
106: /*@
107: SNESAddOptionsChecker - Adds an additional function to check for SNES options.
109: Not Collective
111: Input Parameter:
112: . snescheck - function that checks for options
114: Level: developer
116: .seealso: SNESSetFromOptions()
117: @*/
118: int SNESAddOptionsChecker(int (*snescheck)(SNES))
119: {
121: if (numberofsetfromoptions >= MAXSETFROMOPTIONS) {
122: SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE, "Too many options checkers, only %d allowed", MAXSETFROMOPTIONS);
123: }
125: othersetfromoptions[numberofsetfromoptions++] = snescheck;
126: return(0);
127: }
129: /*@
130: SNESSetFromOptions - Sets various SNES and SLES parameters from user options.
132: Collective on SNES
134: Input Parameter:
135: . snes - the SNES context
137: Options Database Keys:
138: + -snes_type <type> - ls, tr, umls, umtr, test
139: . -snes_stol - convergence tolerance in terms of the norm
140: of the change in the solution between steps
141: . -snes_atol <atol> - absolute tolerance of residual norm
142: . -snes_rtol <rtol> - relative decrease in tolerance norm from initial
143: . -snes_max_it <max_it> - maximum number of iterations
144: . -snes_max_funcs <max_funcs> - maximum number of function evaluations
145: . -snes_max_fail <max_fail> - maximum number of failures
146: . -snes_trtol <trtol> - trust region tolerance
147: . -snes_no_convergence_test - skip convergence test in nonlinear or minimization
148: solver; hence iterations will continue until max_it
149: or some other criterion is reached. Saves expense
150: of convergence test
151: . -snes_monitor - prints residual norm at each iteration
152: . -snes_vecmonitor - plots solution at each iteration
153: . -snes_vecmonitor_update - plots update to solution at each iteration
154: . -snes_xmonitor - plots residual norm at each iteration
155: . -snes_fd - use finite differences to compute Jacobian; very slow, only for testing
156: - -snes_mf_ksp_monitor - if using matrix-free multiply then print h at each KSP iteration
158: Options Database for Eisenstat-Walker method:
159: + -snes_ksp_ew_conv - use Eisenstat-Walker method for determining linear system convergence
160: . -snes_ksp_ew_version ver - version of Eisenstat-Walker method
161: . -snes_ksp_ew_rtol0 <rtol0> - Sets rtol0
162: . -snes_ksp_ew_rtolmax <rtolmax> - Sets rtolmax
163: . -snes_ksp_ew_gamma <gamma> - Sets gamma
164: . -snes_ksp_ew_alpha <alpha> - Sets alpha
165: . -snes_ksp_ew_alpha2 <alpha2> - Sets alpha2
166: - -snes_ksp_ew_threshold <threshold> - Sets threshold
168: Notes:
169: To see all options, run your program with the -help option or consult
170: the users manual.
172: Level: beginner
174: .keywords: SNES, nonlinear, set, options, database
176: .seealso: SNESSetOptionsPrefix()
177: @*/
178: int SNESSetFromOptions(SNES snes)
179: {
180: SLES sles;
181: SNES_KSP_EW_ConvCtx *kctx = (SNES_KSP_EW_ConvCtx *)snes->kspconvctx;
182: PetscTruth flg;
183: int ierr, i;
184: char *deft,type[256];
189: PetscOptionsBegin(snes->comm,snes->prefix,"Nonlinear solver (SNES) options","SNES");
190: if (snes->type_name) {
191: deft = snes->type_name;
192: } else {
193: deft = SNESLS;
194: }
196: if (!SNESRegisterAllCalled) {SNESRegisterAll(PETSC_NULL);}
197: PetscOptionsList("-snes_type","Nonlinear solver method","SNESSetType",SNESList,deft,type,256,&flg);
198: if (flg) {
199: SNESSetType(snes,type);
200: } else if (!snes->type_name) {
201: SNESSetType(snes,deft);
202: }
203: PetscOptionsName("-snes_view","Print detailed information on solver used","SNESView",0);
205: PetscOptionsReal("-snes_stol","Stop if step length less then","SNESSetTolerances",snes->xtol,&snes->xtol,0);
206: PetscOptionsReal("-snes_atol","Stop if function norm less then","SNESSetTolerances",snes->atol,&snes->atol,0);
208: PetscOptionsReal("-snes_rtol","Stop if decrease in function norm less then","SNESSetTolerances",snes->rtol,&snes->rtol,0);
209: PetscOptionsInt("-snes_max_it","Maximum iterations","SNESSetTolerances",snes->max_its,&snes->max_its,PETSC_NULL);
210: PetscOptionsInt("-snes_max_funcs","Maximum function evaluations","SNESSetTolerances",snes->max_funcs,&snes->max_funcs,PETSC_NULL);
211: PetscOptionsInt("-snes_max_fail","Maximum failures","SNESSetTolerances",snes->maxFailures,&snes->maxFailures,PETSC_NULL);
213: PetscOptionsName("-snes_ksp_ew_conv","Use Eisentat-Walker linear system convergence test","SNES_KSP_SetParametersEW",&snes->ksp_ewconv);
215: PetscOptionsInt("-snes_ksp_ew_version","Version 1 or 2","SNES_KSP_SetParametersEW",kctx->version,&kctx->version,0);
216: PetscOptionsReal("-snes_ksp_ew_rtol0","0 <= rtol0 < 1","SNES_KSP_SetParametersEW",kctx->rtol_0,&kctx->rtol_0,0);
217: PetscOptionsReal("-snes_ksp_ew_rtolmax","0 <= rtolmax < 1","SNES_KSP_SetParametersEW",kctx->rtol_max,&kctx->rtol_max,0);
218: PetscOptionsReal("-snes_ksp_ew_gamma","0 <= gamma <= 1","SNES_KSP_SetParametersEW",kctx->gamma,&kctx->gamma,0);
219: PetscOptionsReal("-snes_ksp_ew_alpha","1 < alpha <= 2","SNES_KSP_SetParametersEW",kctx->alpha,&kctx->alpha,0);
220: PetscOptionsReal("-snes_ksp_ew_alpha2","alpha2","SNES_KSP_SetParametersEW",kctx->alpha2,&kctx->alpha2,0);
221: PetscOptionsReal("-snes_ksp_ew_threshold","0 < threshold < 1","SNES_KSP_SetParametersEW",kctx->threshold,&kctx->threshold,0);
223: PetscOptionsName("-snes_no_convergence_test","Don't test for convergence","None",&flg);
224: if (flg) {snes->converged = 0;}
225: PetscOptionsName("-snes_cancelmonitors","Remove all monitors","SNESClearMonitor",&flg);
226: if (flg) {SNESClearMonitor(snes);}
227: PetscOptionsName("-snes_monitor","Monitor norm of function","SNESDefaultMonitor",&flg);
228: if (flg) {SNESSetMonitor(snes,SNESDefaultMonitor,0,0);}
229: PetscOptionsName("-snes_ratiomonitor","Monitor norm of function","SNESSetRatioMonitor",&flg);
230: if (flg) {SNESSetRatioMonitor(snes);}
231: PetscOptionsName("-snes_smonitor","Monitor norm of function (fewer digits)","SNESDefaultSMonitor",&flg);
232: if (flg) {SNESSetMonitor(snes,SNESDefaultSMonitor,0,0);}
233: PetscOptionsName("-snes_vecmonitor","Plot solution at each iteration","SNESVecViewMonitor",&flg);
234: if (flg) {SNESSetMonitor(snes,SNESVecViewMonitor,0,0);}
235: PetscOptionsName("-snes_vecmonitor_update","Plot correction at each iteration","SNESVecViewUpdateMonitor",&flg);
236: if (flg) {SNESSetMonitor(snes,SNESVecViewUpdateMonitor,0,0);}
237: PetscOptionsName("-snes_vecmonitor_residual","Plot residual at each iteration","SNESVecViewResidualMonitor",&flg);
238: if (flg) {SNESSetMonitor(snes,SNESVecViewResidualMonitor,0,0);}
239: PetscOptionsName("-snes_xmonitor","Plot function norm at each iteration","SNESLGMonitor",&flg);
240: if (flg) {SNESSetMonitor(snes,SNESLGMonitor,PETSC_NULL,PETSC_NULL);}
242: PetscOptionsName("-snes_fd","Use finite differences (slow) to compute Jacobian","SNESDefaultComputeJacobian",&flg);
243: if (flg) {
244: SNESSetJacobian(snes,snes->jacobian,snes->jacobian_pre,SNESDefaultComputeJacobian,snes->funP);
245: PetscLogInfo(snes,"SNESSetFromOptions: Setting default finite difference Jacobian matrixn");
246: }
248: for(i = 0; i < numberofsetfromoptions; i++) {
249: (*othersetfromoptions[i])(snes);
250: }
252: if (snes->setfromoptions) {
253: (*snes->setfromoptions)(snes);
254: }
256: PetscOptionsEnd();
258: SNESGetSLES(snes,&sles);
259: SLESSetFromOptions(sles);
261: return(0);
262: }
265: /*@
266: SNESSetApplicationContext - Sets the optional user-defined context for
267: the nonlinear solvers.
269: Collective on SNES
271: Input Parameters:
272: + snes - the SNES context
273: - usrP - optional user context
275: Level: intermediate
277: .keywords: SNES, nonlinear, set, application, context
279: .seealso: SNESGetApplicationContext()
280: @*/
281: int SNESSetApplicationContext(SNES snes,void *usrP)
282: {
285: snes->user = usrP;
286: return(0);
287: }
289: /*@C
290: SNESGetApplicationContext - Gets the user-defined context for the
291: nonlinear solvers.
293: Not Collective
295: Input Parameter:
296: . snes - SNES context
298: Output Parameter:
299: . usrP - user context
301: Level: intermediate
303: .keywords: SNES, nonlinear, get, application, context
305: .seealso: SNESSetApplicationContext()
306: @*/
307: int SNESGetApplicationContext(SNES snes,void **usrP)
308: {
311: *usrP = snes->user;
312: return(0);
313: }
315: /*@
316: SNESGetIterationNumber - Gets the number of nonlinear iterations completed
317: at this time.
319: Not Collective
321: Input Parameter:
322: . snes - SNES context
324: Output Parameter:
325: . iter - iteration number
327: Notes:
328: For example, during the computation of iteration 2 this would return 1.
330: This is useful for using lagged Jacobians (where one does not recompute the
331: Jacobian at each SNES iteration). For example, the code
332: .vb
333: SNESGetIterationNumber(snes,&it);
334: if (!(it % 2)) {
335: [compute Jacobian here]
336: }
337: .ve
338: can be used in your ComputeJacobian() function to cause the Jacobian to be
339: recomputed every second SNES iteration.
341: Level: intermediate
343: .keywords: SNES, nonlinear, get, iteration, number
344: @*/
345: int SNESGetIterationNumber(SNES snes,int* iter)
346: {
350: *iter = snes->iter;
351: return(0);
352: }
354: /*@
355: SNESGetFunctionNorm - Gets the norm of the current function that was set
356: with SNESSSetFunction().
358: Collective on SNES
360: Input Parameter:
361: . snes - SNES context
363: Output Parameter:
364: . fnorm - 2-norm of function
366: Level: intermediate
368: .keywords: SNES, nonlinear, get, function, norm
370: .seealso: SNESGetFunction()
371: @*/
372: int SNESGetFunctionNorm(SNES snes,PetscScalar *fnorm)
373: {
377: *fnorm = snes->norm;
378: return(0);
379: }
381: /*@
382: SNESGetNumberUnsuccessfulSteps - Gets the number of unsuccessful steps
383: attempted by the nonlinear solver.
385: Not Collective
387: Input Parameter:
388: . snes - SNES context
390: Output Parameter:
391: . nfails - number of unsuccessful steps attempted
393: Notes:
394: This counter is reset to zero for each successive call to SNESSolve().
396: Level: intermediate
398: .keywords: SNES, nonlinear, get, number, unsuccessful, steps
399: @*/
400: int SNESGetNumberUnsuccessfulSteps(SNES snes,int* nfails)
401: {
405: *nfails = snes->numFailures;
406: return(0);
407: }
409: /*@
410: SNESSetMaximumUnsuccessfulSteps - Sets the maximum number of unsuccessful steps
411: attempted by the nonlinear solver before it gives up.
413: Not Collective
415: Input Parameters:
416: + snes - SNES context
417: - maxFails - maximum of unsuccessful steps
419: Level: intermediate
421: .keywords: SNES, nonlinear, set, maximum, unsuccessful, steps
422: @*/
423: int SNESSetMaximumUnsuccessfulSteps(SNES snes, int maxFails)
424: {
427: snes->maxFailures = maxFails;
428: return(0);
429: }
431: /*@
432: SNESGetMaximumUnsuccessfulSteps - Gets the maximum number of unsuccessful steps
433: attempted by the nonlinear solver before it gives up.
435: Not Collective
437: Input Parameter:
438: . snes - SNES context
440: Output Parameter:
441: . maxFails - maximum of unsuccessful steps
443: Level: intermediate
445: .keywords: SNES, nonlinear, get, maximum, unsuccessful, steps
446: @*/
447: int SNESGetMaximumUnsuccessfulSteps(SNES snes, int *maxFails)
448: {
452: *maxFails = snes->maxFailures;
453: return(0);
454: }
456: /*@
457: SNESGetNumberLinearIterations - Gets the total number of linear iterations
458: used by the nonlinear solver.
460: Not Collective
462: Input Parameter:
463: . snes - SNES context
465: Output Parameter:
466: . lits - number of linear iterations
468: Notes:
469: This counter is reset to zero for each successive call to SNESSolve().
471: Level: intermediate
473: .keywords: SNES, nonlinear, get, number, linear, iterations
474: @*/
475: int SNESGetNumberLinearIterations(SNES snes,int* lits)
476: {
480: *lits = snes->linear_its;
481: return(0);
482: }
484: /*@C
485: SNESGetSLES - Returns the SLES context for a SNES solver.
487: Not Collective, but if SNES object is parallel, then SLES object is parallel
489: Input Parameter:
490: . snes - the SNES context
492: Output Parameter:
493: . sles - the SLES context
495: Notes:
496: The user can then directly manipulate the SLES context to set various
497: options, etc. Likewise, the user can then extract and manipulate the
498: KSP and PC contexts as well.
500: Level: beginner
502: .keywords: SNES, nonlinear, get, SLES, context
504: .seealso: SLESGetPC(), SLESGetKSP()
505: @*/
506: int SNESGetSLES(SNES snes,SLES *sles)
507: {
510: *sles = snes->sles;
511: return(0);
512: }
514: static int SNESPublish_Petsc(PetscObject obj)
515: {
516: #if defined(PETSC_HAVE_AMS)
517: SNES v = (SNES) obj;
518: int ierr;
519: #endif
523: #if defined(PETSC_HAVE_AMS)
524: /* if it is already published then return */
525: if (v->amem >=0) return(0);
527: PetscObjectPublishBaseBegin(obj);
528: AMS_Memory_add_field((AMS_Memory)v->amem,"Iteration",&v->iter,1,AMS_INT,AMS_READ,
529: AMS_COMMON,AMS_REDUCT_UNDEF);
530: AMS_Memory_add_field((AMS_Memory)v->amem,"Residual",&v->norm,1,AMS_DOUBLE,AMS_READ,
531: AMS_COMMON,AMS_REDUCT_UNDEF);
532: PetscObjectPublishBaseEnd(obj);
533: #endif
534: return(0);
535: }
537: /* -----------------------------------------------------------*/
538: /*@C
539: SNESCreate - Creates a nonlinear solver context.
541: Collective on MPI_Comm
543: Input Parameters:
544: + comm - MPI communicator
546: Output Parameter:
547: . outsnes - the new SNES context
549: Options Database Keys:
550: + -snes_mf - Activates default matrix-free Jacobian-vector products,
551: and no preconditioning matrix
552: . -snes_mf_operator - Activates default matrix-free Jacobian-vector
553: products, and a user-provided preconditioning matrix
554: as set by SNESSetJacobian()
555: - -snes_fd - Uses (slow!) finite differences to compute Jacobian
557: Level: beginner
559: .keywords: SNES, nonlinear, create, context
561: .seealso: SNESSolve(), SNESDestroy(), SNES
562: @*/
563: int SNESCreate(MPI_Comm comm,SNES *outsnes)
564: {
565: int ierr;
566: SNES snes;
567: SNES_KSP_EW_ConvCtx *kctx;
571: *outsnes = PETSC_NULL;
572: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
573: SNESInitializePackage(PETSC_NULL);
574: #endif
576: PetscHeaderCreate(snes,_p_SNES,int,SNES_COOKIE,0,"SNES",comm,SNESDestroy,SNESView);
577: PetscLogObjectCreate(snes);
578: snes->bops->publish = SNESPublish_Petsc;
579: snes->max_its = 50;
580: snes->max_funcs = 10000;
581: snes->norm = 0.0;
582: snes->rtol = 1.e-8;
583: snes->ttol = 0.0;
584: snes->atol = 1.e-50;
585: snes->xtol = 1.e-8;
586: snes->deltatol = 1.e-12;
587: snes->nfuncs = 0;
588: snes->numFailures = 0;
589: snes->maxFailures = 1;
590: snes->linear_its = 0;
591: snes->numbermonitors = 0;
592: snes->data = 0;
593: snes->view = 0;
594: snes->setupcalled = 0;
595: snes->ksp_ewconv = PETSC_FALSE;
596: snes->vwork = 0;
597: snes->nwork = 0;
598: snes->conv_hist_len = 0;
599: snes->conv_hist_max = 0;
600: snes->conv_hist = PETSC_NULL;
601: snes->conv_hist_its = PETSC_NULL;
602: snes->conv_hist_reset = PETSC_TRUE;
603: snes->reason = SNES_CONVERGED_ITERATING;
605: /* Create context to compute Eisenstat-Walker relative tolerance for KSP */
606: PetscNew(SNES_KSP_EW_ConvCtx,&kctx);
607: PetscLogObjectMemory(snes,sizeof(SNES_KSP_EW_ConvCtx));
608: snes->kspconvctx = (void*)kctx;
609: kctx->version = 2;
610: kctx->rtol_0 = .3; /* Eisenstat and Walker suggest rtol_0=.5, but
611: this was too large for some test cases */
612: kctx->rtol_last = 0;
613: kctx->rtol_max = .9;
614: kctx->gamma = 1.0;
615: kctx->alpha2 = .5*(1.0 + sqrt(5.0));
616: kctx->alpha = kctx->alpha2;
617: kctx->threshold = .1;
618: kctx->lresid_last = 0;
619: kctx->norm_last = 0;
621: SLESCreate(comm,&snes->sles);
622: PetscLogObjectParent(snes,snes->sles)
624: *outsnes = snes;
625: PetscPublishAll(snes);
626: return(0);
627: }
629: /* --------------------------------------------------------------- */
630: /*@C
631: SNESSetFunction - Sets the function evaluation routine and function
632: vector for use by the SNES routines in solving systems of nonlinear
633: equations.
635: Collective on SNES
637: Input Parameters:
638: + snes - the SNES context
639: . func - function evaluation routine
640: . r - vector to store function value
641: - ctx - [optional] user-defined context for private data for the
642: function evaluation routine (may be PETSC_NULL)
644: Calling sequence of func:
645: $ func (SNES snes,Vec x,Vec f,void *ctx);
647: . f - function vector
648: - ctx - optional user-defined function context
650: Notes:
651: The Newton-like methods typically solve linear systems of the form
652: $ f'(x) x = -f(x),
653: where f'(x) denotes the Jacobian matrix and f(x) is the function.
655: Level: beginner
657: .keywords: SNES, nonlinear, set, function
659: .seealso: SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian()
660: @*/
661: int SNESSetFunction(SNES snes,Vec r,int (*func)(SNES,Vec,Vec,void*),void *ctx)
662: {
668: snes->computefunction = func;
669: snes->vec_func = snes->vec_func_always = r;
670: snes->funP = ctx;
671: return(0);
672: }
674: /*@
675: SNESComputeFunction - Calls the function that has been set with
676: SNESSetFunction().
678: Collective on SNES
680: Input Parameters:
681: + snes - the SNES context
682: - x - input vector
684: Output Parameter:
685: . y - function vector, as set by SNESSetFunction()
687: Notes:
688: SNESComputeFunction() is typically used within nonlinear solvers
689: implementations, so most users would not generally call this routine
690: themselves.
692: Level: developer
694: .keywords: SNES, nonlinear, compute, function
696: .seealso: SNESSetFunction(), SNESGetFunction()
697: @*/
698: int SNESComputeFunction(SNES snes,Vec x,Vec y)
699: {
700: int ierr;
709: PetscLogEventBegin(SNES_FunctionEval,snes,x,y,0);
710: PetscStackPush("SNES user function");
711: (*snes->computefunction)(snes,x,y,snes->funP);
712: PetscStackPop;
713: snes->nfuncs++;
714: PetscLogEventEnd(SNES_FunctionEval,snes,x,y,0);
715: return(0);
716: }
718: /*@
719: SNESComputeJacobian - Computes the Jacobian matrix that has been
720: set with SNESSetJacobian().
722: Collective on SNES and Mat
724: Input Parameters:
725: + snes - the SNES context
726: - x - input vector
728: Output Parameters:
729: + A - Jacobian matrix
730: . B - optional preconditioning matrix
731: - flag - flag indicating matrix structure
733: Notes:
734: Most users should not need to explicitly call this routine, as it
735: is used internally within the nonlinear solvers.
737: See SLESSetOperators() for important information about setting the
738: flag parameter.
740: Level: developer
742: .keywords: SNES, compute, Jacobian, matrix
744: .seealso: SNESSetJacobian(), SLESSetOperators()
745: @*/
746: int SNESComputeJacobian(SNES snes,Vec X,Mat *A,Mat *B,MatStructure *flg)
747: {
748: int ierr;
754: if (!snes->computejacobian) return(0);
755: PetscLogEventBegin(SNES_JacobianEval,snes,X,*A,*B);
756: *flg = DIFFERENT_NONZERO_PATTERN;
757: PetscStackPush("SNES user Jacobian function");
758: (*snes->computejacobian)(snes,X,A,B,flg,snes->jacP);
759: PetscStackPop;
760: PetscLogEventEnd(SNES_JacobianEval,snes,X,*A,*B);
761: /* make sure user returned a correct Jacobian and preconditioner */
764: return(0);
765: }
767: /*@C
768: SNESSetJacobian - Sets the function to compute Jacobian as well as the
769: location to store the matrix.
771: Collective on SNES and Mat
773: Input Parameters:
774: + snes - the SNES context
775: . A - Jacobian matrix
776: . B - preconditioner matrix (usually same as the Jacobian)
777: . func - Jacobian evaluation routine
778: - ctx - [optional] user-defined context for private data for the
779: Jacobian evaluation routine (may be PETSC_NULL)
781: Calling sequence of func:
782: $ func (SNES snes,Vec x,Mat *A,Mat *B,int *flag,void *ctx);
784: + x - input vector
785: . A - Jacobian matrix
786: . B - preconditioner matrix, usually the same as A
787: . flag - flag indicating information about the preconditioner matrix
788: structure (same as flag in SLESSetOperators())
789: - ctx - [optional] user-defined Jacobian context
791: Notes:
792: See SLESSetOperators() for important information about setting the flag
793: output parameter in the routine func(). Be sure to read this information!
795: The routine func() takes Mat * as the matrix arguments rather than Mat.
796: This allows the Jacobian evaluation routine to replace A and/or B with a
797: completely new new matrix structure (not just different matrix elements)
798: when appropriate, for instance, if the nonzero structure is changing
799: throughout the global iterations.
801: Level: beginner
803: .keywords: SNES, nonlinear, set, Jacobian, matrix
805: .seealso: SLESSetOperators(), SNESSetFunction()
806: @*/
807: int SNESSetJacobian(SNES snes,Mat A,Mat B,int (*func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
808: {
817: if (func) snes->computejacobian = func;
818: if (ctx) snes->jacP = ctx;
819: if (A) {
820: if (snes->jacobian) {MatDestroy(snes->jacobian);}
821: snes->jacobian = A;
822: ierr = PetscObjectReference((PetscObject)A);
823: }
824: if (B) {
825: if (snes->jacobian_pre) {MatDestroy(snes->jacobian_pre);}
826: snes->jacobian_pre = B;
827: ierr = PetscObjectReference((PetscObject)B);
828: }
829: return(0);
830: }
832: /*@C
833: SNESGetJacobian - Returns the Jacobian matrix and optionally the user
834: provided context for evaluating the Jacobian.
836: Not Collective, but Mat object will be parallel if SNES object is
838: Input Parameter:
839: . snes - the nonlinear solver context
841: Output Parameters:
842: + A - location to stash Jacobian matrix (or PETSC_NULL)
843: . B - location to stash preconditioner matrix (or PETSC_NULL)
844: . ctx - location to stash Jacobian ctx (or PETSC_NULL)
845: - func - location to put Jacobian function (or PETSC_NULL)
847: Level: advanced
849: .seealso: SNESSetJacobian(), SNESComputeJacobian()
850: @*/
851: int SNESGetJacobian(SNES snes,Mat *A,Mat *B,void **ctx,int (**func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*))
852: {
855: if (A) *A = snes->jacobian;
856: if (B) *B = snes->jacobian_pre;
857: if (ctx) *ctx = snes->jacP;
858: if (func) *func = snes->computejacobian;
859: return(0);
860: }
862: /* ----- Routines to initialize and destroy a nonlinear solver ---- */
863: extern int SNESDefaultMatrixFreeCreate2(SNES,Vec,Mat*);
865: /*@
866: SNESSetUp - Sets up the internal data structures for the later use
867: of a nonlinear solver.
869: Collective on SNES
871: Input Parameters:
872: + snes - the SNES context
873: - x - the solution vector
875: Notes:
876: For basic use of the SNES solvers the user need not explicitly call
877: SNESSetUp(), since these actions will automatically occur during
878: the call to SNESSolve(). However, if one wishes to control this
879: phase separately, SNESSetUp() should be called after SNESCreate()
880: and optional routines of the form SNESSetXXX(), but before SNESSolve().
882: Level: advanced
884: .keywords: SNES, nonlinear, setup
886: .seealso: SNESCreate(), SNESSolve(), SNESDestroy()
887: @*/
888: int SNESSetUp(SNES snes,Vec x)
889: {
890: int ierr;
891: PetscTruth flg, iseqtr;
897: snes->vec_sol = snes->vec_sol_always = x;
899: PetscOptionsHasName(snes->prefix,"-snes_mf_operator",&flg);
900: /*
901: This version replaces the user provided Jacobian matrix with a
902: matrix-free version but still employs the user-provided preconditioner matrix
903: */
904: if (flg) {
905: Mat J;
906: MatCreateSNESMF(snes,snes->vec_sol,&J);
907: MatSNESMFSetFromOptions(J);
908: PetscLogInfo(snes,"SNESSetUp: Setting default matrix-free operator routinesn");
909: SNESSetJacobian(snes,J,0,0,0);
910: MatDestroy(J);
911: }
913: #if !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
914: PetscOptionsHasName(snes->prefix,"-snes_mf_operator2",&flg);
915: if (flg) {
916: Mat J;
917: SNESDefaultMatrixFreeCreate2(snes,snes->vec_sol,&J);
918: SNESSetJacobian(snes,J,0,0,0);
919: MatDestroy(J);
920: }
921: #endif
923: PetscOptionsHasName(snes->prefix,"-snes_mf",&flg);
924: /*
925: This version replaces both the user-provided Jacobian and the user-
926: provided preconditioner matrix with the default matrix free version.
927: */
928: if (flg) {
929: Mat J;
930: SLES sles;
931: PC pc;
933: MatCreateSNESMF(snes,snes->vec_sol,&J);
934: MatSNESMFSetFromOptions(J);
935: PetscLogInfo(snes,"SNESSetUp: Setting default matrix-free operator and preconditioner routinesn");
936: SNESSetJacobian(snes,J,J,MatSNESMFComputeJacobian,snes->funP);
937: MatDestroy(J);
939: /* force no preconditioner */
940: SNESGetSLES(snes,&sles);
941: SLESGetPC(sles,&pc);
942: PetscTypeCompare((PetscObject)pc,PCSHELL,&flg);
943: if (!flg) {
944: PCSetType(pc,PCNONE);
945: }
946: }
948: if (!snes->vec_func) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetFunction() first");
949: if (!snes->computefunction) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetFunction() first");
950: if (!snes->jacobian) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetJacobian() first n or use -snes_mf option");
951: if (snes->vec_func == snes->vec_sol) {
952: SETERRQ(PETSC_ERR_ARG_IDN,"Solution vector cannot be function vector");
953: }
955: /* Set the KSP stopping criterion to use the Eisenstat-Walker method */
956: PetscTypeCompare((PetscObject)snes,SNESTR,&iseqtr);
957: if (snes->ksp_ewconv && !iseqtr) {
958: SLES sles; KSP ksp;
959: SNESGetSLES(snes,&sles);
960: SLESGetKSP(sles,&ksp);
961: KSPSetConvergenceTest(ksp,SNES_KSP_EW_Converged_Private,snes);
962: }
964: if (snes->setup) {(*snes->setup)(snes);}
965: snes->setupcalled = 1;
966: return(0);
967: }
969: /*@C
970: SNESDestroy - Destroys the nonlinear solver context that was created
971: with SNESCreate().
973: Collective on SNES
975: Input Parameter:
976: . snes - the SNES context
978: Level: beginner
980: .keywords: SNES, nonlinear, destroy
982: .seealso: SNESCreate(), SNESSolve()
983: @*/
984: int SNESDestroy(SNES snes)
985: {
986: int i,ierr;
990: if (--snes->refct > 0) return(0);
992: /* if memory was published with AMS then destroy it */
993: PetscObjectDepublish(snes);
995: if (snes->destroy) {(*(snes)->destroy)(snes);}
996: if (snes->kspconvctx) {PetscFree(snes->kspconvctx);}
997: if (snes->jacobian) {MatDestroy(snes->jacobian);}
998: if (snes->jacobian_pre) {MatDestroy(snes->jacobian_pre);}
999: SLESDestroy(snes->sles);
1000: if (snes->vwork) {VecDestroyVecs(snes->vwork,snes->nvwork);}
1001: for (i=0; i<snes->numbermonitors; i++) {
1002: if (snes->monitordestroy[i]) {
1003: (*snes->monitordestroy[i])(snes->monitorcontext[i]);
1004: }
1005: }
1006: PetscLogObjectDestroy((PetscObject)snes);
1007: PetscHeaderDestroy((PetscObject)snes);
1008: return(0);
1009: }
1011: /* ----------- Routines to set solver parameters ---------- */
1013: /*@
1014: SNESSetTolerances - Sets various parameters used in convergence tests.
1016: Collective on SNES
1018: Input Parameters:
1019: + snes - the SNES context
1020: . atol - absolute convergence tolerance
1021: . rtol - relative convergence tolerance
1022: . stol - convergence tolerance in terms of the norm
1023: of the change in the solution between steps
1024: . maxit - maximum number of iterations
1025: - maxf - maximum number of function evaluations
1027: Options Database Keys:
1028: + -snes_atol <atol> - Sets atol
1029: . -snes_rtol <rtol> - Sets rtol
1030: . -snes_stol <stol> - Sets stol
1031: . -snes_max_it <maxit> - Sets maxit
1032: - -snes_max_funcs <maxf> - Sets maxf
1034: Notes:
1035: The default maximum number of iterations is 50.
1036: The default maximum number of function evaluations is 1000.
1038: Level: intermediate
1040: .keywords: SNES, nonlinear, set, convergence, tolerances
1042: .seealso: SNESSetTrustRegionTolerance(), SNESSetMinimizationFunctionTolerance()
1043: @*/
1044: int SNESSetTolerances(SNES snes,PetscReal atol,PetscReal rtol,PetscReal stol,int maxit,int maxf)
1045: {
1048: if (atol != PETSC_DEFAULT) snes->atol = atol;
1049: if (rtol != PETSC_DEFAULT) snes->rtol = rtol;
1050: if (stol != PETSC_DEFAULT) snes->xtol = stol;
1051: if (maxit != PETSC_DEFAULT) snes->max_its = maxit;
1052: if (maxf != PETSC_DEFAULT) snes->max_funcs = maxf;
1053: return(0);
1054: }
1056: /*@
1057: SNESGetTolerances - Gets various parameters used in convergence tests.
1059: Not Collective
1061: Input Parameters:
1062: + snes - the SNES context
1063: . atol - absolute convergence tolerance
1064: . rtol - relative convergence tolerance
1065: . stol - convergence tolerance in terms of the norm
1066: of the change in the solution between steps
1067: . maxit - maximum number of iterations
1068: - maxf - maximum number of function evaluations
1070: Notes:
1071: The user can specify PETSC_NULL for any parameter that is not needed.
1073: Level: intermediate
1075: .keywords: SNES, nonlinear, get, convergence, tolerances
1077: .seealso: SNESSetTolerances()
1078: @*/
1079: int SNESGetTolerances(SNES snes,PetscReal *atol,PetscReal *rtol,PetscReal *stol,int *maxit,int *maxf)
1080: {
1083: if (atol) *atol = snes->atol;
1084: if (rtol) *rtol = snes->rtol;
1085: if (stol) *stol = snes->xtol;
1086: if (maxit) *maxit = snes->max_its;
1087: if (maxf) *maxf = snes->max_funcs;
1088: return(0);
1089: }
1091: /*@
1092: SNESSetTrustRegionTolerance - Sets the trust region parameter tolerance.
1094: Collective on SNES
1096: Input Parameters:
1097: + snes - the SNES context
1098: - tol - tolerance
1099:
1100: Options Database Key:
1101: . -snes_trtol <tol> - Sets tol
1103: Level: intermediate
1105: .keywords: SNES, nonlinear, set, trust region, tolerance
1107: .seealso: SNESSetTolerances(), SNESSetMinimizationFunctionTolerance()
1108: @*/
1109: int SNESSetTrustRegionTolerance(SNES snes,PetscReal tol)
1110: {
1113: snes->deltatol = tol;
1114: return(0);
1115: }
1117: /*
1118: Duplicate the lg monitors for SNES from KSP; for some reason with
1119: dynamic libraries things don't work under Sun4 if we just use
1120: macros instead of functions
1121: */
1122: int SNESLGMonitor(SNES snes,int it,PetscReal norm,void *ctx)
1123: {
1128: KSPLGMonitor((KSP)snes,it,norm,ctx);
1129: return(0);
1130: }
1132: int SNESLGMonitorCreate(char *host,char *label,int x,int y,int m,int n,PetscDrawLG *draw)
1133: {
1137: KSPLGMonitorCreate(host,label,x,y,m,n,draw);
1138: return(0);
1139: }
1141: int SNESLGMonitorDestroy(PetscDrawLG draw)
1142: {
1146: KSPLGMonitorDestroy(draw);
1147: return(0);
1148: }
1150: /* ------------ Routines to set performance monitoring options ----------- */
1152: /*@C
1153: SNESSetMonitor - Sets an ADDITIONAL function that is to be used at every
1154: iteration of the nonlinear solver to display the iteration's
1155: progress.
1157: Collective on SNES
1159: Input Parameters:
1160: + snes - the SNES context
1161: . func - monitoring routine
1162: . mctx - [optional] user-defined context for private data for the
1163: monitor routine (use PETSC_NULL if no context is desitre)
1164: - monitordestroy - [optional] routine that frees monitor context
1165: (may be PETSC_NULL)
1167: Calling sequence of func:
1168: $ int func(SNES snes,int its, PetscReal norm,void *mctx)
1170: + snes - the SNES context
1171: . its - iteration number
1172: . norm - 2-norm function value (may be estimated)
1173: - mctx - [optional] monitoring context
1175: Options Database Keys:
1176: + -snes_monitor - sets SNESDefaultMonitor()
1177: . -snes_xmonitor - sets line graph monitor,
1178: uses SNESLGMonitorCreate()
1179: _ -snes_cancelmonitors - cancels all monitors that have
1180: been hardwired into a code by
1181: calls to SNESSetMonitor(), but
1182: does not cancel those set via
1183: the options database.
1185: Notes:
1186: Several different monitoring routines may be set by calling
1187: SNESSetMonitor() multiple times; all will be called in the
1188: order in which they were set.
1190: Level: intermediate
1192: .keywords: SNES, nonlinear, set, monitor
1194: .seealso: SNESDefaultMonitor(), SNESClearMonitor()
1195: @*/
1196: int SNESSetMonitor(SNES snes,int (*func)(SNES,int,PetscReal,void*),void *mctx,int (*monitordestroy)(void *))
1197: {
1200: if (snes->numbermonitors >= MAXSNESMONITORS) {
1201: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1202: }
1204: snes->monitor[snes->numbermonitors] = func;
1205: snes->monitordestroy[snes->numbermonitors] = monitordestroy;
1206: snes->monitorcontext[snes->numbermonitors++] = (void*)mctx;
1207: return(0);
1208: }
1210: /*@C
1211: SNESClearMonitor - Clears all the monitor functions for a SNES object.
1213: Collective on SNES
1215: Input Parameters:
1216: . snes - the SNES context
1218: Options Database:
1219: . -snes_cancelmonitors - cancels all monitors that have been hardwired
1220: into a code by calls to SNESSetMonitor(), but does not cancel those
1221: set via the options database
1223: Notes:
1224: There is no way to clear one specific monitor from a SNES object.
1226: Level: intermediate
1228: .keywords: SNES, nonlinear, set, monitor
1230: .seealso: SNESDefaultMonitor(), SNESSetMonitor()
1231: @*/
1232: int SNESClearMonitor(SNES snes)
1233: {
1236: snes->numbermonitors = 0;
1237: return(0);
1238: }
1240: /*@C
1241: SNESSetConvergenceTest - Sets the function that is to be used
1242: to test for convergence of the nonlinear iterative solution.
1244: Collective on SNES
1246: Input Parameters:
1247: + snes - the SNES context
1248: . func - routine to test for convergence
1249: - cctx - [optional] context for private data for the convergence routine
1250: (may be PETSC_NULL)
1252: Calling sequence of func:
1253: $ int func (SNES snes,PetscReal xnorm,PetscReal gnorm,PetscReal f,SNESConvergedReason *reason,void *cctx)
1255: + snes - the SNES context
1256: . cctx - [optional] convergence context
1257: . reason - reason for convergence/divergence
1258: . xnorm - 2-norm of current iterate
1259: . gnorm - 2-norm of current step
1260: - f - 2-norm of function
1262: Level: advanced
1264: .keywords: SNES, nonlinear, set, convergence, test
1266: .seealso: SNESConverged_LS(), SNESConverged_TR()
1267: @*/
1268: int SNESSetConvergenceTest(SNES snes,int (*func)(SNES,PetscReal,PetscReal,PetscReal,SNESConvergedReason*,void*),void *cctx)
1269: {
1272: (snes)->converged = func;
1273: (snes)->cnvP = cctx;
1274: return(0);
1275: }
1277: /*@C
1278: SNESGetConvergedReason - Gets the reason the SNES iteration was stopped.
1280: Not Collective
1282: Input Parameter:
1283: . snes - the SNES context
1285: Output Parameter:
1286: . reason - negative value indicates diverged, positive value converged, see petscsnes.h or the
1287: manual pages for the individual convergence tests for complete lists
1289: Level: intermediate
1291: Notes: Can only be called after the call the SNESSolve() is complete.
1293: .keywords: SNES, nonlinear, set, convergence, test
1295: .seealso: SNESSetConvergenceTest(), SNESConverged_LS(), SNESConverged_TR(), SNESConvergedReason
1296: @*/
1297: int SNESGetConvergedReason(SNES snes,SNESConvergedReason *reason)
1298: {
1301: *reason = snes->reason;
1302: return(0);
1303: }
1305: /*@
1306: SNESSetConvergenceHistory - Sets the array used to hold the convergence history.
1308: Collective on SNES
1310: Input Parameters:
1311: + snes - iterative context obtained from SNESCreate()
1312: . a - array to hold history
1313: . its - integer array holds the number of linear iterations for each solve.
1314: . na - size of a and its
1315: - reset - PETSC_TRUE indicates each new nonlinear solve resets the history counter to zero,
1316: else it continues storing new values for new nonlinear solves after the old ones
1318: Notes:
1319: If set, this array will contain the function norms computed
1320: at each step.
1322: This routine is useful, e.g., when running a code for purposes
1323: of accurate performance monitoring, when no I/O should be done
1324: during the section of code that is being timed.
1326: Level: intermediate
1328: .keywords: SNES, set, convergence, history
1330: .seealso: SNESGetConvergenceHistory()
1332: @*/
1333: int SNESSetConvergenceHistory(SNES snes,PetscReal *a,int *its,int na,PetscTruth reset)
1334: {
1338: snes->conv_hist = a;
1339: snes->conv_hist_its = its;
1340: snes->conv_hist_max = na;
1341: snes->conv_hist_reset = reset;
1342: return(0);
1343: }
1345: /*@C
1346: SNESGetConvergenceHistory - Gets the array used to hold the convergence history.
1348: Collective on SNES
1350: Input Parameter:
1351: . snes - iterative context obtained from SNESCreate()
1353: Output Parameters:
1354: . a - array to hold history
1355: . its - integer array holds the number of linear iterations (or
1356: negative if not converged) for each solve.
1357: - na - size of a and its
1359: Notes:
1360: The calling sequence for this routine in Fortran is
1361: $ call SNESGetConvergenceHistory(SNES snes, integer na, integer ierr)
1363: This routine is useful, e.g., when running a code for purposes
1364: of accurate performance monitoring, when no I/O should be done
1365: during the section of code that is being timed.
1367: Level: intermediate
1369: .keywords: SNES, get, convergence, history
1371: .seealso: SNESSetConvergencHistory()
1373: @*/
1374: int SNESGetConvergenceHistory(SNES snes,PetscReal **a,int **its,int *na)
1375: {
1378: if (a) *a = snes->conv_hist;
1379: if (its) *its = snes->conv_hist_its;
1380: if (na) *na = snes->conv_hist_len;
1381: return(0);
1382: }
1384: /*@
1385: SNESSetRhsBC - Sets the function which applies boundary conditions
1386: to the Rhs of each system.
1388: Collective on SNES
1390: Input Parameters:
1391: . snes - The nonlinear solver context
1392: . func - The function
1394: Calling sequence of func:
1395: . func (SNES snes, Vec rhs, void *ctx);
1397: . rhs - The current rhs vector
1398: . ctx - The user-context
1400: Level: intermediate
1402: .keywords: SNES, Rhs, boundary conditions
1403: .seealso SNESDefaultRhsBC(), SNESSetSolutionBC(), SNESSetUpdate()
1404: @*/
1405: int SNESSetRhsBC(SNES snes, int (*func)(SNES, Vec, void *))
1406: {
1409: snes->applyrhsbc = func;
1410: return(0);
1411: }
1413: /*@
1414: SNESDefaultRhsBC - The default boundary condition function which does nothing.
1416: Not collective
1418: Input Parameters:
1419: . snes - The nonlinear solver context
1420: . rhs - The Rhs
1421: . ctx - The user-context
1423: Level: beginner
1425: .keywords: SNES, Rhs, boundary conditions
1426: .seealso SNESSetRhsBC(), SNESDefaultSolutionBC(), SNESDefaultUpdate()
1427: @*/
1428: int SNESDefaultRhsBC(SNES snes, Vec rhs, void *ctx)
1429: {
1431: return(0);
1432: }
1434: /*@
1435: SNESSetSolutionBC - Sets the function which applies boundary conditions
1436: to the solution of each system.
1438: Collective on SNES
1440: Input Parameters:
1441: . snes - The nonlinear solver context
1442: . func - The function
1444: Calling sequence of func:
1445: . func (SNES snes, Vec rsol, void *ctx);
1447: . sol - The current solution vector
1448: . ctx - The user-context
1450: Level: intermediate
1452: .keywords: SNES, solution, boundary conditions
1453: .seealso SNESDefaultSolutionBC(), SNESSetRhsBC(), SNESSetUpdate()
1454: @*/
1455: int SNESSetSolutionBC(SNES snes, int (*func)(SNES, Vec, void *))
1456: {
1459: snes->applysolbc = func;
1460: return(0);
1461: }
1463: /*@
1464: SNESDefaultSolutionBC - The default boundary condition function which does nothing.
1466: Not collective
1468: Input Parameters:
1469: . snes - The nonlinear solver context
1470: . sol - The solution
1471: . ctx - The user-context
1473: Level: beginner
1475: .keywords: SNES, solution, boundary conditions
1476: .seealso SNESSetSolutionBC(), SNESDefaultRhsBC(), SNESDefaultUpdate()
1477: @*/
1478: int SNESDefaultSolutionBC(SNES snes, Vec sol, void *ctx)
1479: {
1481: return(0);
1482: }
1484: /*@
1485: SNESSetUpdate - Sets the general-purpose update function called
1486: at the beginning of every step of the iteration.
1488: Collective on SNES
1490: Input Parameters:
1491: . snes - The nonlinear solver context
1492: . func - The function
1494: Calling sequence of func:
1495: . func (TS ts, int step);
1497: . step - The current step of the iteration
1499: Level: intermediate
1501: .keywords: SNES, update
1502: .seealso SNESDefaultUpdate(), SNESSetRhsBC(), SNESSetSolutionBC()
1503: @*/
1504: int SNESSetUpdate(SNES snes, int (*func)(SNES, int))
1505: {
1508: snes->update = func;
1509: return(0);
1510: }
1512: /*@
1513: SNESDefaultUpdate - The default update function which does nothing.
1515: Not collective
1517: Input Parameters:
1518: . snes - The nonlinear solver context
1519: . step - The current step of the iteration
1521: Level: intermediate
1523: .keywords: SNES, update
1524: .seealso SNESSetUpdate(), SNESDefaultRhsBC(), SNESDefaultSolutionBC()
1525: @*/
1526: int SNESDefaultUpdate(SNES snes, int step)
1527: {
1529: return(0);
1530: }
1532: /*
1533: SNESScaleStep_Private - Scales a step so that its length is less than the
1534: positive parameter delta.
1536: Input Parameters:
1537: + snes - the SNES context
1538: . y - approximate solution of linear system
1539: . fnorm - 2-norm of current function
1540: - delta - trust region size
1542: Output Parameters:
1543: + gpnorm - predicted function norm at the new point, assuming local
1544: linearization. The value is zero if the step lies within the trust
1545: region, and exceeds zero otherwise.
1546: - ynorm - 2-norm of the step
1548: Note:
1549: For non-trust region methods such as SNESLS, the parameter delta
1550: is set to be the maximum allowable step size.
1552: .keywords: SNES, nonlinear, scale, step
1553: */
1554: int SNESScaleStep_Private(SNES snes,Vec y,PetscReal *fnorm,PetscReal *delta,PetscReal *gpnorm,PetscReal *ynorm)
1555: {
1556: PetscReal nrm;
1557: PetscScalar cnorm;
1558: int ierr;
1565: VecNorm(y,NORM_2,&nrm);
1566: if (nrm > *delta) {
1567: nrm = *delta/nrm;
1568: *gpnorm = (1.0 - nrm)*(*fnorm);
1569: cnorm = nrm;
1570: VecScale(&cnorm,y);
1571: *ynorm = *delta;
1572: } else {
1573: *gpnorm = 0.0;
1574: *ynorm = nrm;
1575: }
1576: return(0);
1577: }
1579: /*@
1580: SNESSolve - Solves a nonlinear system. Call SNESSolve after calling
1581: SNESCreate() and optional routines of the form SNESSetXXX().
1583: Collective on SNES
1585: Input Parameters:
1586: + snes - the SNES context
1587: - x - the solution vector
1589: Output Parameter:
1590: . its - number of iterations until termination
1592: Notes:
1593: The user should initialize the vector,x, with the initial guess
1594: for the nonlinear solve prior to calling SNESSolve. In particular,
1595: to employ an initial guess of zero, the user should explicitly set
1596: this vector to zero by calling VecSet().
1598: Level: beginner
1600: .keywords: SNES, nonlinear, solve
1602: .seealso: SNESCreate(), SNESDestroy()
1603: @*/
1604: int SNESSolve(SNES snes,Vec x,int *its)
1605: {
1606: int ierr;
1607: PetscTruth flg;
1614: if (!snes->solve) SETERRQ(1,"SNESSetType() or SNESSetFromOptions() must be called before SNESSolve()");
1616: if (!snes->setupcalled) {SNESSetUp(snes,x);}
1617: else {snes->vec_sol = snes->vec_sol_always = x;}
1618: if (snes->conv_hist_reset == PETSC_TRUE) snes->conv_hist_len = 0;
1619: PetscLogEventBegin(SNES_Solve,snes,0,0,0);
1620: snes->nfuncs = 0; snes->linear_its = 0; snes->numFailures = 0;
1621: (*(snes)->solve)(snes,its);
1622: PetscLogEventEnd(SNES_Solve,snes,0,0,0);
1623: PetscOptionsHasName(snes->prefix,"-snes_view",&flg);
1624: if (flg && !PetscPreLoadingOn) { SNESView(snes,PETSC_VIEWER_STDOUT_WORLD); }
1625: return(0);
1626: }
1628: /* --------- Internal routines for SNES Package --------- */
1630: /*@C
1631: SNESSetType - Sets the method for the nonlinear solver.
1633: Collective on SNES
1635: Input Parameters:
1636: + snes - the SNES context
1637: - type - a known method
1639: Options Database Key:
1640: . -snes_type <type> - Sets the method; use -help for a list
1641: of available methods (for instance, ls or tr)
1643: Notes:
1644: See "petsc/include/petscsnes.h" for available methods (for instance)
1645: + SNESLS - Newton's method with line search
1646: (systems of nonlinear equations)
1647: . SNESTR - Newton's method with trust region
1648: (systems of nonlinear equations)
1650: Normally, it is best to use the SNESSetFromOptions() command and then
1651: set the SNES solver type from the options database rather than by using
1652: this routine. Using the options database provides the user with
1653: maximum flexibility in evaluating the many nonlinear solvers.
1654: The SNESSetType() routine is provided for those situations where it
1655: is necessary to set the nonlinear solver independently of the command
1656: line or options database. This might be the case, for example, when
1657: the choice of solver changes during the execution of the program,
1658: and the user's application is taking responsibility for choosing the
1659: appropriate method.
1661: Level: intermediate
1663: .keywords: SNES, set, type
1665: .seealso: SNESType, SNESCreate()
1667: @*/
1668: int SNESSetType(SNES snes,SNESType type)
1669: {
1670: int ierr,(*r)(SNES);
1671: PetscTruth match;
1677: PetscTypeCompare((PetscObject)snes,type,&match);
1678: if (match) return(0);
1680: if (snes->setupcalled) {
1681: ierr = (*(snes)->destroy)(snes);
1682: snes->data = 0;
1683: }
1685: /* Get the function pointers for the iterative method requested */
1686: if (!SNESRegisterAllCalled) {SNESRegisterAll(PETSC_NULL);}
1688: PetscFListFind(snes->comm,SNESList,type,(void (**)(void)) &r);
1690: if (!r) SETERRQ1(1,"Unable to find requested SNES type %s",type);
1692: if (snes->data) {PetscFree(snes->data);}
1693: snes->data = 0;
1694: (*r)(snes);
1695: PetscObjectChangeTypeName((PetscObject)snes,type);
1697: return(0);
1698: }
1701: /* --------------------------------------------------------------------- */
1702: /*@C
1703: SNESRegisterDestroy - Frees the list of nonlinear solvers that were
1704: registered by SNESRegisterDynamic().
1706: Not Collective
1708: Level: advanced
1710: .keywords: SNES, nonlinear, register, destroy
1712: .seealso: SNESRegisterAll(), SNESRegisterAll()
1713: @*/
1714: int SNESRegisterDestroy(void)
1715: {
1719: if (SNESList) {
1720: PetscFListDestroy(&SNESList);
1721: SNESList = 0;
1722: }
1723: SNESRegisterAllCalled = PETSC_FALSE;
1724: return(0);
1725: }
1727: /*@C
1728: SNESGetType - Gets the SNES method type and name (as a string).
1730: Not Collective
1732: Input Parameter:
1733: . snes - nonlinear solver context
1735: Output Parameter:
1736: . type - SNES method (a character string)
1738: Level: intermediate
1740: .keywords: SNES, nonlinear, get, type, name
1741: @*/
1742: int SNESGetType(SNES snes,SNESType *type)
1743: {
1746: *type = snes->type_name;
1747: return(0);
1748: }
1750: /*@C
1751: SNESGetSolution - Returns the vector where the approximate solution is
1752: stored.
1754: Not Collective, but Vec is parallel if SNES is parallel
1756: Input Parameter:
1757: . snes - the SNES context
1759: Output Parameter:
1760: . x - the solution
1762: Level: advanced
1764: .keywords: SNES, nonlinear, get, solution
1766: .seealso: SNESGetFunction(), SNESGetSolutionUpdate()
1767: @*/
1768: int SNESGetSolution(SNES snes,Vec *x)
1769: {
1772: *x = snes->vec_sol_always;
1773: return(0);
1774: }
1776: /*@C
1777: SNESGetSolutionUpdate - Returns the vector where the solution update is
1778: stored.
1780: Not Collective, but Vec is parallel if SNES is parallel
1782: Input Parameter:
1783: . snes - the SNES context
1785: Output Parameter:
1786: . x - the solution update
1788: Level: advanced
1790: .keywords: SNES, nonlinear, get, solution, update
1792: .seealso: SNESGetSolution(), SNESGetFunction
1793: @*/
1794: int SNESGetSolutionUpdate(SNES snes,Vec *x)
1795: {
1798: *x = snes->vec_sol_update_always;
1799: return(0);
1800: }
1802: /*@C
1803: SNESGetFunction - Returns the vector where the function is stored.
1805: Not Collective, but Vec is parallel if SNES is parallel
1807: Input Parameter:
1808: . snes - the SNES context
1810: Output Parameter:
1811: + r - the function (or PETSC_NULL)
1812: . ctx - the function context (or PETSC_NULL)
1813: - func - the function (or PETSC_NULL)
1815: Level: advanced
1817: .keywords: SNES, nonlinear, get, function
1819: .seealso: SNESSetFunction(), SNESGetSolution()
1820: @*/
1821: int SNESGetFunction(SNES snes,Vec *r,void **ctx,int (**func)(SNES,Vec,Vec,void*))
1822: {
1825: if (r) *r = snes->vec_func_always;
1826: if (ctx) *ctx = snes->funP;
1827: if (func) *func = snes->computefunction;
1828: return(0);
1829: }
1831: /*@C
1832: SNESSetOptionsPrefix - Sets the prefix used for searching for all
1833: SNES options in the database.
1835: Collective on SNES
1837: Input Parameter:
1838: + snes - the SNES context
1839: - prefix - the prefix to prepend to all option names
1841: Notes:
1842: A hyphen (-) must NOT be given at the beginning of the prefix name.
1843: The first character of all runtime options is AUTOMATICALLY the hyphen.
1845: Level: advanced
1847: .keywords: SNES, set, options, prefix, database
1849: .seealso: SNESSetFromOptions()
1850: @*/
1851: int SNESSetOptionsPrefix(SNES snes,char *prefix)
1852: {
1857: PetscObjectSetOptionsPrefix((PetscObject)snes,prefix);
1858: SLESSetOptionsPrefix(snes->sles,prefix);
1859: return(0);
1860: }
1862: /*@C
1863: SNESAppendOptionsPrefix - Appends to the prefix used for searching for all
1864: SNES options in the database.
1866: Collective on SNES
1868: Input Parameters:
1869: + snes - the SNES context
1870: - prefix - the prefix to prepend to all option names
1872: Notes:
1873: A hyphen (-) must NOT be given at the beginning of the prefix name.
1874: The first character of all runtime options is AUTOMATICALLY the hyphen.
1876: Level: advanced
1878: .keywords: SNES, append, options, prefix, database
1880: .seealso: SNESGetOptionsPrefix()
1881: @*/
1882: int SNESAppendOptionsPrefix(SNES snes,char *prefix)
1883: {
1885:
1888: PetscObjectAppendOptionsPrefix((PetscObject)snes,prefix);
1889: SLESAppendOptionsPrefix(snes->sles,prefix);
1890: return(0);
1891: }
1893: /*@C
1894: SNESGetOptionsPrefix - Sets the prefix used for searching for all
1895: SNES options in the database.
1897: Not Collective
1899: Input Parameter:
1900: . snes - the SNES context
1902: Output Parameter:
1903: . prefix - pointer to the prefix string used
1905: Notes: On the fortran side, the user should pass in a string 'prifix' of
1906: sufficient length to hold the prefix.
1908: Level: advanced
1910: .keywords: SNES, get, options, prefix, database
1912: .seealso: SNESAppendOptionsPrefix()
1913: @*/
1914: int SNESGetOptionsPrefix(SNES snes,char **prefix)
1915: {
1920: PetscObjectGetOptionsPrefix((PetscObject)snes,prefix);
1921: return(0);
1922: }
1924: /*MC
1925: SNESRegisterDynamic - Adds a method to the nonlinear solver package.
1927: Synopsis:
1928: int SNESRegisterDynamic(char *name_solver,char *path,char *name_create,int (*routine_create)(SNES))
1930: Not collective
1932: Input Parameters:
1933: + name_solver - name of a new user-defined solver
1934: . path - path (either absolute or relative) the library containing this solver
1935: . name_create - name of routine to create method context
1936: - routine_create - routine to create method context
1938: Notes:
1939: SNESRegisterDynamic() may be called multiple times to add several user-defined solvers.
1941: If dynamic libraries are used, then the fourth input argument (routine_create)
1942: is ignored.
1944: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR}, ${BOPT},
1945: and others of the form ${any_environmental_variable} occuring in pathname will be
1946: replaced with appropriate values.
1948: Sample usage:
1949: .vb
1950: SNESRegisterDynamic("my_solver",/home/username/my_lib/lib/libg/solaris/mylib.a,
1951: "MySolverCreate",MySolverCreate);
1952: .ve
1954: Then, your solver can be chosen with the procedural interface via
1955: $ SNESSetType(snes,"my_solver")
1956: or at runtime via the option
1957: $ -snes_type my_solver
1959: Level: advanced
1961: .keywords: SNES, nonlinear, register
1963: .seealso: SNESRegisterAll(), SNESRegisterDestroy()
1964: M*/
1966: int SNESRegister(char *sname,char *path,char *name,int (*function)(SNES))
1967: {
1968: char fullname[256];
1969: int ierr;
1972: PetscFListConcat(path,name,fullname);
1973: PetscFListAdd(&SNESList,sname,fullname,(void (*)(void))function);
1974: return(0);
1975: }