Actual source code: snes.c

  1: /*$Id: snes.c,v 1.229 2001/04/10 19:36:48 bsmith Exp $*/

  3: #include "src/snes/snesimpl.h"      /*I "petscsnes.h"  I*/

  5: PetscTruth SNESRegisterAllCalled = PETSC_FALSE;
  6: PetscFList      SNESList = 0;

  8: /*@C
  9:    SNESView - Prints the SNES data structure.

 11:    Collective on SNES

 13:    Input Parameters:
 14: +  SNES - the SNES context
 15: -  viewer - visualization context

 17:    Options Database Key:
 18: .  -snes_view - Calls SNESView() at end of SNESSolve()

 20:    Notes:
 21:    The available visualization contexts include
 22: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
 23: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
 24:          output where only the first processor opens
 25:          the file.  All other processors send their 
 26:          data to the first processor to print. 

 28:    The user can open an alternative visualization context with
 29:    PetscViewerASCIIOpen() - output to a specified file.

 31:    Level: beginner

 33: .keywords: SNES, view

 35: .seealso: PetscViewerASCIIOpen()
 36: @*/
 37: int SNESView(SNES snes,PetscViewer viewer)
 38: {
 39:   SNES_KSP_EW_ConvCtx *kctx;
 40:   int                 ierr;
 41:   SLES                sles;
 42:   char                *type;
 43:   PetscTruth          isascii,isstring;

 47:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(snes->comm);

 51:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
 52:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
 53:   if (isascii) {
 54:     PetscViewerASCIIPrintf(viewer,"SNES Object:n");
 55:     SNESGetType(snes,&type);
 56:     if (type) {
 57:       PetscViewerASCIIPrintf(viewer,"  type: %sn",type);
 58:     } else {
 59:       PetscViewerASCIIPrintf(viewer,"  type: not set yetn");
 60:     }
 61:     if (snes->view) {
 62:       PetscViewerASCIIPushTab(viewer);
 63:       (*snes->view)(snes,viewer);
 64:       PetscViewerASCIIPopTab(viewer);
 65:     }
 66:     PetscViewerASCIIPrintf(viewer,"  maximum iterations=%d, maximum function evaluations=%dn",snes->max_its,snes->max_funcs);
 67:     PetscViewerASCIIPrintf(viewer,"  tolerances: relative=%g, absolute=%g, truncation=%g, solution=%gn",
 68:                  snes->rtol,snes->atol,snes->trunctol,snes->xtol);
 69:     PetscViewerASCIIPrintf(viewer,"  total number of linear solver iterations=%dn",snes->linear_its);
 70:     PetscViewerASCIIPrintf(viewer,"  total number of function evaluations=%dn",snes->nfuncs);
 71:     if (snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) {
 72:       PetscViewerASCIIPrintf(viewer,"  min function tolerance=%gn",snes->fmin);
 73:     }
 74:     if (snes->ksp_ewconv) {
 75:       kctx = (SNES_KSP_EW_ConvCtx *)snes->kspconvctx;
 76:       if (kctx) {
 77:         PetscViewerASCIIPrintf(viewer,"  Eisenstat-Walker computation of KSP relative tolerance (version %d)n",kctx->version);
 78:         PetscViewerASCIIPrintf(viewer,"    rtol_0=%g, rtol_max=%g, threshold=%gn",kctx->rtol_0,kctx->rtol_max,kctx->threshold);
 79:         PetscViewerASCIIPrintf(viewer,"    gamma=%g, alpha=%g, alpha2=%gn",kctx->gamma,kctx->alpha,kctx->alpha2);
 80:       }
 81:     }
 82:   } else if (isstring) {
 83:     SNESGetType(snes,&type);
 84:     PetscViewerStringSPrintf(viewer," %-3.3s",type);
 85:   }
 86:   SNESGetSLES(snes,&sles);
 87:   PetscViewerASCIIPushTab(viewer);
 88:   SLESView(sles,viewer);
 89:   PetscViewerASCIIPopTab(viewer);
 90:   return(0);
 91: }

 93: /*@
 94:    SNESSetFromOptions - Sets various SNES and SLES parameters from user options.

 96:    Collective on SNES

 98:    Input Parameter:
 99: .  snes - the SNES context

101:    Options Database Keys:
102: +  -snes_type <type> - ls, tr, umls, umtr, test
103: .  -snes_stol - convergence tolerance in terms of the norm
104:                 of the change in the solution between steps
105: .  -snes_atol <atol> - absolute tolerance of residual norm
106: .  -snes_rtol <rtol> - relative decrease in tolerance norm from initial
107: .  -snes_max_it <max_it> - maximum number of iterations
108: .  -snes_max_funcs <max_funcs> - maximum number of function evaluations
109: .  -snes_trtol <trtol> - trust region tolerance
110: .  -snes_no_convergence_test - skip convergence test in nonlinear or minimization 
111:                                solver; hence iterations will continue until max_it
112:                                or some other criterion is reached. Saves expense
113:                                of convergence test
114: .  -snes_monitor - prints residual norm at each iteration 
115: .  -snes_vecmonitor - plots solution at each iteration
116: .  -snes_vecmonitor_update - plots update to solution at each iteration 
117: .  -snes_xmonitor - plots residual norm at each iteration 
118: .  -snes_fd - use finite differences to compute Jacobian; very slow, only for testing
119: -  -snes_mf_ksp_monitor - if using matrix-free multiply then print h at each KSP iteration

121:     Options Database for Eisenstat-Walker method:
122: +  -snes_ksp_eq_conv - use Eisenstat-Walker method for determining linear system convergence
123: .  -snes_ksp_eq_version ver - version of  Eisenstat-Walker method
124: .  -snes_ksp_ew_rtol0 <rtol0> - Sets rtol0
125: .  -snes_ksp_ew_rtolmax <rtolmax> - Sets rtolmax
126: .  -snes_ksp_ew_gamma <gamma> - Sets gamma
127: .  -snes_ksp_ew_alpha <alpha> - Sets alpha
128: .  -snes_ksp_ew_alpha2 <alpha2> - Sets alpha2 
129: -  -snes_ksp_ew_threshold <threshold> - Sets threshold

131:    Notes:
132:    To see all options, run your program with the -help option or consult
133:    the users manual.

135:    Level: beginner

137: .keywords: SNES, nonlinear, set, options, database

139: .seealso: SNESSetOptionsPrefix(), SNESSetTypeFromOptions()
140: @*/
141: int SNESSetFromOptions(SNES snes)
142: {
143:   SLES                sles;
144:   SNES_KSP_EW_ConvCtx *kctx = (SNES_KSP_EW_ConvCtx *)snes->kspconvctx;
145:   PetscTruth          flg;
146:   int                 ierr;
147:   char                *deft,type[256];


152:   PetscOptionsBegin(snes->comm,snes->prefix,"Nonlinear solver (SNES) options","SNES");
153:     if (snes->type_name) {
154:       deft = snes->type_name;
155:     } else {
156:       if (snes->method_class == SNES_NONLINEAR_EQUATIONS) {
157:         deft = SNESEQLS;
158:       } else {
159:         deft = SNESUMTR;
160:       }
161:     }

163:     if (!SNESRegisterAllCalled) {SNESRegisterAll(PETSC_NULL);}
164:     PetscOptionsList("-snes_type","Nonlinear solver method","SNESSetType",SNESList,deft,type,256,&flg);
165:     if (flg) {
166:       SNESSetType(snes,type);
167:     } else if (!snes->type_name) {
168:       SNESSetType(snes,deft);
169:     }

171:     PetscOptionsDouble("-snes_stol","Stop if step length less then","SNESSetTolerances",snes->xtol,&snes->xtol,0);
172:     PetscOptionsDouble("-snes_atol","Stop if function norm less then","SNESSetTolerances",snes->atol,&snes->atol,0);

174:     PetscOptionsDouble("-snes_rtol","Stop if decrease in function norm less then","SNESSetTolerances",snes->rtol,&snes->rtol,0);
175:     PetscOptionsInt("-snes_max_it","Maximum iterations","SNESSetTolerances",snes->max_its,&snes->max_its,PETSC_NULL);
176:     PetscOptionsInt("-snes_max_funcs","Maximum function evaluations","SNESSetTolerances",snes->max_funcs,&snes->max_funcs,PETSC_NULL);
177:     PetscOptionsDouble("-snes_fmin","Minimization function tolerance","SNESSetMinimizationFunctionTolerance",snes->fmin,&snes->fmin,0);

179:     PetscOptionsName("-snes_ksp_ew_conv","Use Eisentat-Walker linear system convergence test","SNES_KSP_SetParametersEW",&snes->ksp_ewconv);

181:     PetscOptionsInt("-snes_ksp_ew_version","Version 1 or 2","SNES_KSP_SetParametersEW",kctx->version,&kctx->version,0);
182:     PetscOptionsDouble("-snes_ksp_ew_rtol0","0 <= rtol0 < 1","SNES_KSP_SetParametersEW",kctx->rtol_0,&kctx->rtol_0,0);
183:     PetscOptionsDouble("-snes_ksp_ew_rtolmax","0 <= rtolmax < 1","SNES_KSP_SetParametersEW",kctx->rtol_max,&kctx->rtol_max,0);
184:     PetscOptionsDouble("-snes_ksp_ew_gamma","0 <= gamma <= 1","SNES_KSP_SetParametersEW",kctx->gamma,&kctx->gamma,0);
185:     PetscOptionsDouble("-snes_ksp_ew_alpha","1 < alpha <= 2","SNES_KSP_SetParametersEW",kctx->alpha,&kctx->alpha,0);
186:     PetscOptionsDouble("-snes_ksp_ew_alpha2","alpha2","SNES_KSP_SetParametersEW",kctx->alpha2,&kctx->alpha2,0);
187:     PetscOptionsDouble("-snes_ksp_ew_threshold","0 < threshold < 1","SNES_KSP_SetParametersEW",kctx->threshold,&kctx->threshold,0);

189:     PetscOptionsName("-snes_no_convergence_test","Don't test for convergence","None",&flg);
190:     if (flg) {snes->converged = 0;}
191:     PetscOptionsName("-snes_cancelmonitors","Remove all monitors","SNESClearMonitor",&flg);
192:     if (flg) {SNESClearMonitor(snes);}
193:     PetscOptionsName("-snes_monitor","Monitor norm of function","SNESDefaultMonitor",&flg);
194:     if (flg) {SNESSetMonitor(snes,SNESDefaultMonitor,0,0);}
195:     PetscOptionsName("-snes_smonitor","Monitor norm of function (fewer digits)","SNESDefaultSMonitor",&flg);
196:     if (flg) {SNESSetMonitor(snes,SNESDefaultSMonitor,0,0);}
197:     PetscOptionsName("-snes_vecmonitor","Plot solution at each iteration","SNESVecViewMonitor",&flg);
198:     if (flg) {SNESSetMonitor(snes,SNESVecViewMonitor,0,0);}
199:     PetscOptionsName("-snes_vecmonitor_update","Plot correction at each iteration","SNESVecViewUpdateMonitor",&flg);
200:     if (flg) {SNESSetMonitor(snes,SNESVecViewUpdateMonitor,0,0);}
201:     PetscOptionsName("-snes_xmonitor","Plot function norm at each iteration","SNESLGMonitor",&flg);
202:     if (flg) {SNESSetMonitor(snes,SNESLGMonitor,PETSC_NULL,PETSC_NULL);}

204:     PetscOptionsName("-snes_fd","Use finite differences (slow) to compute Jacobian","SNESDefaultComputeJacobian",&flg);
205:     if (flg && snes->method_class == SNES_NONLINEAR_EQUATIONS) {
206:       SNESSetJacobian(snes,snes->jacobian,snes->jacobian_pre,SNESDefaultComputeJacobian,snes->funP);
207:       PetscLogInfo(snes,"SNESSetFromOptions: Setting default finite difference Jacobian matrixn");
208:     } else if (flg && snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) {
209:       SNESSetHessian(snes,snes->jacobian,snes->jacobian_pre,SNESDefaultComputeHessian,snes->funP);
210:       PetscLogInfo(snes,"SNESSetFromOptions: Setting default finite difference Hessian matrixn");
211:     }

213:     if (snes->setfromoptions) {
214:       (*snes->setfromoptions)(snes);
215:     }

217:   PetscOptionsEnd();

219:   SNESGetSLES(snes,&sles);
220:   SLESSetFromOptions(sles);

222:   return(0);
223: }


