Actual source code: taosolver_hj.c

petsc-3.10.5 2019-03-28
Report Typos and Errors
  1:  #include <petsc/private/taoimpl.h>

  3: /*@C
  4:    TaoSetHessianRoutine - Sets the function to compute the Hessian as well as the location to store the matrix.

  6:    Logically collective on Tao

  8:    Input Parameters:
  9: +  tao  - the Tao context
 10: .  H    - Matrix used for the hessian
 11: .  Hpre - Matrix that will be used operated on by preconditioner, can be same as H
 12: .  func - Hessian evaluation routine
 13: -  ctx  - [optional] user-defined context for private data for the
 14:          Hessian evaluation routine (may be NULL)

 16:    Calling sequence of func:
 17: $    func(Tao tao,Vec x,Mat H,Mat Hpre,void *ctx);

 19: +  tao  - the Tao  context
 20: .  x    - input vector
 21: .  H    - Hessian matrix
 22: .  Hpre - preconditioner matrix, usually the same as H
 23: -  ctx  - [optional] user-defined Hessian context

 25:    Level: beginner
 26: @*/
 27: PetscErrorCode TaoSetHessianRoutine(Tao tao, Mat H, Mat Hpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, void*), void *ctx)
 28: {

 33:   if (H) {
 36:   }
 37:   if (Hpre) {
 40:   }
 41:   if (ctx) {
 42:     tao->user_hessP = ctx;
 43:   }
 44:   if (func) {
 45:     tao->ops->computehessian = func;
 46:   }
 47:   if (H) {
 48:     PetscObjectReference((PetscObject)H);
 49:     MatDestroy(&tao->hessian);
 50:     tao->hessian = H;
 51:   }
 52:   if (Hpre) {
 53:     PetscObjectReference((PetscObject)Hpre);
 54:     MatDestroy(&tao->hessian_pre);
 55:     tao->hessian_pre = Hpre;
 56:   }
 57:   return(0);
 58: }

 60: PetscErrorCode TaoTestHessian(Tao tao)
 61: {
 62:   Mat               A,B,C,D,hessian;
 63:   Vec               x = tao->solution;
 64:   PetscErrorCode    ierr;
 65:   PetscReal         nrm,gnorm;
 66:   PetscReal         threshold = 1.e-5;
 67:   PetscInt          m,n,M,N;
 68:   PetscBool         complete_print = PETSC_FALSE,test = PETSC_FALSE,flg;
 69:   PetscViewer       viewer,mviewer;
 70:   MPI_Comm          comm;
 71:   PetscInt          tabs;
 72:   static PetscBool  directionsprinted = PETSC_FALSE;
 73:   PetscViewerFormat format;

 76:   PetscObjectOptionsBegin((PetscObject)tao);
 77:   PetscOptionsName("-tao_test_hessian","Compare hand-coded and finite difference Hessians","None",&test);
 78:   PetscOptionsReal("-tao_test_hessian", "Threshold for element difference between hand-coded and finite difference being meaningful","None",threshold,&threshold,NULL);
 79:   PetscOptionsViewer("-tao_test_hessian_view","View difference between hand-coded and finite difference Hessians element entries","None",&mviewer,&format,&complete_print);
 80:   PetscOptionsEnd();
 81:   if (!test) return(0);

 83:   PetscObjectGetComm((PetscObject)tao,&comm);
 84:   PetscViewerASCIIGetStdout(comm,&viewer);
 85:   PetscViewerASCIIGetTab(viewer, &tabs);
 86:   PetscViewerASCIISetTab(viewer, ((PetscObject)tao)->tablevel);
 87:   PetscViewerASCIIPrintf(viewer,"  ---------- Testing Hessian -------------\n");
 88:   if (!complete_print && !directionsprinted) {
 89:     PetscViewerASCIIPrintf(viewer,"  Run with -tao_test_hessian_view and optionally -tao_test_hessian <threshold> to show difference\n");
 90:     PetscViewerASCIIPrintf(viewer,"    of hand-coded and finite difference Hessian entries greater than <threshold>.\n");
 91:   }
 92:   if (!directionsprinted) {
 93:     PetscViewerASCIIPrintf(viewer,"  Testing hand-coded Hessian, if (for double precision runs) ||J - Jfd||_F/||J||_F is\n");
 94:     PetscViewerASCIIPrintf(viewer,"    O(1.e-8), the hand-coded Hessian is probably correct.\n");
 95:     directionsprinted = PETSC_TRUE;
 96:   }
 97:   if (complete_print) {
 98:     PetscViewerPushFormat(mviewer,format);
 99:   }

101:   PetscObjectTypeCompare((PetscObject)tao->hessian,MATMFFD,&flg);
102:   if (!flg) hessian = tao->hessian;
103:   else hessian = tao->hessian_pre;

105:   while (hessian) {
106:     PetscObjectTypeCompareAny((PetscObject)hessian,&flg,MATSEQAIJ,MATMPIAIJ,MATSEQDENSE,MATMPIDENSE,MATSEQBAIJ,MATMPIBAIJ,MATSEQSBAIJ,MATMPIBAIJ,"");
107:     if (flg) {
108:       A    = hessian;
109:       PetscObjectReference((PetscObject)A);
110:     } else {
111:       MatComputeExplicitOperator(hessian,&A);
112:     }

114:     MatCreate(PetscObjectComm((PetscObject)A),&B);
115:     MatGetSize(A,&M,&N);
116:     MatGetLocalSize(A,&m,&n);
117:     MatSetSizes(B,m,n,M,N);
118:     MatSetType(B,((PetscObject)A)->type_name);
119:     MatSetUp(B);
120:     MatSetOption(B,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);

122:     TaoDefaultComputeHessian(tao,x,B,B,NULL);

124:     MatDuplicate(B,MAT_COPY_VALUES,&D);
125:     MatAYPX(D,-1.0,A,DIFFERENT_NONZERO_PATTERN);
126:     MatNorm(D,NORM_FROBENIUS,&nrm);
127:     MatNorm(A,NORM_FROBENIUS,&gnorm);
128:     MatDestroy(&D);
129:     if (!gnorm) gnorm = 1; /* just in case */
130:     PetscViewerASCIIPrintf(viewer,"  ||H - Hfd||_F/||H||_F = %g, ||H - Hfd||_F = %g\n",(double)(nrm/gnorm),(double)nrm);

132:     if (complete_print) {
133:       PetscViewerASCIIPrintf(viewer,"  Hand-coded Hessian ----------\n");
134:       MatView(hessian,mviewer);
135:       PetscViewerASCIIPrintf(viewer,"  Finite difference Hessian ----------\n");
136:       MatView(B,mviewer);
137:     }

139:     if (complete_print) {
140:       PetscInt          Istart, Iend, *ccols, bncols, cncols, j, row;
141:       PetscScalar       *cvals;
142:       const PetscInt    *bcols;
143:       const PetscScalar *bvals;

145:       MatAYPX(B,-1.0,A,DIFFERENT_NONZERO_PATTERN);
146:       MatCreate(PetscObjectComm((PetscObject)A),&C);
147:       MatSetSizes(C,m,n,M,N);
148:       MatSetType(C,((PetscObject)A)->type_name);
149:       MatSetUp(C);
150:       MatSetOption(C,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
151:       MatGetOwnershipRange(B,&Istart,&Iend);

153:       for (row = Istart; row < Iend; row++) {
154:         MatGetRow(B,row,&bncols,&bcols,&bvals);
155:         PetscMalloc2(bncols,&ccols,bncols,&cvals);
156:         for (j = 0, cncols = 0; j < bncols; j++) {
157:           if (PetscAbsScalar(bvals[j]) > threshold) {
158:             ccols[cncols] = bcols[j];
159:             cvals[cncols] = bvals[j];
160:             cncols += 1;
161:           }
162:         }
163:         if (cncols) {
164:           MatSetValues(C,1,&row,cncols,ccols,cvals,INSERT_VALUES);
165:         }
166:         MatRestoreRow(B,row,&bncols,&bcols,&bvals);
167:         PetscFree2(ccols,cvals);
168:       }
169:       MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);
170:       MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);
171:       PetscViewerASCIIPrintf(viewer,"  Finite-difference minus hand-coded Hessian with tolerance %g ----------\n",(double)threshold);
172:       MatView(C,mviewer);
173:       MatDestroy(&C);
174:     }
175:     MatDestroy(&A);
176:     MatDestroy(&B);

178:     if (hessian != tao->hessian_pre) {
179:       hessian = tao->hessian_pre;
180:       PetscViewerASCIIPrintf(viewer,"  ---------- Testing Hessian for preconditioner -------------\n");
181:     } else hessian = NULL;
182:   }
183:   if (complete_print) {
184:     PetscViewerPopFormat(mviewer);
185:   }
186:   PetscViewerASCIISetTab(viewer,tabs);
187:   return(0);
188: }

