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