226: /*@
227:    SNESSetApplicationContext - Sets the optional user-defined context for 
228:    the nonlinear solvers.  

230:    Collective on SNES

232:    Input Parameters:
233: +  snes - the SNES context
234: -  usrP - optional user context

236:    Level: intermediate

238: .keywords: SNES, nonlinear, set, application, context

240: .seealso: SNESGetApplicationContext()
241: @*/
242: int SNESSetApplicationContext(SNES snes,void *usrP)
243: {
246:   snes->user                = usrP;
247:   return(0);
248: }

250: /*@C
251:    SNESGetApplicationContext - Gets the user-defined context for the 
252:    nonlinear solvers.  

254:    Not Collective

256:    Input Parameter:
257: .  snes - SNES context

259:    Output Parameter:
260: .  usrP - user context

262:    Level: intermediate

264: .keywords: SNES, nonlinear, get, application, context

266: .seealso: SNESSetApplicationContext()
267: @*/
268: int SNESGetApplicationContext(SNES snes,void **usrP)
269: {
272:   *usrP = snes->user;
273:   return(0);
274: }

276: /*@
277:    SNESGetIterationNumber - Gets the number of nonlinear iterations completed
278:    at this time.

280:    Not Collective

282:    Input Parameter:
283: .  snes - SNES context

285:    Output Parameter:
286: .  iter - iteration number

288:    Notes:
289:    For example, during the computation of iteration 2 this would return 1.

291:    This is useful for using lagged Jacobians (where one does not recompute the 
292:    Jacobian at each SNES iteration). For example, the code
293: .vb
294:       SNESGetIterationNumber(snes,&it);
295:       if (!(it % 2)) {
296:         [compute Jacobian here]
297:       }
298: .ve
299:    can be used in your ComputeJacobian() function to cause the Jacobian to be
300:    recomputed every second SNES iteration.

302:    Level: intermediate

304: .keywords: SNES, nonlinear, get, iteration, number
305: @*/
306: int SNESGetIterationNumber(SNES snes,int* iter)
307: {
311:   *iter = snes->iter;
312:   return(0);
313: }

315: /*@
316:    SNESGetFunctionNorm - Gets the norm of the current function that was set
317:    with SNESSSetFunction().

319:    Collective on SNES

321:    Input Parameter:
322: .  snes - SNES context

324:    Output Parameter:
325: .  fnorm - 2-norm of function

327:    Note:
328:    SNESGetFunctionNorm() is valid for SNES_NONLINEAR_EQUATIONS methods only.
329:    A related routine for SNES_UNCONSTRAINED_MINIMIZATION methods is
330:    SNESGetGradientNorm().

332:    Level: intermediate

334: .keywords: SNES, nonlinear, get, function, norm

336: .seealso: SNESGetFunction()
337: @*/
338: int SNESGetFunctionNorm(SNES snes,Scalar *fnorm)
339: {
343:   if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
344:     SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"For SNES_NONLINEAR_EQUATIONS only");
345:   }
346:   *fnorm = snes->norm;
347:   return(0);
348: }

350: /*@
351:    SNESGetGradientNorm - Gets the norm of the current gradient that was set
352:    with SNESSSetGradient().

354:    Collective on SNES

356:    Input Parameter:
357: .  snes - SNES context

359:    Output Parameter:
360: .  fnorm - 2-norm of gradient

362:    Note:
363:    SNESGetGradientNorm() is valid for SNES_UNCONSTRAINED_MINIMIZATION 
364:    methods only.  A related routine for SNES_NONLINEAR_EQUATIONS methods
365:    is SNESGetFunctionNorm().

367:    Level: intermediate

369: .keywords: SNES, nonlinear, get, gradient, norm

371: .seelso: SNESSetGradient()
372: @*/
373: int SNESGetGradientNorm(SNES snes,Scalar *gnorm)
374: {
378:   if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
379:     SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"For SNES_UNCONSTRAINED_MINIMIZATION only");
380:   }
381:   *gnorm = snes->norm;
382:   return(0);
383: }

385: /*@
386:    SNESGetNumberUnsuccessfulSteps - Gets the number of unsuccessful steps
387:    attempted by the nonlinear solver.

389:    Not Collective

391:    Input Parameter:
392: .  snes - SNES context

394:    Output Parameter:
395: .  nfails - number of unsuccessful steps attempted

397:    Notes:
398:    This counter is reset to zero for each successive call to SNESSolve().

400:    Level: intermediate

402: .keywords: SNES, nonlinear, get, number, unsuccessful, steps
403: @*/
404: int SNESGetNumberUnsuccessfulSteps(SNES snes,int* nfails)
405: {
409:   *nfails = snes->nfailures;
410:   return(0);
411: }

413: /*@
414:    SNESGetNumberLinearIterations - Gets the total number of linear iterations
415:    used by the nonlinear solver.

417:    Not Collective

419:    Input Parameter:
420: .  snes - SNES context

422:    Output Parameter:
423: .  lits - number of linear iterations

425:    Notes:
426:    This counter is reset to zero for each successive call to SNESSolve().

428:    Level: intermediate

430: .keywords: SNES, nonlinear, get, number, linear, iterations
431: @*/
432: int SNESGetNumberLinearIterations(SNES snes,int* lits)
433: {
437:   *lits = snes->linear_its;
438:   return(0);
439: }

441: /*@C
442:    SNESGetSLES - Returns the SLES context for a SNES solver.

444:    Not Collective, but if SNES object is parallel, then SLES object is parallel

446:    Input Parameter:
447: .  snes - the SNES context

449:    Output Parameter:
450: .  sles - the SLES context

452:    Notes:
453:    The user can then directly manipulate the SLES context to set various
454:    options, etc.  Likewise, the user can then extract and manipulate the 
455:    KSP and PC contexts as well.

457:    Level: beginner

459: .keywords: SNES, nonlinear, get, SLES, context

461: .seealso: SLESGetPC(), SLESGetKSP()
462: @*/
463: int SNESGetSLES(SNES snes,SLES *sles)
464: {
467:   *sles = snes->sles;
468:   return(0);
469: }

471: static int SNESPublish_Petsc(PetscObject obj)
472: {
473: #if defined(PETSC_HAVE_AMS)
474:   SNES          v = (SNES) obj;
475:   int          ierr;
476: #endif


480: #if defined(PETSC_HAVE_AMS)
481:   /* if it is already published then return */
482:   if (v->amem >=0) return(0);

484:   PetscObjectPublishBaseBegin(obj);
485:   AMS_Memory_add_field((AMS_Memory)v->amem,"Iteration",&v->iter,1,AMS_INT,AMS_READ,
486:                                 AMS_COMMON,AMS_REDUCT_UNDEF);
487:   AMS_Memory_add_field((AMS_Memory)v->amem,"Residual",&v->norm,1,AMS_DOUBLE,AMS_READ,
488:                                 AMS_COMMON,AMS_REDUCT_UNDEF);
489:   PetscObjectPublishBaseEnd(obj);
490: #endif
491:   return(0);
492: }