190: /*@C
191:    TaoComputeHessian - Computes the Hessian matrix that has been
192:    set with TaoSetHessianRoutine().

194:    Collective on Tao

196:    Input Parameters:
197: +  tao - the Tao solver context
198: -  X   - input vector

200:    Output Parameters:
201: +  H    - Hessian matrix
202: -  Hpre - Preconditioning matrix

204:    Options Database Keys:
205: +     -tao_test_hessian - compare the user provided Hessian with one compute via finite differences to check for errors
206: .     -tao_test_hessian <numerical value>  - display entries in the difference between the user provided Hessian and finite difference Hessian that are greater than a certain value to help users detect errors
207: -     -tao_test_hessian_view - display the user provided Hessian, the finite difference Hessian and the difference between them to help users detect the location of errors in the user provided Hessian

209:    Notes:
210:    Most users should not need to explicitly call this routine, as it
211:    is used internally within the minimization solvers.

213:    TaoComputeHessian() is typically used within minimization
214:    implementations, so most users would not generally call this routine
215:    themselves.

217:    Developer Notes:
218:    The Hessian test mechanism follows SNESTestJacobian().

220:    Level: developer

222: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetHessianRoutine()
223: @*/
224: PetscErrorCode TaoComputeHessian(Tao tao, Vec X, Mat H, Mat Hpre)
225: {

232:   if (!tao->ops->computehessian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetHessianRoutine() first");
233:   ++tao->nhess;
234:   VecLockPush(X);
235:   PetscLogEventBegin(TAO_HessianEval,tao,X,H,Hpre);
236:   PetscStackPush("Tao user Hessian function");
237:   (*tao->ops->computehessian)(tao,X,H,Hpre,tao->user_hessP);
238:   PetscStackPop;
239:   PetscLogEventEnd(TAO_HessianEval,tao,X,H,Hpre);
240:   VecLockPop(X);

242:   TaoTestHessian(tao);
243:   return(0);
244: }

246: /*@C
247:    TaoComputeJacobian - Computes the Jacobian matrix that has been
248:    set with TaoSetJacobianRoutine().

250:    Collective on Tao

252:    Input Parameters:
253: +  tao - the Tao solver context
254: -  X   - input vector

256:    Output Parameters:
257: +  J    - Jacobian matrix
258: -  Jpre - Preconditioning matrix

260:    Notes:
261:    Most users should not need to explicitly call this routine, as it
262:    is used internally within the minimization solvers.

264:    TaoComputeJacobian() is typically used within minimization
265:    implementations, so most users would not generally call this routine
266:    themselves.

268:    Level: developer

270: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianRoutine()
271: @*/
272: PetscErrorCode TaoComputeJacobian(Tao tao, Vec X, Mat J, Mat Jpre)
273: {

280:   if (!tao->ops->computejacobian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobian() first");
281:   ++tao->njac;
282:   VecLockPush(X);
283:   PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
284:   PetscStackPush("Tao user Jacobian function");
285:   (*tao->ops->computejacobian)(tao,X,J,Jpre,tao->user_jacP);
286:   PetscStackPop;
287:   PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
288:   VecLockPop(X);
289:   return(0);
290: }

292: /*@C
293:    TaoComputeJacobianState - Computes the Jacobian matrix that has been
294:    set with TaoSetJacobianStateRoutine().

296:    Collective on Tao

298:    Input Parameters:
299: +  tao - the Tao solver context
300: -  X   - input vector

302:    Output Parameters:
303: +  Jpre - Jacobian matrix
304: -  Jinv - Preconditioning matrix

306:    Notes:
307:    Most users should not need to explicitly call this routine, as it
308:    is used internally within the minimization solvers.

310:    TaoComputeJacobianState() is typically used within minimization
311:    implementations, so most users would not generally call this routine
312:    themselves.

314:    Level: developer

316: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
317: @*/
318: PetscErrorCode TaoComputeJacobianState(Tao tao, Vec X, Mat J, Mat Jpre, Mat Jinv)
319: {

326:   if (!tao->ops->computejacobianstate) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianState() first");
327:   ++tao->njac_state;
328:   VecLockPush(X);
329:   PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
330:   PetscStackPush("Tao user Jacobian(state) function");
331:   (*tao->ops->computejacobianstate)(tao,X,J,Jpre,Jinv,tao->user_jac_stateP);
332:   PetscStackPop;
333:   PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
334:   VecLockPop(X);
335:   return(0);
336: }

338: /*@C
339:    TaoComputeJacobianDesign - Computes the Jacobian matrix that has been
340:    set with TaoSetJacobianDesignRoutine().

342:    Collective on Tao

344:    Input Parameters:
345: +  tao - the Tao solver context
346: -  X   - input vector

348:    Output Parameters:
349: .  J - Jacobian matrix

351:    Notes:
352:    Most users should not need to explicitly call this routine, as it
353:    is used internally within the minimization solvers.

355:    TaoComputeJacobianDesign() is typically used within minimization
356:    implementations, so most users would not generally call this routine
357:    themselves.

359:    Level: developer

361: .seealso:  TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianDesignRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
362: @*/
363: PetscErrorCode TaoComputeJacobianDesign(Tao tao, Vec X, Mat J)
364: {

371:   if (!tao->ops->computejacobiandesign) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianDesign() first");
372:   ++tao->njac_design;
373:   VecLockPush(X);
374:   PetscLogEventBegin(TAO_JacobianEval,tao,X,J,NULL);
375:   PetscStackPush("Tao user Jacobian(design) function");
376:   (*tao->ops->computejacobiandesign)(tao,X,J,tao->user_jac_designP);
377:   PetscStackPop;
378:   PetscLogEventEnd(TAO_JacobianEval,tao,X,J,NULL);
379:   VecLockPop(X);
380:   return(0);
381: }

383: /*@C
384:    TaoSetJacobianRoutine - Sets the function to compute the Jacobian as well as the location to store the matrix.

386:    Logically collective on Tao

388:    Input Parameters:
389: +  tao  - the Tao context
390: .  J    - Matrix used for the jacobian
391: .  Jpre - Matrix that will be used operated on by preconditioner, can be same as J
392: .  func - Jacobian evaluation routine
393: -  ctx  - [optional] user-defined context for private data for the
394:           Jacobian evaluation routine (may be NULL)

396:    Calling sequence of func:
397: $    func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);

399: +  tao  - the Tao  context
400: .  x    - input vector
401: .  J    - Jacobian matrix
402: .  Jpre - preconditioning matrix, usually the same as J
403: -  ctx  - [optional] user-defined Jacobian context

405:    Level: intermediate
406: @*/
407: PetscErrorCode TaoSetJacobianRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, void*), void *ctx)
408: {

413:   if (J) {
416:   }
417:   if (Jpre) {
420:   }
421:   if (ctx) {
422:     tao->user_jacP = ctx;
423:   }
424:   if (func) {
425:     tao->ops->computejacobian = func;
426:   }
427:   if (J) {
428:     PetscObjectReference((PetscObject)J);
429:     MatDestroy(&tao->jacobian);
430:     tao->jacobian = J;
431:   }
432:   if (Jpre) {
433:     PetscObjectReference((PetscObject)Jpre);
434:     MatDestroy(&tao->jacobian_pre);
435:     tao->jacobian_pre=Jpre;
436:   }
437:   return(0);
438: }

