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