494: /* -----------------------------------------------------------*/
495: /*@C
496:    SNESCreate - Creates a nonlinear solver context.

498:    Collective on MPI_Comm

500:    Input Parameters:
501: +  comm - MPI communicator
502: -  type - type of method, either 
503:    SNES_NONLINEAR_EQUATIONS (for systems of nonlinear equations) 
504:    or SNES_UNCONSTRAINED_MINIMIZATION (for unconstrained minimization)

506:    Output Parameter:
507: .  outsnes - the new SNES context

509:    Options Database Keys:
510: +   -snes_mf - Activates default matrix-free Jacobian-vector products,
511:                and no preconditioning matrix
512: .   -snes_mf_operator - Activates default matrix-free Jacobian-vector
513:                products, and a user-provided preconditioning matrix
514:                as set by SNESSetJacobian()
515: -   -snes_fd - Uses (slow!) finite differences to compute Jacobian

517:    Level: beginner

519: .keywords: SNES, nonlinear, create, context

521: .seealso: SNESSolve(), SNESDestroy(), SNESProblemType, SNES
522: @*/
523: int SNESCreate(MPI_Comm comm,SNESProblemType type,SNES *outsnes)
524: {
525:   int                 ierr;
526:   SNES                snes;
527:   SNES_KSP_EW_ConvCtx *kctx;

530:   *outsnes = 0;
531:   if (type != SNES_UNCONSTRAINED_MINIMIZATION && type != SNES_NONLINEAR_EQUATIONS){
532:     SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"incorrect method type");
533:   }
534:   PetscHeaderCreate(snes,_p_SNES,int,SNES_COOKIE,0,"SNES",comm,SNESDestroy,SNESView);
535:   PetscLogObjectCreate(snes);
536:   snes->bops->publish     = SNESPublish_Petsc;
537:   snes->max_its           = 50;
538:   snes->max_funcs          = 10000;
539:   snes->norm                  = 0.0;
540:   if (type == SNES_UNCONSTRAINED_MINIMIZATION) {
541:     snes->rtol                  = 1.e-8;
542:     snes->ttol            = 0.0;
543:     snes->atol                  = 1.e-10;
544:   } else {
545:     snes->rtol                  = 1.e-8;
546:     snes->ttol            = 0.0;
547:     snes->atol                  = 1.e-50;
548:   }
549:   snes->xtol                  = 1.e-8;
550:   snes->trunctol          = 1.e-12; /* no longer used */
551:   snes->nfuncs            = 0;
552:   snes->nfailures         = 0;
553:   snes->linear_its        = 0;
554:   snes->numbermonitors    = 0;
555:   snes->data              = 0;
556:   snes->view              = 0;
557:   snes->computeumfunction = 0;
558:   snes->umfunP            = 0;
559:   snes->fc                = 0;
560:   snes->deltatol          = 1.e-12;
561:   snes->fmin              = -1.e30;
562:   snes->method_class      = type;
563:   snes->set_method_called = 0;
564:   snes->setupcalled       = 0;
565:   snes->ksp_ewconv        = PETSC_FALSE;
566:   snes->vwork             = 0;
567:   snes->nwork             = 0;
568:   snes->conv_hist_len     = 0;
569:   snes->conv_hist_max     = 0;
570:   snes->conv_hist         = PETSC_NULL;
571:   snes->conv_hist_its     = PETSC_NULL;
572:   snes->conv_hist_reset   = PETSC_TRUE;
573:   snes->reason            = SNES_CONVERGED_ITERATING;

575:   /* Create context to compute Eisenstat-Walker relative tolerance for KSP */
576:   PetscNew(SNES_KSP_EW_ConvCtx,&kctx);
577:   PetscLogObjectMemory(snes,sizeof(SNES_KSP_EW_ConvCtx));
578:   snes->kspconvctx  = (void*)kctx;
579:   kctx->version     = 2;
580:   kctx->rtol_0      = .3; /* Eisenstat and Walker suggest rtol_0=.5, but 
581:                              this was too large for some test cases */
582:   kctx->rtol_last   = 0;
583:   kctx->rtol_max    = .9;
584:   kctx->gamma       = 1.0;
585:   kctx->alpha2      = .5*(1.0 + sqrt(5.0));
586:   kctx->alpha       = kctx->alpha2;
587:   kctx->threshold   = .1;
588:   kctx->lresid_last = 0;
589:   kctx->norm_last   = 0;

591:   SLESCreate(comm,&snes->sles);
592:   PetscLogObjectParent(snes,snes->sles)

594:   *outsnes = snes;
595:   PetscPublishAll(snes);
596:   return(0);
597: }

599: /* --------------------------------------------------------------- */
600: /*@C
601:    SNESSetFunction - Sets the function evaluation routine and function 
602:    vector for use by the SNES routines in solving systems of nonlinear
603:    equations.

605:    Collective on SNES

607:    Input Parameters:
608: +  snes - the SNES context
609: .  func - function evaluation routine
610: .  r - vector to store function value
611: -  ctx - [optional] user-defined context for private data for the 
612:          function evaluation routine (may be PETSC_NULL)

614:    Calling sequence of func:
615: $    func (SNES snes,Vec x,Vec f,void *ctx);

617: .  f - function vector
618: -  ctx - optional user-defined function context 

620:    Notes:
621:    The Newton-like methods typically solve linear systems of the form
622: $      f'(x) x = -f(x),
623:    where f'(x) denotes the Jacobian matrix and f(x) is the function.

625:    SNESSetFunction() is valid for SNES_NONLINEAR_EQUATIONS methods only.
626:    Analogous routines for SNES_UNCONSTRAINED_MINIMIZATION methods are
627:    SNESSetMinimizationFunction() and SNESSetGradient();

629:    Level: beginner

631: .keywords: SNES, nonlinear, set, function

633: .seealso: SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian()
634: @*/
635: int SNESSetFunction(SNES snes,Vec r,int (*func)(SNES,Vec,Vec,void*),void *ctx)
636: {
641:   if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
642:     SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
643:   }

645:   snes->computefunction     = func;
646:   snes->vec_func            = snes->vec_func_always = r;
647:   snes->funP                = ctx;
648:   return(0);
649: }

651: /*@
652:    SNESComputeFunction - Calls the function that has been set with
653:                          SNESSetFunction().  

655:    Collective on SNES

657:    Input Parameters:
658: +  snes - the SNES context
659: -  x - input vector

661:    Output Parameter:
662: .  y - function vector, as set by SNESSetFunction()

664:    Notes:
665:    SNESComputeFunction() is valid for SNES_NONLINEAR_EQUATIONS methods only.
666:    Analogous routines for SNES_UNCONSTRAINED_MINIMIZATION methods are
667:    SNESComputeMinimizationFunction() and SNESComputeGradient();

669:    SNESComputeFunction() is typically used within nonlinear solvers
670:    implementations, so most users would not generally call this routine
671:    themselves.

673:    Level: developer

675: .keywords: SNES, nonlinear, compute, function

677: .seealso: SNESSetFunction(), SNESGetFunction()
678: @*/
679: int SNESComputeFunction(SNES snes,Vec x,Vec y)
680: {
681:   int    ierr;

689:   if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
690:     SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
691:   }

693:   PetscLogEventBegin(SNES_FunctionEval,snes,x,y,0);
694:   PetscStackPush("SNES user function");
695:   (*snes->computefunction)(snes,x,y,snes->funP);
696:   PetscStackPop;
697:   snes->nfuncs++;
698:   PetscLogEventEnd(SNES_FunctionEval,snes,x,y,0);
699:   return(0);
700: }

702: /*@C
703:    SNESSetMinimizationFunction - Sets the function evaluation routine for 
704:    unconstrained minimization.

706:    Collective on SNES

708:    Input Parameters:
709: +  snes - the SNES context
710: .  func - function evaluation routine
711: -  ctx - [optional] user-defined context for private data for the 
712:          function evaluation routine (may be PETSC_NULL)

714:    Calling sequence of func:
715: $     func (SNES snes,Vec x,PetscReal *f,void *ctx);

717: +  x - input vector
718: .  f - function
719: -  ctx - [optional] user-defined function context 

721:    Level: beginner

723:    Notes:
724:    SNESSetMinimizationFunction() is valid for SNES_UNCONSTRAINED_MINIMIZATION
725:    methods only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is
726:    SNESSetFunction().

728: .keywords: SNES, nonlinear, set, minimization, function

730: .seealso:  SNESGetMinimizationFunction(), SNESComputeMinimizationFunction(),
731:            SNESSetHessian(), SNESSetGradient()
732: @*/
733: int SNESSetMinimizationFunction(SNES snes,int (*func)(SNES,Vec,PetscReal*,void*),void *ctx)
734: {
737:   if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
738:     SETERRQ(PETSC_ERR_ARG_WRONG,"Only for SNES_UNCONSTRAINED_MINIMIZATION");
739:   }
740:   snes->computeumfunction   = func;
741:   snes->umfunP              = ctx;
742:   return(0);
743: }

745: /*@
746:    SNESComputeMinimizationFunction - Computes the function that has been
747:    set with SNESSetMinimizationFunction().

749:    Collective on SNES

751:    Input Parameters:
752: +  snes - the SNES context
753: -  x - input vector

755:    Output Parameter:
756: .  y - function value

758:    Notes:
759:    SNESComputeMinimizationFunction() is valid only for 
760:    SNES_UNCONSTRAINED_MINIMIZATION methods. An analogous routine for 
761:    SNES_NONLINEAR_EQUATIONS methods is SNESComputeFunction().

763:    SNESComputeMinimizationFunction() is typically used within minimization
764:    implementations, so most users would not generally call this routine
765:    themselves.

767:    Level: developer

769: .keywords: SNES, nonlinear, compute, minimization, function

771: .seealso: SNESSetMinimizationFunction(), SNESGetMinimizationFunction(),
772:           SNESComputeGradient(), SNESComputeHessian()
773: @*/
774: int SNESComputeMinimizationFunction(SNES snes,Vec x,PetscReal *y)
775: {
776:   int    ierr;

782:   if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
783:     SETERRQ(PETSC_ERR_ARG_WRONG,"Only for SNES_UNCONSTRAINED_MINIMIZATION");
784:   }

786:   PetscLogEventBegin(SNES_MinimizationFunctionEval,snes,x,y,0);
787:   PetscStackPush("SNES user minimzation function");
788:   (*snes->computeumfunction)(snes,x,y,snes->umfunP);
789:   PetscStackPop;
790:   snes->nfuncs++;
791:   PetscLogEventEnd(SNES_MinimizationFunctionEval,snes,x,y,0);
792:   return(0);
793: }