440: /*@C
441:    TaoSetJacobianStateRoutine - Sets the function to compute the Jacobian
442:    (and its inverse) of the constraint function with respect to the state variables.
443:    Used only for pde-constrained optimization.

445:    Logically collective on Tao

447:    Input Parameters:
448: +  tao  - the Tao context
449: .  J    - Matrix used for the jacobian
450: .  Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.  Only used if Jinv is NULL
451: .  Jinv - [optional] Matrix used to apply the inverse of the state jacobian. Use NULL to default to PETSc KSP solvers to apply the inverse.
452: .  func - Jacobian evaluation routine
453: -  ctx  - [optional] user-defined context for private data for the
454:           Jacobian evaluation routine (may be NULL)

456:    Calling sequence of func:
457: $    func(Tao tao,Vec x,Mat J,Mat Jpre,Mat Jinv,void *ctx);

459: +  tao  - the Tao  context
460: .  x    - input vector
461: .  J    - Jacobian matrix
462: .  Jpre - preconditioner matrix, usually the same as J
463: .  Jinv - inverse of J
464: -  ctx  - [optional] user-defined Jacobian context

466:    Level: intermediate
467: .seealso: TaoComputeJacobianState(), TaoSetJacobianDesignRoutine(), TaoSetStateDesignIS()
468: @*/
469: PetscErrorCode TaoSetJacobianStateRoutine(Tao tao, Mat J, Mat Jpre, Mat Jinv, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, Mat, void*), void *ctx)
470: {

475:   if (J) {
478:   }
479:   if (Jpre) {
482:   }
483:   if (Jinv) {
486:   }
487:   if (ctx) {
488:     tao->user_jac_stateP = ctx;
489:   }
490:   if (func) {
491:     tao->ops->computejacobianstate = func;
492:   }
493:   if (J) {
494:     PetscObjectReference((PetscObject)J);
495:     MatDestroy(&tao->jacobian_state);
496:     tao->jacobian_state = J;
497:   }
498:   if (Jpre) {
499:     PetscObjectReference((PetscObject)Jpre);
500:     MatDestroy(&tao->jacobian_state_pre);
501:     tao->jacobian_state_pre=Jpre;
502:   }
503:   if (Jinv) {
504:     PetscObjectReference((PetscObject)Jinv);
505:     MatDestroy(&tao->jacobian_state_inv);
506:     tao->jacobian_state_inv=Jinv;
507:   }
508:   return(0);
509: }

