1 #define TAO_DLL 2 3 #include <petsc/private/taoimpl.h> /*I "petsctao.h" I*/ 4 5 PetscBool TaoRegisterAllCalled = PETSC_FALSE; 6 PetscFunctionList TaoList = NULL; 7 8 PetscClassId TAO_CLASSID; 9 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; 17 18 const char *TaoSubSetTypes[] = {"subvec","mask","matrixfree","TaoSubSetType","TAO_SUBSET_",NULL}; 19 20 struct _n_TaoMonitorDrawCtx { 21 PetscViewer viewer; 22 PetscInt howoften; /* when > 0 uses iteration % howoften, when negative only final solution plotted */ 23 }; 24 25 /*@ 26 TaoCreate - Creates a TAO solver 27 28 Collective 29 30 Input Parameter: 31 . comm - MPI communicator 32 33 Output Parameter: 34 . newtao - the new Tao context 35 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 48 49 Options Database Keys: 50 . -tao_type - select which method TAO should use 51 52 Level: beginner 53 54 .seealso: TaoSolve(), TaoDestroy() 55 @*/ 56 PetscErrorCode TaoCreate(MPI_Comm comm, Tao *newtao) 57 { 58 PetscErrorCode ierr; 59 Tao tao; 60 61 PetscFunctionBegin; 62 PetscValidPointer(newtao,2); 63 *newtao = NULL; 64 65 ierr = TaoInitializePackage();CHKERRQ(ierr); 66 ierr = TaoLineSearchInitializePackage();CHKERRQ(ierr); 67 68 ierr = PetscHeaderCreate(tao,TAO_CLASSID,"Tao","Optimization solver","Tao",comm,TaoDestroy,TaoView);CHKERRQ(ierr); 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; 89 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; 125 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 tao->crtol = 1e-5; 132 tao->catol = 1e-5; 133 #else 134 tao->gatol = 1e-8; 135 tao->grtol = 1e-8; 136 tao->crtol = 1e-8; 137 tao->catol = 1e-8; 138 #endif 139 tao->gttol = 0.0; 140 tao->steptol = 0.0; 141 tao->trust0 = PETSC_INFINITY; 142 tao->fmin = PETSC_NINFINITY; 143 144 tao->hist_malloc = PETSC_FALSE; 145 tao->hist_reset = PETSC_TRUE; 146 tao->hist_max = 0; 147 tao->hist_len = 0; 148 tao->hist_obj = NULL; 149 tao->hist_resid = NULL; 150 tao->hist_cnorm = NULL; 151 tao->hist_lits = NULL; 152 153 tao->numbermonitors = 0; 154 tao->viewsolution = PETSC_FALSE; 155 tao->viewhessian = PETSC_FALSE; 156 tao->viewgradient = PETSC_FALSE; 157 tao->viewjacobian = PETSC_FALSE; 158 tao->viewconstraints = PETSC_FALSE; 159 160 tao->bounded = PETSC_FALSE; 161 tao->constrained = PETSC_FALSE; 162 tao->eq_constrained = PETSC_FALSE; 163 tao->ineq_constrained = PETSC_FALSE; 164 tao->ineq_doublesided = PETSC_FALSE; 165 166 tao->recycle = PETSC_FALSE; 167 168 tao->header_printed = PETSC_FALSE; 169 170 /* These flags prevents algorithms from overriding user options */ 171 tao->max_it_changed = PETSC_FALSE; 172 tao->max_funcs_changed = PETSC_FALSE; 173 tao->gatol_changed = PETSC_FALSE; 174 tao->grtol_changed = PETSC_FALSE; 175 tao->gttol_changed = PETSC_FALSE; 176 tao->steptol_changed = PETSC_FALSE; 177 tao->trust0_changed = PETSC_FALSE; 178 tao->fmin_changed = PETSC_FALSE; 179 tao->catol_changed = PETSC_FALSE; 180 tao->crtol_changed = PETSC_FALSE; 181 ierr = TaoResetStatistics(tao);CHKERRQ(ierr); 182 *newtao = tao; 183 PetscFunctionReturn(0); 184 } 185 186 /*@ 187 TaoSolve - Solves an optimization problem min F(x) s.t. l <= x <= u 188 189 Collective on Tao 190 191 Input Parameters: 192 . tao - the Tao context 193 194 Notes: 195 The user must set up the Tao with calls to TaoSetInitialVector(), 196 TaoSetObjectiveRoutine(), 197 TaoSetGradientRoutine(), and (if using 2nd order method) TaoSetHessianRoutine(). 198 199 You should call TaoGetConvergedReason() or run with -tao_converged_reason to determine if the optimization algorithm actually succeeded or 200 why it failed. 201 202 Level: beginner 203 204 .seealso: TaoCreate(), TaoSetObjectiveRoutine(), TaoSetGradientRoutine(), TaoSetHessianRoutine(), TaoGetConvergedReason() 205 @*/ 206 PetscErrorCode TaoSolve(Tao tao) 207 { 208 PetscErrorCode ierr; 209 static PetscBool set = PETSC_FALSE; 210 211 PetscFunctionBegin; 212 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 213 ierr = PetscCitationsRegister("@TechReport{tao-user-ref,\n" 214 "title = {Toolkit for Advanced Optimization (TAO) Users Manual},\n" 215 "author = {Todd Munson and Jason Sarich and Stefan Wild and Steve Benson and Lois Curfman McInnes},\n" 216 "Institution = {Argonne National Laboratory},\n" 217 "Year = 2014,\n" 218 "Number = {ANL/MCS-TM-322 - Revision 3.5},\n" 219 "url = {https://www.mcs.anl.gov/research/projects/tao/}\n}\n",&set);CHKERRQ(ierr); 220 tao->header_printed = PETSC_FALSE; 221 ierr = TaoSetUp(tao);CHKERRQ(ierr); 222 ierr = TaoResetStatistics(tao);CHKERRQ(ierr); 223 if (tao->linesearch) { 224 ierr = TaoLineSearchReset(tao->linesearch);CHKERRQ(ierr); 225 } 226 227 ierr = PetscLogEventBegin(TAO_Solve,tao,0,0,0);CHKERRQ(ierr); 228 if (tao->ops->solve) { ierr = (*tao->ops->solve)(tao);CHKERRQ(ierr); } 229 ierr = PetscLogEventEnd(TAO_Solve,tao,0,0,0);CHKERRQ(ierr); 230 231 ierr = VecViewFromOptions(tao->solution,(PetscObject)tao,"-tao_view_solution");CHKERRQ(ierr); 232 233 tao->ntotalits += tao->niter; 234 ierr = TaoViewFromOptions(tao,NULL,"-tao_view");CHKERRQ(ierr); 235 236 if (tao->printreason) { 237 if (tao->reason > 0) { 238 ierr = PetscPrintf(((PetscObject)tao)->comm,"TAO solve converged due to %s iterations %D\n",TaoConvergedReasons[tao->reason],tao->niter);CHKERRQ(ierr); 239 } else { 240 ierr = PetscPrintf(((PetscObject)tao)->comm,"TAO solve did not converge due to %s iteration %D\n",TaoConvergedReasons[tao->reason],tao->niter);CHKERRQ(ierr); 241 } 242 } 243 PetscFunctionReturn(0); 244 } 245 246 /*@ 247 TaoSetUp - Sets up the internal data structures for the later use 248 of a Tao solver 249 250 Collective on tao 251 252 Input Parameters: 253 . tao - the TAO context 254 255 Notes: 256 The user will not need to explicitly call TaoSetUp(), as it will 257 automatically be called in TaoSolve(). However, if the user 258 desires to call it explicitly, it should come after TaoCreate() 259 and any TaoSetSomething() routines, but before TaoSolve(). 260 261 Level: advanced 262 263 .seealso: TaoCreate(), TaoSolve() 264 @*/ 265 PetscErrorCode TaoSetUp(Tao tao) 266 { 267 PetscErrorCode ierr; 268 269 PetscFunctionBegin; 270 PetscValidHeaderSpecific(tao, TAO_CLASSID,1); 271 if (tao->setupcalled) PetscFunctionReturn(0); 272 273 if (!tao->solution) SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetInitialVector"); 274 if (tao->ops->setup) { 275 ierr = (*tao->ops->setup)(tao);CHKERRQ(ierr); 276 } 277 tao->setupcalled = PETSC_TRUE; 278 PetscFunctionReturn(0); 279 } 280 281 /*@C 282 TaoDestroy - Destroys the TAO context that was created with 283 TaoCreate() 284 285 Collective on Tao 286 287 Input Parameter: 288 . tao - the Tao context 289 290 Level: beginner 291 292 .seealso: TaoCreate(), TaoSolve() 293 @*/ 294 PetscErrorCode TaoDestroy(Tao *tao) 295 { 296 PetscErrorCode ierr; 297 298 PetscFunctionBegin; 299 if (!*tao) PetscFunctionReturn(0); 300 PetscValidHeaderSpecific(*tao,TAO_CLASSID,1); 301 if (--((PetscObject)*tao)->refct > 0) {*tao = NULL;PetscFunctionReturn(0);} 302 303 if ((*tao)->ops->destroy) { 304 ierr = (*((*tao))->ops->destroy)(*tao);CHKERRQ(ierr); 305 } 306 ierr = KSPDestroy(&(*tao)->ksp);CHKERRQ(ierr); 307 ierr = TaoLineSearchDestroy(&(*tao)->linesearch);CHKERRQ(ierr); 308 309 if ((*tao)->ops->convergencedestroy) { 310 ierr = (*(*tao)->ops->convergencedestroy)((*tao)->cnvP);CHKERRQ(ierr); 311 if ((*tao)->jacobian_state_inv) { 312 ierr = MatDestroy(&(*tao)->jacobian_state_inv);CHKERRQ(ierr); 313 } 314 } 315 ierr = VecDestroy(&(*tao)->solution);CHKERRQ(ierr); 316 ierr = VecDestroy(&(*tao)->gradient);CHKERRQ(ierr); 317 ierr = VecDestroy(&(*tao)->ls_res);CHKERRQ(ierr); 318 319 if ((*tao)->gradient_norm) { 320 ierr = PetscObjectDereference((PetscObject)(*tao)->gradient_norm);CHKERRQ(ierr); 321 ierr = VecDestroy(&(*tao)->gradient_norm_tmp);CHKERRQ(ierr); 322 } 323 324 ierr = VecDestroy(&(*tao)->XL);CHKERRQ(ierr); 325 ierr = VecDestroy(&(*tao)->XU);CHKERRQ(ierr); 326 ierr = VecDestroy(&(*tao)->IL);CHKERRQ(ierr); 327 ierr = VecDestroy(&(*tao)->IU);CHKERRQ(ierr); 328 ierr = VecDestroy(&(*tao)->DE);CHKERRQ(ierr); 329 ierr = VecDestroy(&(*tao)->DI);CHKERRQ(ierr); 330 ierr = VecDestroy(&(*tao)->constraints_equality);CHKERRQ(ierr); 331 ierr = VecDestroy(&(*tao)->constraints_inequality);CHKERRQ(ierr); 332 ierr = VecDestroy(&(*tao)->stepdirection);CHKERRQ(ierr); 333 ierr = MatDestroy(&(*tao)->hessian_pre);CHKERRQ(ierr); 334 ierr = MatDestroy(&(*tao)->hessian);CHKERRQ(ierr); 335 ierr = MatDestroy(&(*tao)->ls_jac);CHKERRQ(ierr); 336 ierr = MatDestroy(&(*tao)->ls_jac_pre);CHKERRQ(ierr); 337 ierr = MatDestroy(&(*tao)->jacobian_pre);CHKERRQ(ierr); 338 ierr = MatDestroy(&(*tao)->jacobian);CHKERRQ(ierr); 339 ierr = MatDestroy(&(*tao)->jacobian_state_pre);CHKERRQ(ierr); 340 ierr = MatDestroy(&(*tao)->jacobian_state);CHKERRQ(ierr); 341 ierr = MatDestroy(&(*tao)->jacobian_state_inv);CHKERRQ(ierr); 342 ierr = MatDestroy(&(*tao)->jacobian_design);CHKERRQ(ierr); 343 ierr = MatDestroy(&(*tao)->jacobian_equality);CHKERRQ(ierr); 344 ierr = MatDestroy(&(*tao)->jacobian_equality_pre);CHKERRQ(ierr); 345 ierr = MatDestroy(&(*tao)->jacobian_inequality);CHKERRQ(ierr); 346 ierr = MatDestroy(&(*tao)->jacobian_inequality_pre);CHKERRQ(ierr); 347 ierr = ISDestroy(&(*tao)->state_is);CHKERRQ(ierr); 348 ierr = ISDestroy(&(*tao)->design_is);CHKERRQ(ierr); 349 ierr = VecDestroy(&(*tao)->res_weights_v);CHKERRQ(ierr); 350 ierr = TaoCancelMonitors(*tao);CHKERRQ(ierr); 351 if ((*tao)->hist_malloc) { 352 ierr = PetscFree4((*tao)->hist_obj,(*tao)->hist_resid,(*tao)->hist_cnorm,(*tao)->hist_lits);CHKERRQ(ierr); 353 } 354 if ((*tao)->res_weights_n) { 355 ierr = PetscFree((*tao)->res_weights_rows);CHKERRQ(ierr); 356 ierr = PetscFree((*tao)->res_weights_cols);CHKERRQ(ierr); 357 ierr = PetscFree((*tao)->res_weights_w);CHKERRQ(ierr); 358 } 359 ierr = PetscHeaderDestroy(tao);CHKERRQ(ierr); 360 PetscFunctionReturn(0); 361 } 362 363 /*@ 364 TaoSetFromOptions - Sets various Tao parameters from user 365 options. 366 367 Collective on Tao 368 369 Input Parameter: 370 . tao - the Tao solver context 371 372 options Database Keys: 373 + -tao_type <type> - The algorithm that TAO uses (lmvm, nls, etc.) 374 . -tao_gatol <gatol> - absolute error tolerance for ||gradient|| 375 . -tao_grtol <grtol> - relative error tolerance for ||gradient|| 376 . -tao_gttol <gttol> - reduction of ||gradient|| relative to initial gradient 377 . -tao_max_it <max> - sets maximum number of iterations 378 . -tao_max_funcs <max> - sets maximum number of function evaluations 379 . -tao_fmin <fmin> - stop if function value reaches fmin 380 . -tao_steptol <tol> - stop if trust region radius less than <tol> 381 . -tao_trust0 <t> - initial trust region radius 382 . -tao_monitor - prints function value and residual at each iteration 383 . -tao_smonitor - same as tao_monitor, but truncates very small values 384 . -tao_cmonitor - prints function value, residual, and constraint norm at each iteration 385 . -tao_view_solution - prints solution vector at each iteration 386 . -tao_view_ls_residual - prints least-squares residual vector at each iteration 387 . -tao_view_step - prints step direction vector at each iteration 388 . -tao_view_gradient - prints gradient vector at each iteration 389 . -tao_draw_solution - graphically view solution vector at each iteration 390 . -tao_draw_step - graphically view step vector at each iteration 391 . -tao_draw_gradient - graphically view gradient at each iteration 392 . -tao_fd_gradient - use gradient computed with finite differences 393 . -tao_fd_hessian - use hessian computed with finite differences 394 . -tao_mf_hessian - use matrix-free hessian computed with finite differences 395 . -tao_cancelmonitors - cancels all monitors (except those set with command line) 396 . -tao_view - prints information about the Tao after solving 397 - -tao_converged_reason - prints the reason TAO stopped iterating 398 399 Notes: 400 To see all options, run your program with the -help option or consult the 401 user's manual. Should be called after TaoCreate() but before TaoSolve() 402 403 Level: beginner 404 @*/ 405 PetscErrorCode TaoSetFromOptions(Tao tao) 406 { 407 PetscErrorCode ierr; 408 TaoType default_type = TAOLMVM; 409 char type[256], monfilename[PETSC_MAX_PATH_LEN]; 410 PetscViewer monviewer; 411 PetscBool flg; 412 MPI_Comm comm; 413 414 PetscFunctionBegin; 415 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 416 ierr = PetscObjectGetComm((PetscObject)tao,&comm);CHKERRQ(ierr); 417 418 /* So no warnings are given about unused options */ 419 ierr = PetscOptionsHasName(((PetscObject)tao)->options,((PetscObject)tao)->prefix,"-tao_ls_type",&flg);CHKERRQ(ierr); 420 421 ierr = PetscObjectOptionsBegin((PetscObject)tao);CHKERRQ(ierr); 422 { 423 ierr = TaoRegisterAll();CHKERRQ(ierr); 424 if (((PetscObject)tao)->type_name) { 425 default_type = ((PetscObject)tao)->type_name; 426 } 427 /* Check for type from options */ 428 ierr = PetscOptionsFList("-tao_type","Tao Solver type","TaoSetType",TaoList,default_type,type,256,&flg);CHKERRQ(ierr); 429 if (flg) { 430 ierr = TaoSetType(tao,type);CHKERRQ(ierr); 431 } else if (!((PetscObject)tao)->type_name) { 432 ierr = TaoSetType(tao,default_type);CHKERRQ(ierr); 433 } 434 435 ierr = PetscOptionsReal("-tao_catol","Stop if constraints violations within","TaoSetConstraintTolerances",tao->catol,&tao->catol,&flg);CHKERRQ(ierr); 436 if (flg) tao->catol_changed=PETSC_TRUE; 437 ierr = PetscOptionsReal("-tao_crtol","Stop if relative contraint violations within","TaoSetConstraintTolerances",tao->crtol,&tao->crtol,&flg);CHKERRQ(ierr); 438 if (flg) tao->crtol_changed=PETSC_TRUE; 439 ierr = PetscOptionsReal("-tao_gatol","Stop if norm of gradient less than","TaoSetTolerances",tao->gatol,&tao->gatol,&flg);CHKERRQ(ierr); 440 if (flg) tao->gatol_changed=PETSC_TRUE; 441 ierr = PetscOptionsReal("-tao_grtol","Stop if norm of gradient divided by the function value is less than","TaoSetTolerances",tao->grtol,&tao->grtol,&flg);CHKERRQ(ierr); 442 if (flg) tao->grtol_changed=PETSC_TRUE; 443 ierr = 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);CHKERRQ(ierr); 444 if (flg) tao->gttol_changed=PETSC_TRUE; 445 ierr = PetscOptionsInt("-tao_max_it","Stop if iteration number exceeds","TaoSetMaximumIterations",tao->max_it,&tao->max_it,&flg);CHKERRQ(ierr); 446 if (flg) tao->max_it_changed=PETSC_TRUE; 447 ierr = PetscOptionsInt("-tao_max_funcs","Stop if number of function evaluations exceeds","TaoSetMaximumFunctionEvaluations",tao->max_funcs,&tao->max_funcs,&flg);CHKERRQ(ierr); 448 if (flg) tao->max_funcs_changed=PETSC_TRUE; 449 ierr = PetscOptionsReal("-tao_fmin","Stop if function less than","TaoSetFunctionLowerBound",tao->fmin,&tao->fmin,&flg);CHKERRQ(ierr); 450 if (flg) tao->fmin_changed=PETSC_TRUE; 451 ierr = PetscOptionsReal("-tao_steptol","Stop if step size or trust region radius less than","",tao->steptol,&tao->steptol,&flg);CHKERRQ(ierr); 452 if (flg) tao->steptol_changed=PETSC_TRUE; 453 ierr = PetscOptionsReal("-tao_trust0","Initial trust region radius","TaoSetTrustRegionRadius",tao->trust0,&tao->trust0,&flg);CHKERRQ(ierr); 454 if (flg) tao->trust0_changed=PETSC_TRUE; 455 ierr = PetscOptionsString("-tao_view_solution","view solution vector after each evaluation","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);CHKERRQ(ierr); 456 if (flg) { 457 ierr = PetscViewerASCIIOpen(comm,monfilename,&monviewer);CHKERRQ(ierr); 458 ierr = TaoSetMonitor(tao,TaoSolutionMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);CHKERRQ(ierr); 459 } 460 461 ierr = PetscOptionsBool("-tao_converged_reason","Print reason for TAO converged","TaoSolve",tao->printreason,&tao->printreason,NULL);CHKERRQ(ierr); 462 ierr = PetscOptionsString("-tao_view_gradient","view gradient vector after each evaluation","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);CHKERRQ(ierr); 463 if (flg) { 464 ierr = PetscViewerASCIIOpen(comm,monfilename,&monviewer);CHKERRQ(ierr); 465 ierr = TaoSetMonitor(tao,TaoGradientMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);CHKERRQ(ierr); 466 } 467 468 ierr = PetscOptionsString("-tao_view_stepdirection","view step direction vector after each iteration","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);CHKERRQ(ierr); 469 if (flg) { 470 ierr = PetscViewerASCIIOpen(comm,monfilename,&monviewer);CHKERRQ(ierr); 471 ierr = TaoSetMonitor(tao,TaoStepDirectionMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);CHKERRQ(ierr); 472 } 473 474 ierr = PetscOptionsString("-tao_view_residual","view least-squares residual vector after each evaluation","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);CHKERRQ(ierr); 475 if (flg) { 476 ierr = PetscViewerASCIIOpen(comm,monfilename,&monviewer);CHKERRQ(ierr); 477 ierr = TaoSetMonitor(tao,TaoResidualMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);CHKERRQ(ierr); 478 } 479 480 ierr = PetscOptionsString("-tao_monitor","Use the default convergence monitor","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);CHKERRQ(ierr); 481 if (flg) { 482 ierr = PetscViewerASCIIOpen(comm,monfilename,&monviewer);CHKERRQ(ierr); 483 ierr = TaoSetMonitor(tao,TaoMonitorDefault,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);CHKERRQ(ierr); 484 } 485 486 ierr = PetscOptionsString("-tao_gmonitor","Use the convergence monitor with extra globalization info","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);CHKERRQ(ierr); 487 if (flg) { 488 ierr = PetscViewerASCIIOpen(comm,monfilename,&monviewer);CHKERRQ(ierr); 489 ierr = TaoSetMonitor(tao,TaoDefaultGMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);CHKERRQ(ierr); 490 } 491 492 ierr = PetscOptionsString("-tao_smonitor","Use the short convergence monitor","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);CHKERRQ(ierr); 493 if (flg) { 494 ierr = PetscViewerASCIIOpen(comm,monfilename,&monviewer);CHKERRQ(ierr); 495 ierr = TaoSetMonitor(tao,TaoDefaultSMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);CHKERRQ(ierr); 496 } 497 498 ierr = PetscOptionsString("-tao_cmonitor","Use the default convergence monitor with constraint norm","TaoSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);CHKERRQ(ierr); 499 if (flg) { 500 ierr = PetscViewerASCIIOpen(comm,monfilename,&monviewer);CHKERRQ(ierr); 501 ierr = TaoSetMonitor(tao,TaoDefaultCMonitor,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);CHKERRQ(ierr); 502 } 503 504 flg = PETSC_FALSE; 505 ierr = PetscOptionsBool("-tao_cancelmonitors","cancel all monitors and call any registered destroy routines","TaoCancelMonitors",flg,&flg,NULL);CHKERRQ(ierr); 506 if (flg) {ierr = TaoCancelMonitors(tao);CHKERRQ(ierr);} 507 508 flg = PETSC_FALSE; 509 ierr = PetscOptionsBool("-tao_draw_solution","Plot solution vector at each iteration","TaoSetMonitor",flg,&flg,NULL);CHKERRQ(ierr); 510 if (flg) { 511 TaoMonitorDrawCtx drawctx; 512 PetscInt howoften = 1; 513 ierr = TaoMonitorDrawCtxCreate(PetscObjectComm((PetscObject)tao),NULL,NULL,PETSC_DECIDE,PETSC_DECIDE,300,300,howoften,&drawctx);CHKERRQ(ierr); 514 ierr = TaoSetMonitor(tao,TaoDrawSolutionMonitor,drawctx,(PetscErrorCode (*)(void**))TaoMonitorDrawCtxDestroy);CHKERRQ(ierr); 515 } 516 517 flg = PETSC_FALSE; 518 ierr = PetscOptionsBool("-tao_draw_step","plots step direction at each iteration","TaoSetMonitor",flg,&flg,NULL);CHKERRQ(ierr); 519 if (flg) { 520 ierr = TaoSetMonitor(tao,TaoDrawStepMonitor,NULL,NULL);CHKERRQ(ierr); 521 } 522 523 flg = PETSC_FALSE; 524 ierr = PetscOptionsBool("-tao_draw_gradient","plots gradient at each iteration","TaoSetMonitor",flg,&flg,NULL);CHKERRQ(ierr); 525 if (flg) { 526 TaoMonitorDrawCtx drawctx; 527 PetscInt howoften = 1; 528 ierr = TaoMonitorDrawCtxCreate(PetscObjectComm((PetscObject)tao),NULL,NULL,PETSC_DECIDE,PETSC_DECIDE,300,300,howoften,&drawctx);CHKERRQ(ierr); 529 ierr = TaoSetMonitor(tao,TaoDrawGradientMonitor,drawctx,(PetscErrorCode (*)(void**))TaoMonitorDrawCtxDestroy);CHKERRQ(ierr); 530 } 531 flg = PETSC_FALSE; 532 ierr = PetscOptionsBool("-tao_fd_gradient","compute gradient using finite differences","TaoDefaultComputeGradient",flg,&flg,NULL);CHKERRQ(ierr); 533 if (flg) { 534 ierr = TaoSetGradientRoutine(tao,TaoDefaultComputeGradient,NULL);CHKERRQ(ierr); 535 } 536 flg = PETSC_FALSE; 537 ierr = PetscOptionsBool("-tao_fd_hessian","compute hessian using finite differences","TaoDefaultComputeHessian",flg,&flg,NULL);CHKERRQ(ierr); 538 if (flg) { 539 Mat H; 540 541 ierr = MatCreate(PetscObjectComm((PetscObject)tao),&H);CHKERRQ(ierr); 542 ierr = MatSetType(H,MATAIJ);CHKERRQ(ierr); 543 ierr = TaoSetHessianRoutine(tao,H,H,TaoDefaultComputeHessian,NULL);CHKERRQ(ierr); 544 ierr = MatDestroy(&H);CHKERRQ(ierr); 545 } 546 flg = PETSC_FALSE; 547 ierr = PetscOptionsBool("-tao_mf_hessian","compute matrix-free hessian using finite differences","TaoDefaultComputeHessianMFFD",flg,&flg,NULL);CHKERRQ(ierr); 548 if (flg) { 549 Mat H; 550 551 ierr = MatCreate(PetscObjectComm((PetscObject)tao),&H);CHKERRQ(ierr); 552 ierr = TaoSetHessianRoutine(tao,H,H,TaoDefaultComputeHessianMFFD,NULL);CHKERRQ(ierr); 553 ierr = MatDestroy(&H);CHKERRQ(ierr); 554 } 555 flg = PETSC_FALSE; 556 ierr = PetscOptionsBool("-tao_recycle_history","enable recycling/re-using information from the previous TaoSolve() call for some algorithms","TaoSetRecycleHistory",flg,&flg,NULL);CHKERRQ(ierr); 557 if (flg) { 558 ierr = TaoSetRecycleHistory(tao, PETSC_TRUE);CHKERRQ(ierr); 559 } 560 ierr = PetscOptionsEnum("-tao_subset_type","subset type","",TaoSubSetTypes,(PetscEnum)tao->subset_type,(PetscEnum*)&tao->subset_type,NULL);CHKERRQ(ierr); 561 562 if (tao->ops->setfromoptions) { 563 ierr = (*tao->ops->setfromoptions)(PetscOptionsObject,tao);CHKERRQ(ierr); 564 } 565 } 566 ierr = PetscOptionsEnd();CHKERRQ(ierr); 567 PetscFunctionReturn(0); 568 } 569 570 /*@C 571 TaoViewFromOptions - View from Options 572 573 Collective on Tao 574 575 Input Parameters: 576 + A - the Tao context 577 . obj - Optional object 578 - name - command line option 579 580 Level: intermediate 581 .seealso: Tao, TaoView, PetscObjectViewFromOptions(), TaoCreate() 582 @*/ 583 PetscErrorCode TaoViewFromOptions(Tao A,PetscObject obj,const char name[]) 584 { 585 PetscErrorCode ierr; 586 587 PetscFunctionBegin; 588 PetscValidHeaderSpecific(A,TAO_CLASSID,1); 589 ierr = PetscObjectViewFromOptions((PetscObject)A,obj,name);CHKERRQ(ierr); 590 PetscFunctionReturn(0); 591 } 592 593 /*@C 594 TaoView - Prints information about the Tao 595 596 Collective on Tao 597 598 InputParameters: 599 + tao - the Tao context 600 - viewer - visualization context 601 602 Options Database Key: 603 . -tao_view - Calls TaoView() at the end of TaoSolve() 604 605 Notes: 606 The available visualization contexts include 607 + PETSC_VIEWER_STDOUT_SELF - standard output (default) 608 - PETSC_VIEWER_STDOUT_WORLD - synchronized standard 609 output where only the first processor opens 610 the file. All other processors send their 611 data to the first processor to print. 612 613 Level: beginner 614 615 .seealso: PetscViewerASCIIOpen() 616 @*/ 617 PetscErrorCode TaoView(Tao tao, PetscViewer viewer) 618 { 619 PetscErrorCode ierr; 620 PetscBool isascii,isstring; 621 TaoType type; 622 623 PetscFunctionBegin; 624 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 625 if (!viewer) { 626 ierr = PetscViewerASCIIGetStdout(((PetscObject)tao)->comm,&viewer);CHKERRQ(ierr); 627 } 628 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 629 PetscCheckSameComm(tao,1,viewer,2); 630 631 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);CHKERRQ(ierr); 632 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);CHKERRQ(ierr); 633 if (isascii) { 634 ierr = PetscObjectPrintClassNamePrefixType((PetscObject)tao,viewer);CHKERRQ(ierr); 635 636 if (tao->ops->view) { 637 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 638 ierr = (*tao->ops->view)(tao,viewer);CHKERRQ(ierr); 639 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 640 } 641 if (tao->linesearch) { 642 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 643 ierr = TaoLineSearchView(tao->linesearch,viewer);CHKERRQ(ierr); 644 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 645 } 646 if (tao->ksp) { 647 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 648 ierr = KSPView(tao->ksp,viewer);CHKERRQ(ierr); 649 ierr = PetscViewerASCIIPrintf(viewer,"total KSP iterations: %D\n",tao->ksp_tot_its);CHKERRQ(ierr); 650 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 651 } 652 653 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 654 655 if (tao->XL || tao->XU) { 656 ierr = PetscViewerASCIIPrintf(viewer,"Active Set subset type: %s\n",TaoSubSetTypes[tao->subset_type]);CHKERRQ(ierr); 657 } 658 659 ierr = PetscViewerASCIIPrintf(viewer,"convergence tolerances: gatol=%g,",(double)tao->gatol);CHKERRQ(ierr); 660 ierr = PetscViewerASCIIPrintf(viewer," steptol=%g,",(double)tao->steptol);CHKERRQ(ierr); 661 ierr = PetscViewerASCIIPrintf(viewer," gttol=%g\n",(double)tao->gttol);CHKERRQ(ierr); 662 ierr = PetscViewerASCIIPrintf(viewer,"Residual in Function/Gradient:=%g\n",(double)tao->residual);CHKERRQ(ierr); 663 664 if (tao->constrained) { 665 ierr = PetscViewerASCIIPrintf(viewer,"convergence tolerances:");CHKERRQ(ierr); 666 ierr = PetscViewerASCIIPrintf(viewer," catol=%g,",(double)tao->catol);CHKERRQ(ierr); 667 ierr = PetscViewerASCIIPrintf(viewer," crtol=%g\n",(double)tao->crtol);CHKERRQ(ierr); 668 ierr = PetscViewerASCIIPrintf(viewer,"Residual in Constraints:=%g\n",(double)tao->cnorm);CHKERRQ(ierr); 669 } 670 671 if (tao->trust < tao->steptol) { 672 ierr = PetscViewerASCIIPrintf(viewer,"convergence tolerances: steptol=%g\n",(double)tao->steptol);CHKERRQ(ierr); 673 ierr = PetscViewerASCIIPrintf(viewer,"Final trust region radius:=%g\n",(double)tao->trust);CHKERRQ(ierr); 674 } 675 676 if (tao->fmin>-1.e25) { 677 ierr = PetscViewerASCIIPrintf(viewer,"convergence tolerances: function minimum=%g\n",(double)tao->fmin);CHKERRQ(ierr); 678 } 679 ierr = PetscViewerASCIIPrintf(viewer,"Objective value=%g\n",(double)tao->fc);CHKERRQ(ierr); 680 681 ierr = PetscViewerASCIIPrintf(viewer,"total number of iterations=%D, ",tao->niter);CHKERRQ(ierr); 682 ierr = PetscViewerASCIIPrintf(viewer," (max: %D)\n",tao->max_it);CHKERRQ(ierr); 683 684 if (tao->nfuncs>0) { 685 ierr = PetscViewerASCIIPrintf(viewer,"total number of function evaluations=%D,",tao->nfuncs);CHKERRQ(ierr); 686 ierr = PetscViewerASCIIPrintf(viewer," max: %D\n",tao->max_funcs);CHKERRQ(ierr); 687 } 688 if (tao->ngrads>0) { 689 ierr = PetscViewerASCIIPrintf(viewer,"total number of gradient evaluations=%D,",tao->ngrads);CHKERRQ(ierr); 690 ierr = PetscViewerASCIIPrintf(viewer," max: %D\n",tao->max_funcs);CHKERRQ(ierr); 691 } 692 if (tao->nfuncgrads>0) { 693 ierr = PetscViewerASCIIPrintf(viewer,"total number of function/gradient evaluations=%D,",tao->nfuncgrads);CHKERRQ(ierr); 694 ierr = PetscViewerASCIIPrintf(viewer," (max: %D)\n",tao->max_funcs);CHKERRQ(ierr); 695 } 696 if (tao->nhess>0) { 697 ierr = PetscViewerASCIIPrintf(viewer,"total number of Hessian evaluations=%D\n",tao->nhess);CHKERRQ(ierr); 698 } 699 /* if (tao->linear_its>0) { 700 ierr = PetscViewerASCIIPrintf(viewer," total Krylov method iterations=%D\n",tao->linear_its);CHKERRQ(ierr); 701 }*/ 702 if (tao->nconstraints>0) { 703 ierr = PetscViewerASCIIPrintf(viewer,"total number of constraint function evaluations=%D\n",tao->nconstraints);CHKERRQ(ierr); 704 } 705 if (tao->njac>0) { 706 ierr = PetscViewerASCIIPrintf(viewer,"total number of Jacobian evaluations=%D\n",tao->njac);CHKERRQ(ierr); 707 } 708 709 if (tao->reason>0) { 710 ierr = PetscViewerASCIIPrintf(viewer, "Solution converged: ");CHKERRQ(ierr); 711 switch (tao->reason) { 712 case TAO_CONVERGED_GATOL: 713 ierr = PetscViewerASCIIPrintf(viewer," ||g(X)|| <= gatol\n");CHKERRQ(ierr); 714 break; 715 case TAO_CONVERGED_GRTOL: 716 ierr = PetscViewerASCIIPrintf(viewer," ||g(X)||/|f(X)| <= grtol\n");CHKERRQ(ierr); 717 break; 718 case TAO_CONVERGED_GTTOL: 719 ierr = PetscViewerASCIIPrintf(viewer," ||g(X)||/||g(X0)|| <= gttol\n");CHKERRQ(ierr); 720 break; 721 case TAO_CONVERGED_STEPTOL: 722 ierr = PetscViewerASCIIPrintf(viewer," Steptol -- step size small\n");CHKERRQ(ierr); 723 break; 724 case TAO_CONVERGED_MINF: 725 ierr = PetscViewerASCIIPrintf(viewer," Minf -- f < fmin\n");CHKERRQ(ierr); 726 break; 727 case TAO_CONVERGED_USER: 728 ierr = PetscViewerASCIIPrintf(viewer," User Terminated\n");CHKERRQ(ierr); 729 break; 730 default: 731 ierr = PetscViewerASCIIPrintf(viewer,"\n");CHKERRQ(ierr); 732 break; 733 } 734 735 } else { 736 ierr = PetscViewerASCIIPrintf(viewer,"Solver terminated: %d",tao->reason);CHKERRQ(ierr); 737 switch (tao->reason) { 738 case TAO_DIVERGED_MAXITS: 739 ierr = PetscViewerASCIIPrintf(viewer," Maximum Iterations\n");CHKERRQ(ierr); 740 break; 741 case TAO_DIVERGED_NAN: 742 ierr = PetscViewerASCIIPrintf(viewer," NAN or Inf encountered\n");CHKERRQ(ierr); 743 break; 744 case TAO_DIVERGED_MAXFCN: 745 ierr = PetscViewerASCIIPrintf(viewer," Maximum Function Evaluations\n");CHKERRQ(ierr); 746 break; 747 case TAO_DIVERGED_LS_FAILURE: 748 ierr = PetscViewerASCIIPrintf(viewer," Line Search Failure\n");CHKERRQ(ierr); 749 break; 750 case TAO_DIVERGED_TR_REDUCTION: 751 ierr = PetscViewerASCIIPrintf(viewer," Trust Region too small\n");CHKERRQ(ierr); 752 break; 753 case TAO_DIVERGED_USER: 754 ierr = PetscViewerASCIIPrintf(viewer," User Terminated\n");CHKERRQ(ierr); 755 break; 756 default: 757 ierr = PetscViewerASCIIPrintf(viewer,"\n");CHKERRQ(ierr); 758 break; 759 } 760 } 761 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 762 } else if (isstring) { 763 ierr = TaoGetType(tao,&type);CHKERRQ(ierr); 764 ierr = PetscViewerStringSPrintf(viewer," %-3.3s",type);CHKERRQ(ierr); 765 } 766 PetscFunctionReturn(0); 767 } 768 769 /*@ 770 TaoSetRecycleHistory - Sets the boolean flag to enable/disable re-using 771 iterate information from the previous TaoSolve(). This feature is disabled by 772 default. 773 774 For conjugate gradient methods (BNCG), this re-uses the latest search direction 775 from the previous TaoSolve() call when computing the first search direction in a 776 new solution. By default, CG methods set the first search direction to the 777 negative gradient. 778 779 For quasi-Newton family of methods (BQNLS, BQNKLS, BQNKTR, BQNKTL), this re-uses 780 the accumulated quasi-Newton Hessian approximation from the previous TaoSolve() 781 call. By default, QN family of methods reset the initial Hessian approximation to 782 the identity matrix. 783 784 For any other algorithm, this setting has no effect. 785 786 Logically collective on Tao 787 788 Input Parameters: 789 + tao - the Tao context 790 - recycle - boolean flag 791 792 Options Database Keys: 793 . -tao_recycle_history 794 795 Level: intermediate 796 797 .seealso: TaoSetRecycleHistory(), TAOBNCG, TAOBQNLS, TAOBQNKLS, TAOBQNKTR, TAOBQNKTL 798 799 @*/ 800 PetscErrorCode TaoSetRecycleHistory(Tao tao, PetscBool recycle) 801 { 802 PetscFunctionBegin; 803 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 804 tao->recycle = recycle; 805 PetscFunctionReturn(0); 806 } 807 808 /*@ 809 TaoGetRecycleHistory - Retrieve the boolean flag for re-using iterate information 810 from the previous TaoSolve(). This feature is disabled by default. 811 812 Logically collective on Tao 813 814 Input Parameters: 815 , tao - the Tao context 816 817 Output Parameters: 818 , recycle - boolean flag 819 820 Options Database Keys: 821 . -tao_recycle_history 822 823 Level: intermediate 824 825 .seealso: TaoGetRecycleHistory(), TAOBNCG, TAOBQNLS, TAOBQNKLS, TAOBQNKTR, TAOBQNKTL 826 827 @*/ 828 PetscErrorCode TaoGetRecycleHistory(Tao tao, PetscBool *recycle) 829 { 830 PetscFunctionBegin; 831 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 832 *recycle = tao->recycle; 833 PetscFunctionReturn(0); 834 } 835 836 /*@ 837 TaoSetTolerances - Sets parameters used in TAO convergence tests 838 839 Logically collective on Tao 840 841 Input Parameters: 842 + tao - the Tao context 843 . gatol - stop if norm of gradient is less than this 844 . grtol - stop if relative norm of gradient is less than this 845 - gttol - stop if norm of gradient is reduced by this factor 846 847 Options Database Keys: 848 + -tao_gatol <gatol> - Sets gatol 849 . -tao_grtol <grtol> - Sets grtol 850 - -tao_gttol <gttol> - Sets gttol 851 852 Stopping Criteria: 853 $ ||g(X)|| <= gatol 854 $ ||g(X)|| / |f(X)| <= grtol 855 $ ||g(X)|| / ||g(X0)|| <= gttol 856 857 Notes: 858 Use PETSC_DEFAULT to leave one or more tolerances unchanged. 859 860 Level: beginner 861 862 .seealso: TaoGetTolerances() 863 864 @*/ 865 PetscErrorCode TaoSetTolerances(Tao tao, PetscReal gatol, PetscReal grtol, PetscReal gttol) 866 { 867 PetscErrorCode ierr; 868 869 PetscFunctionBegin; 870 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 871 872 if (gatol != PETSC_DEFAULT) { 873 if (gatol<0) { 874 ierr = PetscInfo(tao,"Tried to set negative gatol -- ignored.\n");CHKERRQ(ierr); 875 } else { 876 tao->gatol = PetscMax(0,gatol); 877 tao->gatol_changed=PETSC_TRUE; 878 } 879 } 880 881 if (grtol != PETSC_DEFAULT) { 882 if (grtol<0) { 883 ierr = PetscInfo(tao,"Tried to set negative grtol -- ignored.\n");CHKERRQ(ierr); 884 } else { 885 tao->grtol = PetscMax(0,grtol); 886 tao->grtol_changed=PETSC_TRUE; 887 } 888 } 889 890 if (gttol != PETSC_DEFAULT) { 891 if (gttol<0) { 892 ierr = PetscInfo(tao,"Tried to set negative gttol -- ignored.\n");CHKERRQ(ierr); 893 } else { 894 tao->gttol = PetscMax(0,gttol); 895 tao->gttol_changed=PETSC_TRUE; 896 } 897 } 898 PetscFunctionReturn(0); 899 } 900 901 /*@ 902 TaoSetConstraintTolerances - Sets constraint tolerance parameters used in TAO convergence tests 903 904 Logically collective on Tao 905 906 Input Parameters: 907 + tao - the Tao context 908 . catol - absolute constraint tolerance, constraint norm must be less than catol for used for gatol convergence criteria 909 - crtol - relative contraint tolerance, constraint norm must be less than crtol for used for gatol, gttol convergence criteria 910 911 Options Database Keys: 912 + -tao_catol <catol> - Sets catol 913 - -tao_crtol <crtol> - Sets crtol 914 915 Notes: 916 Use PETSC_DEFAULT to leave any tolerance unchanged. 917 918 Level: intermediate 919 920 .seealso: TaoGetTolerances(), TaoGetConstraintTolerances(), TaoSetTolerances() 921 922 @*/ 923 PetscErrorCode TaoSetConstraintTolerances(Tao tao, PetscReal catol, PetscReal crtol) 924 { 925 PetscErrorCode ierr; 926 927 PetscFunctionBegin; 928 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 929 930 if (catol != PETSC_DEFAULT) { 931 if (catol<0) { 932 ierr = PetscInfo(tao,"Tried to set negative catol -- ignored.\n");CHKERRQ(ierr); 933 } else { 934 tao->catol = PetscMax(0,catol); 935 tao->catol_changed=PETSC_TRUE; 936 } 937 } 938 939 if (crtol != PETSC_DEFAULT) { 940 if (crtol<0) { 941 ierr = PetscInfo(tao,"Tried to set negative crtol -- ignored.\n");CHKERRQ(ierr); 942 } else { 943 tao->crtol = PetscMax(0,crtol); 944 tao->crtol_changed=PETSC_TRUE; 945 } 946 } 947 PetscFunctionReturn(0); 948 } 949 950 /*@ 951 TaoGetConstraintTolerances - Gets constraint tolerance parameters used in TAO convergence tests 952 953 Not ollective 954 955 Input Parameter: 956 . tao - the Tao context 957 958 Output Parameters: 959 + catol - absolute constraint tolerance, constraint norm must be less than catol for used for gatol convergence criteria 960 - crtol - relative contraint tolerance, constraint norm must be less than crtol for used for gatol, gttol convergence criteria 961 962 Level: intermediate 963 964 .seealso: TaoGetTolerances(), TaoSetTolerances(), TaoSetConstraintTolerances() 965 966 @*/ 967 PetscErrorCode TaoGetConstraintTolerances(Tao tao, PetscReal *catol, PetscReal *crtol) 968 { 969 PetscFunctionBegin; 970 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 971 if (catol) *catol = tao->catol; 972 if (crtol) *crtol = tao->crtol; 973 PetscFunctionReturn(0); 974 } 975 976 /*@ 977 TaoSetFunctionLowerBound - Sets a bound on the solution objective value. 978 When an approximate solution with an objective value below this number 979 has been found, the solver will terminate. 980 981 Logically Collective on Tao 982 983 Input Parameters: 984 + tao - the Tao solver context 985 - fmin - the tolerance 986 987 Options Database Keys: 988 . -tao_fmin <fmin> - sets the minimum function value 989 990 Level: intermediate 991 992 .seealso: TaoSetTolerances() 993 @*/ 994 PetscErrorCode TaoSetFunctionLowerBound(Tao tao,PetscReal fmin) 995 { 996 PetscFunctionBegin; 997 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 998 tao->fmin = fmin; 999 tao->fmin_changed=PETSC_TRUE; 1000 PetscFunctionReturn(0); 1001 } 1002 1003 /*@ 1004 TaoGetFunctionLowerBound - Gets the bound on the solution objective value. 1005 When an approximate solution with an objective value below this number 1006 has been found, the solver will terminate. 1007 1008 Not collective on Tao 1009 1010 Input Parameters: 1011 . tao - the Tao solver context 1012 1013 OutputParameters: 1014 . fmin - the minimum function value 1015 1016 Level: intermediate 1017 1018 .seealso: TaoSetFunctionLowerBound() 1019 @*/ 1020 PetscErrorCode TaoGetFunctionLowerBound(Tao tao,PetscReal *fmin) 1021 { 1022 PetscFunctionBegin; 1023 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1024 *fmin = tao->fmin; 1025 PetscFunctionReturn(0); 1026 } 1027 1028 /*@ 1029 TaoSetMaximumFunctionEvaluations - Sets a maximum number of 1030 function evaluations. 1031 1032 Logically Collective on Tao 1033 1034 Input Parameters: 1035 + tao - the Tao solver context 1036 - nfcn - the maximum number of function evaluations (>=0) 1037 1038 Options Database Keys: 1039 . -tao_max_funcs <nfcn> - sets the maximum number of function evaluations 1040 1041 Level: intermediate 1042 1043 .seealso: TaoSetTolerances(), TaoSetMaximumIterations() 1044 @*/ 1045 1046 PetscErrorCode TaoSetMaximumFunctionEvaluations(Tao tao,PetscInt nfcn) 1047 { 1048 PetscFunctionBegin; 1049 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1050 tao->max_funcs = PetscMax(0,nfcn); 1051 tao->max_funcs_changed=PETSC_TRUE; 1052 PetscFunctionReturn(0); 1053 } 1054 1055 /*@ 1056 TaoGetMaximumFunctionEvaluations - Sets a maximum number of 1057 function evaluations. 1058 1059 Not Collective 1060 1061 Input Parameters: 1062 . tao - the Tao solver context 1063 1064 Output Parameters: 1065 . nfcn - the maximum number of function evaluations 1066 1067 Level: intermediate 1068 1069 .seealso: TaoSetMaximumFunctionEvaluations(), TaoGetMaximumIterations() 1070 @*/ 1071 1072 PetscErrorCode TaoGetMaximumFunctionEvaluations(Tao tao,PetscInt *nfcn) 1073 { 1074 PetscFunctionBegin; 1075 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1076 *nfcn = tao->max_funcs; 1077 PetscFunctionReturn(0); 1078 } 1079 1080 /*@ 1081 TaoGetCurrentFunctionEvaluations - Get current number of 1082 function evaluations. 1083 1084 Not Collective 1085 1086 Input Parameters: 1087 . tao - the Tao solver context 1088 1089 Output Parameters: 1090 . nfuncs - the current number of function evaluations 1091 1092 Level: intermediate 1093 1094 .seealso: TaoSetMaximumFunctionEvaluations(), TaoGetMaximumFunctionEvaluations(), TaoGetMaximumIterations() 1095 @*/ 1096 1097 PetscErrorCode TaoGetCurrentFunctionEvaluations(Tao tao,PetscInt *nfuncs) 1098 { 1099 PetscFunctionBegin; 1100 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1101 *nfuncs=PetscMax(tao->nfuncs,tao->nfuncgrads); 1102 PetscFunctionReturn(0); 1103 } 1104 1105 /*@ 1106 TaoSetMaximumIterations - Sets a maximum number of iterates. 1107 1108 Logically Collective on Tao 1109 1110 Input Parameters: 1111 + tao - the Tao solver context 1112 - maxits - the maximum number of iterates (>=0) 1113 1114 Options Database Keys: 1115 . -tao_max_it <its> - sets the maximum number of iterations 1116 1117 Level: intermediate 1118 1119 .seealso: TaoSetTolerances(), TaoSetMaximumFunctionEvaluations() 1120 @*/ 1121 PetscErrorCode TaoSetMaximumIterations(Tao tao,PetscInt maxits) 1122 { 1123 PetscFunctionBegin; 1124 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1125 tao->max_it = PetscMax(0,maxits); 1126 tao->max_it_changed=PETSC_TRUE; 1127 PetscFunctionReturn(0); 1128 } 1129 1130 /*@ 1131 TaoGetMaximumIterations - Sets a maximum number of iterates. 1132 1133 Not Collective 1134 1135 Input Parameters: 1136 . tao - the Tao solver context 1137 1138 Output Parameters: 1139 . maxits - the maximum number of iterates 1140 1141 Level: intermediate 1142 1143 .seealso: TaoSetMaximumIterations(), TaoGetMaximumFunctionEvaluations() 1144 @*/ 1145 PetscErrorCode TaoGetMaximumIterations(Tao tao,PetscInt *maxits) 1146 { 1147 PetscFunctionBegin; 1148 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1149 *maxits = tao->max_it; 1150 PetscFunctionReturn(0); 1151 } 1152 1153 /*@ 1154 TaoSetInitialTrustRegionRadius - Sets the initial trust region radius. 1155 1156 Logically collective on Tao 1157 1158 Input Parameters: 1159 + tao - a TAO optimization solver 1160 - radius - the trust region radius 1161 1162 Level: intermediate 1163 1164 Options Database Key: 1165 . -tao_trust0 <t0> - sets initial trust region radius 1166 1167 .seealso: TaoGetTrustRegionRadius(), TaoSetTrustRegionTolerance() 1168 @*/ 1169 PetscErrorCode TaoSetInitialTrustRegionRadius(Tao tao, PetscReal radius) 1170 { 1171 PetscFunctionBegin; 1172 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1173 tao->trust0 = PetscMax(0.0,radius); 1174 tao->trust0_changed=PETSC_TRUE; 1175 PetscFunctionReturn(0); 1176 } 1177 1178 /*@ 1179 TaoGetInitialTrustRegionRadius - Sets the initial trust region radius. 1180 1181 Not Collective 1182 1183 Input Parameter: 1184 . tao - a TAO optimization solver 1185 1186 Output Parameter: 1187 . radius - the trust region radius 1188 1189 Level: intermediate 1190 1191 .seealso: TaoSetInitialTrustRegionRadius(), TaoGetCurrentTrustRegionRadius() 1192 @*/ 1193 PetscErrorCode TaoGetInitialTrustRegionRadius(Tao tao, PetscReal *radius) 1194 { 1195 PetscFunctionBegin; 1196 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1197 *radius = tao->trust0; 1198 PetscFunctionReturn(0); 1199 } 1200 1201 /*@ 1202 TaoGetCurrentTrustRegionRadius - Gets the current trust region radius. 1203 1204 Not Collective 1205 1206 Input Parameter: 1207 . tao - a TAO optimization solver 1208 1209 Output Parameter: 1210 . radius - the trust region radius 1211 1212 Level: intermediate 1213 1214 .seealso: TaoSetInitialTrustRegionRadius(), TaoGetInitialTrustRegionRadius() 1215 @*/ 1216 PetscErrorCode TaoGetCurrentTrustRegionRadius(Tao tao, PetscReal *radius) 1217 { 1218 PetscFunctionBegin; 1219 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1220 *radius = tao->trust; 1221 PetscFunctionReturn(0); 1222 } 1223 1224 /*@ 1225 TaoGetTolerances - gets the current values of tolerances 1226 1227 Not Collective 1228 1229 Input Parameter: 1230 . tao - the Tao context 1231 1232 Output Parameters: 1233 + gatol - stop if norm of gradient is less than this 1234 . grtol - stop if relative norm of gradient is less than this 1235 - gttol - stop if norm of gradient is reduced by a this factor 1236 1237 Note: NULL can be used as an argument if not all tolerances values are needed 1238 1239 .seealso TaoSetTolerances() 1240 1241 Level: intermediate 1242 @*/ 1243 PetscErrorCode TaoGetTolerances(Tao tao, PetscReal *gatol, PetscReal *grtol, PetscReal *gttol) 1244 { 1245 PetscFunctionBegin; 1246 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1247 if (gatol) *gatol=tao->gatol; 1248 if (grtol) *grtol=tao->grtol; 1249 if (gttol) *gttol=tao->gttol; 1250 PetscFunctionReturn(0); 1251 } 1252 1253 /*@ 1254 TaoGetKSP - Gets the linear solver used by the optimization solver. 1255 Application writers should use TaoGetKSP if they need direct access 1256 to the PETSc KSP object. 1257 1258 Not Collective 1259 1260 Input Parameters: 1261 . tao - the TAO solver 1262 1263 Output Parameters: 1264 . ksp - the KSP linear solver used in the optimization solver 1265 1266 Level: intermediate 1267 1268 @*/ 1269 PetscErrorCode TaoGetKSP(Tao tao, KSP *ksp) 1270 { 1271 PetscFunctionBegin; 1272 *ksp = tao->ksp; 1273 PetscFunctionReturn(0); 1274 } 1275 1276 /*@ 1277 TaoGetLinearSolveIterations - Gets the total number of linear iterations 1278 used by the TAO solver 1279 1280 Not Collective 1281 1282 Input Parameter: 1283 . tao - TAO context 1284 1285 Output Parameter: 1286 . lits - number of linear iterations 1287 1288 Notes: 1289 This counter is reset to zero for each successive call to TaoSolve() 1290 1291 Level: intermediate 1292 1293 .seealso: TaoGetKSP() 1294 @*/ 1295 PetscErrorCode TaoGetLinearSolveIterations(Tao tao,PetscInt *lits) 1296 { 1297 PetscFunctionBegin; 1298 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1299 PetscValidIntPointer(lits,2); 1300 *lits = tao->ksp_tot_its; 1301 PetscFunctionReturn(0); 1302 } 1303 1304 /*@ 1305 TaoGetLineSearch - Gets the line search used by the optimization solver. 1306 Application writers should use TaoGetLineSearch if they need direct access 1307 to the TaoLineSearch object. 1308 1309 Not Collective 1310 1311 Input Parameters: 1312 . tao - the TAO solver 1313 1314 Output Parameters: 1315 . ls - the line search used in the optimization solver 1316 1317 Level: intermediate 1318 1319 @*/ 1320 PetscErrorCode TaoGetLineSearch(Tao tao, TaoLineSearch *ls) 1321 { 1322 PetscFunctionBegin; 1323 *ls = tao->linesearch; 1324 PetscFunctionReturn(0); 1325 } 1326 1327 /*@ 1328 TaoAddLineSearchCounts - Adds the number of function evaluations spent 1329 in the line search to the running total. 1330 1331 Input Parameters: 1332 + tao - the TAO solver 1333 - ls - the line search used in the optimization solver 1334 1335 Level: developer 1336 1337 .seealso: TaoLineSearchApply() 1338 @*/ 1339 PetscErrorCode TaoAddLineSearchCounts(Tao tao) 1340 { 1341 PetscErrorCode ierr; 1342 PetscBool flg; 1343 PetscInt nfeval,ngeval,nfgeval; 1344 1345 PetscFunctionBegin; 1346 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1347 if (tao->linesearch) { 1348 ierr = TaoLineSearchIsUsingTaoRoutines(tao->linesearch,&flg);CHKERRQ(ierr); 1349 if (!flg) { 1350 ierr = TaoLineSearchGetNumberFunctionEvaluations(tao->linesearch,&nfeval,&ngeval,&nfgeval);CHKERRQ(ierr); 1351 tao->nfuncs+=nfeval; 1352 tao->ngrads+=ngeval; 1353 tao->nfuncgrads+=nfgeval; 1354 } 1355 } 1356 PetscFunctionReturn(0); 1357 } 1358 1359 /*@ 1360 TaoGetSolutionVector - Returns the vector with the current TAO solution 1361 1362 Not Collective 1363 1364 Input Parameter: 1365 . tao - the Tao context 1366 1367 Output Parameter: 1368 . X - the current solution 1369 1370 Level: intermediate 1371 1372 Note: The returned vector will be the same object that was passed into TaoSetInitialVector() 1373 @*/ 1374 PetscErrorCode TaoGetSolutionVector(Tao tao, Vec *X) 1375 { 1376 PetscFunctionBegin; 1377 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1378 *X = tao->solution; 1379 PetscFunctionReturn(0); 1380 } 1381 1382 /*@ 1383 TaoGetGradientVector - Returns the vector with the current TAO gradient 1384 1385 Not Collective 1386 1387 Input Parameter: 1388 . tao - the Tao context 1389 1390 Output Parameter: 1391 . G - the current solution 1392 1393 Level: intermediate 1394 @*/ 1395 PetscErrorCode TaoGetGradientVector(Tao tao, Vec *G) 1396 { 1397 PetscFunctionBegin; 1398 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1399 *G = tao->gradient; 1400 PetscFunctionReturn(0); 1401 } 1402 1403 /*@ 1404 TaoResetStatistics - Initialize the statistics used by TAO for all of the solvers. 1405 These statistics include the iteration number, residual norms, and convergence status. 1406 This routine gets called before solving each optimization problem. 1407 1408 Collective on Tao 1409 1410 Input Parameters: 1411 . solver - the Tao context 1412 1413 Level: developer 1414 1415 .seealso: TaoCreate(), TaoSolve() 1416 @*/ 1417 PetscErrorCode TaoResetStatistics(Tao tao) 1418 { 1419 PetscFunctionBegin; 1420 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1421 tao->niter = 0; 1422 tao->nfuncs = 0; 1423 tao->nfuncgrads = 0; 1424 tao->ngrads = 0; 1425 tao->nhess = 0; 1426 tao->njac = 0; 1427 tao->nconstraints = 0; 1428 tao->ksp_its = 0; 1429 tao->ksp_tot_its = 0; 1430 tao->reason = TAO_CONTINUE_ITERATING; 1431 tao->residual = 0.0; 1432 tao->cnorm = 0.0; 1433 tao->step = 0.0; 1434 tao->lsflag = PETSC_FALSE; 1435 if (tao->hist_reset) tao->hist_len=0; 1436 PetscFunctionReturn(0); 1437 } 1438 1439 /*@C 1440 TaoSetUpdate - Sets the general-purpose update function called 1441 at the beginning of every iteration of the nonlinear solve. Specifically 1442 it is called at the top of every iteration, after the new solution and the gradient 1443 is determined, but before the Hessian is computed (if applicable). 1444 1445 Logically Collective on Tao 1446 1447 Input Parameters: 1448 + tao - The tao solver context 1449 - func - The function 1450 1451 Calling sequence of func: 1452 $ func (Tao tao, PetscInt step); 1453 1454 . step - The current step of the iteration 1455 1456 Level: advanced 1457 1458 .seealso TaoSolve() 1459 @*/ 1460 PetscErrorCode TaoSetUpdate(Tao tao, PetscErrorCode (*func)(Tao, PetscInt,void*), void *ctx) 1461 { 1462 PetscFunctionBegin; 1463 PetscValidHeaderSpecific(tao, TAO_CLASSID,1); 1464 tao->ops->update = func; 1465 tao->user_update = ctx; 1466 PetscFunctionReturn(0); 1467 } 1468 1469 /*@C 1470 TaoSetConvergenceTest - Sets the function that is to be used to test 1471 for convergence o fthe iterative minimization solution. The new convergence 1472 testing routine will replace TAO's default convergence test. 1473 1474 Logically Collective on Tao 1475 1476 Input Parameters: 1477 + tao - the Tao object 1478 . conv - the routine to test for convergence 1479 - ctx - [optional] context for private data for the convergence routine 1480 (may be NULL) 1481 1482 Calling sequence of conv: 1483 $ PetscErrorCode conv(Tao tao, void *ctx) 1484 1485 + tao - the Tao object 1486 - ctx - [optional] convergence context 1487 1488 Note: The new convergence testing routine should call TaoSetConvergedReason(). 1489 1490 Level: advanced 1491 1492 .seealso: TaoSetConvergedReason(), TaoGetSolutionStatus(), TaoGetTolerances(), TaoSetMonitor 1493 1494 @*/ 1495 PetscErrorCode TaoSetConvergenceTest(Tao tao, PetscErrorCode (*conv)(Tao,void*), void *ctx) 1496 { 1497 PetscFunctionBegin; 1498 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1499 (tao)->ops->convergencetest = conv; 1500 (tao)->cnvP = ctx; 1501 PetscFunctionReturn(0); 1502 } 1503 1504 /*@C 1505 TaoSetMonitor - Sets an ADDITIONAL function that is to be used at every 1506 iteration of the solver to display the iteration's 1507 progress. 1508 1509 Logically Collective on Tao 1510 1511 Input Parameters: 1512 + tao - the Tao solver context 1513 . mymonitor - monitoring routine 1514 - mctx - [optional] user-defined context for private data for the 1515 monitor routine (may be NULL) 1516 1517 Calling sequence of mymonitor: 1518 $ PetscErrorCode mymonitor(Tao tao,void *mctx) 1519 1520 + tao - the Tao solver context 1521 - mctx - [optional] monitoring context 1522 1523 Options Database Keys: 1524 + -tao_monitor - sets TaoMonitorDefault() 1525 . -tao_smonitor - sets short monitor 1526 . -tao_cmonitor - same as smonitor plus constraint norm 1527 . -tao_view_solution - view solution at each iteration 1528 . -tao_view_gradient - view gradient at each iteration 1529 . -tao_view_ls_residual - view least-squares residual vector at each iteration 1530 - -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. 1531 1532 Notes: 1533 Several different monitoring routines may be set by calling 1534 TaoSetMonitor() multiple times; all will be called in the 1535 order in which they were set. 1536 1537 Fortran Notes: 1538 Only one monitor function may be set 1539 1540 Level: intermediate 1541 1542 .seealso: TaoMonitorDefault(), TaoCancelMonitors(), TaoSetDestroyRoutine() 1543 @*/ 1544 PetscErrorCode TaoSetMonitor(Tao tao, PetscErrorCode (*func)(Tao, void*), void *ctx,PetscErrorCode (*dest)(void**)) 1545 { 1546 PetscErrorCode ierr; 1547 PetscInt i; 1548 PetscBool identical; 1549 1550 PetscFunctionBegin; 1551 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1552 if (tao->numbermonitors >= MAXTAOMONITORS) SETERRQ1(PetscObjectComm((PetscObject)tao),PETSC_ERR_SUP,"Cannot attach another monitor -- max=%d",MAXTAOMONITORS); 1553 1554 for (i=0; i<tao->numbermonitors;i++) { 1555 ierr = PetscMonitorCompare((PetscErrorCode (*)(void))func,ctx,dest,(PetscErrorCode (*)(void))tao->monitor[i],tao->monitorcontext[i],tao->monitordestroy[i],&identical);CHKERRQ(ierr); 1556 if (identical) PetscFunctionReturn(0); 1557 } 1558 tao->monitor[tao->numbermonitors] = func; 1559 tao->monitorcontext[tao->numbermonitors] = (void*)ctx; 1560 tao->monitordestroy[tao->numbermonitors] = dest; 1561 ++tao->numbermonitors; 1562 PetscFunctionReturn(0); 1563 } 1564 1565 /*@ 1566 TaoCancelMonitors - Clears all the monitor functions for a Tao object. 1567 1568 Logically Collective on Tao 1569 1570 Input Parameters: 1571 . tao - the Tao solver context 1572 1573 Options Database: 1574 . -tao_cancelmonitors - cancels all monitors that have been hardwired 1575 into a code by calls to TaoSetMonitor(), but does not cancel those 1576 set via the options database 1577 1578 Notes: 1579 There is no way to clear one specific monitor from a Tao object. 1580 1581 Level: advanced 1582 1583 .seealso: TaoMonitorDefault(), TaoSetMonitor() 1584 @*/ 1585 PetscErrorCode TaoCancelMonitors(Tao tao) 1586 { 1587 PetscInt i; 1588 PetscErrorCode ierr; 1589 1590 PetscFunctionBegin; 1591 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 1592 for (i=0;i<tao->numbermonitors;i++) { 1593 if (tao->monitordestroy[i]) { 1594 ierr = (*tao->monitordestroy[i])(&tao->monitorcontext[i]);CHKERRQ(ierr); 1595 } 1596 } 1597 tao->numbermonitors=0; 1598 PetscFunctionReturn(0); 1599 } 1600 1601 /*@ 1602 TaoMonitorDefault - Default routine for monitoring progress of the 1603 Tao solvers (default). This monitor prints the function value and gradient 1604 norm at each iteration. It can be turned on from the command line using the 1605 -tao_monitor option 1606 1607 Collective on Tao 1608 1609 Input Parameters: 1610 + tao - the Tao context 1611 - ctx - PetscViewer context or NULL 1612 1613 Options Database Keys: 1614 . -tao_monitor 1615 1616 Level: advanced 1617 1618 .seealso: TaoDefaultSMonitor(), TaoSetMonitor() 1619 @*/ 1620 PetscErrorCode TaoMonitorDefault(Tao tao, void *ctx) 1621 { 1622 PetscErrorCode ierr; 1623 PetscInt its, tabs; 1624 PetscReal fct,gnorm; 1625 PetscViewer viewer = (PetscViewer)ctx; 1626 1627 PetscFunctionBegin; 1628 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 1629 its=tao->niter; 1630 fct=tao->fc; 1631 gnorm=tao->residual; 1632 ierr = PetscViewerASCIIGetTab(viewer, &tabs);CHKERRQ(ierr); 1633 ierr = PetscViewerASCIISetTab(viewer, ((PetscObject)tao)->tablevel);CHKERRQ(ierr); 1634 if (its == 0 && ((PetscObject)tao)->prefix && !tao->header_printed) { 1635 ierr = PetscViewerASCIIPrintf(viewer," Iteration information for %s solve.\n",((PetscObject)tao)->prefix);CHKERRQ(ierr); 1636 tao->header_printed = PETSC_TRUE; 1637 } 1638 ierr = PetscViewerASCIIPrintf(viewer,"%3D TAO,",its);CHKERRQ(ierr); 1639 ierr = PetscViewerASCIIPrintf(viewer," Function value: %g,",(double)fct);CHKERRQ(ierr); 1640 if (gnorm >= PETSC_INFINITY) { 1641 ierr = PetscViewerASCIIPrintf(viewer," Residual: Inf \n");CHKERRQ(ierr); 1642 } else { 1643 ierr = PetscViewerASCIIPrintf(viewer," Residual: %g \n",(double)gnorm);CHKERRQ(ierr); 1644 } 1645 ierr = PetscViewerASCIISetTab(viewer, tabs);CHKERRQ(ierr); 1646 PetscFunctionReturn(0); 1647 } 1648 1649 /*@ 1650 TaoDefaultGMonitor - Default routine for monitoring progress of the 1651 Tao solvers (default) with extra detail on the globalization method. 1652 This monitor prints the function value and gradient norm at each 1653 iteration, as well as the step size and trust radius. Note that the 1654 step size and trust radius may be the same for some algorithms. 1655 It can be turned on from the command line using the 1656 -tao_gmonitor option 1657 1658 Collective on Tao 1659 1660 Input Parameters: 1661 + tao - the Tao context 1662 - ctx - PetscViewer context or NULL 1663 1664 Options Database Keys: 1665 . -tao_monitor 1666 1667 Level: advanced 1668 1669 .seealso: TaoDefaultSMonitor(), TaoSetMonitor() 1670 @*/ 1671 PetscErrorCode TaoDefaultGMonitor(Tao tao, void *ctx) 1672 { 1673 PetscErrorCode ierr; 1674 PetscInt its, tabs; 1675 PetscReal fct,gnorm,stp,tr; 1676 PetscViewer viewer = (PetscViewer)ctx; 1677 1678 PetscFunctionBegin; 1679 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 1680 its=tao->niter; 1681 fct=tao->fc; 1682 gnorm=tao->residual; 1683 stp=tao->step; 1684 tr=tao->trust; 1685 ierr = PetscViewerASCIIGetTab(viewer, &tabs);CHKERRQ(ierr); 1686 ierr = PetscViewerASCIISetTab(viewer, ((PetscObject)tao)->tablevel);CHKERRQ(ierr); 1687 if (its == 0 && ((PetscObject)tao)->prefix && !tao->header_printed) { 1688 ierr = PetscViewerASCIIPrintf(viewer," Iteration information for %s solve.\n",((PetscObject)tao)->prefix);CHKERRQ(ierr); 1689 tao->header_printed = PETSC_TRUE; 1690 } 1691 ierr = PetscViewerASCIIPrintf(viewer,"%3D TAO,",its);CHKERRQ(ierr); 1692 ierr = PetscViewerASCIIPrintf(viewer," Function value: %g,",(double)fct);CHKERRQ(ierr); 1693 if (gnorm >= PETSC_INFINITY) { 1694 ierr = PetscViewerASCIIPrintf(viewer," Residual: Inf,");CHKERRQ(ierr); 1695 } else { 1696 ierr = PetscViewerASCIIPrintf(viewer," Residual: %g,",(double)gnorm);CHKERRQ(ierr); 1697 } 1698 ierr = PetscViewerASCIIPrintf(viewer," Step: %g, Trust: %g\n",(double)stp,(double)tr);CHKERRQ(ierr); 1699 ierr = PetscViewerASCIISetTab(viewer, tabs);CHKERRQ(ierr); 1700 PetscFunctionReturn(0); 1701 } 1702 1703 /*@ 1704 TaoDefaultSMonitor - Default routine for monitoring progress of the 1705 solver. Same as TaoMonitorDefault() except 1706 it prints fewer digits of the residual as the residual gets smaller. 1707 This is because the later digits are meaningless and are often 1708 different on different machines; by using this routine different 1709 machines will usually generate the same output. It can be turned on 1710 by using the -tao_smonitor option 1711 1712 Collective on Tao 1713 1714 Input Parameters: 1715 + tao - the Tao context 1716 - ctx - PetscViewer context of type ASCII 1717 1718 Options Database Keys: 1719 . -tao_smonitor 1720 1721 Level: advanced 1722 1723 .seealso: TaoMonitorDefault(), TaoSetMonitor() 1724 @*/ 1725 PetscErrorCode TaoDefaultSMonitor(Tao tao, void *ctx) 1726 { 1727 PetscErrorCode ierr; 1728 PetscInt its, tabs; 1729 PetscReal fct,gnorm; 1730 PetscViewer viewer = (PetscViewer)ctx; 1731 1732 PetscFunctionBegin; 1733 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 1734 its=tao->niter; 1735 fct=tao->fc; 1736 gnorm=tao->residual; 1737 ierr = PetscViewerASCIIGetTab(viewer, &tabs);CHKERRQ(ierr); 1738 ierr = PetscViewerASCIISetTab(viewer, ((PetscObject)tao)->tablevel);CHKERRQ(ierr); 1739 ierr = PetscViewerASCIIPrintf(viewer,"iter = %3D,",its);CHKERRQ(ierr); 1740 ierr = PetscViewerASCIIPrintf(viewer," Function value %g,",(double)fct);CHKERRQ(ierr); 1741 if (gnorm >= PETSC_INFINITY) { 1742 ierr = PetscViewerASCIIPrintf(viewer," Residual: Inf \n");CHKERRQ(ierr); 1743 } else if (gnorm > 1.e-6) { 1744 ierr = PetscViewerASCIIPrintf(viewer," Residual: %g \n",(double)gnorm);CHKERRQ(ierr); 1745 } else if (gnorm > 1.e-11) { 1746 ierr = PetscViewerASCIIPrintf(viewer," Residual: < 1.0e-6 \n");CHKERRQ(ierr); 1747 } else { 1748 ierr = PetscViewerASCIIPrintf(viewer," Residual: < 1.0e-11 \n");CHKERRQ(ierr); 1749 } 1750 ierr = PetscViewerASCIISetTab(viewer, tabs);CHKERRQ(ierr); 1751 PetscFunctionReturn(0); 1752 } 1753 1754 /*@ 1755 TaoDefaultCMonitor - same as TaoMonitorDefault() except 1756 it prints the norm of the constraints function. It can be turned on 1757 from the command line using the -tao_cmonitor option 1758 1759 Collective on Tao 1760 1761 Input Parameters: 1762 + tao - the Tao context 1763 - ctx - PetscViewer context or NULL 1764 1765 Options Database Keys: 1766 . -tao_cmonitor 1767 1768 Level: advanced 1769 1770 .seealso: TaoMonitorDefault(), TaoSetMonitor() 1771 @*/ 1772 PetscErrorCode TaoDefaultCMonitor(Tao tao, void *ctx) 1773 { 1774 PetscErrorCode ierr; 1775 PetscInt its, tabs; 1776 PetscReal fct,gnorm; 1777 PetscViewer viewer = (PetscViewer)ctx; 1778 1779 PetscFunctionBegin; 1780 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 1781 its=tao->niter; 1782 fct=tao->fc; 1783 gnorm=tao->residual; 1784 ierr = PetscViewerASCIIGetTab(viewer, &tabs);CHKERRQ(ierr); 1785 ierr = PetscViewerASCIISetTab(viewer, ((PetscObject)tao)->tablevel);CHKERRQ(ierr); 1786 ierr = PetscViewerASCIIPrintf(viewer,"iter = %D,",its);CHKERRQ(ierr); 1787 ierr = PetscViewerASCIIPrintf(viewer," Function value: %g,",(double)fct);CHKERRQ(ierr); 1788 ierr = PetscViewerASCIIPrintf(viewer," Residual: %g ",(double)gnorm);CHKERRQ(ierr); 1789 ierr = PetscViewerASCIIPrintf(viewer," Constraint: %g \n",(double)tao->cnorm);CHKERRQ(ierr); 1790 ierr = PetscViewerASCIISetTab(viewer, tabs);CHKERRQ(ierr); 1791 PetscFunctionReturn(0); 1792 } 1793 1794 /*@C 1795 TaoSolutionMonitor - Views the solution at each iteration 1796 It can be turned on from the command line using the 1797 -tao_view_solution option 1798 1799 Collective on Tao 1800 1801 Input Parameters: 1802 + tao - the Tao context 1803 - ctx - PetscViewer context or NULL 1804 1805 Options Database Keys: 1806 . -tao_view_solution 1807 1808 Level: advanced 1809 1810 .seealso: TaoDefaultSMonitor(), TaoSetMonitor() 1811 @*/ 1812 PetscErrorCode TaoSolutionMonitor(Tao tao, void *ctx) 1813 { 1814 PetscErrorCode ierr; 1815 PetscViewer viewer = (PetscViewer)ctx; 1816 1817 PetscFunctionBegin; 1818 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 1819 ierr = VecView(tao->solution, viewer);CHKERRQ(ierr); 1820 PetscFunctionReturn(0); 1821 } 1822 1823 /*@C 1824 TaoGradientMonitor - Views the gradient at each iteration 1825 It can be turned on from the command line using the 1826 -tao_view_gradient option 1827 1828 Collective on Tao 1829 1830 Input Parameters: 1831 + tao - the Tao context 1832 - ctx - PetscViewer context or NULL 1833 1834 Options Database Keys: 1835 . -tao_view_gradient 1836 1837 Level: advanced 1838 1839 .seealso: TaoDefaultSMonitor(), TaoSetMonitor() 1840 @*/ 1841 PetscErrorCode TaoGradientMonitor(Tao tao, void *ctx) 1842 { 1843 PetscErrorCode ierr; 1844 PetscViewer viewer = (PetscViewer)ctx; 1845 1846 PetscFunctionBegin; 1847 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 1848 ierr = VecView(tao->gradient, viewer);CHKERRQ(ierr); 1849 PetscFunctionReturn(0); 1850 } 1851 1852 /*@C 1853 TaoStepDirectionMonitor - Views the gradient at each iteration 1854 It can be turned on from the command line using the 1855 -tao_view_gradient option 1856 1857 Collective on Tao 1858 1859 Input Parameters: 1860 + tao - the Tao context 1861 - ctx - PetscViewer context or NULL 1862 1863 Options Database Keys: 1864 . -tao_view_gradient 1865 1866 Level: advanced 1867 1868 .seealso: TaoDefaultSMonitor(), TaoSetMonitor() 1869 @*/ 1870 PetscErrorCode TaoStepDirectionMonitor(Tao tao, void *ctx) 1871 { 1872 PetscErrorCode ierr; 1873 PetscViewer viewer = (PetscViewer)ctx; 1874 1875 PetscFunctionBegin; 1876 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 1877 ierr = VecView(tao->stepdirection, viewer);CHKERRQ(ierr); 1878 PetscFunctionReturn(0); 1879 } 1880 1881 /*@C 1882 TaoDrawSolutionMonitor - Plots the solution at each iteration 1883 It can be turned on from the command line using the 1884 -tao_draw_solution option 1885 1886 Collective on Tao 1887 1888 Input Parameters: 1889 + tao - the Tao context 1890 - ctx - TaoMonitorDraw context 1891 1892 Options Database Keys: 1893 . -tao_draw_solution 1894 1895 Level: advanced 1896 1897 .seealso: TaoSolutionMonitor(), TaoSetMonitor(), TaoDrawGradientMonitor 1898 @*/ 1899 PetscErrorCode TaoDrawSolutionMonitor(Tao tao, void *ctx) 1900 { 1901 PetscErrorCode ierr; 1902 TaoMonitorDrawCtx ictx = (TaoMonitorDrawCtx)ctx; 1903 1904 PetscFunctionBegin; 1905 if (!(((ictx->howoften > 0) && (!(tao->niter % ictx->howoften))) || ((ictx->howoften == -1) && tao->reason))) PetscFunctionReturn(0); 1906 ierr = VecView(tao->solution,ictx->viewer);CHKERRQ(ierr); 1907 PetscFunctionReturn(0); 1908 } 1909 1910 /*@C 1911 TaoDrawGradientMonitor - Plots the gradient at each iteration 1912 It can be turned on from the command line using the 1913 -tao_draw_gradient option 1914 1915 Collective on Tao 1916 1917 Input Parameters: 1918 + tao - the Tao context 1919 - ctx - PetscViewer context 1920 1921 Options Database Keys: 1922 . -tao_draw_gradient 1923 1924 Level: advanced 1925 1926 .seealso: TaoGradientMonitor(), TaoSetMonitor(), TaoDrawSolutionMonitor 1927 @*/ 1928 PetscErrorCode TaoDrawGradientMonitor(Tao tao, void *ctx) 1929 { 1930 PetscErrorCode ierr; 1931 TaoMonitorDrawCtx ictx = (TaoMonitorDrawCtx)ctx; 1932 1933 PetscFunctionBegin; 1934 if (!(((ictx->howoften > 0) && (!(tao->niter % ictx->howoften))) || ((ictx->howoften == -1) && tao->reason))) PetscFunctionReturn(0); 1935 ierr = VecView(tao->gradient,ictx->viewer);CHKERRQ(ierr); 1936 PetscFunctionReturn(0); 1937 } 1938 1939 /*@C 1940 TaoDrawStepMonitor - Plots the step direction at each iteration 1941 It can be turned on from the command line using the 1942 -tao_draw_step option 1943 1944 Collective on Tao 1945 1946 Input Parameters: 1947 + tao - the Tao context 1948 - ctx - PetscViewer context 1949 1950 Options Database Keys: 1951 . -tao_draw_step 1952 1953 Level: advanced 1954 1955 .seealso: TaoSetMonitor(), TaoDrawSolutionMonitor 1956 @*/ 1957 PetscErrorCode TaoDrawStepMonitor(Tao tao, void *ctx) 1958 { 1959 PetscErrorCode ierr; 1960 PetscViewer viewer = (PetscViewer)(ctx); 1961 1962 PetscFunctionBegin; 1963 ierr = VecView(tao->stepdirection, viewer);CHKERRQ(ierr); 1964 PetscFunctionReturn(0); 1965 } 1966 1967 /*@C 1968 TaoResidualMonitor - Views the least-squares residual at each iteration 1969 It can be turned on from the command line using the 1970 -tao_view_ls_residual option 1971 1972 Collective on Tao 1973 1974 Input Parameters: 1975 + tao - the Tao context 1976 - ctx - PetscViewer context or NULL 1977 1978 Options Database Keys: 1979 . -tao_view_ls_residual 1980 1981 Level: advanced 1982 1983 .seealso: TaoDefaultSMonitor(), TaoSetMonitor() 1984 @*/ 1985 PetscErrorCode TaoResidualMonitor(Tao tao, void *ctx) 1986 { 1987 PetscErrorCode ierr; 1988 PetscViewer viewer = (PetscViewer)ctx; 1989 1990 PetscFunctionBegin; 1991 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 1992 ierr = VecView(tao->ls_res,viewer);CHKERRQ(ierr); 1993 PetscFunctionReturn(0); 1994 } 1995 1996 /*@ 1997 TaoDefaultConvergenceTest - Determines whether the solver should continue iterating 1998 or terminate. 1999 2000 Collective on Tao 2001 2002 Input Parameters: 2003 + tao - the Tao context 2004 - dummy - unused dummy context 2005 2006 Output Parameter: 2007 . reason - for terminating 2008 2009 Notes: 2010 This routine checks the residual in the optimality conditions, the 2011 relative residual in the optimity conditions, the number of function 2012 evaluations, and the function value to test convergence. Some 2013 solvers may use different convergence routines. 2014 2015 Level: developer 2016 2017 .seealso: TaoSetTolerances(),TaoGetConvergedReason(),TaoSetConvergedReason() 2018 @*/ 2019 2020 PetscErrorCode TaoDefaultConvergenceTest(Tao tao,void *dummy) 2021 { 2022 PetscInt niter=tao->niter, nfuncs=PetscMax(tao->nfuncs,tao->nfuncgrads); 2023 PetscInt max_funcs=tao->max_funcs; 2024 PetscReal gnorm=tao->residual, gnorm0=tao->gnorm0; 2025 PetscReal f=tao->fc, steptol=tao->steptol,trradius=tao->step; 2026 PetscReal gatol=tao->gatol,grtol=tao->grtol,gttol=tao->gttol; 2027 PetscReal catol=tao->catol,crtol=tao->crtol; 2028 PetscReal fmin=tao->fmin, cnorm=tao->cnorm; 2029 TaoConvergedReason reason=tao->reason; 2030 PetscErrorCode ierr; 2031 2032 PetscFunctionBegin; 2033 PetscValidHeaderSpecific(tao, TAO_CLASSID,1); 2034 if (reason != TAO_CONTINUE_ITERATING) { 2035 PetscFunctionReturn(0); 2036 } 2037 2038 if (PetscIsInfOrNanReal(f)) { 2039 ierr = PetscInfo(tao,"Failed to converged, function value is Inf or NaN\n");CHKERRQ(ierr); 2040 reason = TAO_DIVERGED_NAN; 2041 } else if (f <= fmin && cnorm <=catol) { 2042 ierr = PetscInfo2(tao,"Converged due to function value %g < minimum function value %g\n", (double)f,(double)fmin);CHKERRQ(ierr); 2043 reason = TAO_CONVERGED_MINF; 2044 } else if (gnorm<= gatol && cnorm <=catol) { 2045 ierr = PetscInfo2(tao,"Converged due to residual norm ||g(X)||=%g < %g\n",(double)gnorm,(double)gatol);CHKERRQ(ierr); 2046 reason = TAO_CONVERGED_GATOL; 2047 } else if (f!=0 && PetscAbsReal(gnorm/f) <= grtol && cnorm <= crtol) { 2048 ierr = PetscInfo2(tao,"Converged due to residual ||g(X)||/|f(X)| =%g < %g\n",(double)(gnorm/f),(double)grtol);CHKERRQ(ierr); 2049 reason = TAO_CONVERGED_GRTOL; 2050 } else if (gnorm0 != 0 && ((gttol == 0 && gnorm == 0) || gnorm/gnorm0 < gttol) && cnorm <= crtol) { 2051 ierr = PetscInfo2(tao,"Converged due to relative residual norm ||g(X)||/||g(X0)|| = %g < %g\n",(double)(gnorm/gnorm0),(double)gttol);CHKERRQ(ierr); 2052 reason = TAO_CONVERGED_GTTOL; 2053 } else if (nfuncs > max_funcs) { 2054 ierr = PetscInfo2(tao,"Exceeded maximum number of function evaluations: %D > %D\n", nfuncs,max_funcs);CHKERRQ(ierr); 2055 reason = TAO_DIVERGED_MAXFCN; 2056 } else if (tao->lsflag != 0) { 2057 ierr = PetscInfo(tao,"Tao Line Search failure.\n");CHKERRQ(ierr); 2058 reason = TAO_DIVERGED_LS_FAILURE; 2059 } else if (trradius < steptol && niter > 0) { 2060 ierr = PetscInfo2(tao,"Trust region/step size too small: %g < %g\n", (double)trradius,(double)steptol);CHKERRQ(ierr); 2061 reason = TAO_CONVERGED_STEPTOL; 2062 } else if (niter >= tao->max_it) { 2063 ierr = PetscInfo2(tao,"Exceeded maximum number of iterations: %D > %D\n",niter,tao->max_it);CHKERRQ(ierr); 2064 reason = TAO_DIVERGED_MAXITS; 2065 } else { 2066 reason = TAO_CONTINUE_ITERATING; 2067 } 2068 tao->reason = reason; 2069 PetscFunctionReturn(0); 2070 } 2071 2072 /*@C 2073 TaoSetOptionsPrefix - Sets the prefix used for searching for all 2074 TAO options in the database. 2075 2076 Logically Collective on Tao 2077 2078 Input Parameters: 2079 + tao - the Tao context 2080 - prefix - the prefix string to prepend to all TAO option requests 2081 2082 Notes: 2083 A hyphen (-) must NOT be given at the beginning of the prefix name. 2084 The first character of all runtime options is AUTOMATICALLY the hyphen. 2085 2086 For example, to distinguish between the runtime options for two 2087 different TAO solvers, one could call 2088 .vb 2089 TaoSetOptionsPrefix(tao1,"sys1_") 2090 TaoSetOptionsPrefix(tao2,"sys2_") 2091 .ve 2092 2093 This would enable use of different options for each system, such as 2094 .vb 2095 -sys1_tao_method blmvm -sys1_tao_grtol 1.e-3 2096 -sys2_tao_method lmvm -sys2_tao_grtol 1.e-4 2097 .ve 2098 2099 Level: advanced 2100 2101 .seealso: TaoAppendOptionsPrefix(), TaoGetOptionsPrefix() 2102 @*/ 2103 2104 PetscErrorCode TaoSetOptionsPrefix(Tao tao, const char p[]) 2105 { 2106 PetscErrorCode ierr; 2107 2108 PetscFunctionBegin; 2109 ierr = PetscObjectSetOptionsPrefix((PetscObject)tao,p);CHKERRQ(ierr); 2110 if (tao->linesearch) { 2111 ierr = TaoLineSearchSetOptionsPrefix(tao->linesearch,p);CHKERRQ(ierr); 2112 } 2113 if (tao->ksp) { 2114 ierr = KSPSetOptionsPrefix(tao->ksp,p);CHKERRQ(ierr); 2115 } 2116 PetscFunctionReturn(0); 2117 } 2118 2119 /*@C 2120 TaoAppendOptionsPrefix - Appends to the prefix used for searching for all 2121 TAO options in the database. 2122 2123 Logically Collective on Tao 2124 2125 Input Parameters: 2126 + tao - the Tao solver context 2127 - prefix - the prefix string to prepend to all TAO option requests 2128 2129 Notes: 2130 A hyphen (-) must NOT be given at the beginning of the prefix name. 2131 The first character of all runtime options is AUTOMATICALLY the hyphen. 2132 2133 Level: advanced 2134 2135 .seealso: TaoSetOptionsPrefix(), TaoGetOptionsPrefix() 2136 @*/ 2137 PetscErrorCode TaoAppendOptionsPrefix(Tao tao, const char p[]) 2138 { 2139 PetscErrorCode ierr; 2140 2141 PetscFunctionBegin; 2142 ierr = PetscObjectAppendOptionsPrefix((PetscObject)tao,p);CHKERRQ(ierr); 2143 if (tao->linesearch) { 2144 ierr = TaoLineSearchSetOptionsPrefix(tao->linesearch,p);CHKERRQ(ierr); 2145 } 2146 if (tao->ksp) { 2147 ierr = KSPSetOptionsPrefix(tao->ksp,p);CHKERRQ(ierr); 2148 } 2149 PetscFunctionReturn(0); 2150 } 2151 2152 /*@C 2153 TaoGetOptionsPrefix - Gets the prefix used for searching for all 2154 TAO options in the database 2155 2156 Not Collective 2157 2158 Input Parameters: 2159 . tao - the Tao context 2160 2161 Output Parameters: 2162 . prefix - pointer to the prefix string used is returned 2163 2164 Notes: 2165 On the fortran side, the user should pass in a string 'prefix' of 2166 sufficient length to hold the prefix. 2167 2168 Level: advanced 2169 2170 .seealso: TaoSetOptionsPrefix(), TaoAppendOptionsPrefix() 2171 @*/ 2172 PetscErrorCode TaoGetOptionsPrefix(Tao tao, const char *p[]) 2173 { 2174 return PetscObjectGetOptionsPrefix((PetscObject)tao,p); 2175 } 2176 2177 /*@C 2178 TaoSetType - Sets the method for the unconstrained minimization solver. 2179 2180 Collective on Tao 2181 2182 Input Parameters: 2183 + solver - the Tao solver context 2184 - type - a known method 2185 2186 Options Database Key: 2187 . -tao_type <type> - Sets the method; use -help for a list 2188 of available methods (for instance, "-tao_type lmvm" or "-tao_type tron") 2189 2190 Available methods include: 2191 + nls - Newton's method with line search for unconstrained minimization 2192 . ntr - Newton's method with trust region for unconstrained minimization 2193 . ntl - Newton's method with trust region, line search for unconstrained minimization 2194 . lmvm - Limited memory variable metric method for unconstrained minimization 2195 . cg - Nonlinear conjugate gradient method for unconstrained minimization 2196 . nm - Nelder-Mead algorithm for derivate-free unconstrained minimization 2197 . tron - Newton Trust Region method for bound constrained minimization 2198 . gpcg - Newton Trust Region method for quadratic bound constrained minimization 2199 . blmvm - Limited memory variable metric method for bound constrained minimization 2200 - pounders - Model-based algorithm pounder extended for nonlinear least squares 2201 2202 Level: intermediate 2203 2204 .seealso: TaoCreate(), TaoGetType(), TaoType 2205 2206 @*/ 2207 PetscErrorCode TaoSetType(Tao tao, TaoType type) 2208 { 2209 PetscErrorCode ierr; 2210 PetscErrorCode (*create_xxx)(Tao); 2211 PetscBool issame; 2212 2213 PetscFunctionBegin; 2214 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2215 2216 ierr = PetscObjectTypeCompare((PetscObject)tao,type,&issame);CHKERRQ(ierr); 2217 if (issame) PetscFunctionReturn(0); 2218 2219 ierr = PetscFunctionListFind(TaoList, type, (void(**)(void))&create_xxx);CHKERRQ(ierr); 2220 if (!create_xxx) SETERRQ1(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested Tao type %s",type); 2221 2222 /* Destroy the existing solver information */ 2223 if (tao->ops->destroy) { 2224 ierr = (*tao->ops->destroy)(tao);CHKERRQ(ierr); 2225 } 2226 ierr = KSPDestroy(&tao->ksp);CHKERRQ(ierr); 2227 ierr = TaoLineSearchDestroy(&tao->linesearch);CHKERRQ(ierr); 2228 ierr = VecDestroy(&tao->gradient);CHKERRQ(ierr); 2229 ierr = VecDestroy(&tao->stepdirection);CHKERRQ(ierr); 2230 2231 tao->ops->setup = NULL; 2232 tao->ops->solve = NULL; 2233 tao->ops->view = NULL; 2234 tao->ops->setfromoptions = NULL; 2235 tao->ops->destroy = NULL; 2236 2237 tao->setupcalled = PETSC_FALSE; 2238 2239 ierr = (*create_xxx)(tao);CHKERRQ(ierr); 2240 ierr = PetscObjectChangeTypeName((PetscObject)tao,type);CHKERRQ(ierr); 2241 PetscFunctionReturn(0); 2242 } 2243 2244 /*MC 2245 TaoRegister - Adds a method to the TAO package for unconstrained minimization. 2246 2247 Synopsis: 2248 TaoRegister(char *name_solver,char *path,char *name_Create,PetscErrorCode (*routine_Create)(Tao)) 2249 2250 Not collective 2251 2252 Input Parameters: 2253 + sname - name of a new user-defined solver 2254 - func - routine to Create method context 2255 2256 Notes: 2257 TaoRegister() may be called multiple times to add several user-defined solvers. 2258 2259 Sample usage: 2260 .vb 2261 TaoRegister("my_solver",MySolverCreate); 2262 .ve 2263 2264 Then, your solver can be chosen with the procedural interface via 2265 $ TaoSetType(tao,"my_solver") 2266 or at runtime via the option 2267 $ -tao_type my_solver 2268 2269 Level: advanced 2270 2271 .seealso: TaoRegisterAll(), TaoRegisterDestroy() 2272 M*/ 2273 PetscErrorCode TaoRegister(const char sname[], PetscErrorCode (*func)(Tao)) 2274 { 2275 PetscErrorCode ierr; 2276 2277 PetscFunctionBegin; 2278 ierr = TaoInitializePackage();CHKERRQ(ierr); 2279 ierr = PetscFunctionListAdd(&TaoList,sname, (void (*)(void))func);CHKERRQ(ierr); 2280 PetscFunctionReturn(0); 2281 } 2282 2283 /*@C 2284 TaoRegisterDestroy - Frees the list of minimization solvers that were 2285 registered by TaoRegisterDynamic(). 2286 2287 Not Collective 2288 2289 Level: advanced 2290 2291 .seealso: TaoRegisterAll(), TaoRegister() 2292 @*/ 2293 PetscErrorCode TaoRegisterDestroy(void) 2294 { 2295 PetscErrorCode ierr; 2296 PetscFunctionBegin; 2297 ierr = PetscFunctionListDestroy(&TaoList);CHKERRQ(ierr); 2298 TaoRegisterAllCalled = PETSC_FALSE; 2299 PetscFunctionReturn(0); 2300 } 2301 2302 /*@ 2303 TaoGetIterationNumber - Gets the number of Tao iterations completed 2304 at this time. 2305 2306 Not Collective 2307 2308 Input Parameter: 2309 . tao - Tao context 2310 2311 Output Parameter: 2312 . iter - iteration number 2313 2314 Notes: 2315 For example, during the computation of iteration 2 this would return 1. 2316 2317 Level: intermediate 2318 2319 .seealso: TaoGetLinearSolveIterations(), TaoGetResidualNorm(), TaoGetObjective() 2320 @*/ 2321 PetscErrorCode TaoGetIterationNumber(Tao tao,PetscInt *iter) 2322 { 2323 PetscFunctionBegin; 2324 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2325 PetscValidIntPointer(iter,2); 2326 *iter = tao->niter; 2327 PetscFunctionReturn(0); 2328 } 2329 2330 /*@ 2331 TaoGetObjective - Gets the current value of the objective function 2332 at this time. 2333 2334 Not Collective 2335 2336 Input Parameter: 2337 . tao - Tao context 2338 2339 Output Parameter: 2340 . value - the current value 2341 2342 Level: intermediate 2343 2344 .seealso: TaoGetLinearSolveIterations(), TaoGetIterationNumber(), TaoGetResidualNorm() 2345 @*/ 2346 PetscErrorCode TaoGetObjective(Tao tao,PetscReal *value) 2347 { 2348 PetscFunctionBegin; 2349 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2350 PetscValidRealPointer(value,2); 2351 *value = tao->fc; 2352 PetscFunctionReturn(0); 2353 } 2354 2355 /*@ 2356 TaoGetResidualNorm - Gets the current value of the norm of the residual 2357 at this time. 2358 2359 Not Collective 2360 2361 Input Parameter: 2362 . tao - Tao context 2363 2364 Output Parameter: 2365 . value - the current value 2366 2367 Level: intermediate 2368 2369 Developer Note: This is the 2-norm of the residual, we cannot use TaoGetGradientNorm() because that has 2370 a different meaning. For some reason Tao sometimes calls the gradient the residual. 2371 2372 .seealso: TaoGetLinearSolveIterations(), TaoGetIterationNumber(), TaoGetObjective() 2373 @*/ 2374 PetscErrorCode TaoGetResidualNorm(Tao tao,PetscReal *value) 2375 { 2376 PetscFunctionBegin; 2377 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2378 PetscValidRealPointer(value,2); 2379 *value = tao->residual; 2380 PetscFunctionReturn(0); 2381 } 2382 2383 /*@ 2384 TaoSetIterationNumber - Sets the current iteration number. 2385 2386 Not Collective 2387 2388 Input Parameters: 2389 + tao - Tao context 2390 - iter - iteration number 2391 2392 Level: developer 2393 2394 .seealso: TaoGetLinearSolveIterations() 2395 @*/ 2396 PetscErrorCode TaoSetIterationNumber(Tao tao,PetscInt iter) 2397 { 2398 PetscErrorCode ierr; 2399 2400 PetscFunctionBegin; 2401 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2402 ierr = PetscObjectSAWsTakeAccess((PetscObject)tao);CHKERRQ(ierr); 2403 tao->niter = iter; 2404 ierr = PetscObjectSAWsGrantAccess((PetscObject)tao);CHKERRQ(ierr); 2405 PetscFunctionReturn(0); 2406 } 2407 2408 /*@ 2409 TaoGetTotalIterationNumber - Gets the total number of Tao iterations 2410 completed. This number keeps accumulating if multiple solves 2411 are called with the Tao object. 2412 2413 Not Collective 2414 2415 Input Parameter: 2416 . tao - Tao context 2417 2418 Output Parameter: 2419 . iter - iteration number 2420 2421 Notes: 2422 The total iteration count is updated after each solve, if there is a current 2423 TaoSolve() in progress then those iterations are not yet counted. 2424 2425 Level: intermediate 2426 2427 .seealso: TaoGetLinearSolveIterations() 2428 @*/ 2429 PetscErrorCode TaoGetTotalIterationNumber(Tao tao,PetscInt *iter) 2430 { 2431 PetscFunctionBegin; 2432 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2433 PetscValidIntPointer(iter,2); 2434 *iter = tao->ntotalits; 2435 PetscFunctionReturn(0); 2436 } 2437 2438 /*@ 2439 TaoSetTotalIterationNumber - Sets the current total iteration number. 2440 2441 Not Collective 2442 2443 Input Parameters: 2444 + tao - Tao context 2445 - iter - iteration number 2446 2447 Level: developer 2448 2449 .seealso: TaoGetLinearSolveIterations() 2450 @*/ 2451 PetscErrorCode TaoSetTotalIterationNumber(Tao tao,PetscInt iter) 2452 { 2453 PetscErrorCode ierr; 2454 2455 PetscFunctionBegin; 2456 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2457 ierr = PetscObjectSAWsTakeAccess((PetscObject)tao);CHKERRQ(ierr); 2458 tao->ntotalits = iter; 2459 ierr = PetscObjectSAWsGrantAccess((PetscObject)tao);CHKERRQ(ierr); 2460 PetscFunctionReturn(0); 2461 } 2462 2463 /*@ 2464 TaoSetConvergedReason - Sets the termination flag on a Tao object 2465 2466 Logically Collective on Tao 2467 2468 Input Parameters: 2469 + tao - the Tao context 2470 - reason - one of 2471 $ TAO_CONVERGED_ATOL (2), 2472 $ TAO_CONVERGED_RTOL (3), 2473 $ TAO_CONVERGED_STEPTOL (4), 2474 $ TAO_CONVERGED_MINF (5), 2475 $ TAO_CONVERGED_USER (6), 2476 $ TAO_DIVERGED_MAXITS (-2), 2477 $ TAO_DIVERGED_NAN (-4), 2478 $ TAO_DIVERGED_MAXFCN (-5), 2479 $ TAO_DIVERGED_LS_FAILURE (-6), 2480 $ TAO_DIVERGED_TR_REDUCTION (-7), 2481 $ TAO_DIVERGED_USER (-8), 2482 $ TAO_CONTINUE_ITERATING (0) 2483 2484 Level: intermediate 2485 2486 @*/ 2487 PetscErrorCode TaoSetConvergedReason(Tao tao, TaoConvergedReason reason) 2488 { 2489 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2490 PetscFunctionBegin; 2491 tao->reason = reason; 2492 PetscFunctionReturn(0); 2493 } 2494 2495 /*@ 2496 TaoGetConvergedReason - Gets the reason the Tao iteration was stopped. 2497 2498 Not Collective 2499 2500 Input Parameter: 2501 . tao - the Tao solver context 2502 2503 Output Parameter: 2504 . reason - one of 2505 $ TAO_CONVERGED_GATOL (3) ||g(X)|| < gatol 2506 $ TAO_CONVERGED_GRTOL (4) ||g(X)|| / f(X) < grtol 2507 $ TAO_CONVERGED_GTTOL (5) ||g(X)|| / ||g(X0)|| < gttol 2508 $ TAO_CONVERGED_STEPTOL (6) step size small 2509 $ TAO_CONVERGED_MINF (7) F < F_min 2510 $ TAO_CONVERGED_USER (8) User defined 2511 $ TAO_DIVERGED_MAXITS (-2) its > maxits 2512 $ TAO_DIVERGED_NAN (-4) Numerical problems 2513 $ TAO_DIVERGED_MAXFCN (-5) fevals > max_funcsals 2514 $ TAO_DIVERGED_LS_FAILURE (-6) line search failure 2515 $ TAO_DIVERGED_TR_REDUCTION (-7) trust region failure 2516 $ TAO_DIVERGED_USER(-8) (user defined) 2517 $ TAO_CONTINUE_ITERATING (0) 2518 2519 where 2520 + X - current solution 2521 . X0 - initial guess 2522 . f(X) - current function value 2523 . f(X*) - true solution (estimated) 2524 . g(X) - current gradient 2525 . its - current iterate number 2526 . maxits - maximum number of iterates 2527 . fevals - number of function evaluations 2528 - max_funcsals - maximum number of function evaluations 2529 2530 Level: intermediate 2531 2532 .seealso: TaoSetConvergenceTest(), TaoSetTolerances() 2533 2534 @*/ 2535 PetscErrorCode TaoGetConvergedReason(Tao tao, TaoConvergedReason *reason) 2536 { 2537 PetscFunctionBegin; 2538 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2539 PetscValidPointer(reason,2); 2540 *reason = tao->reason; 2541 PetscFunctionReturn(0); 2542 } 2543 2544 /*@ 2545 TaoGetSolutionStatus - Get the current iterate, objective value, 2546 residual, infeasibility, and termination 2547 2548 Not Collective 2549 2550 Input Parameter: 2551 . tao - the Tao context 2552 2553 Output Parameters: 2554 + iterate - the current iterate number (>=0) 2555 . f - the current function value 2556 . gnorm - the square of the gradient norm, duality gap, or other measure indicating distance from optimality. 2557 . cnorm - the infeasibility of the current solution with regard to the constraints. 2558 . xdiff - the step length or trust region radius of the most recent iterate. 2559 - reason - The termination reason, which can equal TAO_CONTINUE_ITERATING 2560 2561 Level: intermediate 2562 2563 Note: 2564 TAO returns the values set by the solvers in the routine TaoMonitor(). 2565 2566 Note: 2567 If any of the output arguments are set to NULL, no corresponding value will be returned. 2568 2569 .seealso: TaoMonitor(), TaoGetConvergedReason() 2570 @*/ 2571 PetscErrorCode TaoGetSolutionStatus(Tao tao, PetscInt *its, PetscReal *f, PetscReal *gnorm, PetscReal *cnorm, PetscReal *xdiff, TaoConvergedReason *reason) 2572 { 2573 PetscFunctionBegin; 2574 if (its) *its=tao->niter; 2575 if (f) *f=tao->fc; 2576 if (gnorm) *gnorm=tao->residual; 2577 if (cnorm) *cnorm=tao->cnorm; 2578 if (reason) *reason=tao->reason; 2579 if (xdiff) *xdiff=tao->step; 2580 PetscFunctionReturn(0); 2581 } 2582 2583 /*@C 2584 TaoGetType - Gets the current Tao algorithm. 2585 2586 Not Collective 2587 2588 Input Parameter: 2589 . tao - the Tao solver context 2590 2591 Output Parameter: 2592 . type - Tao method 2593 2594 Level: intermediate 2595 2596 @*/ 2597 PetscErrorCode TaoGetType(Tao tao,TaoType *type) 2598 { 2599 PetscFunctionBegin; 2600 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2601 PetscValidPointer(type,2); 2602 *type=((PetscObject)tao)->type_name; 2603 PetscFunctionReturn(0); 2604 } 2605 2606 /*@C 2607 TaoMonitor - Monitor the solver and the current solution. This 2608 routine will record the iteration number and residual statistics, 2609 call any monitors specified by the user, and calls the convergence-check routine. 2610 2611 Input Parameters: 2612 + tao - the Tao context 2613 . its - the current iterate number (>=0) 2614 . f - the current objective function value 2615 . res - the gradient norm, square root of the duality gap, or other measure indicating distince from optimality. This measure will be recorded and 2616 used for some termination tests. 2617 . cnorm - the infeasibility of the current solution with regard to the constraints. 2618 - steplength - multiple of the step direction added to the previous iterate. 2619 2620 Output Parameters: 2621 . reason - The termination reason, which can equal TAO_CONTINUE_ITERATING 2622 2623 Options Database Key: 2624 . -tao_monitor - Use the default monitor, which prints statistics to standard output 2625 2626 .seealso TaoGetConvergedReason(), TaoMonitorDefault(), TaoSetMonitor() 2627 2628 Level: developer 2629 2630 @*/ 2631 PetscErrorCode TaoMonitor(Tao tao, PetscInt its, PetscReal f, PetscReal res, PetscReal cnorm, PetscReal steplength) 2632 { 2633 PetscErrorCode ierr; 2634 PetscInt i; 2635 2636 PetscFunctionBegin; 2637 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2638 tao->fc = f; 2639 tao->residual = res; 2640 tao->cnorm = cnorm; 2641 tao->step = steplength; 2642 if (!its) { 2643 tao->cnorm0 = cnorm; tao->gnorm0 = res; 2644 } 2645 if (PetscIsInfOrNanReal(f) || PetscIsInfOrNanReal(res)) SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_USER, "User provided compute function generated Inf or NaN"); 2646 for (i=0;i<tao->numbermonitors;i++) { 2647 ierr = (*tao->monitor[i])(tao,tao->monitorcontext[i]);CHKERRQ(ierr); 2648 } 2649 PetscFunctionReturn(0); 2650 } 2651 2652 /*@ 2653 TaoSetConvergenceHistory - Sets the array used to hold the convergence history. 2654 2655 Logically Collective on Tao 2656 2657 Input Parameters: 2658 + tao - the Tao solver context 2659 . obj - array to hold objective value history 2660 . resid - array to hold residual history 2661 . cnorm - array to hold constraint violation history 2662 . lits - integer array holds the number of linear iterations for each Tao iteration 2663 . na - size of obj, resid, and cnorm 2664 - reset - PetscTrue indicates each new minimization resets the history counter to zero, 2665 else it continues storing new values for new minimizations after the old ones 2666 2667 Notes: 2668 If set, TAO will fill the given arrays with the indicated 2669 information at each iteration. If 'obj','resid','cnorm','lits' are 2670 *all* NULL then space (using size na, or 1000 if na is PETSC_DECIDE or 2671 PETSC_DEFAULT) is allocated for the history. 2672 If not all are NULL, then only the non-NULL information categories 2673 will be stored, the others will be ignored. 2674 2675 Any convergence information after iteration number 'na' will not be stored. 2676 2677 This routine is useful, e.g., when running a code for purposes 2678 of accurate performance monitoring, when no I/O should be done 2679 during the section of code that is being timed. 2680 2681 Level: intermediate 2682 2683 .seealso: TaoGetConvergenceHistory() 2684 2685 @*/ 2686 PetscErrorCode TaoSetConvergenceHistory(Tao tao, PetscReal obj[], PetscReal resid[], PetscReal cnorm[], PetscInt lits[], PetscInt na,PetscBool reset) 2687 { 2688 PetscErrorCode ierr; 2689 2690 PetscFunctionBegin; 2691 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2692 if (obj) PetscValidRealPointer(obj,2); 2693 if (resid) PetscValidRealPointer(resid,3); 2694 if (cnorm) PetscValidRealPointer(cnorm,4); 2695 if (lits) PetscValidIntPointer(lits,5); 2696 2697 if (na == PETSC_DECIDE || na == PETSC_DEFAULT) na = 1000; 2698 if (!obj && !resid && !cnorm && !lits) { 2699 ierr = PetscCalloc4(na,&obj,na,&resid,na,&cnorm,na,&lits);CHKERRQ(ierr); 2700 tao->hist_malloc = PETSC_TRUE; 2701 } 2702 2703 tao->hist_obj = obj; 2704 tao->hist_resid = resid; 2705 tao->hist_cnorm = cnorm; 2706 tao->hist_lits = lits; 2707 tao->hist_max = na; 2708 tao->hist_reset = reset; 2709 tao->hist_len = 0; 2710 PetscFunctionReturn(0); 2711 } 2712 2713 /*@C 2714 TaoGetConvergenceHistory - Gets the arrays used to hold the convergence history. 2715 2716 Collective on Tao 2717 2718 Input Parameter: 2719 . tao - the Tao context 2720 2721 Output Parameters: 2722 + obj - array used to hold objective value history 2723 . resid - array used to hold residual history 2724 . cnorm - array used to hold constraint violation history 2725 . lits - integer array used to hold linear solver iteration count 2726 - nhist - size of obj, resid, cnorm, and lits 2727 2728 Notes: 2729 This routine must be preceded by calls to TaoSetConvergenceHistory() 2730 and TaoSolve(), otherwise it returns useless information. 2731 2732 The calling sequence for this routine in Fortran is 2733 $ call TaoGetConvergenceHistory(Tao tao, PetscInt nhist, PetscErrorCode ierr) 2734 2735 This routine is useful, e.g., when running a code for purposes 2736 of accurate performance monitoring, when no I/O should be done 2737 during the section of code that is being timed. 2738 2739 Level: advanced 2740 2741 .seealso: TaoSetConvergenceHistory() 2742 2743 @*/ 2744 PetscErrorCode TaoGetConvergenceHistory(Tao tao, PetscReal **obj, PetscReal **resid, PetscReal **cnorm, PetscInt **lits, PetscInt *nhist) 2745 { 2746 PetscFunctionBegin; 2747 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2748 if (obj) *obj = tao->hist_obj; 2749 if (cnorm) *cnorm = tao->hist_cnorm; 2750 if (resid) *resid = tao->hist_resid; 2751 if (nhist) *nhist = tao->hist_len; 2752 PetscFunctionReturn(0); 2753 } 2754 2755 /*@ 2756 TaoSetApplicationContext - Sets the optional user-defined context for 2757 a solver. 2758 2759 Logically Collective on Tao 2760 2761 Input Parameters: 2762 + tao - the Tao context 2763 - usrP - optional user context 2764 2765 Level: intermediate 2766 2767 .seealso: TaoGetApplicationContext(), TaoSetApplicationContext() 2768 @*/ 2769 PetscErrorCode TaoSetApplicationContext(Tao tao,void *usrP) 2770 { 2771 PetscFunctionBegin; 2772 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2773 tao->user = usrP; 2774 PetscFunctionReturn(0); 2775 } 2776 2777 /*@ 2778 TaoGetApplicationContext - Gets the user-defined context for a 2779 TAO solvers. 2780 2781 Not Collective 2782 2783 Input Parameter: 2784 . tao - Tao context 2785 2786 Output Parameter: 2787 . usrP - user context 2788 2789 Level: intermediate 2790 2791 .seealso: TaoSetApplicationContext() 2792 @*/ 2793 PetscErrorCode TaoGetApplicationContext(Tao tao,void *usrP) 2794 { 2795 PetscFunctionBegin; 2796 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2797 *(void**)usrP = tao->user; 2798 PetscFunctionReturn(0); 2799 } 2800 2801 /*@ 2802 TaoSetGradientNorm - Sets the matrix used to define the inner product that measures the size of the gradient. 2803 2804 Collective on tao 2805 2806 Input Parameters: 2807 + tao - the Tao context 2808 - M - gradient norm 2809 2810 Level: beginner 2811 2812 .seealso: TaoGetGradientNorm(), TaoGradientNorm() 2813 @*/ 2814 PetscErrorCode TaoSetGradientNorm(Tao tao, Mat M) 2815 { 2816 PetscErrorCode ierr; 2817 2818 PetscFunctionBegin; 2819 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2820 ierr = PetscObjectReference((PetscObject)M);CHKERRQ(ierr); 2821 ierr = MatDestroy(&tao->gradient_norm);CHKERRQ(ierr); 2822 ierr = VecDestroy(&tao->gradient_norm_tmp);CHKERRQ(ierr); 2823 tao->gradient_norm = M; 2824 ierr = MatCreateVecs(M, NULL, &tao->gradient_norm_tmp);CHKERRQ(ierr); 2825 PetscFunctionReturn(0); 2826 } 2827 2828 /*@ 2829 TaoGetGradientNorm - Returns the matrix used to define the inner product for measuring the size of the gradient. 2830 2831 Not Collective 2832 2833 Input Parameter: 2834 . tao - Tao context 2835 2836 Output Parameter: 2837 . M - gradient norm 2838 2839 Level: beginner 2840 2841 .seealso: TaoSetGradientNorm(), TaoGradientNorm() 2842 @*/ 2843 PetscErrorCode TaoGetGradientNorm(Tao tao, Mat *M) 2844 { 2845 PetscFunctionBegin; 2846 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2847 *M = tao->gradient_norm; 2848 PetscFunctionReturn(0); 2849 } 2850 2851 /*@C 2852 TaoGradientNorm - Compute the norm with respect to the inner product the user has set. 2853 2854 Collective on tao 2855 2856 Input Parameters: 2857 + tao - the Tao context 2858 . gradient - the gradient to be computed 2859 - norm - the norm type 2860 2861 Output Parameter: 2862 . gnorm - the gradient norm 2863 2864 Level: developer 2865 2866 .seealso: TaoSetGradientNorm(), TaoGetGradientNorm() 2867 @*/ 2868 PetscErrorCode TaoGradientNorm(Tao tao, Vec gradient, NormType type, PetscReal *gnorm) 2869 { 2870 PetscErrorCode ierr; 2871 2872 PetscFunctionBegin; 2873 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 2874 PetscValidHeaderSpecific(gradient,VEC_CLASSID,2); 2875 PetscValidLogicalCollectiveEnum(tao,type,3); 2876 PetscValidPointer(gnorm,4); 2877 if (tao->gradient_norm) { 2878 PetscScalar gnorms; 2879 2880 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."); 2881 ierr = MatMult(tao->gradient_norm, gradient, tao->gradient_norm_tmp);CHKERRQ(ierr); 2882 ierr = VecDot(gradient, tao->gradient_norm_tmp, &gnorms);CHKERRQ(ierr); 2883 *gnorm = PetscRealPart(PetscSqrtScalar(gnorms)); 2884 } else { 2885 ierr = VecNorm(gradient, type, gnorm);CHKERRQ(ierr); 2886 } 2887 PetscFunctionReturn(0); 2888 } 2889 2890 /*@C 2891 TaoMonitorDrawCtxCreate - Creates the monitor context for TaoMonitorDrawCtx 2892 2893 Collective on Tao 2894 2895 Output Patameter: 2896 . ctx - the monitor context 2897 2898 Options Database: 2899 . -tao_draw_solution_initial - show initial guess as well as current solution 2900 2901 Level: intermediate 2902 2903 .seealso: TaoMonitorSet(), TaoMonitorDefault(), VecView(), TaoMonitorDrawCtx() 2904 @*/ 2905 PetscErrorCode TaoMonitorDrawCtxCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscInt howoften,TaoMonitorDrawCtx *ctx) 2906 { 2907 PetscErrorCode ierr; 2908 2909 PetscFunctionBegin; 2910 ierr = PetscNew(ctx);CHKERRQ(ierr); 2911 ierr = PetscViewerDrawOpen(comm,host,label,x,y,m,n,&(*ctx)->viewer);CHKERRQ(ierr); 2912 ierr = PetscViewerSetFromOptions((*ctx)->viewer);CHKERRQ(ierr); 2913 (*ctx)->howoften = howoften; 2914 PetscFunctionReturn(0); 2915 } 2916 2917 /*@C 2918 TaoMonitorDrawCtxDestroy - Destroys the monitor context for TaoMonitorDrawSolution() 2919 2920 Collective on Tao 2921 2922 Input Parameters: 2923 . ctx - the monitor context 2924 2925 Level: intermediate 2926 2927 .seealso: TaoMonitorSet(), TaoMonitorDefault(), VecView(), TaoMonitorDrawSolution() 2928 @*/ 2929 PetscErrorCode TaoMonitorDrawCtxDestroy(TaoMonitorDrawCtx *ictx) 2930 { 2931 PetscErrorCode ierr; 2932 2933 PetscFunctionBegin; 2934 ierr = PetscViewerDestroy(&(*ictx)->viewer);CHKERRQ(ierr); 2935 ierr = PetscFree(*ictx);CHKERRQ(ierr); 2936 PetscFunctionReturn(0); 2937 } 2938