795: /*@C
796:    SNESSetGradient - Sets the gradient evaluation routine and gradient
797:    vector for use by the SNES routines.

799:    Collective on SNES

801:    Input Parameters:
802: +  snes - the SNES context
803: .  func - function evaluation routine
804: .  ctx - optional user-defined context for private data for the 
805:          gradient evaluation routine (may be PETSC_NULL)
806: -  r - vector to store gradient value

808:    Calling sequence of func:
809: $     func (SNES, Vec x, Vec g, void *ctx);

811: +  x - input vector
812: .  g - gradient vector
813: -  ctx - optional user-defined gradient context 

815:    Notes:
816:    SNESSetMinimizationFunction() is valid for SNES_UNCONSTRAINED_MINIMIZATION
817:    methods only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is
818:    SNESSetFunction().

820:    Level: beginner

822: .keywords: SNES, nonlinear, set, function

824: .seealso: SNESGetGradient(), SNESComputeGradient(), SNESSetHessian(),
825:           SNESSetMinimizationFunction(),
826: @*/
827: int SNESSetGradient(SNES snes,Vec r,int (*func)(SNES,Vec,Vec,void*),void *ctx)
828: {
833:   if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
834:     SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
835:   }
836:   snes->computefunction     = func;
837:   snes->vec_func            = snes->vec_func_always = r;
838:   snes->funP                = ctx;
839:   return(0);
840: }

842: /*@
843:    SNESComputeGradient - Computes the gradient that has been set with
844:    SNESSetGradient().

846:    Collective on SNES

848:    Input Parameters:
849: +  snes - the SNES context
850: -  x - input vector

852:    Output Parameter:
853: .  y - gradient vector

855:    Notes:
856:    SNESComputeGradient() is valid only for 
857:    SNES_UNCONSTRAINED_MINIMIZATION methods. An analogous routine for 
858:    SNES_NONLINEAR_EQUATIONS methods is SNESComputeFunction().

860:    SNESComputeGradient() is typically used within minimization
861:    implementations, so most users would not generally call this routine
862:    themselves.

864:    Level: developer

866: .keywords: SNES, nonlinear, compute, gradient

868: .seealso:  SNESSetGradient(), SNESGetGradient(), 
869:            SNESComputeMinimizationFunction(), SNESComputeHessian()
870: @*/
871: int SNESComputeGradient(SNES snes,Vec x,Vec y)
872: {
873:   int    ierr;

881:   if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
882:     SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
883:   }
884:   PetscLogEventBegin(SNES_GradientEval,snes,x,y,0);
885:   PetscStackPush("SNES user gradient function");
886:   (*snes->computefunction)(snes,x,y,snes->funP);
887:   PetscStackPop;
888:   PetscLogEventEnd(SNES_GradientEval,snes,x,y,0);
889:   return(0);
890: }

892: /*@
893:    SNESComputeJacobian - Computes the Jacobian matrix that has been
894:    set with SNESSetJacobian().

896:    Collective on SNES and Mat

898:    Input Parameters:
899: +  snes - the SNES context
900: -  x - input vector

902:    Output Parameters:
903: +  A - Jacobian matrix
904: .  B - optional preconditioning matrix
905: -  flag - flag indicating matrix structure

907:    Notes: 
908:    Most users should not need to explicitly call this routine, as it 
909:    is used internally within the nonlinear solvers. 

911:    See SLESSetOperators() for important information about setting the
912:    flag parameter.

914:    SNESComputeJacobian() is valid only for SNES_NONLINEAR_EQUATIONS
915:    methods. An analogous routine for SNES_UNCONSTRAINED_MINIMIZATION 
916:    methods is SNESComputeHessian().

918:    Level: developer

920: .keywords: SNES, compute, Jacobian, matrix

922: .seealso:  SNESSetJacobian(), SLESSetOperators()
923: @*/
924: int SNESComputeJacobian(SNES snes,Vec X,Mat *A,Mat *B,MatStructure *flg)
925: {
926:   int    ierr;

932:   if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
933:     SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
934:   }
935:   if (!snes->computejacobian) return(0);
936:   PetscLogEventBegin(SNES_JacobianEval,snes,X,*A,*B);
937:   *flg = DIFFERENT_NONZERO_PATTERN;
938:   PetscStackPush("SNES user Jacobian function");
939:   (*snes->computejacobian)(snes,X,A,B,flg,snes->jacP);
940:   PetscStackPop;
941:   PetscLogEventEnd(SNES_JacobianEval,snes,X,*A,*B);
942:   /* make sure user returned a correct Jacobian and preconditioner */
945:   return(0);
946: }

948: /*@
949:    SNESComputeHessian - Computes the Hessian matrix that has been
950:    set with SNESSetHessian().

952:    Collective on SNES and Mat

954:    Input Parameters:
955: +  snes - the SNES context
956: -  x - input vector

958:    Output Parameters:
959: +  A - Hessian matrix
960: .  B - optional preconditioning matrix
961: -  flag - flag indicating matrix structure

963:    Notes: 
964:    Most users should not need to explicitly call this routine, as it
965:    is used internally within the nonlinear solvers. 

967:    See SLESSetOperators() for important information about setting the
968:    flag parameter.

970:    SNESComputeHessian() is valid only for 
971:    SNES_UNCONSTRAINED_MINIMIZATION methods. An analogous routine for 
972:    SNES_NONLINEAR_EQUATIONS methods is SNESComputeJacobian().

974:    SNESComputeHessian() is typically used within minimization
975:    implementations, so most users would not generally call this routine
976:    themselves.

978:    Level: developer

980: .keywords: SNES, compute, Hessian, matrix

982: .seealso:  SNESSetHessian(), SLESSetOperators(), SNESComputeGradient(),
983:            SNESComputeMinimizationFunction()
984: @*/
985: int SNESComputeHessian(SNES snes,Vec x,Mat *A,Mat *B,MatStructure *flag)
986: {
987:   int    ierr;

993:   if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
994:     SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
995:   }
996:   if (!snes->computejacobian) return(0);
997:   PetscLogEventBegin(SNES_HessianEval,snes,x,*A,*B);
998:   *flag = DIFFERENT_NONZERO_PATTERN;
999:   PetscStackPush("SNES user Hessian function");
1000:   (*snes->computejacobian)(snes,x,A,B,flag,snes->jacP);
1001:   PetscStackPop;
1002:   PetscLogEventEnd(SNES_HessianEval,snes,x,*A,*B);
1003:   /* make sure user returned a correct Jacobian and preconditioner */
1006:   return(0);
1007: }

1009: /*@C
1010:    SNESSetJacobian - Sets the function to compute Jacobian as well as the
1011:    location to store the matrix.

1013:    Collective on SNES and Mat

1015:    Input Parameters:
1016: +  snes - the SNES context
1017: .  A - Jacobian matrix
1018: .  B - preconditioner matrix (usually same as the Jacobian)
1019: .  func - Jacobian evaluation routine
1020: -  ctx - [optional] user-defined context for private data for the 
1021:          Jacobian evaluation routine (may be PETSC_NULL)

1023:    Calling sequence of func:
1024: $     func (SNES snes,Vec x,Mat *A,Mat *B,int *flag,void *ctx);

1026: +  x - input vector
1027: .  A - Jacobian matrix
1028: .  B - preconditioner matrix, usually the same as A
1029: .  flag - flag indicating information about the preconditioner matrix
1030:    structure (same as flag in SLESSetOperators())
1031: -  ctx - [optional] user-defined Jacobian context

1033:    Notes: 
1034:    See SLESSetOperators() for important information about setting the flag
1035:    output parameter in the routine func().  Be sure to read this information!

1037:    The routine func() takes Mat * as the matrix arguments rather than Mat.  
1038:    This allows the Jacobian evaluation routine to replace A and/or B with a 
1039:    completely new new matrix structure (not just different matrix elements)
1040:    when appropriate, for instance, if the nonzero structure is changing
1041:    throughout the global iterations.

1043:    Level: beginner

1045: .keywords: SNES, nonlinear, set, Jacobian, matrix

1047: .seealso: SLESSetOperators(), SNESSetFunction()
1048: @*/
1049: int SNESSetJacobian(SNES snes,Mat A,Mat B,int (*func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
1050: {
1057:   if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
1058:     SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
1059:   }

1061:   snes->computejacobian = func;
1062:   snes->jacP            = ctx;
1063:   snes->jacobian        = A;
1064:   snes->jacobian_pre    = B;
1065:   return(0);
1066: }

1068: /*@C
1069:    SNESGetJacobian - Returns the Jacobian matrix and optionally the user 
1070:    provided context for evaluating the Jacobian.

1072:    Not Collective, but Mat object will be parallel if SNES object is

1074:    Input Parameter:
1075: .  snes - the nonlinear solver context

1077:    Output Parameters:
1078: +  A - location to stash Jacobian matrix (or PETSC_NULL)
1079: .  B - location to stash preconditioner matrix (or PETSC_NULL)
1080: .  ctx - location to stash Jacobian ctx (or PETSC_NULL)
1081: -  func - location to put Jacobian function (or PETSC_NULL)

1083:    Level: advanced

1085: .seealso: SNESSetJacobian(), SNESComputeJacobian()
1086: @*/
1087: int SNESGetJacobian(SNES snes,Mat *A,Mat *B,void **ctx,int (**func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*))
1088: {
1091:   if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
1092:     SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
1093:   }
1094:   if (A)    *A    = snes->jacobian;
1095:   if (B)    *B    = snes->jacobian_pre;
1096:   if (ctx)  *ctx  = snes->jacP;
1097:   if (func) *func = snes->computejacobian;
1098:   return(0);
1099: }

