1 #include <petsc/private/taoimpl.h> /*I "petsctao.h" I*/ 2 3 /*@ 4 TaoSetInitialVector - Sets the initial guess for the solve 5 6 Logically collective on Tao 7 8 Input Parameters: 9 + tao - the Tao context 10 - x0 - the initial guess 11 12 Level: beginner 13 .seealso: TaoCreate(), TaoSolve() 14 @*/ 15 16 PetscErrorCode TaoSetInitialVector(Tao tao, Vec x0) 17 { 18 PetscErrorCode ierr; 19 20 PetscFunctionBegin; 21 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 22 if (x0) { 23 PetscValidHeaderSpecific(x0,VEC_CLASSID,2); 24 PetscObjectReference((PetscObject)x0); 25 } 26 ierr = VecDestroy(&tao->solution);CHKERRQ(ierr); 27 tao->solution = x0; 28 PetscFunctionReturn(0); 29 } 30 31 PetscErrorCode TaoTestGradient(Tao tao,Vec g1) 32 { 33 Vec x = tao->solution,g2,g3; 34 PetscBool complete_print = PETSC_FALSE,test = PETSC_FALSE; 35 PetscReal hcnorm,fdnorm,hcmax,fdmax,diffmax,diffnorm; 36 PetscScalar dot; 37 MPI_Comm comm; 38 PetscViewer viewer; 39 PetscInt tabs; 40 static PetscBool directionsprinted = PETSC_FALSE; 41 PetscErrorCode ierr; 42 43 PetscFunctionBegin; 44 ierr = PetscObjectOptionsBegin((PetscObject)tao);CHKERRQ(ierr); 45 ierr = PetscOptionsName("-tao_test_gradient","Compare hand-coded and finite difference Gradients","None",&test);CHKERRQ(ierr); 46 ierr = PetscOptionsBool("-tao_test_gradient_view","View difference between hand-coded and finite difference Gradients element entries","None",complete_print,&complete_print,NULL);CHKERRQ(ierr); 47 ierr = PetscOptionsEnd();CHKERRQ(ierr); 48 if (!test) PetscFunctionReturn(0); 49 50 ierr = PetscObjectGetComm((PetscObject)tao,&comm);CHKERRQ(ierr); 51 ierr = PetscViewerASCIIGetStdout(comm,&viewer);CHKERRQ(ierr); 52 ierr = PetscViewerASCIIGetTab(viewer, &tabs);CHKERRQ(ierr); 53 ierr = PetscViewerASCIISetTab(viewer, ((PetscObject)tao)->tablevel);CHKERRQ(ierr); 54 ierr = PetscViewerASCIIPrintf(viewer," ---------- Testing Gradient -------------\n");CHKERRQ(ierr); 55 if (!complete_print && !directionsprinted) { 56 ierr = PetscViewerASCIIPrintf(viewer," Run with -tao_test_gradient_view and optionally -tao_test_gradient <threshold> to show difference\n");CHKERRQ(ierr); 57 ierr = PetscViewerASCIIPrintf(viewer," of hand-coded and finite difference gradient entries greater than <threshold>.\n");CHKERRQ(ierr); 58 } 59 if (!directionsprinted) { 60 ierr = PetscViewerASCIIPrintf(viewer," Testing hand-coded Gradient, if (for double precision runs) ||G - Gfd||_F/||G||_F is\n");CHKERRQ(ierr); 61 ierr = PetscViewerASCIIPrintf(viewer," O(1.e-8), the hand-coded Hessian is probably correct.\n");CHKERRQ(ierr); 62 directionsprinted = PETSC_TRUE; 63 } 64 65 ierr = VecDuplicate(x,&g2);CHKERRQ(ierr); 66 ierr = VecDuplicate(x,&g3);CHKERRQ(ierr); 67 68 /* Compute finite difference gradient, assume the gradient is already computed by TaoComputeGradient() and put into g1 */ 69 ierr = TaoDefaultComputeGradient(tao,x,g2,NULL);CHKERRQ(ierr); 70 71 ierr = VecNorm(g2,NORM_2,&fdnorm);CHKERRQ(ierr); 72 ierr = VecNorm(g1,NORM_2,&hcnorm);CHKERRQ(ierr); 73 ierr = VecNorm(g2,NORM_INFINITY,&fdmax);CHKERRQ(ierr); 74 ierr = VecNorm(g1,NORM_INFINITY,&hcmax);CHKERRQ(ierr); 75 ierr = VecDot(g1,g2,&dot);CHKERRQ(ierr); 76 ierr = VecCopy(g1,g3);CHKERRQ(ierr); 77 ierr = VecAXPY(g3,-1.0,g2);CHKERRQ(ierr); 78 ierr = VecNorm(g3,NORM_2,&diffnorm);CHKERRQ(ierr); 79 ierr = VecNorm(g3,NORM_INFINITY,&diffmax);CHKERRQ(ierr); 80 ierr = PetscViewerASCIIPrintf(viewer," ||Gfd|| %g, ||G|| = %g, angle cosine = (Gfd'G)/||Gfd||||G|| = %g\n", (double)fdnorm, (double)hcnorm, (double)(PetscRealPart(dot)/(fdnorm*hcnorm)));CHKERRQ(ierr); 81 ierr = PetscViewerASCIIPrintf(viewer," 2-norm ||G - Gfd||/||G|| = %g, ||G - Gfd|| = %g\n",(double)(diffnorm/PetscMax(hcnorm,fdnorm)),(double)diffnorm);CHKERRQ(ierr); 82 ierr = PetscViewerASCIIPrintf(viewer," max-norm ||G - Gfd||/||G|| = %g, ||G - Gfd|| = %g\n",(double)(diffmax/PetscMax(hcmax,fdmax)),(double)diffmax);CHKERRQ(ierr); 83 84 if (complete_print) { 85 ierr = PetscViewerASCIIPrintf(viewer," Hand-coded gradient ----------\n");CHKERRQ(ierr); 86 ierr = VecView(g1,viewer);CHKERRQ(ierr); 87 ierr = PetscViewerASCIIPrintf(viewer," Finite difference gradient ----------\n");CHKERRQ(ierr); 88 ierr = VecView(g2,viewer);CHKERRQ(ierr); 89 ierr = PetscViewerASCIIPrintf(viewer," Hand-coded minus finite-difference gradient ----------\n");CHKERRQ(ierr); 90 ierr = VecView(g3,viewer);CHKERRQ(ierr); 91 } 92 ierr = VecDestroy(&g2);CHKERRQ(ierr); 93 ierr = VecDestroy(&g3);CHKERRQ(ierr); 94 PetscFunctionReturn(0); 95 } 96 97 /*@ 98 TaoComputeGradient - Computes the gradient of the objective function 99 100 Collective on Tao 101 102 Input Parameters: 103 + tao - the Tao context 104 - X - input vector 105 106 Output Parameter: 107 . G - gradient vector 108 109 Options Database Keys: 110 + -tao_test_gradient - compare the user provided gradient with one compute via finite differences to check for errors 111 - -tao_test_gradient_display - display the user provided gradient, the finite difference gradient and the difference between them to help users detect the location of errors in the user provided gradient 112 113 Notes: 114 TaoComputeGradient() is typically used within minimization implementations, 115 so most users would not generally call this routine themselves. 116 117 Level: advanced 118 119 .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetGradientRoutine() 120 @*/ 121 PetscErrorCode TaoComputeGradient(Tao tao, Vec X, Vec G) 122 { 123 PetscErrorCode ierr; 124 PetscReal dummy; 125 126 PetscFunctionBegin; 127 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 128 PetscValidHeaderSpecific(X,VEC_CLASSID,2); 129 PetscValidHeaderSpecific(G,VEC_CLASSID,2); 130 PetscCheckSameComm(tao,1,X,2); 131 PetscCheckSameComm(tao,1,G,3); 132 ierr = VecLockPush(X);CHKERRQ(ierr); 133 if (tao->ops->computegradient) { 134 ierr = PetscLogEventBegin(Tao_GradientEval,tao,X,G,NULL);CHKERRQ(ierr); 135 PetscStackPush("Tao user gradient evaluation routine"); 136 ierr = (*tao->ops->computegradient)(tao,X,G,tao->user_gradP);CHKERRQ(ierr); 137 PetscStackPop; 138 ierr = PetscLogEventEnd(Tao_GradientEval,tao,X,G,NULL);CHKERRQ(ierr); 139 tao->ngrads++; 140 } else if (tao->ops->computeobjectiveandgradient) { 141 ierr = PetscLogEventBegin(Tao_ObjGradientEval,tao,X,G,NULL);CHKERRQ(ierr); 142 PetscStackPush("Tao user objective/gradient evaluation routine"); 143 ierr = (*tao->ops->computeobjectiveandgradient)(tao,X,&dummy,G,tao->user_objgradP);CHKERRQ(ierr); 144 PetscStackPop; 145 ierr = PetscLogEventEnd(Tao_ObjGradientEval,tao,X,G,NULL);CHKERRQ(ierr); 146 tao->nfuncgrads++; 147 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"TaoSetGradientRoutine() has not been called"); 148 ierr = VecLockPop(X);CHKERRQ(ierr); 149 150 ierr = TaoTestGradient(tao,G);CHKERRQ(ierr); 151 PetscFunctionReturn(0); 152 } 153 154 /*@ 155 TaoComputeObjective - Computes the objective function value at a given point 156 157 Collective on Tao 158 159 Input Parameters: 160 + tao - the Tao context 161 - X - input vector 162 163 Output Parameter: 164 . f - Objective value at X 165 166 Notes: 167 TaoComputeObjective() is typically used within minimization implementations, 168 so most users would not generally call this routine themselves. 169 170 Level: advanced 171 172 .seealso: TaoComputeGradient(), TaoComputeObjectiveAndGradient(), TaoSetObjectiveRoutine() 173 @*/ 174 PetscErrorCode TaoComputeObjective(Tao tao, Vec X, PetscReal *f) 175 { 176 PetscErrorCode ierr; 177 Vec temp; 178 179 PetscFunctionBegin; 180 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 181 PetscValidHeaderSpecific(X,VEC_CLASSID,2); 182 PetscCheckSameComm(tao,1,X,2); 183 ierr = VecLockPush(X);CHKERRQ(ierr); 184 if (tao->ops->computeobjective) { 185 ierr = PetscLogEventBegin(Tao_ObjectiveEval,tao,X,NULL,NULL);CHKERRQ(ierr); 186 PetscStackPush("Tao user objective evaluation routine"); 187 ierr = (*tao->ops->computeobjective)(tao,X,f,tao->user_objP);CHKERRQ(ierr); 188 PetscStackPop; 189 ierr = PetscLogEventEnd(Tao_ObjectiveEval,tao,X,NULL,NULL);CHKERRQ(ierr); 190 tao->nfuncs++; 191 } else if (tao->ops->computeobjectiveandgradient) { 192 ierr = PetscInfo(tao,"Duplicating variable vector in order to call func/grad routine\n");CHKERRQ(ierr); 193 ierr = VecDuplicate(X,&temp);CHKERRQ(ierr); 194 ierr = PetscLogEventBegin(Tao_ObjGradientEval,tao,X,NULL,NULL);CHKERRQ(ierr); 195 PetscStackPush("Tao user objective/gradient evaluation routine"); 196 ierr = (*tao->ops->computeobjectiveandgradient)(tao,X,f,temp,tao->user_objgradP);CHKERRQ(ierr); 197 PetscStackPop; 198 ierr = PetscLogEventEnd(Tao_ObjGradientEval,tao,X,NULL,NULL);CHKERRQ(ierr); 199 ierr = VecDestroy(&temp);CHKERRQ(ierr); 200 tao->nfuncgrads++; 201 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"TaoSetObjectiveRoutine() has not been called"); 202 ierr = PetscInfo1(tao,"TAO Function evaluation: %20.19e\n",(double)(*f));CHKERRQ(ierr); 203 ierr = VecLockPop(X);CHKERRQ(ierr); 204 PetscFunctionReturn(0); 205 } 206 207 /*@ 208 TaoComputeObjectiveAndGradient - Computes the objective function value at a given point 209 210 Collective on Tao 211 212 Input Parameters: 213 + tao - the Tao context 214 - X - input vector 215 216 Output Parameter: 217 + f - Objective value at X 218 - g - Gradient vector at X 219 220 Notes: 221 TaoComputeObjectiveAndGradient() is typically used within minimization implementations, 222 so most users would not generally call this routine themselves. 223 224 Level: advanced 225 226 .seealso: TaoComputeGradient(), TaoComputeObjectiveAndGradient(), TaoSetObjectiveRoutine() 227 @*/ 228 PetscErrorCode TaoComputeObjectiveAndGradient(Tao tao, Vec X, PetscReal *f, Vec G) 229 { 230 PetscErrorCode ierr; 231 232 PetscFunctionBegin; 233 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 234 PetscValidHeaderSpecific(X,VEC_CLASSID,2); 235 PetscValidHeaderSpecific(G,VEC_CLASSID,4); 236 PetscCheckSameComm(tao,1,X,2); 237 PetscCheckSameComm(tao,1,G,4); 238 ierr = VecLockPush(X);CHKERRQ(ierr); 239 if (tao->ops->computeobjectiveandgradient) { 240 ierr = PetscLogEventBegin(Tao_ObjGradientEval,tao,X,G,NULL);CHKERRQ(ierr); 241 if (tao->ops->computegradient == TaoDefaultComputeGradient) { 242 ierr = TaoComputeObjective(tao,X,f);CHKERRQ(ierr); 243 ierr = TaoDefaultComputeGradient(tao,X,G,NULL);CHKERRQ(ierr); 244 } else { 245 PetscStackPush("Tao user objective/gradient evaluation routine"); 246 ierr = (*tao->ops->computeobjectiveandgradient)(tao,X,f,G,tao->user_objgradP);CHKERRQ(ierr); 247 PetscStackPop; 248 } 249 ierr = PetscLogEventEnd(Tao_ObjGradientEval,tao,X,G,NULL);CHKERRQ(ierr); 250 tao->nfuncgrads++; 251 } else if (tao->ops->computeobjective && tao->ops->computegradient) { 252 ierr = PetscLogEventBegin(Tao_ObjectiveEval,tao,X,NULL,NULL);CHKERRQ(ierr); 253 PetscStackPush("Tao user objective evaluation routine"); 254 ierr = (*tao->ops->computeobjective)(tao,X,f,tao->user_objP);CHKERRQ(ierr); 255 PetscStackPop; 256 ierr = PetscLogEventEnd(Tao_ObjectiveEval,tao,X,NULL,NULL);CHKERRQ(ierr); 257 tao->nfuncs++; 258 ierr = PetscLogEventBegin(Tao_GradientEval,tao,X,G,NULL);CHKERRQ(ierr); 259 PetscStackPush("Tao user gradient evaluation routine"); 260 ierr = (*tao->ops->computegradient)(tao,X,G,tao->user_gradP);CHKERRQ(ierr); 261 PetscStackPop; 262 ierr = PetscLogEventEnd(Tao_GradientEval,tao,X,G,NULL);CHKERRQ(ierr); 263 tao->ngrads++; 264 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"TaoSetObjectiveRoutine() or TaoSetGradientRoutine() not set"); 265 ierr = PetscInfo1(tao,"TAO Function evaluation: %20.19e\n",(double)(*f));CHKERRQ(ierr); 266 ierr = VecLockPop(X);CHKERRQ(ierr); 267 PetscFunctionReturn(0); 268 } 269 270 /*@C 271 TaoSetObjectiveRoutine - Sets the function evaluation routine for minimization 272 273 Logically collective on Tao 274 275 Input Parameter: 276 + tao - the Tao context 277 . func - the objective function 278 - ctx - [optional] user-defined context for private data for the function evaluation 279 routine (may be NULL) 280 281 Calling sequence of func: 282 $ func (Tao tao, Vec x, PetscReal *f, void *ctx); 283 284 + x - input vector 285 . f - function value 286 - ctx - [optional] user-defined function context 287 288 Level: beginner 289 290 .seealso: TaoSetGradientRoutine(), TaoSetHessianRoutine() TaoSetObjectiveAndGradientRoutine() 291 @*/ 292 PetscErrorCode TaoSetObjectiveRoutine(Tao tao, PetscErrorCode (*func)(Tao, Vec, PetscReal*,void*),void *ctx) 293 { 294 PetscFunctionBegin; 295 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 296 tao->user_objP = ctx; 297 tao->ops->computeobjective = func; 298 PetscFunctionReturn(0); 299 } 300 301 /*@C 302 TaoSetSeparableObjectiveRoutine - Sets the function evaluation routine for least-square applications 303 304 Logically collective on Tao 305 306 Input Parameter: 307 + tao - the Tao context 308 . func - the objective function evaluation routine 309 - ctx - [optional] user-defined context for private data for the function evaluation 310 routine (may be NULL) 311 312 Calling sequence of func: 313 $ func (Tao tao, Vec x, Vec f, void *ctx); 314 315 + x - input vector 316 . f - function value vector 317 - ctx - [optional] user-defined function context 318 319 Level: beginner 320 321 .seealso: TaoSetObjectiveRoutine(), TaoSetJacobianRoutine() 322 @*/ 323 PetscErrorCode TaoSetSeparableObjectiveRoutine(Tao tao, Vec sepobj, PetscErrorCode (*func)(Tao, Vec, Vec, void*),void *ctx) 324 { 325 PetscFunctionBegin; 326 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 327 PetscValidHeaderSpecific(sepobj, VEC_CLASSID,2); 328 tao->user_sepobjP = ctx; 329 tao->sep_objective = sepobj; 330 tao->ops->computeseparableobjective = func; 331 PetscFunctionReturn(0); 332 } 333 334 /*@ 335 TaoSetSeparableObjectiveWeights - Give weights for the separable objective values. A vector can be used if only diagonal terms are used, otherwise a matrix can be give. If this function is not used, or if sigma_v and sigma_w are both NULL, then the default identity matrix will be used for weights. 336 337 Collective on Tao 338 339 Input Parameters: 340 + tao - the Tao context 341 . sigma_v - vector of weights (diagonal terms only) 342 . n - the number of weights (if using off-diagonal) 343 . rows - index list of rows for sigma_w 344 . cols - index list of columns for sigma_w 345 - vals - array of weights 346 347 348 349 Note: Either sigma_v or sigma_w (or both) should be NULL 350 351 Level: intermediate 352 353 .seealso: TaoSetSeparableObjectiveRoutine() 354 @*/ 355 PetscErrorCode TaoSetSeparableObjectiveWeights(Tao tao, Vec sigma_v, PetscInt n, PetscInt *rows, PetscInt *cols, PetscReal *vals) 356 { 357 PetscErrorCode ierr; 358 PetscInt i; 359 PetscFunctionBegin; 360 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 361 ierr = VecDestroy(&tao->sep_weights_v);CHKERRQ(ierr); 362 tao->sep_weights_v=sigma_v; 363 if (sigma_v) { 364 ierr = PetscObjectReference((PetscObject)sigma_v);CHKERRQ(ierr); 365 } 366 if (vals) { 367 if (tao->sep_weights_n) { 368 ierr = PetscFree(tao->sep_weights_rows);CHKERRQ(ierr); 369 ierr = PetscFree(tao->sep_weights_cols);CHKERRQ(ierr); 370 ierr = PetscFree(tao->sep_weights_w);CHKERRQ(ierr); 371 } 372 ierr = PetscMalloc1(n,&tao->sep_weights_rows);CHKERRQ(ierr); 373 ierr = PetscMalloc1(n,&tao->sep_weights_cols);CHKERRQ(ierr); 374 ierr = PetscMalloc1(n,&tao->sep_weights_w);CHKERRQ(ierr); 375 tao->sep_weights_n=n; 376 for (i=0;i<n;i++) { 377 tao->sep_weights_rows[i]=rows[i]; 378 tao->sep_weights_cols[i]=cols[i]; 379 tao->sep_weights_w[i]=vals[i]; 380 } 381 } else { 382 tao->sep_weights_n=0; 383 tao->sep_weights_rows=0; 384 tao->sep_weights_cols=0; 385 } 386 PetscFunctionReturn(0); 387 } 388 /*@ 389 TaoComputeSeparableObjective - Computes a separable objective function vector at a given point (for least-square applications) 390 391 Collective on Tao 392 393 Input Parameters: 394 + tao - the Tao context 395 - X - input vector 396 397 Output Parameter: 398 . f - Objective vector at X 399 400 Notes: 401 TaoComputeSeparableObjective() is typically used within minimization implementations, 402 so most users would not generally call this routine themselves. 403 404 Level: advanced 405 406 .seealso: TaoSetSeparableObjectiveRoutine() 407 @*/ 408 PetscErrorCode TaoComputeSeparableObjective(Tao tao, Vec X, Vec F) 409 { 410 PetscErrorCode ierr; 411 412 PetscFunctionBegin; 413 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 414 PetscValidHeaderSpecific(X,VEC_CLASSID,2); 415 PetscValidHeaderSpecific(F,VEC_CLASSID,3); 416 PetscCheckSameComm(tao,1,X,2); 417 PetscCheckSameComm(tao,1,F,3); 418 if (tao->ops->computeseparableobjective) { 419 ierr = PetscLogEventBegin(Tao_ObjectiveEval,tao,X,NULL,NULL);CHKERRQ(ierr); 420 PetscStackPush("Tao user separable objective evaluation routine"); 421 ierr = (*tao->ops->computeseparableobjective)(tao,X,F,tao->user_sepobjP);CHKERRQ(ierr); 422 PetscStackPop; 423 ierr = PetscLogEventEnd(Tao_ObjectiveEval,tao,X,NULL,NULL);CHKERRQ(ierr); 424 tao->nfuncs++; 425 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"TaoSetSeparableObjectiveRoutine() has not been called"); 426 ierr = PetscInfo(tao,"TAO separable function evaluation.\n");CHKERRQ(ierr); 427 PetscFunctionReturn(0); 428 } 429 430 /*@C 431 TaoSetGradientRoutine - Sets the gradient evaluation routine for minimization 432 433 Logically collective on Tao 434 435 Input Parameter: 436 + tao - the Tao context 437 . func - the gradient function 438 - ctx - [optional] user-defined context for private data for the gradient evaluation 439 routine (may be NULL) 440 441 Calling sequence of func: 442 $ func (Tao tao, Vec x, Vec g, void *ctx); 443 444 + x - input vector 445 . g - gradient value (output) 446 - ctx - [optional] user-defined function context 447 448 Level: beginner 449 450 .seealso: TaoSetObjectiveRoutine(), TaoSetHessianRoutine() TaoSetObjectiveAndGradientRoutine() 451 @*/ 452 PetscErrorCode TaoSetGradientRoutine(Tao tao, PetscErrorCode (*func)(Tao, Vec, Vec, void*),void *ctx) 453 { 454 PetscFunctionBegin; 455 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 456 tao->user_gradP = ctx; 457 tao->ops->computegradient = func; 458 PetscFunctionReturn(0); 459 } 460 461 /*@C 462 TaoSetObjectiveAndGradientRoutine - Sets a combined objective function and gradient evaluation routine for minimization 463 464 Logically collective on Tao 465 466 Input Parameter: 467 + tao - the Tao context 468 . func - the gradient function 469 - ctx - [optional] user-defined context for private data for the gradient evaluation 470 routine (may be NULL) 471 472 Calling sequence of func: 473 $ func (Tao tao, Vec x, PetscReal *f, Vec g, void *ctx); 474 475 + x - input vector 476 . f - objective value (output) 477 . g - gradient value (output) 478 - ctx - [optional] user-defined function context 479 480 Level: beginner 481 482 .seealso: TaoSetObjectiveRoutine(), TaoSetHessianRoutine() TaoSetObjectiveAndGradientRoutine() 483 @*/ 484 PetscErrorCode TaoSetObjectiveAndGradientRoutine(Tao tao, PetscErrorCode (*func)(Tao, Vec, PetscReal *, Vec, void*), void *ctx) 485 { 486 PetscFunctionBegin; 487 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 488 tao->user_objgradP = ctx; 489 tao->ops->computeobjectiveandgradient = func; 490 PetscFunctionReturn(0); 491 } 492 493 /*@ 494 TaoIsObjectiveDefined -- Checks to see if the user has 495 declared an objective-only routine. Useful for determining when 496 it is appropriate to call TaoComputeObjective() or 497 TaoComputeObjectiveAndGradient() 498 499 Collective on Tao 500 501 Input Parameter: 502 + tao - the Tao context 503 - ctx - PETSC_TRUE if objective function routine is set by user, 504 PETSC_FALSE otherwise 505 Level: developer 506 507 .seealso: TaoSetObjectiveRoutine(), TaoIsGradientDefined(), TaoIsObjectiveAndGradientDefined() 508 @*/ 509 PetscErrorCode TaoIsObjectiveDefined(Tao tao, PetscBool *flg) 510 { 511 PetscFunctionBegin; 512 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 513 if (tao->ops->computeobjective == 0) *flg = PETSC_FALSE; 514 else *flg = PETSC_TRUE; 515 PetscFunctionReturn(0); 516 } 517 518 /*@ 519 TaoIsGradientDefined -- Checks to see if the user has 520 declared an objective-only routine. Useful for determining when 521 it is appropriate to call TaoComputeGradient() or 522 TaoComputeGradientAndGradient() 523 524 Not Collective 525 526 Input Parameter: 527 + tao - the Tao context 528 - ctx - PETSC_TRUE if gradient routine is set by user, PETSC_FALSE otherwise 529 Level: developer 530 531 .seealso: TaoSetGradientRoutine(), TaoIsObjectiveDefined(), TaoIsObjectiveAndGradientDefined() 532 @*/ 533 PetscErrorCode TaoIsGradientDefined(Tao tao, PetscBool *flg) 534 { 535 PetscFunctionBegin; 536 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 537 if (tao->ops->computegradient == 0) *flg = PETSC_FALSE; 538 else *flg = PETSC_TRUE; 539 PetscFunctionReturn(0); 540 } 541 542 /*@ 543 TaoIsObjectiveAndGradientDefined -- Checks to see if the user has 544 declared a joint objective/gradient routine. Useful for determining when 545 it is appropriate to call TaoComputeObjective() or 546 TaoComputeObjectiveAndGradient() 547 548 Not Collective 549 550 Input Parameter: 551 + tao - the Tao context 552 - ctx - PETSC_TRUE if objective/gradient routine is set by user, PETSC_FALSE otherwise 553 Level: developer 554 555 .seealso: TaoSetObjectiveAndGradientRoutine(), TaoIsObjectiveDefined(), TaoIsGradientDefined() 556 @*/ 557 PetscErrorCode TaoIsObjectiveAndGradientDefined(Tao tao, PetscBool *flg) 558 { 559 PetscFunctionBegin; 560 PetscValidHeaderSpecific(tao,TAO_CLASSID,1); 561 if (tao->ops->computeobjectiveandgradient == 0) *flg = PETSC_FALSE; 562 else *flg = PETSC_TRUE; 563 PetscFunctionReturn(0); 564 } 565