511: /*@C
512:    TaoSetJacobianDesignRoutine - Sets the function to compute the Jacobian of
513:    the constraint function with respect to the design variables.  Used only for
514:    pde-constrained optimization.

516:    Logically collective on Tao

518:    Input Parameters:
519: +  tao  - the Tao context
520: .  J    - Matrix used for the jacobian
521: .  func - Jacobian evaluation routine
522: -  ctx  - [optional] user-defined context for private data for the
523:           Jacobian evaluation routine (may be NULL)

525:    Calling sequence of func:
526: $    func(Tao tao,Vec x,Mat J,void *ctx);

528: +  tao - the Tao  context
529: .  x   - input vector
530: .  J   - Jacobian matrix
531: -  ctx - [optional] user-defined Jacobian context

533:    Level: intermediate

535: .seealso: TaoComputeJacobianDesign(), TaoSetJacobianStateRoutine(), TaoSetStateDesignIS()
536: @*/
537: PetscErrorCode TaoSetJacobianDesignRoutine(Tao tao, Mat J, PetscErrorCode (*func)(Tao, Vec, Mat, void*), void *ctx)
538: {

543:   if (J) {
546:   }
547:   if (ctx) {
548:     tao->user_jac_designP = ctx;
549:   }
550:   if (func) {
551:     tao->ops->computejacobiandesign = func;
552:   }
553:   if (J) {
554:     PetscObjectReference((PetscObject)J);
555:     MatDestroy(&tao->jacobian_design);
556:     tao->jacobian_design = J;
557:   }
558:   return(0);
559: }