1101: /*@C
1102:    SNESSetHessian - Sets the function to compute Hessian as well as the
1103:    location to store the matrix.

1105:    Collective on SNES and Mat

1107:    Input Parameters:
1108: +  snes - the SNES context
1109: .  A - Hessian matrix
1110: .  B - preconditioner matrix (usually same as the Hessian)
1111: .  func - Jacobian evaluation routine
1112: -  ctx - [optional] user-defined context for private data for the 
1113:          Hessian evaluation routine (may be PETSC_NULL)

1115:    Calling sequence of func:
1116: $    func (SNES snes,Vec x,Mat *A,Mat *B,int *flag,void *ctx);

1118: +  x - input vector
1119: .  A - Hessian matrix
1120: .  B - preconditioner matrix, usually the same as A
1121: .  flag - flag indicating information about the preconditioner matrix
1122:    structure (same as flag in SLESSetOperators())
1123: -  ctx - [optional] user-defined Hessian context

1125:    Notes: 
1126:    See SLESSetOperators() for important information about setting the flag
1127:    output parameter in the routine func().  Be sure to read this information!

1129:    The function func() takes Mat * as the matrix arguments rather than Mat.  
1130:    This allows the Hessian evaluation routine to replace A and/or B with a 
1131:    completely new new matrix structure (not just different matrix elements)
1132:    when appropriate, for instance, if the nonzero structure is changing
1133:    throughout the global iterations.

1135:    Level: beginner

1137: .keywords: SNES, nonlinear, set, Hessian, matrix

1139: .seealso: SNESSetMinimizationFunction(), SNESSetGradient(), SLESSetOperators()
1140: @*/
1141: int SNESSetHessian(SNES snes,Mat A,Mat B,int (*func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
1142: {
1149:   if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
1150:     SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
1151:   }
1152:   snes->computejacobian = func;
1153:   snes->jacP            = ctx;
1154:   snes->jacobian        = A;
1155:   snes->jacobian_pre    = B;
1156:   return(0);
1157: }

1159: /*@
1160:    SNESGetHessian - Returns the Hessian matrix and optionally the user 
1161:    provided context for evaluating the Hessian.

1163:    Not Collective, but Mat object is parallel if SNES object is parallel

1165:    Input Parameter:
1166: .  snes - the nonlinear solver context

1168:    Output Parameters:
1169: +  A - location to stash Hessian matrix (or PETSC_NULL)
1170: .  B - location to stash preconditioner matrix (or PETSC_NULL)
1171: -  ctx - location to stash Hessian ctx (or PETSC_NULL)

1173:    Level: advanced

1175: .seealso: SNESSetHessian(), SNESComputeHessian()

1177: .keywords: SNES, get, Hessian
1178: @*/
1179: int SNESGetHessian(SNES snes,Mat *A,Mat *B,void **ctx)
1180: {
1183:   if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION){
1184:     SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
1185:   }
1186:   if (A)   *A = snes->jacobian;
1187:   if (B)   *B = snes->jacobian_pre;
1188:   if (ctx) *ctx = snes->jacP;
1189:   return(0);
1190: }

1192: /* ----- Routines to initialize and destroy a nonlinear solver ---- */

1194: /*@
1195:    SNESSetUp - Sets up the internal data structures for the later use
1196:    of a nonlinear solver.

1198:    Collective on SNES

1200:    Input Parameters:
1201: +  snes - the SNES context
1202: -  x - the solution vector

1204:    Notes:
1205:    For basic use of the SNES solvers the user need not explicitly call
1206:    SNESSetUp(), since these actions will automatically occur during
1207:    the call to SNESSolve().  However, if one wishes to control this
1208:    phase separately, SNESSetUp() should be called after SNESCreate()
1209:    and optional routines of the form SNESSetXXX(), but before SNESSolve().  

1211:    Level: advanced

1213: .keywords: SNES, nonlinear, setup

1215: .seealso: SNESCreate(), SNESSolve(), SNESDestroy()
1216: @*/
1217: int SNESSetUp(SNES snes,Vec x)
1218: {
1219:   int        ierr;
1220:   PetscTruth flg;

1226:   snes->vec_sol = snes->vec_sol_always = x;

1228:   PetscOptionsHasName(snes->prefix,"-snes_mf_operator",&flg);
1229:   /*
1230:       This version replaces the user provided Jacobian matrix with a
1231:       matrix-free version but still employs the user-provided preconditioner matrix
1232:   */
1233:   if (flg) {
1234:     Mat J;
1235:     MatCreateSNESMF(snes,snes->vec_sol,&J);
1236:     PetscLogObjectParent(snes,J);
1237:     snes->mfshell  = J;
1238:     snes->jacobian = J;
1239:     PetscLogInfo(snes,"SNESSetUp: Setting default matrix-free operator routinesn");
1240:     MatSNESMFSetFromOptions(J);
1241:   }
1242:   PetscOptionsHasName(snes->prefix,"-snes_mf",&flg);
1243:   /*
1244:       This version replaces both the user-provided Jacobian and the user-
1245:       provided preconditioner matrix with the default matrix free version.
1246:    */
1247:   if (flg) {
1248:     Mat  J;
1249:     SLES sles;
1250:     PC   pc;

1252:     MatCreateSNESMF(snes,snes->vec_sol,&J);
1253:     PetscLogObjectParent(snes,J);
1254:     snes->mfshell = J;
1255:     PetscLogInfo(snes,"SNESSetUp: Setting default matrix-free operator and preconditioner routinesn");
1256:     if (snes->method_class == SNES_NONLINEAR_EQUATIONS) {
1257:       SNESSetJacobian(snes,J,J,MatSNESMFFormJacobian,snes->funP);
1258:     } else if (snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) {
1259:       SNESSetHessian(snes,J,J,MatSNESMFFormJacobian,snes->funP);
1260:     } else {
1261:       SETERRQ(PETSC_ERR_SUP,"Method class doesn't support matrix-free option");
1262:     }
1263:     MatSNESMFSetFromOptions(J);
1264:     /* force no preconditioner */
1265:     SNESGetSLES(snes,&sles);
1266:     SLESGetPC(sles,&pc);
1267:     PCSetType(pc,PCNONE);
1268:   }

1270:   if ((snes->method_class == SNES_NONLINEAR_EQUATIONS)) {
1271:     PetscTruth iseqtr;

1273:     if (!snes->vec_func) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetFunction() first");
1274:     if (!snes->computefunction) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetFunction() first");
1275:     if (!snes->jacobian) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetJacobian() first n or use -snes_mf option");
1276:     if (snes->vec_func == snes->vec_sol) {
1277:       SETERRQ(PETSC_ERR_ARG_IDN,"Solution vector cannot be function vector");
1278:     }

1280:     /* Set the KSP stopping criterion to use the Eisenstat-Walker method */
1281:     PetscTypeCompare((PetscObject)snes,SNESEQTR,&iseqtr);
1282:     if (snes->ksp_ewconv && !iseqtr) {
1283:       SLES sles; KSP ksp;
1284:       SNESGetSLES(snes,&sles);
1285:       SLESGetKSP(sles,&ksp);
1286:       KSPSetConvergenceTest(ksp,SNES_KSP_EW_Converged_Private,snes);
1287:     }
1288:   } else if ((snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION)) {
1289:     if (!snes->vec_func) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetGradient() first");
1290:     if (!snes->computefunction) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetGradient() first");
1291:     if (!snes->computeumfunction) {
1292:       SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetMinimizationFunction() first");
1293:     }
1294:     if (!snes->jacobian) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetHessian()");
1295:   } else {
1296:     SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Unknown method class");
1297:   }
1298:   if (snes->setup) {(*snes->setup)(snes);}
1299:   snes->setupcalled = 1;
1300:   return(0);
1301: }

1303: /*@C
1304:    SNESDestroy - Destroys the nonlinear solver context that was created
1305:    with SNESCreate().

1307:    Collective on SNES

1309:    Input Parameter:
1310: .  snes - the SNES context

1312:    Level: beginner

1314: .keywords: SNES, nonlinear, destroy

1316: .seealso: SNESCreate(), SNESSolve()
1317: @*/
1318: int SNESDestroy(SNES snes)
1319: {
1320:   int i,ierr;

1324:   if (--snes->refct > 0) return(0);

1326:   /* if memory was published with AMS then destroy it */
1327:   PetscObjectDepublish(snes);

1329:   if (snes->destroy) {(*(snes)->destroy)(snes);}
1330:   if (snes->kspconvctx) {PetscFree(snes->kspconvctx);}
1331:   if (snes->mfshell) {MatDestroy(snes->mfshell);}
1332:   SLESDestroy(snes->sles);
1333:   if (snes->vwork) {VecDestroyVecs(snes->vwork,snes->nvwork);}
1334:   for (i=0; i<snes->numbermonitors; i++) {
1335:     if (snes->monitordestroy[i]) {
1336:       (*snes->monitordestroy[i])(snes->monitorcontext[i]);
1337:     }
1338:   }
1339:   PetscLogObjectDestroy((PetscObject)snes);
1340:   PetscHeaderDestroy((PetscObject)snes);
1341:   return(0);
1342: }

1344: /* ----------- Routines to set solver parameters ---------- */

1346: /*@
1347:    SNESSetTolerances - Sets various parameters used in convergence tests.

1349:    Collective on SNES

1351:    Input Parameters:
1352: +  snes - the SNES context
1353: .  atol - absolute convergence tolerance
1354: .  rtol - relative convergence tolerance
1355: .  stol -  convergence tolerance in terms of the norm
1356:            of the change in the solution between steps
1357: .  maxit - maximum number of iterations
1358: -  maxf - maximum number of function evaluations

1360:    Options Database Keys: 
1361: +    -snes_atol <atol> - Sets atol
1362: .    -snes_rtol <rtol> - Sets rtol
1363: .    -snes_stol <stol> - Sets stol
1364: .    -snes_max_it <maxit> - Sets maxit
1365: -    -snes_max_funcs <maxf> - Sets maxf

1367:    Notes:
1368:    The default maximum number of iterations is 50.
1369:    The default maximum number of function evaluations is 1000.

1371:    Level: intermediate

1373: .keywords: SNES, nonlinear, set, convergence, tolerances

1375: .seealso: SNESSetTrustRegionTolerance(), SNESSetMinimizationFunctionTolerance()
1376: @*/
1377: int SNESSetTolerances(SNES snes,PetscReal atol,PetscReal rtol,PetscReal stol,int maxit,int maxf)
1378: {
1381:   if (atol != PETSC_DEFAULT)  snes->atol      = atol;
1382:   if (rtol != PETSC_DEFAULT)  snes->rtol      = rtol;
1383:   if (stol != PETSC_DEFAULT)  snes->xtol      = stol;
1384:   if (maxit != PETSC_DEFAULT) snes->max_its   = maxit;
1385:   if (maxf != PETSC_DEFAULT)  snes->max_funcs = maxf;
1386:   return(0);
1387: }

