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) { 17 delta = 2.0 / gnorm2; 18 } else { 19 delta = 2.0 * PetscAbsScalar(bnk->f) / gnorm2; 20 } 21 ierr = MatSymBrdnSetDelta(bqnk->B, delta);CHKERRQ(ierr); 22 ierr = MatLMVMUpdate(bqnk->B, tao->solution, bnk->unprojected_gradient);CHKERRQ(ierr); 23 PetscFunctionReturn(0); 24 } 25 26 static PetscErrorCode TaoBQNLSComputeStep(Tao tao, PetscBool shift, KSPConvergedReason *ksp_reason, PetscInt *step_type) 27 { 28 TAO_BNK *bnk = (TAO_BNK *)tao->data; 29 TAO_BQNK *bqnk = (TAO_BQNK*)bnk->ctx; 30 PetscErrorCode ierr; 31 PetscInt nupdates; 32 33 PetscFunctionBegin; 34 ierr = MatSolve(bqnk->B, tao->gradient, tao->stepdirection);CHKERRQ(ierr); 35 ierr = VecScale(tao->stepdirection, -1.0);CHKERRQ(ierr); 36 ierr = TaoBNKBoundStep(tao, bnk->as_type, tao->stepdirection);CHKERRQ(ierr); 37 *ksp_reason = KSP_CONVERGED_ATOL; 38 ierr = MatLMVMGetUpdateCount(bqnk->B, &nupdates);CHKERRQ(ierr); 39 if (nupdates == 0) { 40 *step_type = BNK_SCALED_GRADIENT; 41 } else { 42 *step_type = BNK_BFGS; 43 } 44 PetscFunctionReturn(0); 45 } 46 47 static PetscErrorCode TaoSetFromOptions_BQNLS(PetscOptionItems *PetscOptionsObject,Tao tao) 48 { 49 TAO_BNK *bnk = (TAO_BNK *)tao->data; 50 TAO_BQNK *bqnk = (TAO_BQNK*)bnk->ctx; 51 PetscErrorCode ierr; 52 KSPType ksp_type; 53 PetscBool is_spd; 54 55 PetscFunctionBegin; 56 ierr = PetscOptionsHead(PetscOptionsObject,"Quasi-Newton-Krylov method for bound constrained optimization");CHKERRQ(ierr); 57 ierr = PetscOptionsEList("-tao_bqnls_as_type", "active set estimation method", "", BNK_AS, BNK_AS_TYPES, BNK_AS[bnk->as_type], &bnk->as_type, 0);CHKERRQ(ierr); 58 ierr = PetscOptionsReal("-tao_bqnls_epsilon", "(developer) tolerance used when computing actual and predicted reduction", "", bnk->epsilon, &bnk->epsilon,NULL);CHKERRQ(ierr); 59 ierr = PetscOptionsReal("-tao_bqnls_as_tol", "(developer) initial tolerance used when estimating actively bounded variables", "", bnk->as_tol, &bnk->as_tol,NULL);CHKERRQ(ierr); 60 ierr = PetscOptionsReal("-tao_bqnls_as_step", "(developer) step length used when estimating actively bounded variables", "", bnk->as_step, &bnk->as_step,NULL);CHKERRQ(ierr); 61 ierr = PetscOptionsInt("-tao_bqnls_max_cg_its", "number of BNCG iterations to take for each Newton step", "", bnk->max_cg_its, &bnk->max_cg_its,NULL);CHKERRQ(ierr); 62 ierr = PetscOptionsTail();CHKERRQ(ierr); 63 ierr = TaoSetFromOptions(bnk->bncg);CHKERRQ(ierr); 64 ierr = TaoLineSearchSetFromOptions(tao->linesearch);CHKERRQ(ierr); 65 ierr = KSPSetFromOptions(tao->ksp);CHKERRQ(ierr); 66 ierr = KSPGetType(tao->ksp,&ksp_type);CHKERRQ(ierr); 67 bnk->is_nash = bnk->is_gltr = bnk->is_stcg = PETSC_FALSE; 68 ierr = MatSetFromOptions(bqnk->B);CHKERRQ(ierr); 69 ierr = MatGetOption(bqnk->B, MAT_SPD, &is_spd);CHKERRQ(ierr); 70 if (!is_spd) SETERRQ(PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "LMVM matrix must be symmetric positive-definite"); 71 PetscFunctionReturn(0); 72 } 73 74 PETSC_EXTERN PetscErrorCode TaoCreate_BQNLS(Tao tao) 75 { 76 TAO_BNK *bnk; 77 TAO_BQNK *bqnk; 78 PetscErrorCode ierr; 79 80 PetscFunctionBegin; 81 ierr = TaoCreate_BQNK(tao);CHKERRQ(ierr); 82 ierr = KSPSetOptionsPrefix(tao->ksp, "unused");CHKERRQ(ierr); 83 tao->ops->solve = TaoSolve_BNLS; 84 tao->ops->setfromoptions = TaoSetFromOptions_BQNLS; 85 86 bnk = (TAO_BNK*)tao->data; 87 bnk->update_type = BNK_UPDATE_STEP; 88 bnk->computehessian = TaoBQNLSComputeHessian; 89 bnk->computestep = TaoBQNLSComputeStep; 90 91 bqnk = (TAO_BQNK*)bnk->ctx; 92 ierr = MatSetOptionsPrefix(bqnk->B, "tao_bqnls_");CHKERRQ(ierr); 93 ierr = MatSetType(bqnk->B, MATLMVMBFGS);CHKERRQ(ierr); 94 PetscFunctionReturn(0); 95 } 96