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