561: /*@
562:    TaoSetStateDesignIS - Indicate to the Tao which variables in the
563:    solution vector are state variables and which are design.  Only applies to
564:    pde-constrained optimization.

566:    Logically Collective on Tao

568:    Input Parameters:
569: +  tao  - The Tao context
570: .  s_is - the index set corresponding to the state variables
571: -  d_is - the index set corresponding to the design variables

573:    Level: intermediate

575: .seealso: TaoSetJacobianStateRoutine(), TaoSetJacobianDesignRoutine()
576: @*/
577: PetscErrorCode TaoSetStateDesignIS(Tao tao, IS s_is, IS d_is)
578: {

582:   PetscObjectReference((PetscObject)s_is);
583:   ISDestroy(&tao->state_is);
584:   tao->state_is = s_is;
585:   PetscObjectReference((PetscObject)(d_is));
586:   ISDestroy(&tao->design_is);
587:   tao->design_is = d_is;
588:   return(0);
589: }

591: /*@C
592:    TaoComputeJacobianEquality - Computes the Jacobian matrix that has been
593:    set with TaoSetJacobianEqualityRoutine().

595:    Collective on Tao

597:    Input Parameters:
598: +  tao - the Tao solver context
599: -  X   - input vector

601:    Output Parameters:
602: +  J    - Jacobian matrix
603: -  Jpre - Preconditioning matrix

605:    Notes:
606:    Most users should not need to explicitly call this routine, as it
607:    is used internally within the minimization solvers.

609:    Level: developer

611: .seealso:  TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
612: @*/
613: PetscErrorCode TaoComputeJacobianEquality(Tao tao, Vec X, Mat J, Mat Jpre)
614: {

621:   if (!tao->ops->computejacobianequality) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianEquality() first");
622:   ++tao->njac_equality;
623:   VecLockPush(X);
624:   PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
625:   PetscStackPush("Tao user Jacobian(equality) function");
626:   (*tao->ops->computejacobianequality)(tao,X,J,Jpre,tao->user_jac_equalityP);
627:   PetscStackPop;
628:   PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
629:   VecLockPop(X);
630:   return(0);
631: }

