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