Actual source code: taosolver.c

petsc-master 2020-07-09
Report Typos and Errors
  1: #define TAO_DLL

  3:  #include <petsc/private/taoimpl.h>

  5: PetscBool TaoRegisterAllCalled = PETSC_FALSE;
  6: PetscFunctionList TaoList = NULL;

  8: PetscClassId TAO_CLASSID;

 10: PetscLogEvent TAO_Solve;
 11: PetscLogEvent TAO_ObjectiveEval;
 12: PetscLogEvent TAO_GradientEval;
 13: PetscLogEvent TAO_ObjGradEval;
 14: PetscLogEvent TAO_HessianEval;
 15: PetscLogEvent TAO_JacobianEval;
 16: PetscLogEvent TAO_ConstraintsEval;

 18: const char *TaoSubSetTypes[] = {"subvec","mask","matrixfree","TaoSubSetType","TAO_SUBSET_",NULL};

 20: struct _n_TaoMonitorDrawCtx {
 21:   PetscViewer viewer;
 22:   PetscInt    howoften;  /* when > 0 uses iteration % howoften, when negative only final solution plotted */
 23: };

 25: /*@
 26:   TaoCreate - Creates a TAO solver

 28:   Collective

 30:   Input Parameter:
 31: . comm - MPI communicator

 33:   Output Parameter:
 34: . newtao - the new Tao context

 36:   Available methods include:
 37: +    nls - Newton's method with line search for unconstrained minimization
 38: .    ntr - Newton's method with trust region for unconstrained minimization
 39: .    ntl - Newton's method with trust region, line search for unconstrained minimization
 40: .    lmvm - Limited memory variable metric method for unconstrained minimization
 41: .    cg - Nonlinear conjugate gradient method for unconstrained minimization
 42: .    nm - Nelder-Mead algorithm for derivate-free unconstrained minimization
 43: .    tron - Newton Trust Region method for bound constrained minimization
 44: .    gpcg - Newton Trust Region method for quadratic bound constrained minimization
 45: .    blmvm - Limited memory variable metric method for bound constrained minimization
 46: .    lcl - Linearly constrained Lagrangian method for pde-constrained minimization
 47: -    pounders - Model-based algorithm for nonlinear least squares

 49:    Options Database Keys:
 50: .   -tao_type - select which method TAO should use

 52:    Level: beginner

 54: .seealso: TaoSolve(), TaoDestroy()
 55: @*/
 56: PetscErrorCode TaoCreate(MPI_Comm comm, Tao *newtao)
 57: {
 59:   Tao            tao;

 63:   *newtao = NULL;

 65:   TaoInitializePackage();
 66:   TaoLineSearchInitializePackage();

 68:   PetscHeaderCreate(tao,TAO_CLASSID,"Tao","Optimization solver","Tao",comm,TaoDestroy,TaoView);
 69:   tao->ops->computeobjective = NULL;
 70:   tao->ops->computeobjectiveandgradient = NULL;
 71:   tao->ops->computegradient = NULL;
 72:   tao->ops->computehessian = NULL;
 73:   tao->ops->computeresidual = NULL;
 74:   tao->ops->computeresidualjacobian = NULL;
 75:   tao->ops->computeconstraints = NULL;
 76:   tao->ops->computejacobian = NULL;
 77:   tao->ops->computejacobianequality = NULL;
 78:   tao->ops->computejacobianinequality = NULL;
 79:   tao->ops->computeequalityconstraints = NULL;
 80:   tao->ops->computeinequalityconstraints = NULL;
 81:   tao->ops->convergencetest = TaoDefaultConvergenceTest;
 82:   tao->ops->convergencedestroy = NULL;
 83:   tao->ops->computedual = NULL;
 84:   tao->ops->setup = NULL;
 85:   tao->ops->solve = NULL;
 86:   tao->ops->view = NULL;
 87:   tao->ops->setfromoptions = NULL;
 88:   tao->ops->destroy = NULL;

 90:   tao->solution=NULL;
 91:   tao->gradient=NULL;
 92:   tao->ls_res = NULL;
 93:   tao->ls_jac = NULL;
 94:   tao->constraints=NULL;
 95:   tao->constraints_equality=NULL;
 96:   tao->constraints_inequality=NULL;
 97:   tao->res_weights_v=NULL;
 98:   tao->res_weights_w=NULL;
 99:   tao->stepdirection=NULL;
100:   tao->niter=0;
101:   tao->ntotalits=0;
102:   tao->XL = NULL;
103:   tao->XU = NULL;
104:   tao->IL = NULL;
105:   tao->IU = NULL;
106:   tao->DI = NULL;
107:   tao->DE = NULL;
108:   tao->gradient_norm = NULL;
109:   tao->gradient_norm_tmp = NULL;
110:   tao->hessian = NULL;
111:   tao->hessian_pre = NULL;
112:   tao->jacobian = NULL;
113:   tao->jacobian_pre = NULL;
114:   tao->jacobian_state = NULL;
115:   tao->jacobian_state_pre = NULL;
116:   tao->jacobian_state_inv = NULL;
117:   tao->jacobian_design = NULL;
118:   tao->jacobian_design_pre = NULL;
119:   tao->jacobian_equality = NULL;
120:   tao->jacobian_equality_pre = NULL;
121:   tao->jacobian_inequality = NULL;
122:   tao->jacobian_inequality_pre = NULL;
123:   tao->state_is = NULL;
124:   tao->design_is = NULL;

126:   tao->max_it     = 10000;
127:   tao->max_funcs   = 10000;
128: #if defined(PETSC_USE_REAL_SINGLE)
129:   tao->gatol       = 1e-5;
130:   tao->grtol       = 1e-5;
131: #else
132:   tao->gatol       = 1e-8;
133:   tao->grtol       = 1e-8;
134: #endif
135:   tao->crtol       = 0.0;
136:   tao->catol       = 0.0;
137:   tao->gttol       = 0.0;
138:   tao->steptol     = 0.0;
139:   tao->trust0      = PETSC_INFINITY;
140:   tao->fmin        = PETSC_NINFINITY;
141:   tao->hist_malloc = PETSC_FALSE;
142:   tao->hist_reset = PETSC_TRUE;
143:   tao->hist_max = 0;
144:   tao->hist_len = 0;
145:   tao->hist_obj = NULL;
146:   tao->hist_resid = NULL;
147:   tao->hist_cnorm = NULL;
148:   tao->hist_lits = NULL;

150:   tao->numbermonitors=0;
151:   tao->viewsolution=PETSC_FALSE;
152:   tao->viewhessian=PETSC_FALSE;
153:   tao->viewgradient=PETSC_FALSE;
154:   tao->viewjacobian=PETSC_FALSE;
155:   tao->viewconstraints = PETSC_FALSE;

157:   tao->bounded = PETSC_FALSE;

159:   tao->header_printed = PETSC_FALSE;

161:   /* These flags prevents algorithms from overriding user options */
162:   tao->max_it_changed   =PETSC_FALSE;
163:   tao->max_funcs_changed=PETSC_FALSE;
164:   tao->gatol_changed    =PETSC_FALSE;
165:   tao->grtol_changed    =PETSC_FALSE;
166:   tao->gttol_changed    =PETSC_FALSE;
167:   tao->steptol_changed  =PETSC_FALSE;
168:   tao->trust0_changed   =PETSC_FALSE;
169:   tao->fmin_changed     =PETSC_FALSE;
170:   tao->catol_changed    =PETSC_FALSE;
171:   tao->crtol_changed    =PETSC_FALSE;
172:   TaoResetStatistics(tao);
173:   *newtao = tao;
174:   return(0);
175: }

177: /*@
178:   TaoSolve - Solves an optimization problem min F(x) s.t. l <= x <= u

180:   Collective on Tao

182:   Input Parameters:
183: . tao - the Tao context

185:   Notes:
186:   The user must set up the Tao with calls to TaoSetInitialVector(),
187:   TaoSetObjectiveRoutine(),
188:   TaoSetGradientRoutine(), and (if using 2nd order method) TaoSetHessianRoutine().

190:   You should call TaoGetConvergedReason() or run with -tao_converged_reason to determine if the optimization algorithm actually succeeded or
191:   why it failed.

193:   Level: beginner

195: .seealso: TaoCreate(), TaoSetObjectiveRoutine(), TaoSetGradientRoutine(), TaoSetHessianRoutine(), TaoGetConvergedReason()
196:  @*/
197: PetscErrorCode TaoSolve(Tao tao)
198: {
199:   PetscErrorCode   ierr;
200:   static PetscBool set = PETSC_FALSE;

204:   PetscCitationsRegister("@TechReport{tao-user-ref,\n"
205:                                 "title   = {Toolkit for Advanced Optimization (TAO) Users Manual},\n"
206:                                 "author  = {Todd Munson and Jason Sarich and Stefan Wild and Steve Benson and Lois Curfman McInnes},\n"
207:                                 "Institution = {Argonne National Laboratory},\n"
208:                                 "Year   = 2014,\n"
209:                                 "Number = {ANL/MCS-TM-322 - Revision 3.5},\n"
210:                                 "url    = {https://www.mcs.anl.gov/research/projects/tao/}\n}\n",&set);
211:   tao->header_printed = PETSC_FALSE;
212:   TaoSetUp(tao);
213:   TaoResetStatistics(tao);
214:   if (tao->linesearch) {
215:     TaoLineSearchReset(tao->linesearch);
216:   }

218:   PetscLogEventBegin(TAO_Solve,tao,0,0,0);
219:   if (tao->ops->solve){ (*tao->ops->solve)(tao); }
220:   PetscLogEventEnd(TAO_Solve,tao,0,0,0);

222:   VecViewFromOptions(tao->solution,(PetscObject)tao,"-tao_view_solution");

224:   tao->ntotalits += tao->niter;
225:   TaoViewFromOptions(tao,NULL,"-tao_view");

227:   if (tao->printreason) {
228:     if (tao->reason > 0) {
229:       PetscPrintf(((PetscObject)tao)->comm,"TAO solve converged due to %s iterations %D\n",TaoConvergedReasons[tao->reason],tao->niter);
230:     } else {
231:       PetscPrintf(((PetscObject)tao)->comm,"TAO solve did not converge due to %s iteration %D\n",TaoConvergedReasons[tao->reason],tao->niter);
232:     }
233:   }
234:   return(0);
235: }

237: /*@
238:   TaoSetUp - Sets up the internal data structures for the later use
239:   of a Tao solver

241:   Collective on tao

243:   Input Parameters:
244: . tao - the TAO context

246:   Notes:
247:   The user will not need to explicitly call TaoSetUp(), as it will
248:   automatically be called in TaoSolve().  However, if the user
249:   desires to call it explicitly, it should come after TaoCreate()
250:   and any TaoSetSomething() routines, but before TaoSolve().

252:   Level: advanced

254: .seealso: TaoCreate(), TaoSolve()
255: @*/
256: PetscErrorCode TaoSetUp(Tao tao)
257: {

262:   if (tao->setupcalled) return(0);

264:   if (!tao->solution) SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetInitialVector");
265:   if (tao->ops->setup) {
266:     (*tao->ops->setup)(tao);
267:   }
268:   tao->setupcalled = PETSC_TRUE;
269:   return(0);
270: }

272: /*@
273:   TaoDestroy - Destroys the TAO context that was created with
274:   TaoCreate()

276:   Collective on Tao

278:   Input Parameter:
279: . tao - the Tao context

281:   Level: beginner

283: .seealso: TaoCreate(), TaoSolve()
284: @*/
285: PetscErrorCode TaoDestroy(Tao *tao)
286: {

290:   if (!*tao) return(0);
292:   if (--((PetscObject)*tao)->refct > 0) {*tao = NULL;return(0);}

294:   if ((*tao)->ops->destroy) {
295:     (*((*tao))->ops->destroy)(*tao);
296:   }
297:   KSPDestroy(&(*tao)->ksp);
298:   TaoLineSearchDestroy(&(*tao)->linesearch);

300:   if ((*tao)->ops->convergencedestroy) {
301:     (*(*tao)->ops->convergencedestroy)((*tao)->cnvP);
302:     if ((*tao)->jacobian_state_inv) {
303:       MatDestroy(&(*tao)->jacobian_state_inv);
304:     }
305:   }
306:   VecDestroy(&(*tao)->solution);
307:   VecDestroy(&(*tao)->gradient);
308:   VecDestroy(&(*tao)->ls_res);

310:   if ((*tao)->gradient_norm) {
311:     PetscObjectDereference((PetscObject)(*tao)->gradient_norm);
312:     VecDestroy(&(*tao)->gradient_norm_tmp);
313:   }

315:   VecDestroy(&(*tao)->XL);
316:   VecDestroy(&(*tao)->XU);
317:   VecDestroy(&(*tao)->IL);
318:   VecDestroy(&(*tao)->IU);
319:   VecDestroy(&(*tao)->DE);
320:   VecDestroy(&(*tao)->DI);
321:   VecDestroy(&(*tao)->constraints_equality);
322:   VecDestroy(&(*tao)->constraints_inequality);
323:   VecDestroy(&(*tao)->stepdirection);
324:   MatDestroy(&(*tao)->hessian_pre);
325:   MatDestroy(&(*tao)->hessian);
326:   MatDestroy(&(*tao)->ls_jac);
327:   MatDestroy(&(*tao)->ls_jac_pre);
328:   MatDestroy(&(*tao)->jacobian_pre);
329:   MatDestroy(&(*tao)->jacobian);
330:   MatDestroy(&(*tao)->jacobian_state_pre);
331:   MatDestroy(&(*tao)->jacobian_state);
332:   MatDestroy(&(*tao)->jacobian_state_inv);
333:   MatDestroy(&(*tao)->jacobian_design);
334:   MatDestroy(&(*tao)->jacobian_equality);
335:   MatDestroy(&(*tao)->jacobian_equality_pre);
336:   MatDestroy(&(*tao)->jacobian_inequality);
337:   MatDestroy(&(*tao)->jacobian_inequality_pre);
338:   ISDestroy(&(*tao)->state_is);
339:   ISDestroy(&(*tao)->design_is);
340:   VecDestroy(&(*tao)->res_weights_v);
341:   TaoCancelMonitors(*tao);
342:   if ((*tao)->hist_malloc) {
343:     PetscFree4((*tao)->hist_obj,(*tao)->hist_resid,(*tao)->hist_cnorm,(*tao)->hist_lits);
344:   }
345:   if ((*tao)->res_weights_n) {
346:     PetscFree((*tao)->res_weights_rows);
347:     PetscFree((*tao)->res_weights_cols);
348:     PetscFree((*tao)->res_weights_w);
349:   }
350:   PetscHeaderDestroy(tao);
351:   return(0);
352: }

