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