1389: /*@
1390:    SNESGetTolerances - Gets various parameters used in convergence tests.

1392:    Not Collective

1394:    Input Parameters:
1395: +  snes - the SNES context
1396: .  atol - absolute convergence tolerance
1397: .  rtol - relative convergence tolerance
1398: .  stol -  convergence tolerance in terms of the norm
1399:            of the change in the solution between steps
1400: .  maxit - maximum number of iterations
1401: -  maxf - maximum number of function evaluations

1403:    Notes:
1404:    The user can specify PETSC_NULL for any parameter that is not needed.

1406:    Level: intermediate

1408: .keywords: SNES, nonlinear, get, convergence, tolerances

1410: .seealso: SNESSetTolerances()
1411: @*/
1412: int SNESGetTolerances(SNES snes,PetscReal *atol,PetscReal *rtol,PetscReal *stol,int *maxit,int *maxf)
1413: {
1416:   if (atol)  *atol  = snes->atol;
1417:   if (rtol)  *rtol  = snes->rtol;
1418:   if (stol)  *stol  = snes->xtol;
1419:   if (maxit) *maxit = snes->max_its;
1420:   if (maxf)  *maxf  = snes->max_funcs;
1421:   return(0);
1422: }

1424: /*@
1425:    SNESSetTrustRegionTolerance - Sets the trust region parameter tolerance.  

1427:    Collective on SNES

1429:    Input Parameters:
1430: +  snes - the SNES context
1431: -  tol - tolerance
1432:    
1433:    Options Database Key: 
1434: .  -snes_trtol <tol> - Sets tol

1436:    Level: intermediate

1438: .keywords: SNES, nonlinear, set, trust region, tolerance

1440: .seealso: SNESSetTolerances(), SNESSetMinimizationFunctionTolerance()
1441: @*/
1442: int SNESSetTrustRegionTolerance(SNES snes,PetscReal tol)
1443: {
1446:   snes->deltatol = tol;
1447:   return(0);
1448: }

1450: /*@
1451:    SNESSetMinimizationFunctionTolerance - Sets the minimum allowable function tolerance
1452:    for unconstrained minimization solvers.
1453:    
1454:    Collective on SNES

1456:    Input Parameters:
1457: +  snes - the SNES context
1458: -  ftol - minimum function tolerance

1460:    Options Database Key: 
1461: .  -snes_fmin <ftol> - Sets ftol

1463:    Note:
1464:    SNESSetMinimizationFunctionTolerance() is valid for SNES_UNCONSTRAINED_MINIMIZATION
1465:    methods only.

1467:    Level: intermediate

1469: .keywords: SNES, nonlinear, set, minimum, convergence, function, tolerance

1471: .seealso: SNESSetTolerances(), SNESSetTrustRegionTolerance()
1472: @*/
1473: int SNESSetMinimizationFunctionTolerance(SNES snes,PetscReal ftol)
1474: {
1477:   snes->fmin = ftol;
1478:   return(0);
1479: }
1480: /* 
1481:    Duplicate the lg monitors for SNES from KSP; for some reason with 
1482:    dynamic libraries things don't work under Sun4 if we just use 
1483:    macros instead of functions
1484: */
1485: int SNESLGMonitor(SNES snes,int it,PetscReal norm,void *ctx)
1486: {

1491:   KSPLGMonitor((KSP)snes,it,norm,ctx);
1492:   return(0);
1493: }

1495: int SNESLGMonitorCreate(char *host,char *label,int x,int y,int m,int n,PetscDrawLG *draw)
1496: {

1500:   KSPLGMonitorCreate(host,label,x,y,m,n,draw);
1501:   return(0);
1502: }

1504: int SNESLGMonitorDestroy(PetscDrawLG draw)
1505: {

1509:   KSPLGMonitorDestroy(draw);
1510:   return(0);
1511: }

1513: /* ------------ Routines to set performance monitoring options ----------- */

1515: /*@C
1516:    SNESSetMonitor - Sets an ADDITIONAL function that is to be used at every
1517:    iteration of the nonlinear solver to display the iteration's 
1518:    progress.   

1520:    Collective on SNES

1522:    Input Parameters:
1523: +  snes - the SNES context
1524: .  func - monitoring routine
1525: .  mctx - [optional] user-defined context for private data for the 
1526:           monitor routine (use PETSC_NULL if no context is desitre)
1527: -  monitordestroy - [optional] routine that frees monitor context
1528:           (may be PETSC_NULL)

1530:    Calling sequence of func:
1531: $     int func(SNES snes,int its, PetscReal norm,void *mctx)

1533: +    snes - the SNES context
1534: .    its - iteration number
1535: .    norm - 2-norm function value (may be estimated)
1536:             for SNES_NONLINEAR_EQUATIONS methods
1537: .    norm - 2-norm gradient value (may be estimated)
1538:             for SNES_UNCONSTRAINED_MINIMIZATION methods
1539: -    mctx - [optional] monitoring context

1541:    Options Database Keys:
1542: +    -snes_monitor        - sets SNESDefaultMonitor()
1543: .    -snes_xmonitor       - sets line graph monitor,
1544:                             uses SNESLGMonitorCreate()
1545: _    -snes_cancelmonitors - cancels all monitors that have
1546:                             been hardwired into a code by 
1547:                             calls to SNESSetMonitor(), but
1548:                             does not cancel those set via
1549:                             the options database.

1551:    Notes: 
1552:    Several different monitoring routines may be set by calling
1553:    SNESSetMonitor() multiple times; all will be called in the 
1554:    order in which they were set.

1556:    Level: intermediate

1558: .keywords: SNES, nonlinear, set, monitor

1560: .seealso: SNESDefaultMonitor(), SNESClearMonitor()
1561: @*/
1562: int SNESSetMonitor(SNES snes,int (*func)(SNES,int,PetscReal,void*),void *mctx,int (*monitordestroy)(void *))
1563: {
1566:   if (snes->numbermonitors >= MAXSNESMONITORS) {
1567:     SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1568:   }

1570:   snes->monitor[snes->numbermonitors]           = func;
1571:   snes->monitordestroy[snes->numbermonitors]    = monitordestroy;
1572:   snes->monitorcontext[snes->numbermonitors++]  = (void*)mctx;
1573:   return(0);
1574: }

1576: /*@C
1577:    SNESClearMonitor - Clears all the monitor functions for a SNES object.

1579:    Collective on SNES

1581:    Input Parameters:
1582: .  snes - the SNES context

1584:    Options Database:
1585: .  -snes_cancelmonitors - cancels all monitors that have been hardwired
1586:     into a code by calls to SNESSetMonitor(), but does not cancel those 
1587:     set via the options database

1589:    Notes: 
1590:    There is no way to clear one specific monitor from a SNES object.

1592:    Level: intermediate

1594: .keywords: SNES, nonlinear, set, monitor

1596: .seealso: SNESDefaultMonitor(), SNESSetMonitor()
1597: @*/
1598: int SNESClearMonitor(SNES snes)
1599: {
1602:   snes->numbermonitors = 0;
1603:   return(0);
1604: }

1606: /*@C
1607:    SNESSetConvergenceTest - Sets the function that is to be used 
1608:    to test for convergence of the nonlinear iterative solution.   

1610:    Collective on SNES

1612:    Input Parameters:
1613: +  snes - the SNES context
1614: .  func - routine to test for convergence
1615: -  cctx - [optional] context for private data for the convergence routine 
1616:           (may be PETSC_NULL)

1618:    Calling sequence of func:
1619: $     int func (SNES snes,PetscReal xnorm,PetscReal gnorm,PetscReal f,SNESConvergedReason *reason,void *cctx)

1621: +    snes - the SNES context
1622: .    cctx - [optional] convergence context
1623: .    reason - reason for convergence/divergence
1624: .    xnorm - 2-norm of current iterate
1625: .    gnorm - 2-norm of current step (SNES_NONLINEAR_EQUATIONS methods)
1626: .    f - 2-norm of function (SNES_NONLINEAR_EQUATIONS methods)
1627: .    gnorm - 2-norm of current gradient (SNES_UNCONSTRAINED_MINIMIZATION methods)
1628: -    f - function value (SNES_UNCONSTRAINED_MINIMIZATION methods)

1630:    Level: advanced

1632: .keywords: SNES, nonlinear, set, convergence, test

1634: .seealso: SNESConverged_EQ_LS(), SNESConverged_EQ_TR(), 
1635:           SNESConverged_UM_LS(), SNESConverged_UM_TR()
1636: @*/
1637: int SNESSetConvergenceTest(SNES snes,int (*func)(SNES,PetscReal,PetscReal,PetscReal,SNESConvergedReason*,void*),void *cctx)
1638: {
1641:   (snes)->converged = func;
1642:   (snes)->cnvP      = cctx;
1643:   return(0);
1644: }

