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