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: }