1646: /*@C
1647:    SNESGetConvergedReason - Gets the reason the SNES iteration was stopped.

1649:    Not Collective

1651:    Input Parameter:
1652: .  snes - the SNES context

1654:    Output Parameter:
1655: .  reason - negative value indicates diverged, positive value converged, see petscsnes.h or the 
1656:             manual pages for the individual convergence tests for complete lists

1658:    Level: intermediate

1660:    Notes: Can only be called after the call the SNESSolve() is complete.

1662: .keywords: SNES, nonlinear, set, convergence, test

1664: .seealso: SNESSetConvergenceTest(), SNESConverged_EQ_LS(), SNESConverged_EQ_TR(), 
1665:           SNESConverged_UM_LS(), SNESConverged_UM_TR(), SNESConvergedReason
1666: @*/
1667: int SNESGetConvergedReason(SNES snes,SNESConvergedReason *reason)
1668: {
1671:   *reason = snes->reason;
1672:   return(0);
1673: }

1675: /*@
1676:    SNESSetConvergenceHistory - Sets the array used to hold the convergence history.

1678:    Collective on SNES

1680:    Input Parameters:
1681: +  snes - iterative context obtained from SNESCreate()
1682: .  a   - array to hold history
1683: .  its - integer array holds the number of linear iterations for each solve.
1684: .  na  - size of a and its
1685: -  reset - PETSC_TRUE indicates each new nonlinear solve resets the history counter to zero,
1686:            else it continues storing new values for new nonlinear solves after the old ones

1688:    Notes:
1689:    If set, this array will contain the function norms (for
1690:    SNES_NONLINEAR_EQUATIONS methods) or gradient norms
1691:    (for SNES_UNCONSTRAINED_MINIMIZATION methods) computed
1692:    at each step.

1694:    This routine is useful, e.g., when running a code for purposes
1695:    of accurate performance monitoring, when no I/O should be done
1696:    during the section of code that is being timed.

1698:    Level: intermediate

1700: .keywords: SNES, set, convergence, history

1702: .seealso: SNESGetConvergenceHistory()

1704: @*/
1705: int SNESSetConvergenceHistory(SNES snes,PetscReal *a,int *its,int na,PetscTruth reset)
1706: {
1710:   snes->conv_hist       = a;
1711:   snes->conv_hist_its   = its;
1712:   snes->conv_hist_max   = na;
1713:   snes->conv_hist_reset = reset;
1714:   return(0);
1715: }

1717: /*@C
1718:    SNESGetConvergenceHistory - Gets the array used to hold the convergence history.

1720:    Collective on SNES

1722:    Input Parameter:
1723: .  snes - iterative context obtained from SNESCreate()

1725:    Output Parameters:
1726: .  a   - array to hold history
1727: .  its - integer array holds the number of linear iterations (or
1728:          negative if not converged) for each solve.
1729: -  na  - size of a and its

1731:    Notes:
1732:     The calling sequence for this routine in Fortran is
1733: $   call SNESGetConvergenceHistory(SNES snes, integer na, integer ierr)

1735:    This routine is useful, e.g., when running a code for purposes
1736:    of accurate performance monitoring, when no I/O should be done
1737:    during the section of code that is being timed.

1739:    Level: intermediate

1741: .keywords: SNES, get, convergence, history

1743: .seealso: SNESSetConvergencHistory()

1745: @*/
1746: int SNESGetConvergenceHistory(SNES snes,PetscReal **a,int **its,int *na)
1747: {
1750:   if (a)   *a   = snes->conv_hist;
1751:   if (its) *its = snes->conv_hist_its;
1752:   if (na) *na   = snes->conv_hist_len;
1753:   return(0);
1754: }

1756: /*
1757:    SNESScaleStep_Private - Scales a step so that its length is less than the
1758:    positive parameter delta.

1760:     Input Parameters:
1761: +   snes - the SNES context
1762: .   y - approximate solution of linear system
1763: .   fnorm - 2-norm of current function
1764: -   delta - trust region size

1766:     Output Parameters:
1767: +   gpnorm - predicted function norm at the new point, assuming local 
1768:     linearization.  The value is zero if the step lies within the trust 
1769:     region, and exceeds zero otherwise.
1770: -   ynorm - 2-norm of the step

1772:     Note:
1773:     For non-trust region methods such as SNESEQLS, the parameter delta 
1774:     is set to be the maximum allowable step size.  

1776: .keywords: SNES, nonlinear, scale, step
1777: */
1778: int SNESScaleStep_Private(SNES snes,Vec y,PetscReal *fnorm,PetscReal *delta,
1779:                           PetscReal *gpnorm,PetscReal *ynorm)
1780: {
1781:   PetscReal norm;
1782:   Scalar cnorm;
1783:   int    ierr;


1790:   VecNorm(y,NORM_2,&norm);
1791:   if (norm > *delta) {
1792:      norm = *delta/norm;
1793:      *gpnorm = (1.0 - norm)*(*fnorm);
1794:      cnorm = norm;
1795:      VecScale(&cnorm,y);
1796:      *ynorm = *delta;
1797:   } else {
1798:      *gpnorm = 0.0;
1799:      *ynorm = norm;
1800:   }
1801:   return(0);
1802: }

1804: /*@
1805:    SNESSolve - Solves a nonlinear system.  Call SNESSolve after calling 
1806:    SNESCreate() and optional routines of the form SNESSetXXX().

1808:    Collective on SNES

1810:    Input Parameters:
1811: +  snes - the SNES context
1812: -  x - the solution vector

1814:    Output Parameter:
1815: .  its - number of iterations until termination

1817:    Notes:
1818:    The user should initialize the vector,x, with the initial guess
1819:    for the nonlinear solve prior to calling SNESSolve.  In particular,
1820:    to employ an initial guess of zero, the user should explicitly set
1821:    this vector to zero by calling VecSet().

1823:    Level: beginner

1825: .keywords: SNES, nonlinear, solve

1827: .seealso: SNESCreate(), SNESDestroy()
1828: @*/
1829: int SNESSolve(SNES snes,Vec x,int *its)
1830: {
1831:   int        ierr;
1832:   PetscTruth flg;

1839:   if (!snes->solve) SETERRQ(1,"SNESSetType() or SNESSetFromOptions() must be called before SNESSolve()");

1841:   if (!snes->setupcalled) {SNESSetUp(snes,x);}
1842:   else {snes->vec_sol = snes->vec_sol_always = x;}
1843:   if (snes->conv_hist_reset == PETSC_TRUE) snes->conv_hist_len = 0;
1844:   PetscLogEventBegin(SNES_Solve,snes,0,0,0);
1845:   snes->nfuncs = 0; snes->linear_its = 0; snes->nfailures = 0;
1846:   (*(snes)->solve)(snes,its);
1847:   PetscLogEventEnd(SNES_Solve,snes,0,0,0);
1848:   PetscOptionsHasName(snes->prefix,"-snes_view",&flg);
1849:   if (flg && !PetscPreLoadingOn) { SNESView(snes,PETSC_VIEWER_STDOUT_WORLD); }
1850:   return(0);
1851: }

1853: /* --------- Internal routines for SNES Package --------- */

1855: /*@C
1856:    SNESSetType - Sets the method for the nonlinear solver.  

1858:    Collective on SNES

1860:    Input Parameters:
1861: +  snes - the SNES context
1862: -  type - a known method

1864:    Options Database Key:
1865: .  -snes_type <type> - Sets the method; use -help for a list
1866:    of available methods (for instance, ls or tr)

1868:    Notes:
1869:    See "petsc/include/petscsnes.h" for available methods (for instance)
1870: +    SNESEQLS - Newton's method with line search
1871:      (systems of nonlinear equations)
1872: .    SNESEQTR - Newton's method with trust region
1873:      (systems of nonlinear equations)
1874: .    SNESUMTR - Newton's method with trust region 
1875:      (unconstrained minimization)
1876: -    SNESUMLS - Newton's method with line search
1877:      (unconstrained minimization)

1879:   Normally, it is best to use the SNESSetFromOptions() command and then
1880:   set the SNES solver type from the options database rather than by using
1881:   this routine.  Using the options database provides the user with
1882:   maximum flexibility in evaluating the many nonlinear solvers.
1883:   The SNESSetType() routine is provided for those situations where it
1884:   is necessary to set the nonlinear solver independently of the command
1885:   line or options database.  This might be the case, for example, when
1886:   the choice of solver changes during the execution of the program,
1887:   and the user's application is taking responsibility for choosing the
1888:   appropriate method.

1890:   Level: intermediate

1892: .keywords: SNES, set, type

1894: .seealso: SNESType, SNESCreate()

1896: @*/
1897: int SNESSetType(SNES snes,SNESType type)
1898: {
1899:   int        ierr,(*r)(SNES);
1900:   PetscTruth match;


1906:   PetscTypeCompare((PetscObject)snes,type,&match);
1907:   if (match) return(0);

1909:   if (snes->setupcalled) {
1910:     ierr       = (*(snes)->destroy)(snes);
1911:     snes->data = 0;
1912:   }

1914:   /* Get the function pointers for the iterative method requested */
1915:   if (!SNESRegisterAllCalled) {SNESRegisterAll(PETSC_NULL);}

1917:    PetscFListFind(snes->comm,SNESList,type,(void (**)(void)) &r);

1919:   if (!r) SETERRQ1(1,"Unable to find requested SNES type %s",type);

1921:   if (snes->data) {PetscFree(snes->data);}
1922:   snes->data = 0;
1923:   (*r)(snes);

1925:   PetscObjectChangeTypeName((PetscObject)snes,type);
1926:   snes->set_method_called = 1;

1928:   return(0);
1929: }


1932: /* --------------------------------------------------------------------- */
1933: /*@C
1934:    SNESRegisterDestroy - Frees the list of nonlinear solvers that were
1935:    registered by SNESRegisterDynamic().

1937:    Not Collective

1939:    Level: advanced

1941: .keywords: SNES, nonlinear, register, destroy

1943: .seealso: SNESRegisterAll(), SNESRegisterAll()
1944: @*/
1945: int SNESRegisterDestroy(void)
1946: {

1950:   if (SNESList) {
1951:     PetscFListDestroy(&SNESList);
1952:     SNESList = 0;
1953:   }
1954:   SNESRegisterAllCalled = PETSC_FALSE;
1955:   return(0);
1956: }

