1 #include <../src/tao/bound/impls/bqnk/bqnk.h> /*I "petsctao.h" I*/ /*I "petscmat.h" I*/ 2 #include <petscksp.h> 3 4 static const char *BQNK_INIT[64] = {"constant", "direction"}; 5 static const char *BNK_UPDATE[64] = {"step", "reduction", "interpolation"}; 6 static const char *BNK_AS[64] = {"none", "bertsekas"}; 7 8 static PetscErrorCode TaoBQNKComputeHessian(Tao tao) 9 { 10 TAO_BNK *bnk = (TAO_BNK *)tao->data; 11 TAO_BQNK *bqnk = (TAO_BQNK*)bnk->ctx; 12 PetscErrorCode ierr; 13 PetscReal gnorm2, delta; 14 15 PetscFunctionBegin; 16 /* Alias the LMVM matrix into the TAO hessian */ 17 if (tao->hessian) { 18 ierr = MatDestroy(&tao->hessian);CHKERRQ(ierr); 19 } 20 if (tao->hessian_pre) { 21 ierr = MatDestroy(&tao->hessian_pre);CHKERRQ(ierr); 22 } 23 ierr = PetscObjectReference((PetscObject)bqnk->B);CHKERRQ(ierr); 24 tao->hessian = bqnk->B; 25 ierr = PetscObjectReference((PetscObject)bqnk->B);CHKERRQ(ierr); 26 tao->hessian_pre = bqnk->B; 27 /* Update the Hessian with the latest solution */ 28 if (bqnk->is_spd) { 29 gnorm2 = bnk->gnorm*bnk->gnorm; 30 if (gnorm2 == 0.0) gnorm2 = PETSC_MACHINE_EPSILON; 31 if (bnk->f == 0.0) { 32 delta = 2.0 / gnorm2; 33 } else { 34 delta = 2.0 * PetscAbsScalar(bnk->f) / gnorm2; 35 } 36 ierr = MatLMVMSymBroydenSetDelta(bqnk->B, delta);CHKERRQ(ierr); 37 } 38 ierr = MatLMVMUpdate(tao->hessian, tao->solution, bnk->unprojected_gradient);CHKERRQ(ierr); 39 ierr = MatLMVMResetShift(tao->hessian);CHKERRQ(ierr); 40 /* Prepare the reduced sub-matrices for the inactive set */ 41 ierr = MatDestroy(&bnk->H_inactive);CHKERRQ(ierr); 42 if (bnk->active_idx) { 43 ierr = MatCreateSubMatrixVirtual(tao->hessian, bnk->inactive_idx, bnk->inactive_idx, &bnk->H_inactive);CHKERRQ(ierr); 44 ierr = PCLMVMSetIS(bqnk->pc, bnk->inactive_idx);CHKERRQ(ierr); 45 } else { 46 ierr = PetscObjectReference((PetscObject)tao->hessian);CHKERRQ(ierr); 47 bnk->H_inactive = tao->hessian; 48 ierr = PCLMVMClearIS(bqnk->pc);CHKERRQ(ierr); 49 } 50 ierr = MatDestroy(&bnk->Hpre_inactive);CHKERRQ(ierr); 51 ierr = PetscObjectReference((PetscObject)bnk->H_inactive);CHKERRQ(ierr); 52 bnk->Hpre_inactive = bnk->H_inactive; 53 PetscFunctionReturn(0); 54 } 55 56 static PetscErrorCode TaoBQNKComputeStep(Tao tao, PetscBool shift, KSPConvergedReason *ksp_reason, PetscInt *step_type) 57 { 58 TAO_BNK *bnk = (TAO_BNK *)tao->data; 59 TAO_BQNK *bqnk = (TAO_BQNK*)bnk->ctx; 60 PetscErrorCode ierr; 61 62 PetscFunctionBegin; 63 ierr = TaoBNKComputeStep(tao, shift, ksp_reason, step_type);CHKERRQ(ierr); 64 if (*ksp_reason < 0) { 65 /* Krylov solver failed to converge so reset the LMVM matrix */ 66 ierr = MatLMVMReset(bqnk->B, PETSC_FALSE);CHKERRQ(ierr); 67 ierr = MatLMVMUpdate(bqnk->B, tao->solution, bnk->unprojected_gradient);CHKERRQ(ierr); 68 } 69 PetscFunctionReturn(0); 70 } 71 72 PetscErrorCode TaoSolve_BQNK(Tao tao) 73 { 74 TAO_BNK *bnk = (TAO_BNK *)tao->data; 75 TAO_BQNK *bqnk = (TAO_BQNK*)bnk->ctx; 76 Mat_LMVM *lmvm = (Mat_LMVM*)bqnk->B->data; 77 Mat_LMVM *J0; 78 Mat_SymBrdn *diag_ctx; 79 PetscBool flg = PETSC_FALSE; 80 PetscErrorCode ierr; 81 82 PetscFunctionBegin; 83 if (!tao->recycle) { 84 ierr = MatLMVMReset(bqnk->B, PETSC_FALSE);CHKERRQ(ierr); 85 lmvm->nresets = 0; 86 if (lmvm->J0) { 87 ierr = PetscObjectBaseTypeCompare((PetscObject)lmvm->J0, MATLMVM, &flg);CHKERRQ(ierr); 88 if (flg) { 89 J0 = (Mat_LMVM*)lmvm->J0->data; 90 J0->nresets = 0; 91 } 92 } 93 flg = PETSC_FALSE; 94 ierr = PetscObjectTypeCompareAny((PetscObject)bqnk->B, &flg, MATLMVMSYMBROYDEN, MATLMVMSYMBADBROYDEN, MATLMVMBFGS, MATLMVMDFP, ""); 95 if (flg) { 96 diag_ctx = (Mat_SymBrdn*)lmvm->ctx; 97 J0 = (Mat_LMVM*)diag_ctx->D->data; 98 J0->nresets = 0; 99 } 100 } 101 ierr = (*bqnk->solve)(tao);CHKERRQ(ierr); 102 PetscFunctionReturn(0); 103 } 104 105 PetscErrorCode TaoSetUp_BQNK(Tao tao) 106 { 107 TAO_BNK *bnk = (TAO_BNK *)tao->data; 108 TAO_BQNK *bqnk = (TAO_BQNK*)bnk->ctx; 109 PetscErrorCode ierr; 110 PetscInt n, N; 111 PetscBool is_lmvm, is_sym, is_spd; 112 113 PetscFunctionBegin; 114 ierr = TaoSetUp_BNK(tao);CHKERRQ(ierr); 115 ierr = VecGetLocalSize(tao->solution,&n);CHKERRQ(ierr); 116 ierr = VecGetSize(tao->solution,&N);CHKERRQ(ierr); 117 ierr = MatSetSizes(bqnk->B, n, n, N, N);CHKERRQ(ierr); 118 ierr = MatLMVMAllocate(bqnk->B,tao->solution,bnk->unprojected_gradient);CHKERRQ(ierr); 119 ierr = PetscObjectBaseTypeCompare((PetscObject)bqnk->B, MATLMVM, &is_lmvm);CHKERRQ(ierr); 120 if (!is_lmvm) SETERRQ(PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "Matrix must be an LMVM-type"); 121 ierr = MatGetOption(bqnk->B, MAT_SYMMETRIC, &is_sym);CHKERRQ(ierr); 122 if (!is_sym) SETERRQ(PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "LMVM matrix must be symmetric"); 123 ierr = MatGetOption(bqnk->B, MAT_SPD, &is_spd);CHKERRQ(ierr); 124 ierr = KSPGetPC(tao->ksp, &bqnk->pc);CHKERRQ(ierr); 125 ierr = PCSetType(bqnk->pc, PCLMVM);CHKERRQ(ierr); 126 ierr = PCLMVMSetMatLMVM(bqnk->pc, bqnk->B);CHKERRQ(ierr); 127 PetscFunctionReturn(0); 128 } 129 130 static PetscErrorCode TaoSetFromOptions_BQNK(PetscOptionItems *PetscOptionsObject,Tao tao) 131 { 132 TAO_BNK *bnk = (TAO_BNK *)tao->data; 133 TAO_BQNK *bqnk = (TAO_BQNK*)bnk->ctx; 134 PetscErrorCode ierr; 135 KSPType ksp_type; 136 137 PetscFunctionBegin; 138 ierr = PetscOptionsHead(PetscOptionsObject,"Quasi-Newton-Krylov method for bound constrained optimization");CHKERRQ(ierr); 139 ierr = PetscOptionsEList("-tao_bqnk_init_type", "radius initialization type", "", BQNK_INIT, BQNK_INIT_TYPES, BQNK_INIT[bnk->init_type], &bnk->init_type, NULL);CHKERRQ(ierr); 140 ierr = PetscOptionsEList("-tao_bqnk_update_type", "radius update type", "", BNK_UPDATE, BNK_UPDATE_TYPES, BNK_UPDATE[bnk->update_type], &bnk->update_type, NULL);CHKERRQ(ierr); 141 ierr = PetscOptionsEList("-tao_bqnk_as_type", "active set estimation method", "", BNK_AS, BNK_AS_TYPES, BNK_AS[bnk->as_type], &bnk->as_type, NULL);CHKERRQ(ierr); 142 ierr = PetscOptionsReal("-tao_bqnk_sval", "(developer) Hessian perturbation starting value", "", bnk->sval, &bnk->sval,NULL);CHKERRQ(ierr); 143 ierr = PetscOptionsReal("-tao_bqnk_imin", "(developer) minimum initial Hessian perturbation", "", bnk->imin, &bnk->imin,NULL);CHKERRQ(ierr); 144 ierr = PetscOptionsReal("-tao_bqnk_imax", "(developer) maximum initial Hessian perturbation", "", bnk->imax, &bnk->imax,NULL);CHKERRQ(ierr); 145 ierr = PetscOptionsReal("-tao_bqnk_imfac", "(developer) initial merit factor for Hessian perturbation", "", bnk->imfac, &bnk->imfac,NULL);CHKERRQ(ierr); 146 ierr = PetscOptionsReal("-tao_bqnk_pmin", "(developer) minimum Hessian perturbation", "", bnk->pmin, &bnk->pmin,NULL);CHKERRQ(ierr); 147 ierr = PetscOptionsReal("-tao_bqnk_pmax", "(developer) maximum Hessian perturbation", "", bnk->pmax, &bnk->pmax,NULL);CHKERRQ(ierr); 148 ierr = PetscOptionsReal("-tao_bqnk_pgfac", "(developer) Hessian perturbation growth factor", "", bnk->pgfac, &bnk->pgfac,NULL);CHKERRQ(ierr); 149 ierr = PetscOptionsReal("-tao_bqnk_psfac", "(developer) Hessian perturbation shrink factor", "", bnk->psfac, &bnk->psfac,NULL);CHKERRQ(ierr); 150 ierr = PetscOptionsReal("-tao_bqnk_pmgfac", "(developer) merit growth factor for Hessian perturbation", "", bnk->pmgfac, &bnk->pmgfac,NULL);CHKERRQ(ierr); 151 ierr = PetscOptionsReal("-tao_bqnk_pmsfac", "(developer) merit shrink factor for Hessian perturbation", "", bnk->pmsfac, &bnk->pmsfac,NULL);CHKERRQ(ierr); 152 ierr = PetscOptionsReal("-tao_bqnk_eta1", "(developer) threshold for rejecting step (-tao_bqnk_update_type reduction)", "", bnk->eta1, &bnk->eta1,NULL);CHKERRQ(ierr); 153 ierr = PetscOptionsReal("-tao_bqnk_eta2", "(developer) threshold for accepting marginal step (-tao_bqnk_update_type reduction)", "", bnk->eta2, &bnk->eta2,NULL);CHKERRQ(ierr); 154 ierr = PetscOptionsReal("-tao_bqnk_eta3", "(developer) threshold for accepting reasonable step (-tao_bqnk_update_type reduction)", "", bnk->eta3, &bnk->eta3,NULL);CHKERRQ(ierr); 155 ierr = PetscOptionsReal("-tao_bqnk_eta4", "(developer) threshold for accepting good step (-tao_bqnk_update_type reduction)", "", bnk->eta4, &bnk->eta4,NULL);CHKERRQ(ierr); 156 ierr = PetscOptionsReal("-tao_bqnk_alpha1", "(developer) radius reduction factor for rejected step (-tao_bqnk_update_type reduction)", "", bnk->alpha1, &bnk->alpha1,NULL);CHKERRQ(ierr); 157 ierr = PetscOptionsReal("-tao_bqnk_alpha2", "(developer) radius reduction factor for marginally accepted bad step (-tao_bqnk_update_type reduction)", "", bnk->alpha2, &bnk->alpha2,NULL);CHKERRQ(ierr); 158 ierr = PetscOptionsReal("-tao_bqnk_alpha3", "(developer) radius increase factor for reasonable accepted step (-tao_bqnk_update_type reduction)", "", bnk->alpha3, &bnk->alpha3,NULL);CHKERRQ(ierr); 159 ierr = PetscOptionsReal("-tao_bqnk_alpha4", "(developer) radius increase factor for good accepted step (-tao_bqnk_update_type reduction)", "", bnk->alpha4, &bnk->alpha4,NULL);CHKERRQ(ierr); 160 ierr = PetscOptionsReal("-tao_bqnk_alpha5", "(developer) radius increase factor for very good accepted step (-tao_bqnk_update_type reduction)", "", bnk->alpha5, &bnk->alpha5,NULL);CHKERRQ(ierr); 161 ierr = PetscOptionsReal("-tao_bqnk_nu1", "(developer) threshold for small line-search step length (-tao_bqnk_update_type step)", "", bnk->nu1, &bnk->nu1,NULL);CHKERRQ(ierr); 162 ierr = PetscOptionsReal("-tao_bqnk_nu2", "(developer) threshold for reasonable line-search step length (-tao_bqnk_update_type step)", "", bnk->nu2, &bnk->nu2,NULL);CHKERRQ(ierr); 163 ierr = PetscOptionsReal("-tao_bqnk_nu3", "(developer) threshold for large line-search step length (-tao_bqnk_update_type step)", "", bnk->nu3, &bnk->nu3,NULL);CHKERRQ(ierr); 164 ierr = PetscOptionsReal("-tao_bqnk_nu4", "(developer) threshold for very large line-search step length (-tao_bqnk_update_type step)", "", bnk->nu4, &bnk->nu4,NULL);CHKERRQ(ierr); 165 ierr = PetscOptionsReal("-tao_bqnk_omega1", "(developer) radius reduction factor for very small line-search step length (-tao_bqnk_update_type step)", "", bnk->omega1, &bnk->omega1,NULL);CHKERRQ(ierr); 166 ierr = PetscOptionsReal("-tao_bqnk_omega2", "(developer) radius reduction factor for small line-search step length (-tao_bqnk_update_type step)", "", bnk->omega2, &bnk->omega2,NULL);CHKERRQ(ierr); 167 ierr = PetscOptionsReal("-tao_bqnk_omega3", "(developer) radius factor for decent line-search step length (-tao_bqnk_update_type step)", "", bnk->omega3, &bnk->omega3,NULL);CHKERRQ(ierr); 168 ierr = PetscOptionsReal("-tao_bqnk_omega4", "(developer) radius increase factor for large line-search step length (-tao_bqnk_update_type step)", "", bnk->omega4, &bnk->omega4,NULL);CHKERRQ(ierr); 169 ierr = PetscOptionsReal("-tao_bqnk_omega5", "(developer) radius increase factor for very large line-search step length (-tao_bqnk_update_type step)", "", bnk->omega5, &bnk->omega5,NULL);CHKERRQ(ierr); 170 ierr = PetscOptionsReal("-tao_bqnk_mu1", "(developer) threshold for accepting very good step (-tao_bqnk_update_type interpolation)", "", bnk->mu1, &bnk->mu1,NULL);CHKERRQ(ierr); 171 ierr = PetscOptionsReal("-tao_bqnk_mu2", "(developer) threshold for accepting good step (-tao_bqnk_update_type interpolation)", "", bnk->mu2, &bnk->mu2,NULL);CHKERRQ(ierr); 172 ierr = PetscOptionsReal("-tao_bqnk_gamma1", "(developer) radius reduction factor for rejected very bad step (-tao_bqnk_update_type interpolation)", "", bnk->gamma1, &bnk->gamma1,NULL);CHKERRQ(ierr); 173 ierr = PetscOptionsReal("-tao_bqnk_gamma2", "(developer) radius reduction factor for rejected bad step (-tao_bqnk_update_type interpolation)", "", bnk->gamma2, &bnk->gamma2,NULL);CHKERRQ(ierr); 174 ierr = PetscOptionsReal("-tao_bqnk_gamma3", "(developer) radius increase factor for accepted good step (-tao_bqnk_update_type interpolation)", "", bnk->gamma3, &bnk->gamma3,NULL);CHKERRQ(ierr); 175 ierr = PetscOptionsReal("-tao_bqnk_gamma4", "(developer) radius increase factor for accepted very good step (-tao_bqnk_update_type interpolation)", "", bnk->gamma4, &bnk->gamma4,NULL);CHKERRQ(ierr); 176 ierr = PetscOptionsReal("-tao_bqnk_theta", "(developer) trust region interpolation factor (-tao_bqnk_update_type interpolation)", "", bnk->theta, &bnk->theta,NULL);CHKERRQ(ierr); 177 ierr = PetscOptionsReal("-tao_bqnk_min_radius", "(developer) lower bound on initial radius", "", bnk->min_radius, &bnk->min_radius,NULL);CHKERRQ(ierr); 178 ierr = PetscOptionsReal("-tao_bqnk_max_radius", "(developer) upper bound on radius", "", bnk->max_radius, &bnk->max_radius,NULL);CHKERRQ(ierr); 179 ierr = PetscOptionsReal("-tao_bqnk_epsilon", "(developer) tolerance used when computing actual and predicted reduction", "", bnk->epsilon, &bnk->epsilon,NULL);CHKERRQ(ierr); 180 ierr = PetscOptionsReal("-tao_bqnk_as_tol", "(developer) initial tolerance used when estimating actively bounded variables", "", bnk->as_tol, &bnk->as_tol,NULL);CHKERRQ(ierr); 181 ierr = PetscOptionsReal("-tao_bqnk_as_step", "(developer) step length used when estimating actively bounded variables", "", bnk->as_step, &bnk->as_step,NULL);CHKERRQ(ierr); 182 ierr = PetscOptionsInt("-tao_bqnk_max_cg_its", "number of BNCG iterations to take for each Newton step", "", bnk->max_cg_its, &bnk->max_cg_its,NULL);CHKERRQ(ierr); 183 ierr = PetscOptionsTail();CHKERRQ(ierr); 184 ierr = TaoSetFromOptions(bnk->bncg);CHKERRQ(ierr); 185 ierr = TaoLineSearchSetFromOptions(tao->linesearch);CHKERRQ(ierr); 186 ierr = KSPSetFromOptions(tao->ksp);CHKERRQ(ierr); 187 ierr = KSPGetType(tao->ksp,&ksp_type);CHKERRQ(ierr); 188 ierr = PetscStrcmp(ksp_type,KSPNASH,&bnk->is_nash);CHKERRQ(ierr); 189 ierr = PetscStrcmp(ksp_type,KSPSTCG,&bnk->is_stcg);CHKERRQ(ierr); 190 ierr = PetscStrcmp(ksp_type,KSPGLTR,&bnk->is_gltr);CHKERRQ(ierr); 191 if (bnk->init_type == BNK_INIT_INTERPOLATION) bnk->init_type = BNK_INIT_DIRECTION; 192 ierr = MatSetFromOptions(bqnk->B);CHKERRQ(ierr); 193 ierr = MatGetOption(bqnk->B, MAT_SPD, &bqnk->is_spd);CHKERRQ(ierr); 194 PetscFunctionReturn(0); 195 } 196 197 static PetscErrorCode TaoView_BQNK(Tao tao, PetscViewer viewer) 198 { 199 TAO_BNK *bnk = (TAO_BNK*)tao->data; 200 TAO_BQNK *bqnk = (TAO_BQNK*)bnk->ctx; 201 PetscErrorCode ierr; 202 PetscBool isascii; 203 204 PetscFunctionBegin; 205 ierr = TaoView_BNK(tao, viewer);CHKERRQ(ierr); 206 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);CHKERRQ(ierr); 207 if (isascii) { 208 ierr = PetscViewerPushFormat(viewer, PETSC_VIEWER_ASCII_INFO);CHKERRQ(ierr); 209 ierr = MatView(bqnk->B, viewer);CHKERRQ(ierr); 210 ierr = PetscViewerPopFormat(viewer);CHKERRQ(ierr); 211 } 212 PetscFunctionReturn(0); 213 } 214 215 static PetscErrorCode TaoDestroy_BQNK(Tao tao) 216 { 217 TAO_BNK *bnk = (TAO_BNK*)tao->data; 218 TAO_BQNK *bqnk = (TAO_BQNK*)bnk->ctx; 219 PetscErrorCode ierr; 220 221 PetscFunctionBegin; 222 ierr = MatDestroy(&bnk->Hpre_inactive);CHKERRQ(ierr); 223 ierr = MatDestroy(&bnk->H_inactive);CHKERRQ(ierr); 224 ierr = MatDestroy(&bqnk->B);CHKERRQ(ierr); 225 ierr = PetscFree(bnk->ctx);CHKERRQ(ierr); 226 ierr = TaoDestroy_BNK(tao);CHKERRQ(ierr); 227 PetscFunctionReturn(0); 228 } 229 230 PETSC_INTERN PetscErrorCode TaoCreate_BQNK(Tao tao) 231 { 232 TAO_BNK *bnk; 233 TAO_BQNK *bqnk; 234 PetscErrorCode ierr; 235 236 PetscFunctionBegin; 237 ierr = TaoCreate_BNK(tao);CHKERRQ(ierr); 238 ierr = KSPSetOptionsPrefix(tao->ksp, "tao_bqnk_");CHKERRQ(ierr); 239 tao->ops->solve = TaoSolve_BQNK; 240 tao->ops->setfromoptions = TaoSetFromOptions_BQNK; 241 tao->ops->destroy = TaoDestroy_BQNK; 242 tao->ops->view = TaoView_BQNK; 243 tao->ops->setup = TaoSetUp_BQNK; 244 245 bnk = (TAO_BNK *)tao->data; 246 bnk->computehessian = TaoBQNKComputeHessian; 247 bnk->computestep = TaoBQNKComputeStep; 248 bnk->init_type = BNK_INIT_DIRECTION; 249 250 ierr = PetscNewLog(tao,&bqnk);CHKERRQ(ierr); 251 bnk->ctx = (void*)bqnk; 252 bqnk->is_spd = PETSC_TRUE; 253 254 ierr = MatCreate(PetscObjectComm((PetscObject)tao), &bqnk->B);CHKERRQ(ierr); 255 ierr = PetscObjectIncrementTabLevel((PetscObject)bqnk->B, (PetscObject)tao, 1);CHKERRQ(ierr); 256 ierr = MatSetOptionsPrefix(bqnk->B, "tao_bqnk_");CHKERRQ(ierr); 257 ierr = MatSetType(bqnk->B, MATLMVMSR1);CHKERRQ(ierr); 258 PetscFunctionReturn(0); 259 } 260 261 /*@ 262 TaoGetLMVMMatrix - Returns a pointer to the internal LMVM matrix. Valid 263 only for quasi-Newton family of methods. 264 265 Input Parameters: 266 . tao - Tao solver context 267 268 Output Parameters: 269 . B - LMVM matrix 270 271 Level: advanced 272 273 .seealso: TAOBQNLS, TAOBQNKLS, TAOBQNKTL, TAOBQNKTR, MATLMVM, TaoSetLMVMMatrix() 274 @*/ 275 PetscErrorCode TaoGetLMVMMatrix(Tao tao, Mat *B) 276 { 277 TAO_BNK *bnk = (TAO_BNK*)tao->data; 278 TAO_BQNK *bqnk = (TAO_BQNK*)bnk->ctx; 279 PetscErrorCode ierr; 280 PetscBool flg = PETSC_FALSE; 281 282 PetscFunctionBegin; 283 ierr = PetscObjectTypeCompareAny((PetscObject)tao, &flg, TAOBQNLS, TAOBQNKLS, TAOBQNKTR, TAOBQNKTL, "");CHKERRQ(ierr); 284 if (!flg) SETERRQ(PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "LMVM Matrix only exists for quasi-Newton algorithms"); 285 *B = bqnk->B; 286 PetscFunctionReturn(0); 287 } 288 289 /*@ 290 TaoSetLMVMMatrix - Sets an external LMVM matrix into the Tao solver. Valid 291 only for quasi-Newton family of methods. 292 293 QN family of methods create their own LMVM matrices and users who wish to 294 manipulate this matrix should use TaoGetLMVMMatrix() instead. 295 296 Input Parameters: 297 + tao - Tao solver context 298 - B - LMVM matrix 299 300 Level: advanced 301 302 .seealso: TAOBQNLS, TAOBQNKLS, TAOBQNKTL, TAOBQNKTR, MATLMVM, TaoGetLMVMMatrix() 303 @*/ 304 PetscErrorCode TaoSetLMVMMatrix(Tao tao, Mat B) 305 { 306 TAO_BNK *bnk = (TAO_BNK*)tao->data; 307 TAO_BQNK *bqnk = (TAO_BQNK*)bnk->ctx; 308 PetscErrorCode ierr; 309 PetscBool flg = PETSC_FALSE; 310 311 PetscFunctionBegin; 312 ierr = PetscObjectTypeCompareAny((PetscObject)tao, &flg, TAOBQNLS, TAOBQNKLS, TAOBQNKTR, TAOBQNKTL, "");CHKERRQ(ierr); 313 if (!flg) SETERRQ(PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "LMVM Matrix only exists for quasi-Newton algorithms"); 314 ierr = PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &flg);CHKERRQ(ierr); 315 if (!flg) SETERRQ(PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "Given matrix is not an LMVM matrix"); 316 if (bqnk->B) { 317 ierr = MatDestroy(&bqnk->B);CHKERRQ(ierr); 318 } 319 ierr = PetscObjectReference((PetscObject)B);CHKERRQ(ierr); 320 bqnk->B = B; 321 PetscFunctionReturn(0); 322 } 323