354: /*@
355:   TaoSetFromOptions - Sets various Tao parameters from user
356:   options.

358:   Collective on Tao

360:   Input Paremeter:
361: . tao - the Tao solver context

363:   options Database Keys:
364: + -tao_type <type> - The algorithm that TAO uses (lmvm, nls, etc.)
365: . -tao_gatol <gatol> - absolute error tolerance for ||gradient||
366: . -tao_grtol <grtol> - relative error tolerance for ||gradient||
367: . -tao_gttol <gttol> - reduction of ||gradient|| relative to initial gradient
368: . -tao_max_it <max> - sets maximum number of iterations
369: . -tao_max_funcs <max> - sets maximum number of function evaluations
370: . -tao_fmin <fmin> - stop if function value reaches fmin
371: . -tao_steptol <tol> - stop if trust region radius less than <tol>
372: . -tao_trust0 <t> - initial trust region radius
373: . -tao_monitor - prints function value and residual at each iteration
374: . -tao_smonitor - same as tao_monitor, but truncates very small values
375: . -tao_cmonitor - prints function value, residual, and constraint norm at each iteration
376: . -tao_view_solution - prints solution vector at each iteration
377: . -tao_view_ls_residual - prints least-squares residual vector at each iteration
378: . -tao_view_step - prints step direction vector at each iteration
379: . -tao_view_gradient - prints gradient vector at each iteration
380: . -tao_draw_solution - graphically view solution vector at each iteration
381: . -tao_draw_step - graphically view step vector at each iteration
382: . -tao_draw_gradient - graphically view gradient at each iteration
383: . -tao_fd_gradient - use gradient computed with finite differences
384: . -tao_fd_hessian - use hessian computed with finite differences
385: . -tao_mf_hessian - use matrix-free hessian computed with finite differences
386: . -tao_cancelmonitors - cancels all monitors (except those set with command line)
387: . -tao_view - prints information about the Tao after solving
388: - -tao_converged_reason - prints the reason TAO stopped iterating

390:   Notes:
391:   To see all options, run your program with the -help option or consult the
392:   user's manual. Should be called after TaoCreate() but before TaoSolve()

394:   Level: beginner
395: @*/
396: PetscErrorCode TaoSetFromOptions(Tao tao)
397: {
399:   TaoType        default_type = TAOLMVM;
400:   char           type[256], monfilename[PETSC_MAX_PATH_LEN];
401:   PetscViewer    monviewer;
402:   PetscBool      flg;
403:   MPI_Comm       comm;

407:   PetscObjectGetComm((PetscObject)tao,&comm);

409:   /* So no warnings are given about unused options */
410:   PetscOptionsHasName(((PetscObject)tao)->options,((PetscObject)tao)->prefix,"-tao_ls_type",&flg);

412:   PetscObjectOptionsBegin((PetscObject)tao);
413:   {
414:     TaoRegisterAll();
415:     if (((PetscObject)tao)->type_name) {
416:       default_type = ((PetscObject)tao)->type_name;
417:     }
418:     /* Check for type from options */
419:     PetscOptionsFList("-tao_type","Tao Solver type","TaoSetType",TaoList,default_type,type,256,&flg);
420:     if (flg) {
421:       TaoSetType(tao,type);
422:     } else if (!((PetscObject)tao)->type_name) {
423:       TaoSetType(tao,default_type);
424:     }

426:     PetscOptionsReal("-tao_catol","Stop if constraints violations within","TaoSetConstraintTolerances",tao->catol,&tao->catol,&flg);
427:     if (flg) tao->catol_changed=PETSC_TRUE;
428:     PetscOptionsReal("-tao_crtol","Stop if relative contraint violations within","TaoSetConstraintTolerances",tao->crtol,&tao->crtol,&flg);
429:     if (flg) tao->crtol_changed=PETSC_TRUE;
430:     PetscOptionsReal("-tao_gatol","Stop if norm of gradient less than","TaoSetTolerances",tao->gatol,&tao->gatol,&flg);
431:     if (flg) tao->gatol_changed=PETSC_TRUE;
432:     PetscOptionsReal("-tao_grtol","Stop if norm of gradient divided by the function value is less than","TaoSetTolerances",tao->grtol,&tao->grtol,&flg);
433:     if (flg) tao->grtol_changed=PETSC_TRUE;
434:     PetscOptionsReal("-tao_gttol","Stop if the norm of the gradient is less than the norm of the initial gradient times tol","TaoSetTolerances",tao->gttol,&tao->gttol,&flg);
435:     if (flg) tao->gttol_changed=PETSC_TRUE;
436:     PetscOptionsInt("-tao_max_it","Stop if iteration number exceeds","TaoSetMaximumIterations",tao->max_it,&tao->max_it,&flg);
437:     if (flg) tao->max_it_changed=PETSC_TRUE;
438:     PetscOptionsInt("-tao_max_funcs","Stop if number of function evaluations exceeds","TaoSetMaximumFunctionEvaluations",tao->max_funcs,&tao->max_funcs,&flg);
439:     if (flg) tao->max_funcs_changed=PETSC_TRUE;
440:     PetscOptionsReal("-tao_fmin","Stop if function less than","TaoSetFunctionLowerBound",tao->fmin,&tao->fmin,&flg);
441:     if (flg) tao->fmin_changed=PETSC_TRUE;
442:     PetscOptionsReal("-tao_steptol","Stop if step size or trust region radius less than","",tao->steptol,&tao->steptol,&flg);
443:     if (flg) tao->steptol_changed=PETSC_TRUE;
444:     PetscOptionsReal("-tao_trust0","Initial trust region radius","TaoSetTrustRegionRadius",tao->trust0,&tao->trust0,&flg);
445:     if (flg) tao->trust0_changed=PETSC_TRUE;
446:     PetscOptionsString("-tao_view_solution","view solution vector after each evaluation","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);
447:     if (flg) {
448:       PetscViewerASCIIOpen(comm,monfilename,&monviewer);
449:       TaoSetMonitor(tao,TaoSolutionMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
450:     }

452:     PetscOptionsBool("-tao_converged_reason","Print reason for TAO converged","TaoSolve",tao->printreason,&tao->printreason,NULL);
453:     PetscOptionsString("-tao_view_gradient","view gradient vector after each evaluation","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);
454:     if (flg) {
455:       PetscViewerASCIIOpen(comm,monfilename,&monviewer);
456:       TaoSetMonitor(tao,TaoGradientMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
457:     }

459:     PetscOptionsString("-tao_view_stepdirection","view step direction vector after each iteration","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);
460:     if (flg) {
461:       PetscViewerASCIIOpen(comm,monfilename,&monviewer);
462:       TaoSetMonitor(tao,TaoStepDirectionMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
463:     }

465:     PetscOptionsString("-tao_view_residual","view least-squares residual vector after each evaluation","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);
466:     if (flg) {
467:       PetscViewerASCIIOpen(comm,monfilename,&monviewer);
468:       TaoSetMonitor(tao,TaoResidualMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
469:     }

471:     PetscOptionsString("-tao_monitor","Use the default convergence monitor","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);
472:     if (flg) {
473:       PetscViewerASCIIOpen(comm,monfilename,&monviewer);
474:       TaoSetMonitor(tao,TaoMonitorDefault,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
475:     }

477:     PetscOptionsString("-tao_gmonitor","Use the convergence monitor with extra globalization info","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);
478:     if (flg) {
479:       PetscViewerASCIIOpen(comm,monfilename,&monviewer);
480:       TaoSetMonitor(tao,TaoDefaultGMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
481:     }

483:     PetscOptionsString("-tao_smonitor","Use the short convergence monitor","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);
484:     if (flg) {
485:       PetscViewerASCIIOpen(comm,monfilename,&monviewer);
486:       TaoSetMonitor(tao,TaoDefaultSMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
487:     }

489:     PetscOptionsString("-tao_cmonitor","Use the default convergence monitor with constraint norm","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);
490:     if (flg) {
491:       PetscViewerASCIIOpen(comm,monfilename,&monviewer);
492:       TaoSetMonitor(tao,TaoDefaultCMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
493:     }


496:     flg = PETSC_FALSE;
497:     PetscOptionsBool("-tao_cancelmonitors","cancel all monitors and call any registered destroy routines","TaoCancelMonitors",flg,&flg,NULL);
498:     if (flg) {TaoCancelMonitors(tao);}

500:     flg = PETSC_FALSE;
501:     PetscOptionsBool("-tao_draw_solution","Plot solution vector at each iteration","TaoSetMonitor",flg,&flg,NULL);
502:     if (flg) {
503:       TaoMonitorDrawCtx drawctx;
504:       PetscInt          howoften = 1;
505:       TaoMonitorDrawCtxCreate(PetscObjectComm((PetscObject)tao),NULL,NULL,PETSC_DECIDE,PETSC_DECIDE,300,300,howoften,&drawctx);
506:       TaoSetMonitor(tao,TaoDrawSolutionMonitor,drawctx,(PetscErrorCode (*)(void**))TaoMonitorDrawCtxDestroy);
507:     }

509:     flg = PETSC_FALSE;
510:     PetscOptionsBool("-tao_draw_step","plots step direction at each iteration","TaoSetMonitor",flg,&flg,NULL);
511:     if (flg) {
512:       TaoSetMonitor(tao,TaoDrawStepMonitor,NULL,NULL);
513:     }

515:     flg = PETSC_FALSE;
516:     PetscOptionsBool("-tao_draw_gradient","plots gradient at each iteration","TaoSetMonitor",flg,&flg,NULL);
517:     if (flg) {
518:       TaoMonitorDrawCtx drawctx;
519:       PetscInt          howoften = 1;
520:       TaoMonitorDrawCtxCreate(PetscObjectComm((PetscObject)tao),NULL,NULL,PETSC_DECIDE,PETSC_DECIDE,300,300,howoften,&drawctx);
521:       TaoSetMonitor(tao,TaoDrawGradientMonitor,drawctx,(PetscErrorCode (*)(void**))TaoMonitorDrawCtxDestroy);
522:     }
523:     flg = PETSC_FALSE;
524:     PetscOptionsBool("-tao_fd_gradient","compute gradient using finite differences","TaoDefaultComputeGradient",flg,&flg,NULL);
525:     if (flg) {
526:       TaoSetGradientRoutine(tao,TaoDefaultComputeGradient,NULL);
527:     }
528:     flg = PETSC_FALSE;
529:     PetscOptionsBool("-tao_fd_hessian","compute hessian using finite differences","TaoDefaultComputeHessian",flg,&flg,NULL);
530:     if (flg) {
531:       Mat H;

533:       MatCreate(PetscObjectComm((PetscObject)tao),&H);
534:       MatSetType(H,MATAIJ);
535:       TaoSetHessianRoutine(tao,H,H,TaoDefaultComputeHessian,NULL);
536:       MatDestroy(&H);
537:     }
538:     flg = PETSC_FALSE;
539:     PetscOptionsBool("-tao_mf_hessian","compute matrix-free hessian using finite differences","TaoDefaultComputeHessianMFFD",flg,&flg,NULL);
540:     if (flg) {
541:       Mat H;

543:       MatCreate(PetscObjectComm((PetscObject)tao),&H);
544:       TaoSetHessianRoutine(tao,H,H,TaoDefaultComputeHessianMFFD,NULL);
545:       MatDestroy(&H);
546:     }
547:     PetscOptionsEnum("-tao_subset_type","subset type","",TaoSubSetTypes,(PetscEnum)tao->subset_type,(PetscEnum*)&tao->subset_type,NULL);

549:     if (tao->ops->setfromoptions) {
550:       (*tao->ops->setfromoptions)(PetscOptionsObject,tao);
551:     }
552:   }
553:   PetscOptionsEnd();
554:   return(0);
555: }

557: /*@C
558:    TaoViewFromOptions - View from Options

560:    Collective on Tao

562:    Input Parameters:
563: +  A - the  Tao context
564: .  obj - Optional object
565: -  name - command line option

567:    Level: intermediate
568: .seealso:  Tao, TaoView, PetscObjectViewFromOptions(), TaoCreate()
569: @*/
570: PetscErrorCode  TaoViewFromOptions(Tao A,PetscObject obj,const char name[])
571: {

576:   PetscObjectViewFromOptions((PetscObject)A,obj,name);
577:   return(0);
578: }

580: /*@C
581:   TaoView - Prints information about the Tao

583:   Collective on Tao

585:   InputParameters:
586: + tao - the Tao context
587: - viewer - visualization context

589:   Options Database Key:
590: . -tao_view - Calls TaoView() at the end of TaoSolve()

592:   Notes:
593:   The available visualization contexts include
594: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
595: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
596:          output where only the first processor opens
597:          the file.  All other processors send their
598:          data to the first processor to print.

600:   Level: beginner

602: .seealso: PetscViewerASCIIOpen()
603: @*/
604: PetscErrorCode TaoView(Tao tao, PetscViewer viewer)
605: {
606:   PetscErrorCode      ierr;
607:   PetscBool           isascii,isstring;
608:   TaoType             type;

612:   if (!viewer) {
613:     PetscViewerASCIIGetStdout(((PetscObject)tao)->comm,&viewer);
614:   }

618:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
619:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);
620:   if (isascii) {
621:     PetscObjectPrintClassNamePrefixType((PetscObject)tao,viewer);

623:     if (tao->ops->view) {
624:       PetscViewerASCIIPushTab(viewer);
625:       (*tao->ops->view)(tao,viewer);
626:       PetscViewerASCIIPopTab(viewer);
627:     }
628:     if (tao->linesearch) {
629:       PetscViewerASCIIPushTab(viewer);
630:       TaoLineSearchView(tao->linesearch,viewer);
631:       PetscViewerASCIIPopTab(viewer);
632:     }
633:     if (tao->ksp) {
634:       PetscViewerASCIIPushTab(viewer);
635:       KSPView(tao->ksp,viewer);
636:       PetscViewerASCIIPrintf(viewer,"total KSP iterations: %D\n",tao->ksp_tot_its);
637:       PetscViewerASCIIPopTab(viewer);
638:     }

640:     PetscViewerASCIIPushTab(viewer);

642:     if (tao->XL || tao->XU) {
643:       PetscViewerASCIIPrintf(viewer,"Active Set subset type: %s\n",TaoSubSetTypes[tao->subset_type]);
644:     }

646:     PetscViewerASCIIPrintf(viewer,"convergence tolerances: gatol=%g,",(double)tao->gatol);
647:     PetscViewerASCIIPrintf(viewer," steptol=%g,",(double)tao->steptol);
648:     PetscViewerASCIIPrintf(viewer," gttol=%g\n",(double)tao->gttol);
649:     PetscViewerASCIIPrintf(viewer,"Residual in Function/Gradient:=%g\n",(double)tao->residual);

651:     if (tao->cnorm>0 || tao->catol>0 || tao->crtol>0){
652:       ierr=PetscViewerASCIIPrintf(viewer,"convergence tolerances:");
653:       ierr=PetscViewerASCIIPrintf(viewer," catol=%g,",(double)tao->catol);
654:       ierr=PetscViewerASCIIPrintf(viewer," crtol=%g\n",(double)tao->crtol);
655:       PetscViewerASCIIPrintf(viewer,"Residual in Constraints:=%g\n",(double)tao->cnorm);
656:     }

658:     if (tao->trust < tao->steptol){
659:       ierr=PetscViewerASCIIPrintf(viewer,"convergence tolerances: steptol=%g\n",(double)tao->steptol);
660:       ierr=PetscViewerASCIIPrintf(viewer,"Final trust region radius:=%g\n",(double)tao->trust);
661:     }

663:     if (tao->fmin>-1.e25){
664:       ierr=PetscViewerASCIIPrintf(viewer,"convergence tolerances: function minimum=%g\n",(double)tao->fmin);
665:     }
666:     PetscViewerASCIIPrintf(viewer,"Objective value=%g\n",(double)tao->fc);

668:     PetscViewerASCIIPrintf(viewer,"total number of iterations=%D,          ",tao->niter);
669:     PetscViewerASCIIPrintf(viewer,"              (max: %D)\n",tao->max_it);

671:     if (tao->nfuncs>0){
672:       PetscViewerASCIIPrintf(viewer,"total number of function evaluations=%D,",tao->nfuncs);
673:       PetscViewerASCIIPrintf(viewer,"                max: %D\n",tao->max_funcs);
674:     }
675:     if (tao->ngrads>0){
676:       PetscViewerASCIIPrintf(viewer,"total number of gradient evaluations=%D,",tao->ngrads);
677:       PetscViewerASCIIPrintf(viewer,"                max: %D\n",tao->max_funcs);
678:     }
679:     if (tao->nfuncgrads>0){
680:       PetscViewerASCIIPrintf(viewer,"total number of function/gradient evaluations=%D,",tao->nfuncgrads);
681:       PetscViewerASCIIPrintf(viewer,"    (max: %D)\n",tao->max_funcs);
682:     }
683:     if (tao->nhess>0){
684:       PetscViewerASCIIPrintf(viewer,"total number of Hessian evaluations=%D\n",tao->nhess);
685:     }
686:     /*  if (tao->linear_its>0){
687:      PetscViewerASCIIPrintf(viewer,"  total Krylov method iterations=%D\n",tao->linear_its);
688:      }*/
689:     if (tao->nconstraints>0){
690:       PetscViewerASCIIPrintf(viewer,"total number of constraint function evaluations=%D\n",tao->nconstraints);
691:     }
692:     if (tao->njac>0){
693:       PetscViewerASCIIPrintf(viewer,"total number of Jacobian evaluations=%D\n",tao->njac);
694:     }

696:     if (tao->reason>0){
697:       PetscViewerASCIIPrintf(viewer,    "Solution converged: ");
698:       switch (tao->reason) {
699:       case TAO_CONVERGED_GATOL:
700:         PetscViewerASCIIPrintf(viewer," ||g(X)|| <= gatol\n");
701:         break;
702:       case TAO_CONVERGED_GRTOL:
703:         PetscViewerASCIIPrintf(viewer," ||g(X)||/|f(X)| <= grtol\n");
704:         break;
705:       case TAO_CONVERGED_GTTOL:
706:         PetscViewerASCIIPrintf(viewer," ||g(X)||/||g(X0)|| <= gttol\n");
707:         break;
708:       case TAO_CONVERGED_STEPTOL:
709:         PetscViewerASCIIPrintf(viewer," Steptol -- step size small\n");
710:         break;
711:       case TAO_CONVERGED_MINF:
712:         PetscViewerASCIIPrintf(viewer," Minf --  f < fmin\n");
713:         break;
714:       case TAO_CONVERGED_USER:
715:         PetscViewerASCIIPrintf(viewer," User Terminated\n");
716:         break;
717:       default:
718:         PetscViewerASCIIPrintf(viewer,"\n");
719:         break;
720:       }

722:     } else {
723:       PetscViewerASCIIPrintf(viewer,"Solver terminated: %d",tao->reason);
724:       switch (tao->reason) {
725:       case TAO_DIVERGED_MAXITS:
726:         PetscViewerASCIIPrintf(viewer," Maximum Iterations\n");
727:         break;
728:       case TAO_DIVERGED_NAN:
729:         PetscViewerASCIIPrintf(viewer," NAN or Inf encountered\n");
730:         break;
731:       case TAO_DIVERGED_MAXFCN:
732:         PetscViewerASCIIPrintf(viewer," Maximum Function Evaluations\n");
733:         break;
734:       case TAO_DIVERGED_LS_FAILURE:
735:         PetscViewerASCIIPrintf(viewer," Line Search Failure\n");
736:         break;
737:       case TAO_DIVERGED_TR_REDUCTION:
738:         PetscViewerASCIIPrintf(viewer," Trust Region too small\n");
739:         break;
740:       case TAO_DIVERGED_USER:
741:         PetscViewerASCIIPrintf(viewer," User Terminated\n");
742:         break;
743:       default:
744:         PetscViewerASCIIPrintf(viewer,"\n");
745:         break;
746:       }
747:     }
748:     PetscViewerASCIIPopTab(viewer);
749:   } else if (isstring) {
750:     TaoGetType(tao,&type);
751:     PetscViewerStringSPrintf(viewer," %-3.3s",type);
752:   }
753:   return(0);
754: }

756: /*@
757:   TaoSetTolerances - Sets parameters used in TAO convergence tests

759:   Logically collective on Tao

761:   Input Parameters:
762: + tao - the Tao context
763: . gatol - stop if norm of gradient is less than this
764: . grtol - stop if relative norm of gradient is less than this
765: - gttol - stop if norm of gradient is reduced by this factor

767:   Options Database Keys:
768: + -tao_gatol <gatol> - Sets gatol
769: . -tao_grtol <grtol> - Sets grtol
770: - -tao_gttol <gttol> - Sets gttol

772:   Stopping Criteria:
773: $ ||g(X)||                            <= gatol
774: $ ||g(X)|| / |f(X)|                   <= grtol
775: $ ||g(X)|| / ||g(X0)||                <= gttol

777:   Notes:
778:   Use PETSC_DEFAULT to leave one or more tolerances unchanged.

780:   Level: beginner

782: .seealso: TaoGetTolerances()

784: @*/
785: PetscErrorCode TaoSetTolerances(Tao tao, PetscReal gatol, PetscReal grtol, PetscReal gttol)
786: {


792:   if (gatol != PETSC_DEFAULT) {
793:     if (gatol<0) {
794:       PetscInfo(tao,"Tried to set negative gatol -- ignored.\n");
795:     } else {
796:       tao->gatol = PetscMax(0,gatol);
797:       tao->gatol_changed=PETSC_TRUE;
798:     }
799:   }

801:   if (grtol != PETSC_DEFAULT) {
802:     if (grtol<0) {
803:       PetscInfo(tao,"Tried to set negative grtol -- ignored.\n");
804:     } else {
805:       tao->grtol = PetscMax(0,grtol);
806:       tao->grtol_changed=PETSC_TRUE;
807:     }
808:   }

810:   if (gttol != PETSC_DEFAULT) {
811:     if (gttol<0) {
812:       PetscInfo(tao,"Tried to set negative gttol -- ignored.\n");
813:     } else {
814:       tao->gttol = PetscMax(0,gttol);
815:       tao->gttol_changed=PETSC_TRUE;
816:     }
817:   }
818:   return(0);
819: }

821: /*@
822:   TaoSetConstraintTolerances - Sets constraint tolerance parameters used in TAO  convergence tests

824:   Logically collective on Tao

826:   Input Parameters:
827: + tao - the Tao context
828: . catol - absolute constraint tolerance, constraint norm must be less than catol for used for gatol convergence criteria
829: - crtol - relative contraint tolerance, constraint norm must be less than crtol for used for gatol, gttol convergence criteria

831:   Options Database Keys:
832: + -tao_catol <catol> - Sets catol
833: - -tao_crtol <crtol> - Sets crtol

835:   Notes:
836:   Use PETSC_DEFAULT to leave any tolerance unchanged.

838:   Level: intermediate

840: .seealso: TaoGetTolerances(), TaoGetConstraintTolerances(), TaoSetTolerances()

842: @*/
843: PetscErrorCode TaoSetConstraintTolerances(Tao tao, PetscReal catol, PetscReal crtol)
844: {


850:   if (catol != PETSC_DEFAULT) {
851:     if (catol<0) {
852:       PetscInfo(tao,"Tried to set negative catol -- ignored.\n");
853:     } else {
854:       tao->catol = PetscMax(0,catol);
855:       tao->catol_changed=PETSC_TRUE;
856:     }
857:   }

859:   if (crtol != PETSC_DEFAULT) {
860:     if (crtol<0) {
861:       PetscInfo(tao,"Tried to set negative crtol -- ignored.\n");
862:     } else {
863:       tao->crtol = PetscMax(0,crtol);
864:       tao->crtol_changed=PETSC_TRUE;
865:     }
866:   }
867:   return(0);
868: }

870: /*@
871:   TaoGetConstraintTolerances - Gets constraint tolerance parameters used in TAO  convergence tests

873:   Not ollective

875:   Input Parameter:
876: . tao - the Tao context

878:   Output Parameter:
879: + catol - absolute constraint tolerance, constraint norm must be less than catol for used for gatol convergence criteria
880: - crtol - relative contraint tolerance, constraint norm must be less than crtol for used for gatol, gttol convergence criteria

882:   Level: intermediate

884: .seealso: TaoGetTolerances(), TaoSetTolerances(), TaoSetConstraintTolerances()

886: @*/
887: PetscErrorCode TaoGetConstraintTolerances(Tao tao, PetscReal *catol, PetscReal *crtol)
888: {
891:   if (catol) *catol = tao->catol;
892:   if (crtol) *crtol = tao->crtol;
893:   return(0);
894: }

896: /*@
897:    TaoSetFunctionLowerBound - Sets a bound on the solution objective value.
898:    When an approximate solution with an objective value below this number
899:    has been found, the solver will terminate.

901:    Logically Collective on Tao

903:    Input Parameters:
904: +  tao - the Tao solver context
905: -  fmin - the tolerance

907:    Options Database Keys:
908: .    -tao_fmin <fmin> - sets the minimum function value

910:    Level: intermediate

912: .seealso: TaoSetTolerances()
913: @*/
914: PetscErrorCode TaoSetFunctionLowerBound(Tao tao,PetscReal fmin)
915: {
918:   tao->fmin = fmin;
919:   tao->fmin_changed=PETSC_TRUE;
920:   return(0);
921: }

923: /*@
924:    TaoGetFunctionLowerBound - Gets the bound on the solution objective value.
925:    When an approximate solution with an objective value below this number
926:    has been found, the solver will terminate.

928:    Not collective on Tao

930:    Input Parameters:
931: .  tao - the Tao solver context

933:    OutputParameters:
934: .  fmin - the minimum function value

936:    Level: intermediate

938: .seealso: TaoSetFunctionLowerBound()
939: @*/
940: PetscErrorCode TaoGetFunctionLowerBound(Tao tao,PetscReal *fmin)
941: {
944:   *fmin = tao->fmin;
945:   return(0);
946: }

948: /*@
949:    TaoSetMaximumFunctionEvaluations - Sets a maximum number of
950:    function evaluations.

952:    Logically Collective on Tao

954:    Input Parameters:
955: +  tao - the Tao solver context
956: -  nfcn - the maximum number of function evaluations (>=0)

958:    Options Database Keys:
959: .    -tao_max_funcs <nfcn> - sets the maximum number of function evaluations

961:    Level: intermediate

963: .seealso: TaoSetTolerances(), TaoSetMaximumIterations()
964: @*/

966: PetscErrorCode TaoSetMaximumFunctionEvaluations(Tao tao,PetscInt nfcn)
967: {
970:   tao->max_funcs = PetscMax(0,nfcn);
971:   tao->max_funcs_changed=PETSC_TRUE;
972:   return(0);
973: }

975: /*@
976:    TaoGetMaximumFunctionEvaluations - Sets a maximum number of
977:    function evaluations.

979:    Not Collective

981:    Input Parameters:
982: .  tao - the Tao solver context

984:    Output Parameters:
985: .  nfcn - the maximum number of function evaluations

987:    Level: intermediate

989: .seealso: TaoSetMaximumFunctionEvaluations(), TaoGetMaximumIterations()
990: @*/

992: PetscErrorCode TaoGetMaximumFunctionEvaluations(Tao tao,PetscInt *nfcn)
993: {
996:   *nfcn = tao->max_funcs;
997:   return(0);
998: }

1000: /*@
1001:    TaoGetCurrentFunctionEvaluations - Get current number of
1002:    function evaluations.

1004:    Not Collective

1006:    Input Parameters:
1007: .  tao - the Tao solver context

1009:    Output Parameters:
1010: .  nfuncs - the current number of function evaluations

1012:    Level: intermediate

1014: .seealso: TaoSetMaximumFunctionEvaluations(), TaoGetMaximumFunctionEvaluations(), TaoGetMaximumIterations()
1015: @*/

1017: PetscErrorCode TaoGetCurrentFunctionEvaluations(Tao tao,PetscInt *nfuncs)
1018: {
1021:   *nfuncs=PetscMax(tao->nfuncs,tao->nfuncgrads);
1022:   return(0);
1023: }

1025: /*@
1026:    TaoSetMaximumIterations - Sets a maximum number of iterates.

1028:    Logically Collective on Tao

1030:    Input Parameters:
1031: +  tao - the Tao solver context
1032: -  maxits - the maximum number of iterates (>=0)

1034:    Options Database Keys:
1035: .    -tao_max_it <its> - sets the maximum number of iterations

1037:    Level: intermediate

1039: .seealso: TaoSetTolerances(), TaoSetMaximumFunctionEvaluations()
1040: @*/
1041: PetscErrorCode TaoSetMaximumIterations(Tao tao,PetscInt maxits)
1042: {
1045:   tao->max_it = PetscMax(0,maxits);
1046:   tao->max_it_changed=PETSC_TRUE;
1047:   return(0);
1048: }

1050: /*@
1051:    TaoGetMaximumIterations - Sets a maximum number of iterates.

1053:    Not Collective

1055:    Input Parameters:
1056: .  tao - the Tao solver context

1058:    Output Parameters:
1059: .  maxits - the maximum number of iterates

1061:    Level: intermediate

1063: .seealso: TaoSetMaximumIterations(), TaoGetMaximumFunctionEvaluations()
1064: @*/
1065: PetscErrorCode TaoGetMaximumIterations(Tao tao,PetscInt *maxits)
1066: {
1069:   *maxits = tao->max_it;
1070:   return(0);
1071: }

1073: /*@
1074:    TaoSetInitialTrustRegionRadius - Sets the initial trust region radius.

1076:    Logically collective on Tao

1078:    Input Parameter:
1079: +  tao - a TAO optimization solver
1080: -  radius - the trust region radius

1082:    Level: intermediate

1084:    Options Database Key:
1085: .  -tao_trust0 <t0> - sets initial trust region radius

1087: .seealso: TaoGetTrustRegionRadius(), TaoSetTrustRegionTolerance()
1088: @*/
1089: PetscErrorCode TaoSetInitialTrustRegionRadius(Tao tao, PetscReal radius)
1090: {
1093:   tao->trust0 = PetscMax(0.0,radius);
1094:   tao->trust0_changed=PETSC_TRUE;
1095:   return(0);
1096: }

1098: /*@
1099:    TaoGetInitialTrustRegionRadius - Sets the initial trust region radius.

1101:    Not Collective

1103:    Input Parameter:
1104: .  tao - a TAO optimization solver

1106:    Output Parameter:
1107: .  radius - the trust region radius

1109:    Level: intermediate

1111: .seealso: TaoSetInitialTrustRegionRadius(), TaoGetCurrentTrustRegionRadius()
1112: @*/
1113: PetscErrorCode TaoGetInitialTrustRegionRadius(Tao tao, PetscReal *radius)
1114: {
1117:   *radius = tao->trust0;
1118:   return(0);
1119: }

1121: /*@
1122:    TaoGetCurrentTrustRegionRadius - Gets the current trust region radius.

1124:    Not Collective

1126:    Input Parameter:
1127: .  tao - a TAO optimization solver

1129:    Output Parameter:
1130: .  radius - the trust region radius

1132:    Level: intermediate

1134: .seealso: TaoSetInitialTrustRegionRadius(), TaoGetInitialTrustRegionRadius()
1135: @*/
1136: PetscErrorCode TaoGetCurrentTrustRegionRadius(Tao tao, PetscReal *radius)
1137: {
1140:   *radius = tao->trust;
1141:   return(0);
1142: }

1144: /*@
1145:   TaoGetTolerances - gets the current values of tolerances

1147:   Not Collective

1149:   Input Parameters:
1150: . tao - the Tao context

1152:   Output Parameters:
1153: + gatol - stop if norm of gradient is less than this
1154: . grtol - stop if relative norm of gradient is less than this
1155: - gttol - stop if norm of gradient is reduced by a this factor

1157:   Note: NULL can be used as an argument if not all tolerances values are needed

1159: .seealso TaoSetTolerances()

1161:   Level: intermediate
1162: @*/
1163: PetscErrorCode TaoGetTolerances(Tao tao, PetscReal *gatol, PetscReal *grtol, PetscReal *gttol)
1164: {
1167:   if (gatol) *gatol=tao->gatol;
1168:   if (grtol) *grtol=tao->grtol;
1169:   if (gttol) *gttol=tao->gttol;
1170:   return(0);
1171: }

1173: /*@
1174:   TaoGetKSP - Gets the linear solver used by the optimization solver.
1175:   Application writers should use TaoGetKSP if they need direct access
1176:   to the PETSc KSP object.

1178:   Not Collective

1180:    Input Parameters:
1181: .  tao - the TAO solver

1183:    Output Parameters:
1184: .  ksp - the KSP linear solver used in the optimization solver

1186:    Level: intermediate

1188: @*/
1189: PetscErrorCode TaoGetKSP(Tao tao, KSP *ksp)
1190: {
1192:   *ksp = tao->ksp;
1193:   return(0);
1194: }

1196: /*@
1197:    TaoGetLinearSolveIterations - Gets the total number of linear iterations
1198:    used by the TAO solver

1200:    Not Collective

1202:    Input Parameter:
1203: .  tao - TAO context

1205:    Output Parameter:
1206: .  lits - number of linear iterations

1208:    Notes:
1209:    This counter is reset to zero for each successive call to TaoSolve()

1211:    Level: intermediate

1213: .seealso:  TaoGetKSP()
1214: @*/
1215: PetscErrorCode  TaoGetLinearSolveIterations(Tao tao,PetscInt *lits)
1216: {
1220:   *lits = tao->ksp_tot_its;
1221:   return(0);
1222: }

1224: /*@
1225:   TaoGetLineSearch - Gets the line search used by the optimization solver.
1226:   Application writers should use TaoGetLineSearch if they need direct access
1227:   to the TaoLineSearch object.

1229:   Not Collective

1231:    Input Parameters:
1232: .  tao - the TAO solver

1234:    Output Parameters:
1235: .  ls - the line search used in the optimization solver

1237:    Level: intermediate

1239: @*/
1240: PetscErrorCode TaoGetLineSearch(Tao tao, TaoLineSearch *ls)
1241: {
1243:   *ls = tao->linesearch;
1244:   return(0);
1245: }

1247: /*@
1248:   TaoAddLineSearchCounts - Adds the number of function evaluations spent
1249:   in the line search to the running total.

1251:    Input Parameters:
1252: +  tao - the TAO solver
1253: -  ls - the line search used in the optimization solver

1255:    Level: developer

1257: .seealso: TaoLineSearchApply()
1258: @*/
1259: PetscErrorCode TaoAddLineSearchCounts(Tao tao)
1260: {
1262:   PetscBool      flg;
1263:   PetscInt       nfeval,ngeval,nfgeval;

1267:   if (tao->linesearch) {
1268:     TaoLineSearchIsUsingTaoRoutines(tao->linesearch,&flg);
1269:     if (!flg) {
1270:       TaoLineSearchGetNumberFunctionEvaluations(tao->linesearch,&nfeval,&ngeval,&nfgeval);
1271:       tao->nfuncs+=nfeval;
1272:       tao->ngrads+=ngeval;
1273:       tao->nfuncgrads+=nfgeval;
1274:     }
1275:   }
1276:   return(0);
1277: }

1279: /*@
1280:   TaoGetSolutionVector - Returns the vector with the current TAO solution

1282:   Not Collective

1284:   Input Parameter:
1285: . tao - the Tao context

1287:   Output Parameter:
1288: . X - the current solution

1290:   Level: intermediate

1292:   Note:  The returned vector will be the same object that was passed into TaoSetInitialVector()
1293: @*/
1294: PetscErrorCode TaoGetSolutionVector(Tao tao, Vec *X)
1295: {
1298:   *X = tao->solution;
1299:   return(0);
1300: }

1302: /*@
1303:   TaoGetGradientVector - Returns the vector with the current TAO gradient

1305:   Not Collective

1307:   Input Parameter:
1308: . tao - the Tao context

1310:   Output Parameter:
1311: . G - the current solution

1313:   Level: intermediate
1314: @*/
1315: PetscErrorCode TaoGetGradientVector(Tao tao, Vec *G)
1316: {
1319:   *G = tao->gradient;
1320:   return(0);
1321: }

1323: /*@
1324:    TaoResetStatistics - Initialize the statistics used by TAO for all of the solvers.
1325:    These statistics include the iteration number, residual norms, and convergence status.
1326:    This routine gets called before solving each optimization problem.

1328:    Collective on Tao

1330:    Input Parameters:
1331: .  solver - the Tao context

1333:    Level: developer

1335: .seealso: TaoCreate(), TaoSolve()
1336: @*/
1337: PetscErrorCode TaoResetStatistics(Tao tao)
1338: {
1341:   tao->niter        = 0;
1342:   tao->nfuncs       = 0;
1343:   tao->nfuncgrads   = 0;
1344:   tao->ngrads       = 0;
1345:   tao->nhess        = 0;
1346:   tao->njac         = 0;
1347:   tao->nconstraints = 0;
1348:   tao->ksp_its      = 0;
1349:   tao->ksp_tot_its  = 0;
1350:   tao->reason       = TAO_CONTINUE_ITERATING;
1351:   tao->residual     = 0.0;
1352:   tao->cnorm        = 0.0;
1353:   tao->step         = 0.0;
1354:   tao->lsflag       = PETSC_FALSE;
1355:   if (tao->hist_reset) tao->hist_len=0;
1356:   return(0);
1357: }

1359: /*@C
1360:   TaoSetUpdate - Sets the general-purpose update function called
1361:   at the beginning of every iteration of the nonlinear solve. Specifically
1362:   it is called at the top of every iteration, after the new solution and the gradient
1363:   is determined, but before the Hessian is computed (if applicable).

1365:   Logically Collective on Tao

1367:   Input Parameters:
1368: + tao - The tao solver context
1369: - func - The function

1371:   Calling sequence of func:
1372: $ func (Tao tao, PetscInt step);

1374: . step - The current step of the iteration

1376:   Level: advanced

1378: .seealso TaoSolve()
1379: @*/
1380: PetscErrorCode  TaoSetUpdate(Tao tao, PetscErrorCode (*func)(Tao, PetscInt,void*), void *ctx)
1381: {
1384:   tao->ops->update = func;
1385:   tao->user_update = ctx;
1386:   return(0);
1387: }

1389: /*@C
1390:   TaoSetConvergenceTest - Sets the function that is to be used to test
1391:   for convergence o fthe iterative minimization solution.  The new convergence
1392:   testing routine will replace TAO's default convergence test.

1394:   Logically Collective on Tao

1396:   Input Parameters:
1397: + tao - the Tao object
1398: . conv - the routine to test for convergence
1399: - ctx - [optional] context for private data for the convergence routine
1400:         (may be NULL)

1402:   Calling sequence of conv:
1403: $   PetscErrorCode conv(Tao tao, void *ctx)

1405: + tao - the Tao object
1406: - ctx - [optional] convergence context

1408:   Note: The new convergence testing routine should call TaoSetConvergedReason().

1410:   Level: advanced

1412: .seealso: TaoSetConvergedReason(), TaoGetSolutionStatus(), TaoGetTolerances(), TaoSetMonitor

1414: @*/
1415: PetscErrorCode TaoSetConvergenceTest(Tao tao, PetscErrorCode (*conv)(Tao,void*), void *ctx)
1416: {
1419:   (tao)->ops->convergencetest = conv;
1420:   (tao)->cnvP = ctx;
1421:   return(0);
1422: }

1424: /*@C
1425:    TaoSetMonitor - Sets an ADDITIONAL function that is to be used at every
1426:    iteration of the solver to display the iteration's
1427:    progress.

1429:    Logically Collective on Tao

1431:    Input Parameters:
1432: +  tao - the Tao solver context
1433: .  mymonitor - monitoring routine
1434: -  mctx - [optional] user-defined context for private data for the
1435:           monitor routine (may be NULL)

1437:    Calling sequence of mymonitor:
1438: $     PetscErrorCode mymonitor(Tao tao,void *mctx)

1440: +    tao - the Tao solver context
1441: -    mctx - [optional] monitoring context


1444:    Options Database Keys:
1445: +    -tao_monitor        - sets TaoMonitorDefault()
1446: .    -tao_smonitor       - sets short monitor
1447: .    -tao_cmonitor       - same as smonitor plus constraint norm
1448: .    -tao_view_solution   - view solution at each iteration
1449: .    -tao_view_gradient   - view gradient at each iteration
1450: .    -tao_view_ls_residual - view least-squares residual vector at each iteration
1451: -    -tao_cancelmonitors - cancels all monitors that have been hardwired into a code by calls to TaoSetMonitor(), but does not cancel those set via the options database.


1454:    Notes:
1455:    Several different monitoring routines may be set by calling
1456:    TaoSetMonitor() multiple times; all will be called in the
1457:    order in which they were set.

1459:    Fortran Notes:
1460:     Only one monitor function may be set

1462:    Level: intermediate

1464: .seealso: TaoMonitorDefault(), TaoCancelMonitors(),  TaoSetDestroyRoutine()
1465: @*/
1466: PetscErrorCode TaoSetMonitor(Tao tao, PetscErrorCode (*func)(Tao, void*), void *ctx,PetscErrorCode (*dest)(void**))
1467: {
1469:   PetscInt       i;
1470:   PetscBool      identical;

1474:   if (tao->numbermonitors >= MAXTAOMONITORS) SETERRQ1(PetscObjectComm((PetscObject)tao),1,"Cannot attach another monitor -- max=",MAXTAOMONITORS);

1476:   for (i=0; i<tao->numbermonitors;i++) {
1477:     PetscMonitorCompare((PetscErrorCode (*)(void))func,ctx,dest,(PetscErrorCode (*)(void))tao->monitor[i],tao->monitorcontext[i],tao->monitordestroy[i],&identical);
1478:     if (identical) return(0);
1479:   }
1480:   tao->monitor[tao->numbermonitors] = func;
1481:   tao->monitorcontext[tao->numbermonitors] = (void*)ctx;
1482:   tao->monitordestroy[tao->numbermonitors] = dest;
1483:   ++tao->numbermonitors;
1484:   return(0);
1485: }

1487: /*@
1488:    TaoCancelMonitors - Clears all the monitor functions for a Tao object.

1490:    Logically Collective on Tao

1492:    Input Parameters:
1493: .  tao - the Tao solver context

1495:    Options Database:
1496: .  -tao_cancelmonitors - cancels all monitors that have been hardwired
1497:     into a code by calls to TaoSetMonitor(), but does not cancel those
1498:     set via the options database

1500:    Notes:
1501:    There is no way to clear one specific monitor from a Tao object.

1503:    Level: advanced

1505: .seealso: TaoMonitorDefault(), TaoSetMonitor()
1506: @*/
1507: PetscErrorCode TaoCancelMonitors(Tao tao)
1508: {
1509:   PetscInt       i;

1514:   for (i=0;i<tao->numbermonitors;i++) {
1515:     if (tao->monitordestroy[i]) {
1516:       (*tao->monitordestroy[i])(&tao->monitorcontext[i]);
1517:     }
1518:   }
1519:   tao->numbermonitors=0;
1520:   return(0);
1521: }

1523: /*@
1524:    TaoMonitorDefault - Default routine for monitoring progress of the
1525:    Tao solvers (default).  This monitor prints the function value and gradient
1526:    norm at each iteration.  It can be turned on from the command line using the
1527:    -tao_monitor option

1529:    Collective on Tao

1531:    Input Parameters:
1532: +  tao - the Tao context
1533: -  ctx - PetscViewer context or NULL

1535:    Options Database Keys:
1536: .  -tao_monitor

1538:    Level: advanced

1540: .seealso: TaoDefaultSMonitor(), TaoSetMonitor()
1541: @*/
1542: PetscErrorCode TaoMonitorDefault(Tao tao, void *ctx)
1543: {
1545:   PetscInt       its, tabs;
1546:   PetscReal      fct,gnorm;
1547:   PetscViewer    viewer = (PetscViewer)ctx;

1551:   its=tao->niter;
1552:   fct=tao->fc;
1553:   gnorm=tao->residual;
1554:   PetscViewerASCIIGetTab(viewer, &tabs);
1555:   PetscViewerASCIISetTab(viewer, ((PetscObject)tao)->tablevel);
1556:   if (its == 0 && ((PetscObject)tao)->prefix && !tao->header_printed) {
1557:      PetscViewerASCIIPrintf(viewer,"  Iteration information for %s solve.\n",((PetscObject)tao)->prefix);
1558:      tao->header_printed = PETSC_TRUE;
1559:    }
1560:   ierr=PetscViewerASCIIPrintf(viewer,"%3D TAO,",its);
1561:   ierr=PetscViewerASCIIPrintf(viewer,"  Function value: %g,",(double)fct);
1562:   if (gnorm >= PETSC_INFINITY) {
1563:     ierr=PetscViewerASCIIPrintf(viewer,"  Residual: Inf \n");
1564:   } else {
1565:     ierr=PetscViewerASCIIPrintf(viewer,"  Residual: %g \n",(double)gnorm);
1566:   }
1567:   PetscViewerASCIISetTab(viewer, tabs);
1568:   return(0);
1569: }

1571: /*@
1572:    TaoDefaultGMonitor - Default routine for monitoring progress of the
1573:    Tao solvers (default) with extra detail on the globalization method.
1574:    This monitor prints the function value and gradient norm at each
1575:    iteration, as well as the step size and trust radius. Note that the
1576:    step size and trust radius may be the same for some algorithms.
1577:    It can be turned on from the command line using the
1578:    -tao_gmonitor option

1580:    Collective on Tao

1582:    Input Parameters:
1583: +  tao - the Tao context
1584: -  ctx - PetscViewer context or NULL

1586:    Options Database Keys:
1587: .  -tao_monitor

1589:    Level: advanced

1591: .seealso: TaoDefaultSMonitor(), TaoSetMonitor()
1592: @*/
1593: PetscErrorCode TaoDefaultGMonitor(Tao tao, void *ctx)
1594: {
1596:   PetscInt       its, tabs;
1597:   PetscReal      fct,gnorm,stp,tr;
1598:   PetscViewer    viewer = (PetscViewer)ctx;

1602:   its=tao->niter;
1603:   fct=tao->fc;
1604:   gnorm=tao->residual;
1605:   stp=tao->step;
1606:   tr=tao->trust;
1607:   PetscViewerASCIIGetTab(viewer, &tabs);
1608:   PetscViewerASCIISetTab(viewer, ((PetscObject)tao)->tablevel);
1609:   if (its == 0 && ((PetscObject)tao)->prefix && !tao->header_printed) {
1610:      PetscViewerASCIIPrintf(viewer,"  Iteration information for %s solve.\n",((PetscObject)tao)->prefix);
1611:      tao->header_printed = PETSC_TRUE;
1612:    }
1613:   ierr=PetscViewerASCIIPrintf(viewer,"%3D TAO,",its);
1614:   ierr=PetscViewerASCIIPrintf(viewer,"  Function value: %g,",(double)fct);
1615:   if (gnorm >= PETSC_INFINITY) {
1616:     ierr=PetscViewerASCIIPrintf(viewer,"  Residual: Inf,");
1617:   } else {
1618:     ierr=PetscViewerASCIIPrintf(viewer,"  Residual: %g,",(double)gnorm);
1619:   }
1620:   PetscViewerASCIIPrintf(viewer,"  Step: %g,  Trust: %g\n",(double)stp,(double)tr);
1621:   PetscViewerASCIISetTab(viewer, tabs);
1622:   return(0);
1623: }

1625: /*@
1626:    TaoDefaultSMonitor - Default routine for monitoring progress of the
1627:    solver. Same as TaoMonitorDefault() except
1628:    it prints fewer digits of the residual as the residual gets smaller.
1629:    This is because the later digits are meaningless and are often
1630:    different on different machines; by using this routine different
1631:    machines will usually generate the same output. It can be turned on
1632:    by using the -tao_smonitor option

1634:    Collective on Tao

1636:    Input Parameters:
1637: +  tao - the Tao context
1638: -  ctx - PetscViewer context of type ASCII

1640:    Options Database Keys:
1641: .  -tao_smonitor

1643:    Level: advanced

1645: .seealso: TaoMonitorDefault(), TaoSetMonitor()
1646: @*/
1647: PetscErrorCode TaoDefaultSMonitor(Tao tao, void *ctx)
1648: {
1650:   PetscInt       its, tabs;
1651:   PetscReal      fct,gnorm;
1652:   PetscViewer    viewer = (PetscViewer)ctx;

1656:   its=tao->niter;
1657:   fct=tao->fc;
1658:   gnorm=tao->residual;
1659:   PetscViewerASCIIGetTab(viewer, &tabs);
1660:   PetscViewerASCIISetTab(viewer, ((PetscObject)tao)->tablevel);
1661:   ierr=PetscViewerASCIIPrintf(viewer,"iter = %3D,",its);
1662:   ierr=PetscViewerASCIIPrintf(viewer," Function value %g,",(double)fct);
1663:   if (gnorm >= PETSC_INFINITY) {
1664:     ierr=PetscViewerASCIIPrintf(viewer," Residual: Inf \n");
1665:   } else if (gnorm > 1.e-6) {
1666:     ierr=PetscViewerASCIIPrintf(viewer," Residual: %g \n",(double)gnorm);
1667:   } else if (gnorm > 1.e-11) {
1668:     ierr=PetscViewerASCIIPrintf(viewer," Residual: < 1.0e-6 \n");
1669:   } else {
1670:     ierr=PetscViewerASCIIPrintf(viewer," Residual: < 1.0e-11 \n");
1671:   }
1672:   PetscViewerASCIISetTab(viewer, tabs);
1673:   return(0);
1674: }

1676: /*@
1677:    TaoDefaultCMonitor - same as TaoMonitorDefault() except
1678:    it prints the norm of the constraints function. It can be turned on
1679:    from the command line using the -tao_cmonitor option

1681:    Collective on Tao

1683:    Input Parameters:
1684: +  tao - the Tao context
1685: -  ctx - PetscViewer context or NULL

1687:    Options Database Keys:
1688: .  -tao_cmonitor

1690:    Level: advanced

1692: .seealso: TaoMonitorDefault(), TaoSetMonitor()
1693: @*/
1694: PetscErrorCode TaoDefaultCMonitor(Tao tao, void *ctx)
1695: {
1697:   PetscInt       its, tabs;
1698:   PetscReal      fct,gnorm;
1699:   PetscViewer    viewer = (PetscViewer)ctx;

1703:   its=tao->niter;
1704:   fct=tao->fc;
1705:   gnorm=tao->residual;
1706:   PetscViewerASCIIGetTab(viewer, &tabs);
1707:   PetscViewerASCIISetTab(viewer, ((PetscObject)tao)->tablevel);
1708:   ierr=PetscViewerASCIIPrintf(viewer,"iter = %D,",its);
1709:   ierr=PetscViewerASCIIPrintf(viewer," Function value: %g,",(double)fct);
1710:   ierr=PetscViewerASCIIPrintf(viewer,"  Residual: %g ",(double)gnorm);
1711:   PetscViewerASCIIPrintf(viewer,"  Constraint: %g \n",(double)tao->cnorm);
1712:   PetscViewerASCIISetTab(viewer, tabs);
1713:   return(0);
1714: }

1716: /*@C
1717:    TaoSolutionMonitor - Views the solution at each iteration
1718:    It can be turned on from the command line using the
1719:    -tao_view_solution option

1721:    Collective on Tao

1723:    Input Parameters:
1724: +  tao - the Tao context
1725: -  ctx - PetscViewer context or NULL

1727:    Options Database Keys:
1728: .  -tao_view_solution

1730:    Level: advanced

1732: .seealso: TaoDefaultSMonitor(), TaoSetMonitor()
1733: @*/
1734: PetscErrorCode TaoSolutionMonitor(Tao tao, void *ctx)
1735: {
1737:   PetscViewer    viewer  = (PetscViewer)ctx;

1741:   VecView(tao->solution, viewer);
1742:   return(0);
1743: }

1745: /*@C
1746:    TaoGradientMonitor - Views the gradient at each iteration
1747:    It can be turned on from the command line using the
1748:    -tao_view_gradient option

1750:    Collective on Tao

1752:    Input Parameters:
1753: +  tao - the Tao context
1754: -  ctx - PetscViewer context or NULL

1756:    Options Database Keys:
1757: .  -tao_view_gradient

1759:    Level: advanced

1761: .seealso: TaoDefaultSMonitor(), TaoSetMonitor()
1762: @*/
1763: PetscErrorCode TaoGradientMonitor(Tao tao, void *ctx)
1764: {
1766:   PetscViewer    viewer = (PetscViewer)ctx;

1770:   VecView(tao->gradient, viewer);
1771:   return(0);
1772: }

1774: /*@C
1775:    TaoStepDirectionMonitor - Views the gradient at each iteration
1776:    It can be turned on from the command line using the
1777:    -tao_view_gradient option

1779:    Collective on Tao

1781:    Input Parameters:
1782: +  tao - the Tao context
1783: -  ctx - PetscViewer context or NULL

1785:    Options Database Keys:
1786: .  -tao_view_gradient

1788:    Level: advanced

1790: .seealso: TaoDefaultSMonitor(), TaoSetMonitor()
1791: @*/
1792: PetscErrorCode TaoStepDirectionMonitor(Tao tao, void *ctx)
1793: {
1795:   PetscViewer    viewer = (PetscViewer)ctx;

1799:   VecView(tao->stepdirection, viewer);
1800:   return(0);
1801: }

1803: /*@C
1804:    TaoDrawSolutionMonitor - Plots the solution at each iteration
1805:    It can be turned on from the command line using the
1806:    -tao_draw_solution option

1808:    Collective on Tao

1810:    Input Parameters:
1811: +  tao - the Tao context
1812: -  ctx - TaoMonitorDraw context

1814:    Options Database Keys:
1815: .  -tao_draw_solution

1817:    Level: advanced

1819: .seealso: TaoSolutionMonitor(), TaoSetMonitor(), TaoDrawGradientMonitor
1820: @*/
1821: PetscErrorCode TaoDrawSolutionMonitor(Tao tao, void *ctx)
1822: {
1823:   PetscErrorCode    ierr;
1824:   TaoMonitorDrawCtx ictx = (TaoMonitorDrawCtx)ctx;

1827:   if (!(((ictx->howoften > 0) && (!(tao->niter % ictx->howoften))) || ((ictx->howoften == -1) && tao->reason))) return(0);
1828:   VecView(tao->solution,ictx->viewer);
1829:   return(0);
1830: }

1832: /*@C
1833:    TaoDrawGradientMonitor - Plots the gradient at each iteration
1834:    It can be turned on from the command line using the
1835:    -tao_draw_gradient option

1837:    Collective on Tao

1839:    Input Parameters:
1840: +  tao - the Tao context
1841: -  ctx - PetscViewer context

1843:    Options Database Keys:
1844: .  -tao_draw_gradient

1846:    Level: advanced

1848: .seealso: TaoGradientMonitor(), TaoSetMonitor(), TaoDrawSolutionMonitor
1849: @*/
1850: PetscErrorCode TaoDrawGradientMonitor(Tao tao, void *ctx)
1851: {
1852:   PetscErrorCode    ierr;
1853:   TaoMonitorDrawCtx ictx = (TaoMonitorDrawCtx)ctx;

1856:   if (!(((ictx->howoften > 0) && (!(tao->niter % ictx->howoften))) || ((ictx->howoften == -1) && tao->reason))) return(0);
1857:   VecView(tao->gradient,ictx->viewer);
1858:   return(0);
1859: }

1861: /*@C
1862:    TaoDrawStepMonitor - Plots the step direction at each iteration
1863:    It can be turned on from the command line using the
1864:    -tao_draw_step option

1866:    Collective on Tao

1868:    Input Parameters:
1869: +  tao - the Tao context
1870: -  ctx - PetscViewer context

1872:    Options Database Keys:
1873: .  -tao_draw_step

1875:    Level: advanced

1877: .seealso: TaoSetMonitor(), TaoDrawSolutionMonitor
1878: @*/
1879: PetscErrorCode TaoDrawStepMonitor(Tao tao, void *ctx)
1880: {
1882:   PetscViewer    viewer = (PetscViewer)(ctx);

1885:   VecView(tao->stepdirection, viewer);
1886:   return(0);
1887: }

1889: /*@C
1890:    TaoResidualMonitor - Views the least-squares residual at each iteration
1891:    It can be turned on from the command line using the
1892:    -tao_view_ls_residual option

1894:    Collective on Tao

1896:    Input Parameters:
1897: +  tao - the Tao context
1898: -  ctx - PetscViewer context or NULL

1900:    Options Database Keys:
1901: .  -tao_view_ls_residual

1903:    Level: advanced

1905: .seealso: TaoDefaultSMonitor(), TaoSetMonitor()
1906: @*/
1907: PetscErrorCode TaoResidualMonitor(Tao tao, void *ctx)
1908: {
1910:   PetscViewer    viewer  = (PetscViewer)ctx;

1914:   VecView(tao->ls_res,viewer);
1915:   return(0);
1916: }

1918: /*@
1919:    TaoDefaultConvergenceTest - Determines whether the solver should continue iterating
1920:    or terminate.

1922:    Collective on Tao

1924:    Input Parameters:
1925: +  tao - the Tao context
1926: -  dummy - unused dummy context

1928:    Output Parameter:
1929: .  reason - for terminating

1931:    Notes:
1932:    This routine checks the residual in the optimality conditions, the
1933:    relative residual in the optimity conditions, the number of function
1934:    evaluations, and the function value to test convergence.  Some
1935:    solvers may use different convergence routines.

1937:    Level: developer

1939: .seealso: TaoSetTolerances(),TaoGetConvergedReason(),TaoSetConvergedReason()
1940: @*/

1942: PetscErrorCode TaoDefaultConvergenceTest(Tao tao,void *dummy)
1943: {
1944:   PetscInt           niter=tao->niter, nfuncs=PetscMax(tao->nfuncs,tao->nfuncgrads);
1945:   PetscInt           max_funcs=tao->max_funcs;
1946:   PetscReal          gnorm=tao->residual, gnorm0=tao->gnorm0;
1947:   PetscReal          f=tao->fc, steptol=tao->steptol,trradius=tao->step;
1948:   PetscReal          gatol=tao->gatol,grtol=tao->grtol,gttol=tao->gttol;
1949:   PetscReal          catol=tao->catol,crtol=tao->crtol;
1950:   PetscReal          fmin=tao->fmin, cnorm=tao->cnorm;
1951:   TaoConvergedReason reason=tao->reason;
1952:   PetscErrorCode     ierr;

1956:   if (reason != TAO_CONTINUE_ITERATING) {
1957:     return(0);
1958:   }

1960:   if (PetscIsInfOrNanReal(f)) {
1961:     PetscInfo(tao,"Failed to converged, function value is Inf or NaN\n");
1962:     reason = TAO_DIVERGED_NAN;
1963:   } else if (f <= fmin && cnorm <=catol) {
1964:     PetscInfo2(tao,"Converged due to function value %g < minimum function value %g\n", (double)f,(double)fmin);
1965:     reason = TAO_CONVERGED_MINF;
1966:   } else if (gnorm<= gatol && cnorm <=catol) {
1967:     PetscInfo2(tao,"Converged due to residual norm ||g(X)||=%g < %g\n",(double)gnorm,(double)gatol);
1968:     reason = TAO_CONVERGED_GATOL;
1969:   } else if ( f!=0 && PetscAbsReal(gnorm/f) <= grtol && cnorm <= crtol) {
1970:     PetscInfo2(tao,"Converged due to residual ||g(X)||/|f(X)| =%g < %g\n",(double)(gnorm/f),(double)grtol);
1971:     reason = TAO_CONVERGED_GRTOL;
1972:   } else if (gnorm0 != 0 && ((gttol == 0 && gnorm == 0) || gnorm/gnorm0 < gttol) && cnorm <= crtol) {
1973:     PetscInfo2(tao,"Converged due to relative residual norm ||g(X)||/||g(X0)|| = %g < %g\n",(double)(gnorm/gnorm0),(double)gttol);
1974:     reason = TAO_CONVERGED_GTTOL;
1975:   } else if (nfuncs > max_funcs){
1976:     PetscInfo2(tao,"Exceeded maximum number of function evaluations: %D > %D\n", nfuncs,max_funcs);
1977:     reason = TAO_DIVERGED_MAXFCN;
1978:   } else if ( tao->lsflag != 0 ){
1979:     PetscInfo(tao,"Tao Line Search failure.\n");
1980:     reason = TAO_DIVERGED_LS_FAILURE;
1981:   } else if (trradius < steptol && niter > 0){
1982:     PetscInfo2(tao,"Trust region/step size too small: %g < %g\n", (double)trradius,(double)steptol);
1983:     reason = TAO_CONVERGED_STEPTOL;
1984:   } else if (niter >= tao->max_it) {
1985:     PetscInfo2(tao,"Exceeded maximum number of iterations: %D > %D\n",niter,tao->max_it);
1986:     reason = TAO_DIVERGED_MAXITS;
1987:   } else {
1988:     reason = TAO_CONTINUE_ITERATING;
1989:   }
1990:   tao->reason = reason;
1991:   return(0);
1992: }

1994: /*@C
1995:    TaoSetOptionsPrefix - Sets the prefix used for searching for all
1996:    TAO options in the database.


1999:    Logically Collective on Tao

2001:    Input Parameters:
2002: +  tao - the Tao context
2003: -  prefix - the prefix string to prepend to all TAO option requests

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

2009:    For example, to distinguish between the runtime options for two
2010:    different TAO solvers, one could call
2011: .vb
2012:       TaoSetOptionsPrefix(tao1,"sys1_")
2013:       TaoSetOptionsPrefix(tao2,"sys2_")
2014: .ve

2016:    This would enable use of different options for each system, such as
2017: .vb
2018:       -sys1_tao_method blmvm -sys1_tao_gtol 1.e-3
2019:       -sys2_tao_method lmvm  -sys2_tao_gtol 1.e-4
2020: .ve


2023:    Level: advanced

2025: .seealso: TaoAppendOptionsPrefix(), TaoGetOptionsPrefix()
2026: @*/

2028: PetscErrorCode TaoSetOptionsPrefix(Tao tao, const char p[])
2029: {

2033:   PetscObjectSetOptionsPrefix((PetscObject)tao,p);
2034:   if (tao->linesearch) {
2035:     TaoLineSearchSetOptionsPrefix(tao->linesearch,p);
2036:   }
2037:   if (tao->ksp) {
2038:     KSPSetOptionsPrefix(tao->ksp,p);
2039:   }
2040:   return(0);
2041: }

2043: /*@C
2044:    TaoAppendOptionsPrefix - Appends to the prefix used for searching for all
2045:    TAO options in the database.


2048:    Logically Collective on Tao

2050:    Input Parameters:
2051: +  tao - the Tao solver context
2052: -  prefix - the prefix string to prepend to all TAO option requests

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


2059:    Level: advanced

2061: .seealso: TaoSetOptionsPrefix(), TaoGetOptionsPrefix()
2062: @*/
2063: PetscErrorCode TaoAppendOptionsPrefix(Tao tao, const char p[])
2064: {

2068:   PetscObjectAppendOptionsPrefix((PetscObject)tao,p);
2069:   if (tao->linesearch) {
2070:     TaoLineSearchSetOptionsPrefix(tao->linesearch,p);
2071:   }
2072:   if (tao->ksp) {
2073:     KSPSetOptionsPrefix(tao->ksp,p);
2074:   }
2075:   return(0);
2076: }

2078: /*@C
2079:   TaoGetOptionsPrefix - Gets the prefix used for searching for all
2080:   TAO options in the database

2082:   Not Collective

2084:   Input Parameters:
2085: . tao - the Tao context

2087:   Output Parameters:
2088: . prefix - pointer to the prefix string used is returned

2090:   Notes:
2091:     On the fortran side, the user should pass in a string 'prefix' of
2092:   sufficient length to hold the prefix.

2094:   Level: advanced

2096: .seealso: TaoSetOptionsPrefix(), TaoAppendOptionsPrefix()
2097: @*/
2098: PetscErrorCode TaoGetOptionsPrefix(Tao tao, const char *p[])
2099: {
2100:    return PetscObjectGetOptionsPrefix((PetscObject)tao,p);
2101: }

2103: /*@C
2104:    TaoSetType - Sets the method for the unconstrained minimization solver.

2106:    Collective on Tao

2108:    Input Parameters:
2109: +  solver - the Tao solver context
2110: -  type - a known method

2112:    Options Database Key:
2113: .  -tao_type <type> - Sets the method; use -help for a list
2114:    of available methods (for instance, "-tao_type lmvm" or "-tao_type tron")

2116:    Available methods include:
2117: +    nls - Newton's method with line search for unconstrained minimization
2118: .    ntr - Newton's method with trust region for unconstrained minimization
2119: .    ntl - Newton's method with trust region, line search for unconstrained minimization
2120: .    lmvm - Limited memory variable metric method for unconstrained minimization
2121: .    cg - Nonlinear conjugate gradient method for unconstrained minimization
2122: .    nm - Nelder-Mead algorithm for derivate-free unconstrained minimization
2123: .    tron - Newton Trust Region method for bound constrained minimization
2124: .    gpcg - Newton Trust Region method for quadratic bound constrained minimization
2125: .    blmvm - Limited memory variable metric method for bound constrained minimization
2126: -    pounders - Model-based algorithm pounder extended for nonlinear least squares

2128:   Level: intermediate

2130: .seealso: TaoCreate(), TaoGetType(), TaoType

2132: @*/
2133: PetscErrorCode TaoSetType(Tao tao, TaoType type)
2134: {
2136:   PetscErrorCode (*create_xxx)(Tao);
2137:   PetscBool      issame;


2142:   PetscObjectTypeCompare((PetscObject)tao,type,&issame);
2143:   if (issame) return(0);

2145:   PetscFunctionListFind(TaoList, type, (void(**)(void))&create_xxx);
2146:   if (!create_xxx) SETERRQ1(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested Tao type %s",type);

2148:   /* Destroy the existing solver information */
2149:   if (tao->ops->destroy) {
2150:     (*tao->ops->destroy)(tao);
2151:   }
2152:   KSPDestroy(&tao->ksp);
2153:   TaoLineSearchDestroy(&tao->linesearch);
2154:   VecDestroy(&tao->gradient);
2155:   VecDestroy(&tao->stepdirection);

2157:   tao->ops->setup = NULL;
2158:   tao->ops->solve = NULL;
2159:   tao->ops->view  = NULL;
2160:   tao->ops->setfromoptions = NULL;
2161:   tao->ops->destroy = NULL;

2163:   tao->setupcalled = PETSC_FALSE;

2165:   (*create_xxx)(tao);
2166:   PetscObjectChangeTypeName((PetscObject)tao,type);
2167:   return(0);
2168: }

2170: /*MC
2171:    TaoRegister - Adds a method to the TAO package for unconstrained minimization.

2173:    Synopsis:
2174:    TaoRegister(char *name_solver,char *path,char *name_Create,PetscErrorCode (*routine_Create)(Tao))

2176:    Not collective

2178:    Input Parameters:
2179: +  sname - name of a new user-defined solver
2180: -  func - routine to Create method context

2182:    Notes:
2183:    TaoRegister() may be called multiple times to add several user-defined solvers.

2185:    Sample usage:
2186: .vb
2187:    TaoRegister("my_solver",MySolverCreate);
2188: .ve

2190:    Then, your solver can be chosen with the procedural interface via
2191: $     TaoSetType(tao,"my_solver")
2192:    or at runtime via the option
2193: $     -tao_type my_solver

2195:    Level: advanced

2197: .seealso: TaoRegisterAll(), TaoRegisterDestroy()
2198: M*/
2199: PetscErrorCode TaoRegister(const char sname[], PetscErrorCode (*func)(Tao))
2200: {

2204:   TaoInitializePackage();
2205:   PetscFunctionListAdd(&TaoList,sname, (void (*)(void))func);
2206:   return(0);
2207: }

2209: /*@C
2210:    TaoRegisterDestroy - Frees the list of minimization solvers that were
2211:    registered by TaoRegisterDynamic().

2213:    Not Collective

2215:    Level: advanced

2217: .seealso: TaoRegisterAll(), TaoRegister()
2218: @*/
2219: PetscErrorCode TaoRegisterDestroy(void)
2220: {
2223:   PetscFunctionListDestroy(&TaoList);
2224:   TaoRegisterAllCalled = PETSC_FALSE;
2225:   return(0);
2226: }

2228: /*@
2229:    TaoGetIterationNumber - Gets the number of Tao iterations completed
2230:    at this time.

2232:    Not Collective

2234:    Input Parameter:
2235: .  tao - Tao context

2237:    Output Parameter:
2238: .  iter - iteration number

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


2244:    Level: intermediate

2246: .seealso:   TaoGetLinearSolveIterations(), TaoGetResidualNorm(), TaoGetObjective()
2247: @*/
2248: PetscErrorCode  TaoGetIterationNumber(Tao tao,PetscInt *iter)
2249: {
2253:   *iter = tao->niter;
2254:   return(0);
2255: }

2257: /*@
2258:    TaoGetObjective - Gets the current value of the objective function
2259:    at this time.

2261:    Not Collective

2263:    Input Parameter:
2264: .  tao - Tao context

2266:    Output Parameter:
2267: .  value - the current value

2269:    Level: intermediate

2271: .seealso:   TaoGetLinearSolveIterations(), TaoGetIterationNumber(), TaoGetResidualNorm()
2272: @*/
2273: PetscErrorCode  TaoGetObjective(Tao tao,PetscReal *value)
2274: {
2278:   *value = tao->fc;
2279:   return(0);
2280: }

2282: /*@
2283:    TaoGetResidualNorm - Gets the current value of the norm of the residual
2284:    at this time.

2286:    Not Collective

2288:    Input Parameter:
2289: .  tao - Tao context

2291:    Output Parameter:
2292: .  value - the current value

2294:    Level: intermediate

2296:    Developer Note: This is the 2-norm of the residual, we cannot use TaoGetGradientNorm() because that has
2297:                    a different meaning. For some reason Tao sometimes calls the gradient the residual.

2299: .seealso:   TaoGetLinearSolveIterations(), TaoGetIterationNumber(), TaoGetObjective()
2300: @*/
2301: PetscErrorCode  TaoGetResidualNorm(Tao tao,PetscReal *value)
2302: {
2306:   *value = tao->residual;
2307:   return(0);
2308: }

2310: /*@
2311:    TaoSetIterationNumber - Sets the current iteration number.

2313:    Not Collective

2315:    Input Parameter:
2316: +  tao - Tao context
2317: -  iter - iteration number

2319:    Level: developer

2321: .seealso:   TaoGetLinearSolveIterations()
2322: @*/
2323: PetscErrorCode  TaoSetIterationNumber(Tao tao,PetscInt iter)
2324: {

2329:   PetscObjectSAWsTakeAccess((PetscObject)tao);
2330:   tao->niter = iter;
2331:   PetscObjectSAWsGrantAccess((PetscObject)tao);
2332:   return(0);
2333: }

2335: /*@
2336:    TaoGetTotalIterationNumber - Gets the total number of Tao iterations
2337:    completed. This number keeps accumulating if multiple solves
2338:    are called with the Tao object.

2340:    Not Collective

2342:    Input Parameter:
2343: .  tao - Tao context

2345:    Output Parameter:
2346: .  iter - iteration number

2348:    Notes:
2349:    The total iteration count is updated after each solve, if there is a current
2350:    TaoSolve() in progress then those iterations are not yet counted.

2352:    Level: intermediate

2354: .seealso:   TaoGetLinearSolveIterations()
2355: @*/
2356: PetscErrorCode  TaoGetTotalIterationNumber(Tao tao,PetscInt *iter)
2357: {
2361:   *iter = tao->ntotalits;
2362:   return(0);
2363: }

2365: /*@
2366:    TaoSetTotalIterationNumber - Sets the current total iteration number.

2368:    Not Collective

2370:    Input Parameter:
2371: +  tao - Tao context
2372: -  iter - iteration number

2374:    Level: developer

2376: .seealso:   TaoGetLinearSolveIterations()
2377: @*/
2378: PetscErrorCode  TaoSetTotalIterationNumber(Tao tao,PetscInt iter)
2379: {

2384:   PetscObjectSAWsTakeAccess((PetscObject)tao);
2385:   tao->ntotalits = iter;
2386:   PetscObjectSAWsGrantAccess((PetscObject)tao);
2387:   return(0);
2388: }

2390: /*@
2391:   TaoSetConvergedReason - Sets the termination flag on a Tao object

2393:   Logically Collective on Tao

2395:   Input Parameters:
2396: + tao - the Tao context
2397: - reason - one of
2398: $     TAO_CONVERGED_ATOL (2),
2399: $     TAO_CONVERGED_RTOL (3),
2400: $     TAO_CONVERGED_STEPTOL (4),
2401: $     TAO_CONVERGED_MINF (5),
2402: $     TAO_CONVERGED_USER (6),
2403: $     TAO_DIVERGED_MAXITS (-2),
2404: $     TAO_DIVERGED_NAN (-4),
2405: $     TAO_DIVERGED_MAXFCN (-5),
2406: $     TAO_DIVERGED_LS_FAILURE (-6),
2407: $     TAO_DIVERGED_TR_REDUCTION (-7),
2408: $     TAO_DIVERGED_USER (-8),
2409: $     TAO_CONTINUE_ITERATING (0)

2411:    Level: intermediate

2413: @*/
2414: PetscErrorCode TaoSetConvergedReason(Tao tao, TaoConvergedReason reason)
2415: {
2418:   tao->reason = reason;
2419:   return(0);
2420: }

2422: /*@
2423:    TaoGetConvergedReason - Gets the reason the Tao iteration was stopped.

2425:    Not Collective

2427:    Input Parameter:
2428: .  tao - the Tao solver context

2430:    Output Parameter:
2431: .  reason - one of
2432: $  TAO_CONVERGED_GATOL (3)           ||g(X)|| < gatol
2433: $  TAO_CONVERGED_GRTOL (4)           ||g(X)|| / f(X)  < grtol
2434: $  TAO_CONVERGED_GTTOL (5)           ||g(X)|| / ||g(X0)|| < gttol
2435: $  TAO_CONVERGED_STEPTOL (6)         step size small
2436: $  TAO_CONVERGED_MINF (7)            F < F_min
2437: $  TAO_CONVERGED_USER (8)            User defined
2438: $  TAO_DIVERGED_MAXITS (-2)          its > maxits
2439: $  TAO_DIVERGED_NAN (-4)             Numerical problems
2440: $  TAO_DIVERGED_MAXFCN (-5)          fevals > max_funcsals
2441: $  TAO_DIVERGED_LS_FAILURE (-6)      line search failure
2442: $  TAO_DIVERGED_TR_REDUCTION (-7)    trust region failure
2443: $  TAO_DIVERGED_USER(-8)             (user defined)
2444:  $  TAO_CONTINUE_ITERATING (0)

2446:    where
2447: +  X - current solution
2448: .  X0 - initial guess
2449: .  f(X) - current function value
2450: .  f(X*) - true solution (estimated)
2451: .  g(X) - current gradient
2452: .  its - current iterate number
2453: .  maxits - maximum number of iterates
2454: .  fevals - number of function evaluations
2455: -  max_funcsals - maximum number of function evaluations

2457:    Level: intermediate

2459: .seealso: TaoSetConvergenceTest(), TaoSetTolerances()

2461: @*/
2462: PetscErrorCode TaoGetConvergedReason(Tao tao, TaoConvergedReason *reason)
2463: {
2467:   *reason = tao->reason;
2468:   return(0);
2469: }

2471: /*@
2472:   TaoGetSolutionStatus - Get the current iterate, objective value,
2473:   residual, infeasibility, and termination

2475:   Not Collective

2477:    Input Parameters:
2478: .  tao - the Tao context

2480:    Output Parameters:
2481: +  iterate - the current iterate number (>=0)
2482: .  f - the current function value
2483: .  gnorm - the square of the gradient norm, duality gap, or other measure indicating distance from optimality.
2484: .  cnorm - the infeasibility of the current solution with regard to the constraints.
2485: .  xdiff - the step length or trust region radius of the most recent iterate.
2486: -  reason - The termination reason, which can equal TAO_CONTINUE_ITERATING

2488:    Level: intermediate

2490:    Note:
2491:    TAO returns the values set by the solvers in the routine TaoMonitor().

2493:    Note:
2494:    If any of the output arguments are set to NULL, no corresponding value will be returned.

2496: .seealso: TaoMonitor(), TaoGetConvergedReason()
2497: @*/
2498: PetscErrorCode TaoGetSolutionStatus(Tao tao, PetscInt *its, PetscReal *f, PetscReal *gnorm, PetscReal *cnorm, PetscReal *xdiff, TaoConvergedReason *reason)
2499: {
2501:   if (its) *its=tao->niter;
2502:   if (f) *f=tao->fc;
2503:   if (gnorm) *gnorm=tao->residual;
2504:   if (cnorm) *cnorm=tao->cnorm;
2505:   if (reason) *reason=tao->reason;
2506:   if (xdiff) *xdiff=tao->step;
2507:   return(0);
2508: }

2510: /*@C
2511:    TaoGetType - Gets the current Tao algorithm.

2513:    Not Collective

2515:    Input Parameter:
2516: .  tao - the Tao solver context

2518:    Output Parameter:
2519: .  type - Tao method

2521:    Level: intermediate

2523: @*/
2524: PetscErrorCode TaoGetType(Tao tao,TaoType *type)
2525: {
2529:   *type=((PetscObject)tao)->type_name;
2530:   return(0);
2531: }

2533: /*@C
2534:   TaoMonitor - Monitor the solver and the current solution.  This
2535:   routine will record the iteration number and residual statistics,
2536:   call any monitors specified by the user, and calls the convergence-check routine.

2538:    Input Parameters:
2539: +  tao - the Tao context
2540: .  its - the current iterate number (>=0)
2541: .  f - the current objective function value
2542: .  res - the gradient norm, square root of the duality gap, or other measure indicating distince from optimality.  This measure will be recorded and
2543:           used for some termination tests.
2544: .  cnorm - the infeasibility of the current solution with regard to the constraints.
2545: -  steplength - multiple of the step direction added to the previous iterate.

2547:    Output Parameters:
2548: .  reason - The termination reason, which can equal TAO_CONTINUE_ITERATING

2550:    Options Database Key:
2551: .  -tao_monitor - Use the default monitor, which prints statistics to standard output

2553: .seealso TaoGetConvergedReason(), TaoMonitorDefault(), TaoSetMonitor()

2555:    Level: developer

2557: @*/
2558: PetscErrorCode TaoMonitor(Tao tao, PetscInt its, PetscReal f, PetscReal res, PetscReal cnorm, PetscReal steplength)
2559: {
2561:   PetscInt       i;

2565:   tao->fc = f;
2566:   tao->residual = res;
2567:   tao->cnorm = cnorm;
2568:   tao->step = steplength;
2569:   if (!its) {
2570:     tao->cnorm0 = cnorm; tao->gnorm0 = res;
2571:   }
2572:   if (PetscIsInfOrNanReal(f) || PetscIsInfOrNanReal(res)) SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_USER, "User provided compute function generated Inf or NaN");
2573:   for (i=0;i<tao->numbermonitors;i++) {
2574:     (*tao->monitor[i])(tao,tao->monitorcontext[i]);
2575:   }
2576:   return(0);
2577: }

2579: /*@
2580:    TaoSetConvergenceHistory - Sets the array used to hold the convergence history.

2582:    Logically Collective on Tao

2584:    Input Parameters:
2585: +  tao - the Tao solver context
2586: .  obj   - array to hold objective value history
2587: .  resid - array to hold residual history
2588: .  cnorm - array to hold constraint violation history
2589: .  lits - integer array holds the number of linear iterations for each Tao iteration
2590: .  na  - size of obj, resid, and cnorm
2591: -  reset - PetscTrue indicates each new minimization resets the history counter to zero,
2592:            else it continues storing new values for new minimizations after the old ones

2594:    Notes:
2595:    If set, TAO will fill the given arrays with the indicated
2596:    information at each iteration.  If 'obj','resid','cnorm','lits' are
2597:    *all* NULL then space (using size na, or 1000 if na is PETSC_DECIDE or
2598:    PETSC_DEFAULT) is allocated for the history.
2599:    If not all are NULL, then only the non-NULL information categories
2600:    will be stored, the others will be ignored.

2602:    Any convergence information after iteration number 'na' will not be stored.

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

2608:    Level: intermediate

2610: .seealso: TaoGetConvergenceHistory()

2612: @*/
2613: PetscErrorCode TaoSetConvergenceHistory(Tao tao, PetscReal obj[], PetscReal resid[], PetscReal cnorm[], PetscInt lits[], PetscInt na,PetscBool reset)
2614: {


2624:   if (na == PETSC_DECIDE || na == PETSC_DEFAULT) na = 1000;
2625:   if (!obj && !resid && !cnorm && !lits) {
2626:     PetscCalloc4(na,&obj,na,&resid,na,&cnorm,na,&lits);
2627:     tao->hist_malloc = PETSC_TRUE;
2628:   }

2630:   tao->hist_obj = obj;
2631:   tao->hist_resid = resid;
2632:   tao->hist_cnorm = cnorm;
2633:   tao->hist_lits = lits;
2634:   tao->hist_max   = na;
2635:   tao->hist_reset = reset;
2636:   tao->hist_len = 0;
2637:   return(0);
2638: }

2640: /*@C
2641:    TaoGetConvergenceHistory - Gets the arrays used to hold the convergence history.

2643:    Collective on Tao

2645:    Input Parameter:
2646: .  tao - the Tao context

2648:    Output Parameters:
2649: +  obj   - array used to hold objective value history
2650: .  resid - array used to hold residual history
2651: .  cnorm - array used to hold constraint violation history
2652: .  lits  - integer array used to hold linear solver iteration count
2653: -  nhist  - size of obj, resid, cnorm, and lits

2655:    Notes:
2656:     This routine must be preceded by calls to TaoSetConvergenceHistory()
2657:     and TaoSolve(), otherwise it returns useless information.

2659:     The calling sequence for this routine in Fortran is
2660: $   call TaoGetConvergenceHistory(Tao tao, PetscInt nhist, PetscErrorCode ierr)

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

2666:    Level: advanced

2668: .seealso: TaoSetConvergenceHistory()

2670: @*/
2671: PetscErrorCode TaoGetConvergenceHistory(Tao tao, PetscReal **obj, PetscReal **resid, PetscReal **cnorm, PetscInt **lits, PetscInt *nhist)
2672: {
2675:   if (obj)   *obj   = tao->hist_obj;
2676:   if (cnorm) *cnorm = tao->hist_cnorm;
2677:   if (resid) *resid = tao->hist_resid;
2678:   if (nhist) *nhist = tao->hist_len;
2679:   return(0);
2680: }

2682: /*@
2683:    TaoSetApplicationContext - Sets the optional user-defined context for
2684:    a solver.

2686:    Logically Collective on Tao

2688:    Input Parameters:
2689: +  tao  - the Tao context
2690: -  usrP - optional user context

2692:    Level: intermediate

2694: .seealso: TaoGetApplicationContext(), TaoSetApplicationContext()
2695: @*/
2696: PetscErrorCode  TaoSetApplicationContext(Tao tao,void *usrP)
2697: {
2700:   tao->user = usrP;
2701:   return(0);
2702: }

2704: /*@
2705:    TaoGetApplicationContext - Gets the user-defined context for a
2706:    TAO solvers.

2708:    Not Collective

2710:    Input Parameter:
2711: .  tao  - Tao context

2713:    Output Parameter:
2714: .  usrP - user context

2716:    Level: intermediate

2718: .seealso: TaoSetApplicationContext()
2719: @*/
2720: PetscErrorCode  TaoGetApplicationContext(Tao tao,void *usrP)
2721: {
2724:   *(void**)usrP = tao->user;
2725:   return(0);
2726: }

2728: /*@
2729:    TaoSetGradientNorm - Sets the matrix used to define the inner product that measures the size of the gradient.

2731:    Collective on tao

2733:    Input Parameters:
2734: +  tao  - the Tao context
2735: -  M    - gradient norm

2737:    Level: beginner

2739: .seealso: TaoGetGradientNorm(), TaoGradientNorm()
2740: @*/
2741: PetscErrorCode  TaoSetGradientNorm(Tao tao, Mat M)
2742: {

2747:   PetscObjectReference((PetscObject)M);
2748:   MatDestroy(&tao->gradient_norm);
2749:   VecDestroy(&tao->gradient_norm_tmp);
2750:   tao->gradient_norm = M;
2751:   MatCreateVecs(M, NULL, &tao->gradient_norm_tmp);
2752:   return(0);
2753: }

2755: /*@
2756:    TaoGetGradientNorm - Returns the matrix used to define the inner product for measuring the size of the gradient.

2758:    Not Collective

2760:    Input Parameter:
2761: .  tao  - Tao context

2763:    Output Parameter:
2764: .  M - gradient norm

2766:    Level: beginner

2768: .seealso: TaoSetGradientNorm(), TaoGradientNorm()
2769: @*/
2770: PetscErrorCode  TaoGetGradientNorm(Tao tao, Mat *M)
2771: {
2774:   *M = tao->gradient_norm;
2775:   return(0);
2776: }

2778: /*c
2779:    TaoGradientNorm - Compute the norm with respect to the inner product the user has set.

2781:    Collective on tao

2783:    Input Parameter:
2784: .  tao      - the Tao context
2785: .  gradient - the gradient to be computed
2786: .  norm     - the norm type

2788:    Output Parameter:
2789: .  gnorm    - the gradient norm

2791:    Level: developer

2793: .seealso: TaoSetGradientNorm(), TaoGetGradientNorm()
2794: @*/
2795: PetscErrorCode  TaoGradientNorm(Tao tao, Vec gradient, NormType type, PetscReal *gnorm)
2796: {

2804:   if (tao->gradient_norm) {
2805:     PetscScalar gnorms;

2807:     if (type != NORM_2) SETERRQ(PetscObjectComm((PetscObject)gradient), PETSC_ERR_ARG_WRONG, "Norm type must be NORM_2 if an inner product for the gradient norm is set.");
2808:     MatMult(tao->gradient_norm, gradient, tao->gradient_norm_tmp);
2809:     VecDot(gradient, tao->gradient_norm_tmp, &gnorms);
2810:     *gnorm = PetscRealPart(PetscSqrtScalar(gnorms));
2811:   } else {
2812:     VecNorm(gradient, type, gnorm);
2813:   }
2814:   return(0);
2815: }

2817: /*@C
2818:    TaoMonitorDrawCtxCreate - Creates the monitor context for TaoMonitorDrawCtx

2820:    Collective on Tao

2822:    Output Patameter:
2823: .    ctx - the monitor context

2825:    Options Database:
2826: .   -tao_draw_solution_initial - show initial guess as well as current solution

2828:    Level: intermediate

2830: .seealso: TaoMonitorSet(), TaoMonitorDefault(), VecView(), TaoMonitorDrawCtx()
2831: @*/
2832: PetscErrorCode  TaoMonitorDrawCtxCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscInt howoften,TaoMonitorDrawCtx *ctx)
2833: {
2834:   PetscErrorCode   ierr;

2837:   PetscNew(ctx);
2838:   PetscViewerDrawOpen(comm,host,label,x,y,m,n,&(*ctx)->viewer);
2839:   PetscViewerSetFromOptions((*ctx)->viewer);
2840:   (*ctx)->howoften = howoften;
2841:   return(0);
2842: }

2844: /*@C
2845:    TaoMonitorDrawCtxDestroy - Destroys the monitor context for TaoMonitorDrawSolution()

2847:    Collective on Tao

2849:    Input Parameters:
2850: .    ctx - the monitor context

2852:    Level: intermediate

2854: .seealso: TaoMonitorSet(), TaoMonitorDefault(), VecView(), TaoMonitorDrawSolution()
2855: @*/
2856: PetscErrorCode  TaoMonitorDrawCtxDestroy(TaoMonitorDrawCtx *ictx)
2857: {

2861:   PetscViewerDestroy(&(*ictx)->viewer);
2862:   PetscFree(*ictx);
2863:   return(0);
2864: }