633: /*@C
634:    TaoComputeJacobianInequality - Computes the Jacobian matrix that has been
635:    set with TaoSetJacobianInequalityRoutine().

637:    Collective on Tao

639:    Input Parameters:
640: +  tao - the Tao solver context
641: -  X   - input vector

643:    Output Parameters:
644: +  J    - Jacobian matrix
645: -  Jpre - Preconditioning matrix

647:    Notes:
648:    Most users should not need to explicitly call this routine, as it
649:    is used internally within the minimization solvers.

651:    Level: developer

653: .seealso:  TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
654: @*/
655: PetscErrorCode TaoComputeJacobianInequality(Tao tao, Vec X, Mat J, Mat Jpre)
656: {

663:   if (!tao->ops->computejacobianinequality) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianInequality() first");
664:   ++tao->njac_inequality;
665:   VecLockPush(X);
666:   PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
667:   PetscStackPush("Tao user Jacobian(inequality) function");
668:   (*tao->ops->computejacobianinequality)(tao,X,J,Jpre,tao->user_jac_inequalityP);
669:   PetscStackPop;
670:   PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
671:   VecLockPop(X);
672:   return(0);
673: }

675: /*@C
676:    TaoSetJacobianEqualityRoutine - Sets the function to compute the Jacobian
677:    (and its inverse) of the constraint function with respect to the equality variables.
678:    Used only for pde-constrained optimization.

680:    Logically collective on Tao

682:    Input Parameters:
683: +  tao  - the Tao context
684: .  J    - Matrix used for the jacobian
685: .  Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.
686: .  func - Jacobian evaluation routine
687: -  ctx  - [optional] user-defined context for private data for the
688:           Jacobian evaluation routine (may be NULL)

690:    Calling sequence of func:
691: $    func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);

693: +  tao  - the Tao  context
694: .  x    - input vector
695: .  J    - Jacobian matrix
696: .  Jpre - preconditioner matrix, usually the same as J
697: -  ctx  - [optional] user-defined Jacobian context

699:    Level: intermediate

701: .seealso: TaoComputeJacobianEquality(), TaoSetJacobianDesignRoutine(), TaoSetEqualityDesignIS()
702: @*/
703: PetscErrorCode TaoSetJacobianEqualityRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat,void*), void *ctx)
704: {

709:   if (J) {
712:   }
713:   if (Jpre) {
716:   }
717:   if (ctx) {
718:     tao->user_jac_equalityP = ctx;
719:   }
720:   if (func) {
721:     tao->ops->computejacobianequality = func;
722:   }
723:   if (J) {
724:     PetscObjectReference((PetscObject)J);
725:     MatDestroy(&tao->jacobian_equality);
726:     tao->jacobian_equality = J;
727:   }
728:   if (Jpre) {
729:     PetscObjectReference((PetscObject)Jpre);
730:     MatDestroy(&tao->jacobian_equality_pre);
731:     tao->jacobian_equality_pre=Jpre;
732:   }
733:   return(0);
734: }