1958: /*@C
1959:    SNESGetType - Gets the SNES method type and name (as a string).

1961:    Not Collective

1963:    Input Parameter:
1964: .  snes - nonlinear solver context

1966:    Output Parameter:
1967: .  type - SNES method (a charactor string)

1969:    Level: intermediate

1971: .keywords: SNES, nonlinear, get, type, name
1972: @*/
1973: int SNESGetType(SNES snes,SNESType *type)
1974: {
1977:   *type = snes->type_name;
1978:   return(0);
1979: }

1981: /*@C
1982:    SNESGetSolution - Returns the vector where the approximate solution is
1983:    stored.

1985:    Not Collective, but Vec is parallel if SNES is parallel

1987:    Input Parameter:
1988: .  snes - the SNES context

1990:    Output Parameter:
1991: .  x - the solution

1993:    Level: advanced

1995: .keywords: SNES, nonlinear, get, solution

1997: .seealso: SNESGetFunction(), SNESGetGradient(), SNESGetSolutionUpdate()
1998: @*/
1999: int SNESGetSolution(SNES snes,Vec *x)
2000: {
2003:   *x = snes->vec_sol_always;
2004:   return(0);
2005: }

2007: /*@C
2008:    SNESGetSolutionUpdate - Returns the vector where the solution update is
2009:    stored. 

2011:    Not Collective, but Vec is parallel if SNES is parallel

2013:    Input Parameter:
2014: .  snes - the SNES context

2016:    Output Parameter:
2017: .  x - the solution update

2019:    Level: advanced

2021: .keywords: SNES, nonlinear, get, solution, update

2023: .seealso: SNESGetSolution(), SNESGetFunction
2024: @*/
2025: int SNESGetSolutionUpdate(SNES snes,Vec *x)
2026: {
2029:   *x = snes->vec_sol_update_always;
2030:   return(0);
2031: }

2033: /*@C
2034:    SNESGetFunction - Returns the vector where the function is stored.

2036:    Not Collective, but Vec is parallel if SNES is parallel

2038:    Input Parameter:
2039: .  snes - the SNES context

2041:    Output Parameter:
2042: +  r - the function (or PETSC_NULL)
2043: .  ctx - the function context (or PETSC_NULL)
2044: -  func - the function (or PETSC_NULL)

2046:    Notes:
2047:    SNESGetFunction() is valid for SNES_NONLINEAR_EQUATIONS methods only
2048:    Analogous routines for SNES_UNCONSTRAINED_MINIMIZATION methods are
2049:    SNESGetMinimizationFunction() and SNESGetGradient();

2051:    Level: advanced

2053: .keywords: SNES, nonlinear, get, function

2055: .seealso: SNESSetFunction(), SNESGetSolution(), SNESGetMinimizationFunction(),
2056:           SNESGetGradient()

2058: @*/
2059: int SNESGetFunction(SNES snes,Vec *r,void **ctx,int (**func)(SNES,Vec,Vec,void*))
2060: {
2063:   if (snes->method_class != SNES_NONLINEAR_EQUATIONS) {
2064:     SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_NONLINEAR_EQUATIONS only");
2065:   }
2066:   if (r)    *r    = snes->vec_func_always;
2067:   if (ctx)  *ctx  = snes->funP;
2068:   if (func) *func = snes->computefunction;
2069:   return(0);
2070: }

2072: /*@C
2073:    SNESGetGradient - Returns the vector where the gradient is stored.

2075:    Not Collective, but Vec is parallel if SNES is parallel

2077:    Input Parameter:
2078: .  snes - the SNES context

2080:    Output Parameter:
2081: +  r - the gradient (or PETSC_NULL)
2082: -  ctx - the gradient context (or PETSC_NULL)

2084:    Notes:
2085:    SNESGetGradient() is valid for SNES_UNCONSTRAINED_MINIMIZATION methods 
2086:    only.  An analogous routine for SNES_NONLINEAR_EQUATIONS methods is
2087:    SNESGetFunction().

2089:    Level: advanced

2091: .keywords: SNES, nonlinear, get, gradient

2093: .seealso: SNESGetMinimizationFunction(), SNESGetSolution(), SNESGetFunction(),
2094:           SNESSetGradient(), SNESSetFunction()

2096: @*/
2097: int SNESGetGradient(SNES snes,Vec *r,void **ctx)
2098: {
2101:   if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
2102:     SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
2103:   }
2104:   if (r)   *r = snes->vec_func_always;
2105:   if (ctx) *ctx = snes->funP;
2106:   return(0);
2107: }

2109: /*@C
2110:    SNESGetMinimizationFunction - Returns the scalar function value for 
2111:    unconstrained minimization problems.

2113:    Not Collective

2115:    Input Parameter:
2116: .  snes - the SNES context

2118:    Output Parameter:
2119: +  r - the function (or PETSC_NULL)
2120: -  ctx - the function context (or PETSC_NULL)

2122:    Notes:
2123:    SNESGetMinimizationFunction() is valid for SNES_UNCONSTRAINED_MINIMIZATION
2124:    methods only.  An analogous routine for SNES_NONLINEAR_EQUATIONS methods is
2125:    SNESGetFunction().

2127:    Level: advanced

2129: .keywords: SNES, nonlinear, get, function

2131: .seealso: SNESGetGradient(), SNESGetSolution(), SNESGetFunction(), SNESSetFunction()

2133: @*/
2134: int SNESGetMinimizationFunction(SNES snes,PetscReal *r,void **ctx)
2135: {
2139:   if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) {
2140:     SETERRQ(PETSC_ERR_ARG_WRONG,"For SNES_UNCONSTRAINED_MINIMIZATION only");
2141:   }
2142:   if (r)   *r = snes->fc;
2143:   if (ctx) *ctx = snes->umfunP;
2144:   return(0);
2145: }

2147: /*@C
2148:    SNESSetOptionsPrefix - Sets the prefix used for searching for all 
2149:    SNES options in the database.

2151:    Collective on SNES

2153:    Input Parameter:
2154: +  snes - the SNES context
2155: -  prefix - the prefix to prepend to all option names

2157:    Notes:
2158:    A hyphen (-) must NOT be given at the beginning of the prefix name.
2159:    The first character of all runtime options is AUTOMATICALLY the hyphen.

2161:    Level: advanced

2163: .keywords: SNES, set, options, prefix, database

2165: .seealso: SNESSetFromOptions()
2166: @*/
2167: int SNESSetOptionsPrefix(SNES snes,char *prefix)
2168: {

2173:   PetscObjectSetOptionsPrefix((PetscObject)snes,prefix);
2174:   SLESSetOptionsPrefix(snes->sles,prefix);
2175:   return(0);
2176: }

2178: /*@C
2179:    SNESAppendOptionsPrefix - Appends to the prefix used for searching for all 
2180:    SNES options in the database.

2182:    Collective on SNES

2184:    Input Parameters:
2185: +  snes - the SNES context
2186: -  prefix - the prefix to prepend to all option names

2188:    Notes:
2189:    A hyphen (-) must NOT be given at the beginning of the prefix name.
2190:    The first character of all runtime options is AUTOMATICALLY the hyphen.

2192:    Level: advanced

2194: .keywords: SNES, append, options, prefix, database

2196: .seealso: SNESGetOptionsPrefix()
2197: @*/
2198: int SNESAppendOptionsPrefix(SNES snes,char *prefix)
2199: {
2201: 
2204:   PetscObjectAppendOptionsPrefix((PetscObject)snes,prefix);
2205:   SLESAppendOptionsPrefix(snes->sles,prefix);
2206:   return(0);
2207: }

2209: /*@C
2210:    SNESGetOptionsPrefix - Sets the prefix used for searching for all 
2211:    SNES options in the database.

2213:    Not Collective

2215:    Input Parameter:
2216: .  snes - the SNES context

2218:    Output Parameter:
2219: .  prefix - pointer to the prefix string used

2221:    Notes: On the fortran side, the user should pass in a string 'prifix' of
2222:    sufficient length to hold the prefix.

2224:    Level: advanced

2226: .keywords: SNES, get, options, prefix, database

2228: .seealso: SNESAppendOptionsPrefix()
2229: @*/
2230: int SNESGetOptionsPrefix(SNES snes,char **prefix)
2231: {

2236:   PetscObjectGetOptionsPrefix((PetscObject)snes,prefix);
2237:   return(0);
2238: }

2240: /*MC
2241:    SNESRegisterDynamic - Adds a method to the nonlinear solver package.

2243:    Synopsis:
2244:    SNESRegisterDynamic(char *name_solver,char *path,char *name_create,int (*routine_create)(SNES))

2246:    Not collective

2248:    Input Parameters:
2249: +  name_solver - name of a new user-defined solver
2250: .  path - path (either absolute or relative) the library containing this solver
2251: .  name_create - name of routine to create method context
2252: -  routine_create - routine to create method context

2254:    Notes:
2255:    SNESRegisterDynamic() may be called multiple times to add several user-defined solvers.

2257:    If dynamic libraries are used, then the fourth input argument (routine_create)
2258:    is ignored.

2260:    Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR}, ${BOPT},
2261:    and others of the form ${any_environmental_variable} occuring in pathname will be 
2262:    replaced with appropriate values.

2264:    Sample usage:
2265: .vb
2266:    SNESRegisterDynamic("my_solver",/home/username/my_lib/lib/libg/solaris/mylib.a,
2267:                 "MySolverCreate",MySolverCreate);
2268: .ve

2270:    Then, your solver can be chosen with the procedural interface via
2271: $     SNESSetType(snes,"my_solver")
2272:    or at runtime via the option
2273: $     -snes_type my_solver

2275:    Level: advanced

2277: .keywords: SNES, nonlinear, register

2279: .seealso: SNESRegisterAll(), SNESRegisterDestroy()
2280: M*/

2282: int SNESRegister(char *sname,char *path,char *name,int (*function)(SNES))
2283: {
2284:   char fullname[256];
2285:   int  ierr;

2288:   PetscFListConcat(path,name,fullname);
2289:   PetscFListAdd(&SNESList,sname,fullname,(void (*)())function);
2290:   return(0);
2291: }