1 #include <../src/tao/bound/impls/bqnk/bqnk.h>
2
3 static const char *BNK_AS[64] = {"none", "bertsekas"};
4
TaoBQNLSComputeHessian(Tao tao)5 static PetscErrorCode TaoBQNLSComputeHessian(Tao tao)
6 {
7 TAO_BNK *bnk = (TAO_BNK *)tao->data;
8 TAO_BQNK *bqnk = (TAO_BQNK *)bnk->ctx;
9 PetscReal gnorm2, delta;
10
11 PetscFunctionBegin;
12 /* Compute the initial scaling and update the approximation */
13 gnorm2 = bnk->gnorm * bnk->gnorm;
14 if (gnorm2 == 0.0) gnorm2 = PETSC_MACHINE_EPSILON;
15 if (bnk->f == 0.0) delta = 2.0 / gnorm2;
16 else delta = 2.0 * PetscAbsScalar(bnk->f) / gnorm2;
17 PetscCall(MatLMVMSymBroydenSetDelta(bqnk->B, delta));
18 PetscCall(MatLMVMUpdate(bqnk->B, tao->solution, bnk->unprojected_gradient));
19 PetscFunctionReturn(PETSC_SUCCESS);
20 }
21
TaoBQNLSComputeStep(Tao tao,PetscBool shift,KSPConvergedReason * ksp_reason,PetscInt * step_type)22 static PetscErrorCode TaoBQNLSComputeStep(Tao tao, PetscBool shift, KSPConvergedReason *ksp_reason, PetscInt *step_type)
23 {
24 TAO_BNK *bnk = (TAO_BNK *)tao->data;
25 TAO_BQNK *bqnk = (TAO_BQNK *)bnk->ctx;
26 PetscInt nupdates;
27
28 PetscFunctionBegin;
29 PetscCall(MatSolve(bqnk->B, tao->gradient, tao->stepdirection));
30 PetscCall(VecScale(tao->stepdirection, -1.0));
31 PetscCall(TaoBNKBoundStep(tao, bnk->as_type, tao->stepdirection));
32 *ksp_reason = KSP_CONVERGED_ATOL;
33 PetscCall(MatLMVMGetUpdateCount(bqnk->B, &nupdates));
34 if (nupdates == 0) *step_type = BNK_SCALED_GRADIENT;
35 else *step_type = BNK_BFGS;
36 PetscFunctionReturn(PETSC_SUCCESS);
37 }
38
TaoSetFromOptions_BQNLS(Tao tao,PetscOptionItems PetscOptionsObject)39 static PetscErrorCode TaoSetFromOptions_BQNLS(Tao tao, PetscOptionItems PetscOptionsObject)
40 {
41 TAO_BNK *bnk = (TAO_BNK *)tao->data;
42 TAO_BQNK *bqnk = (TAO_BQNK *)bnk->ctx;
43 PetscBool is_set, is_spd;
44
45 PetscFunctionBegin;
46 PetscOptionsHeadBegin(PetscOptionsObject, "Quasi-Newton-Krylov method for bound constrained optimization");
47 PetscCall(PetscOptionsEList("-tao_bnk_as_type", "active set estimation method", "", BNK_AS, BNK_AS_TYPES, BNK_AS[bnk->as_type], &bnk->as_type, NULL));
48 PetscCall(PetscOptionsReal("-tao_bnk_epsilon", "(developer) tolerance used when computing actual and predicted reduction", "", bnk->epsilon, &bnk->epsilon, NULL));
49 PetscCall(PetscOptionsReal("-tao_bnk_as_tol", "(developer) initial tolerance used when estimating actively bounded variables", "", bnk->as_tol, &bnk->as_tol, NULL));
50 PetscCall(PetscOptionsReal("-tao_bnk_as_step", "(developer) step length used when estimating actively bounded variables", "", bnk->as_step, &bnk->as_step, NULL));
51 PetscCall(PetscOptionsInt("-tao_bnk_max_cg_its", "number of BNCG iterations to take for each Newton step", "", bnk->max_cg_its, &bnk->max_cg_its, NULL));
52 PetscOptionsHeadEnd();
53
54 PetscCall(TaoSetOptionsPrefix(bnk->bncg, ((PetscObject)tao)->prefix));
55 PetscCall(TaoAppendOptionsPrefix(bnk->bncg, "tao_bnk_"));
56 PetscCall(TaoSetFromOptions(bnk->bncg));
57
58 PetscCall(MatSetOptionsPrefix(bqnk->B, ((PetscObject)tao)->prefix));
59 PetscCall(MatAppendOptionsPrefix(bqnk->B, "tao_bqnls_"));
60 PetscCall(MatSetFromOptions(bqnk->B));
61 PetscCall(MatIsSPDKnown(bqnk->B, &is_set, &is_spd));
62 PetscCheck(is_set && is_spd, PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "LMVM matrix must be symmetric positive-definite");
63 PetscFunctionReturn(PETSC_SUCCESS);
64 }
65
66 /*MC
67 TAOBQNLS - Bounded Quasi-Newton Line Search method for nonlinear minimization with bound
68 constraints. This method approximates the action of the inverse-Hessian with a
69 limited memory quasi-Newton formula. The quasi-Newton matrix and its options are
70 accessible via the prefix `-tao_bqnls_`
71
72 Options Database Keys:
73 + -tao_bnk_max_cg_its - maximum number of bounded conjugate-gradient iterations taken in each Newton loop
74 . -tao_bnk_as_type - active-set estimation method ("none", "bertsekas")
75 . -tao_bnk_epsilon - (developer) tolerance for small pred/actual ratios that trigger automatic step acceptance
76 . -tao_bnk_as_tol - (developer) initial tolerance used in estimating bounded active variables (-as_type bertsekas)
77 - -tao_bnk_as_step - (developer) trial step length used in estimating bounded active variables (-as_type bertsekas)
78
79 Level: beginner
80
81 .seealso: `TAOBNK`
82 M*/
TaoCreate_BQNLS(Tao tao)83 PETSC_EXTERN PetscErrorCode TaoCreate_BQNLS(Tao tao)
84 {
85 TAO_BNK *bnk;
86 TAO_BQNK *bqnk;
87
88 PetscFunctionBegin;
89 PetscCall(TaoCreate_BQNK(tao));
90 tao->ops->setfromoptions = TaoSetFromOptions_BQNLS;
91
92 bnk = (TAO_BNK *)tao->data;
93 bnk->update_type = BNK_UPDATE_STEP;
94 bnk->computehessian = TaoBQNLSComputeHessian;
95 bnk->computestep = TaoBQNLSComputeStep;
96
97 bqnk = (TAO_BQNK *)bnk->ctx;
98 bqnk->solve = TaoSolve_BNLS;
99 PetscCall(MatSetType(bqnk->B, MATLMVMBFGS));
100 PetscFunctionReturn(PETSC_SUCCESS);
101 }
102