736: /*@C
737:    TaoSetJacobianInequalityRoutine - Sets the function to compute the Jacobian
738:    (and its inverse) of the constraint function with respect to the inequality variables.
739:    Used only for pde-constrained optimization.

741:    Logically collective on Tao

743:    Input Parameters:
744: +  tao  - the Tao context
745: .  J    - Matrix used for the jacobian
746: .  Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.
747: .  func - Jacobian evaluation routine
748: -  ctx  - [optional] user-defined context for private data for the
749:           Jacobian evaluation routine (may be NULL)

751:    Calling sequence of func:
752: $    func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);

754: +  tao  - the Tao  context
755: .  x    - input vector
756: .  J    - Jacobian matrix
757: .  Jpre - preconditioner matrix, usually the same as J
758: -  ctx  - [optional] user-defined Jacobian context

760:    Level: intermediate

762: .seealso: TaoComputeJacobianInequality(), TaoSetJacobianDesignRoutine(), TaoSetInequalityDesignIS()
763: @*/
764: PetscErrorCode TaoSetJacobianInequalityRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat,void*), void *ctx)
765: {

770:   if (J) {
773:   }
774:   if (Jpre) {
777:   }
778:   if (ctx) {
779:     tao->user_jac_inequalityP = ctx;
780:   }
781:   if (func) {
782:     tao->ops->computejacobianinequality = func;
783:   }
784:   if (J) {
785:     PetscObjectReference((PetscObject)J);
786:     MatDestroy(&tao->jacobian_inequality);
787:     tao->jacobian_inequality = J;
788:   }
789:   if (Jpre) {
790:     PetscObjectReference((PetscObject)Jpre);
791:     MatDestroy(&tao->jacobian_inequality_pre);
792:     tao->jacobian_inequality_pre=Jpre;
793:   }
794:   return(0);
795: }