xref: /petsc/src/tao/interface/taosolver_hj.c (revision 050fc7a34c6d42d9b367a0bc7bbb4624be070c15)
1 #include <petsc-private/taosolverimpl.h> /*I "taosolver.h" I*/
2 
3 #undef __FUNCT__
4 #define __FUNCT__ "TaoSetHessianRoutine"
5 /*@C
6    TaoSetHessianRoutine - Sets the function to compute the Hessian as well as the location to store the matrix.
7 
8    Logically collective on TaoSolver
9 
10    Input Parameters:
11 +  tao - the TaoSolver context
12 .  H - Matrix used for the hessian
13 .  Hpre - Matrix that will be used operated on by preconditioner, can be same as H
14 .  hess - Hessian evaluation routine
15 -  ctx - [optional] user-defined context for private data for the
16          Hessian evaluation routine (may be NULL)
17 
18    Calling sequence of hess:
19 $    hess (TaoSolver tao,Vec x,Mat *H,Mat *Hpre,MatStructure *flag,void *ctx);
20 
21 +  tao - the TaoSolver  context
22 .  x - input vector
23 .  H - Hessian matrix
24 .  Hpre - preconditioner matrix, usually the same as H
25 .  flag - flag indicating information about the preconditioner matrix
26    structure (see below)
27 -  ctx - [optional] user-defined Hessian context
28 
29 
30    Notes:
31 
32    The function hess() takes Mat * as the matrix arguments rather than Mat.
33    This allows the Hessian evaluation routine to replace A and/or B with a
34    completely new new matrix structure (not just different matrix elements)
35    when appropriate, for instance, if the nonzero structure is changing
36    throughout the global iterations.
37 
38    The flag can be used to eliminate unnecessary work in the preconditioner
39    during the repeated solution of linear systems of the same size.  The
40    available options are
41 $    SAME_PRECONDITIONER -
42 $      Hpre is identical during successive linear solves.
43 $      This option is intended for folks who are using
44 $      different Amat and Pmat matrices and want to reuse the
45 $      same preconditioner matrix.  For example, this option
46 $      saves work by not recomputing incomplete factorization
47 $      for ILU/ICC preconditioners.
48 $    SAME_NONZERO_PATTERN -
49 $      Hpre has the same nonzero structure during
50 $      successive linear solves.
51 $    DIFFERENT_NONZERO_PATTERN -
52 $      Hpre does not have the same nonzero structure.
53 
54    Caution:
55    If you specify SAME_NONZERO_PATTERN, the software believes your assertion
56    and does not check the structure of the matrix.  If you erroneously
57    claim that the structure is the same when it actually is not, the new
58    preconditioner will not function correctly.  Thus, use this optimization
59    feature carefully!
60 
61    If in doubt about whether your preconditioner matrix has changed
62    structure or not, use the flag DIFFERENT_NONZERO_PATTERN.
63 
64    Level: beginner
65 
66 @*/
67 PetscErrorCode TaoSetHessianRoutine(TaoSolver tao, Mat H, Mat Hpre, PetscErrorCode (*func)(TaoSolver, Vec, Mat*, Mat *, MatStructure *, void*), void *ctx)
68 {
69   PetscErrorCode ierr;
70 
71   PetscFunctionBegin;
72   PetscValidHeaderSpecific(tao,TAOSOLVER_CLASSID,1);
73   if (H) {
74     PetscValidHeaderSpecific(H,MAT_CLASSID,2);
75     PetscCheckSameComm(tao,1,H,2);
76   }
77   if (Hpre) {
78     PetscValidHeaderSpecific(Hpre,MAT_CLASSID,3);
79     PetscCheckSameComm(tao,1,Hpre,3);
80   }
81   if (ctx) {
82     tao->user_hessP = ctx;
83   }
84   if (func) {
85     tao->ops->computehessian = func;
86   }
87   if (H) {
88     ierr = PetscObjectReference((PetscObject)H);CHKERRQ(ierr);
89     ierr = MatDestroy(&tao->hessian);CHKERRQ(ierr);
90     tao->hessian = H;
91   }
92   if (Hpre) {
93     ierr = PetscObjectReference((PetscObject)Hpre);CHKERRQ(ierr);
94     ierr = MatDestroy(&tao->hessian_pre);CHKERRQ(ierr);
95     tao->hessian_pre = Hpre;
96   }
97   PetscFunctionReturn(0);
98 }
99 
100 #undef __FUNCT__
101 #define __FUNCT__ "TaoComputeHessian"
102 /*@C
103    TaoComputeHessian - Computes the Hessian matrix that has been
104    set with TaoSetHessianRoutine().
105 
106    Collective on TaoSolver
107 
108    Input Parameters:
109 +  solver - the TaoSolver solver context
110 -  xx - input vector
111 
112    Output Parameters:
113 +  H - Hessian matrix
114 .  Hpre - Preconditioning matrix
115 -  flag - flag indicating matrix structure (SAME_NONZERO_PATTERN, DIFFERENT_NONZERO_PATTERN, or SAME_PRECONDITIONER)
116 
117    Notes:
118    Most users should not need to explicitly call this routine, as it
119    is used internally within the minimization solvers.
120 
121    TaoComputeHessian() is typically used within minimization
122    implementations, so most users would not generally call this routine
123    themselves.
124 
125    Level: developer
126 
127 .seealso:  TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetHessian()
128 
129 @*/
130 PetscErrorCode TaoComputeHessian(TaoSolver tao, Vec X, Mat *H, Mat *Hpre, MatStructure *flg)
131 {
132   PetscErrorCode ierr;
133 
134   PetscFunctionBegin;
135   PetscValidHeaderSpecific(tao,TAOSOLVER_CLASSID,1);
136   PetscValidHeaderSpecific(X, VEC_CLASSID,2);
137   PetscValidPointer(flg,5);
138   PetscCheckSameComm(tao,1,X,2);
139 
140   if (!tao->ops->computehessian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetHessian() first");
141   *flg = DIFFERENT_NONZERO_PATTERN;
142   ++tao->nhess;
143   ierr = PetscLogEventBegin(TaoSolver_HessianEval,tao,X,*H,*Hpre);CHKERRQ(ierr);
144   PetscStackPush("TaoSolver user Hessian function");
145   CHKMEMQ;
146   ierr = (*tao->ops->computehessian)(tao,X,H,Hpre,flg,tao->user_hessP);CHKERRQ(ierr);
147   CHKMEMQ;
148   PetscStackPop;
149   ierr = PetscLogEventEnd(TaoSolver_HessianEval,tao,X,*H,*Hpre);CHKERRQ(ierr);
150   PetscFunctionReturn(0);
151 }
152 
153 #undef __FUNCT__
154 #define __FUNCT__ "TaoComputeJacobian"
155 /*@C
156    TaoComputeJacobian - Computes the Jacobian matrix that has been
157    set with TaoSetJacobianRoutine().
158 
159    Collective on TaoSolver
160 
161    Input Parameters:
162 +  solver - the TaoSolver solver context
163 -  xx - input vector
164 
165    Output Parameters:
166 +  H - Jacobian matrix
167 .  Hpre - Preconditioning matrix
168 -  flag - flag indicating matrix structure (SAME_NONZERO_PATTERN, DIFFERENT_NONZERO_PATTERN, or SAME_PRECONDITIONER)
169 
170    Notes:
171    Most users should not need to explicitly call this routine, as it
172    is used internally within the minimization solvers.
173 
174    TaoComputeJacobian() is typically used within minimization
175    implementations, so most users would not generally call this routine
176    themselves.
177 
178    Level: developer
179 
180 .seealso:  TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobian()
181 
182 @*/
183 PetscErrorCode TaoComputeJacobian(TaoSolver tao, Vec X, Mat *J, Mat *Jpre, MatStructure *flg)
184 {
185   PetscErrorCode ierr;
186 
187   PetscFunctionBegin;
188   PetscValidHeaderSpecific(tao,TAOSOLVER_CLASSID,1);
189   PetscValidHeaderSpecific(X, VEC_CLASSID,2);
190   PetscValidPointer(flg,5);
191   PetscCheckSameComm(tao,1,X,2);
192 
193   if (!tao->ops->computejacobian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobian() first");
194   *flg = DIFFERENT_NONZERO_PATTERN;
195   ++tao->njac;
196   ierr = PetscLogEventBegin(TaoSolver_JacobianEval,tao,X,*J,*Jpre);CHKERRQ(ierr);
197   PetscStackPush("TaoSolver user Jacobian function");
198   CHKMEMQ;
199   ierr = (*tao->ops->computejacobian)(tao,X,J,Jpre,flg,tao->user_jacP);CHKERRQ(ierr);
200   CHKMEMQ;
201   PetscStackPop;
202   ierr = PetscLogEventEnd(TaoSolver_JacobianEval,tao,X,*J,*Jpre);CHKERRQ(ierr);
203   PetscFunctionReturn(0);
204 }
205 
206 #undef __FUNCT__
207 #define __FUNCT__ "TaoComputeJacobianState"
208 /*@C
209    TaoComputeJacobianState - Computes the Jacobian matrix that has been
210    set with TaoSetJacobianStateRoutine().
211 
212    Collective on TaoSolver
213 
214    Input Parameters:
215 +  solver - the TaoSolver solver context
216 -  xx - input vector
217 
218    Output Parameters:
219 +  H - Jacobian matrix
220 .  Hpre - Preconditioning matrix
221 -  flag - flag indicating matrix structure (SAME_NONZERO_PATTERN, DIFFERENT_NONZERO_PATTERN, or SAME_PRECONDITIONER)
222 
223    Notes:
224    Most users should not need to explicitly call this routine, as it
225    is used internally within the minimization solvers.
226 
227    TaoComputeJacobianState() is typically used within minimization
228    implementations, so most users would not generally call this routine
229    themselves.
230 
231    Level: developer
232 
233 .seealso:  TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
234 
235 @*/
236 PetscErrorCode TaoComputeJacobianState(TaoSolver tao, Vec X, Mat *J, Mat *Jpre, Mat *Jinv, MatStructure *flg)
237 {
238   PetscErrorCode ierr;
239 
240   PetscFunctionBegin;
241   PetscValidHeaderSpecific(tao,TAOSOLVER_CLASSID,1);
242   PetscValidHeaderSpecific(X, VEC_CLASSID,2);
243   PetscValidPointer(flg,5);
244   PetscCheckSameComm(tao,1,X,2);
245 
246   if (!tao->ops->computejacobianstate) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianState() first");
247   *flg = DIFFERENT_NONZERO_PATTERN;
248   ++tao->njac_state;
249   ierr = PetscLogEventBegin(TaoSolver_JacobianEval,tao,X,*J,*Jpre);CHKERRQ(ierr);
250   PetscStackPush("TaoSolver user Jacobian(state) function");
251   CHKMEMQ;
252   ierr = (*tao->ops->computejacobianstate)(tao,X,J,Jpre,Jinv,flg,tao->user_jac_stateP);CHKERRQ(ierr);
253   CHKMEMQ;
254   PetscStackPop;
255   ierr = PetscLogEventEnd(TaoSolver_JacobianEval,tao,X,*J,*Jpre);CHKERRQ(ierr);
256   PetscFunctionReturn(0);
257 }
258 
259 #undef __FUNCT__
260 #define __FUNCT__ "TaoComputeJacobianDesign"
261 /*@C
262    TaoComputeJacobianDesign - Computes the Jacobian matrix that has been
263    set with TaoSetJacobianDesignRoutine().
264 
265    Collective on TaoSolver
266 
267    Input Parameters:
268 +  solver - the TaoSolver solver context
269 -  xx - input vector
270 
271    Output Parameters:
272 .  H - Jacobian matrix
273 
274    Notes:
275    Most users should not need to explicitly call this routine, as it
276    is used internally within the minimization solvers.
277 
278    TaoComputeJacobianDesign() is typically used within minimization
279    implementations, so most users would not generally call this routine
280    themselves.
281 
282    Level: developer
283 
284 .seealso:  TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianDesignRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
285 
286 @*/
287 PetscErrorCode TaoComputeJacobianDesign(TaoSolver tao, Vec X, Mat *J)
288 {
289   PetscErrorCode ierr;
290 
291   PetscFunctionBegin;
292   PetscValidHeaderSpecific(tao,TAOSOLVER_CLASSID,1);
293   PetscValidHeaderSpecific(X, VEC_CLASSID,2);
294   PetscCheckSameComm(tao,1,X,2);
295 
296   if (!tao->ops->computejacobiandesign) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianDesign() first");
297   ++tao->njac_design;
298   ierr = PetscLogEventBegin(TaoSolver_JacobianEval,tao,X,*J,NULL);CHKERRQ(ierr);
299   PetscStackPush("TaoSolver user Jacobian(design) function");
300   CHKMEMQ;
301   ierr = (*tao->ops->computejacobiandesign)(tao,X,J,tao->user_jac_designP);CHKERRQ(ierr);
302   CHKMEMQ;
303   PetscStackPop;
304   ierr = PetscLogEventEnd(TaoSolver_JacobianEval,tao,X,*J,NULL);CHKERRQ(ierr);
305   PetscFunctionReturn(0);
306 }
307 
308 #undef __FUNCT__
309 #define __FUNCT__ "TaoSetJacobianRoutine"
310 /*@C
311    TaoSetJacobianRoutine - Sets the function to compute the Jacobian as well as the location to store the matrix.
312 
313    Logically collective on TaoSolver
314 
315    Input Parameters:
316 +  tao - the TaoSolver context
317 .  J - Matrix used for the jacobian
318 .  Jpre - Matrix that will be used operated on by preconditioner, can be same as J
319 .  jac - Jacobian evaluation routine
320 -  ctx - [optional] user-defined context for private data for the
321          Jacobian evaluation routine (may be NULL)
322 
323    Calling sequence of jac:
324 $    jac (TaoSolver tao,Vec x,Mat *J,Mat *Jpre,MatStructure *flag,void *ctx);
325 
326 +  tao - the TaoSolver  context
327 .  x - input vector
328 .  J - Jacobian matrix
329 .  Jpre - preconditioner matrix, usually the same as J
330 .  flag - flag indicating information about the preconditioner matrix
331    structure (see below)
332 -  ctx - [optional] user-defined Jacobian context
333 
334    Notes:
335 
336    The function jac() takes Mat * as the matrix arguments rather than Mat.
337    This allows the Jacobian evaluation routine to replace A and/or B with a
338    completely new new matrix structure (not just different matrix elements)
339    when appropriate, for instance, if the nonzero structure is changing
340    throughout the global iterations.
341 
342    The flag can be used to eliminate unnecessary work in the preconditioner
343    during the repeated solution of linear systems of the same size.  The
344    available options are
345 $    SAME_PRECONDITIONER -
346 $      Jpre is identical during successive linear solves.
347 $      This option is intended for folks who are using
348 $      different Amat and Pmat matrices and want to reuse the
349 $      same preconditioner matrix.  For example, this option
350 $      saves work by not recomputing incomplete factorization
351 $      for ILU/ICC preconditioners.
352 $    SAME_NONZERO_PATTERN -
353 $      Jpre has the same nonzero structure during
354 $      successive linear solves.
355 $    DIFFERENT_NONZERO_PATTERN -
356 $      Jpre does not have the same nonzero structure.
357 
358    Caution:
359    If you specify SAME_NONZERO_PATTERN, the software believes your assertion
360    and does not check the structure of the matrix.  If you erroneously
361    claim that the structure is the same when it actually is not, the new
362    preconditioner will not function correctly.  Thus, use this optimization
363    feature carefully!
364 
365    If in doubt about whether your preconditioner matrix has changed
366    structure or not, use the flag DIFFERENT_NONZERO_PATTERN.
367 
368    Level: intermediate
369 
370 @*/
371 PetscErrorCode TaoSetJacobianRoutine(TaoSolver tao, Mat J, Mat Jpre, PetscErrorCode (*func)(TaoSolver, Vec, Mat*, Mat *, MatStructure *, void*), void *ctx)
372 {
373   PetscErrorCode ierr;
374   PetscFunctionBegin;
375   PetscValidHeaderSpecific(tao,TAOSOLVER_CLASSID,1);
376   if (J) {
377     PetscValidHeaderSpecific(J,MAT_CLASSID,2);
378     PetscCheckSameComm(tao,1,J,2);
379   }
380   if (Jpre) {
381     PetscValidHeaderSpecific(Jpre,MAT_CLASSID,3);
382     PetscCheckSameComm(tao,1,Jpre,3);
383   }
384   if (ctx) {
385     tao->user_jacP = ctx;
386   }
387   if (func) {
388     tao->ops->computejacobian = func;
389   }
390   if (J) {
391     ierr = PetscObjectReference((PetscObject)J);CHKERRQ(ierr);
392     ierr = MatDestroy(&tao->jacobian);CHKERRQ(ierr);
393     tao->jacobian = J;
394   }
395   if (Jpre) {
396     ierr = PetscObjectReference((PetscObject)Jpre);CHKERRQ(ierr);
397     ierr = MatDestroy(&tao->jacobian_pre);CHKERRQ(ierr);
398     tao->jacobian_pre=Jpre;
399   }
400   PetscFunctionReturn(0);
401 }
402 
403 #undef __FUNCT__
404 #define __FUNCT__ "TaoSetJacobianStateRoutine"
405 /*@C
406    TaoSetJacobianStateRoutine - Sets the function to compute the Jacobian
407    (and its inverse) of the constraint function with respect to the state variables.
408    Used only for pde-constrained optimization.
409 
410    Logically collective on TaoSolver
411 
412    Input Parameters:
413 +  tao - the TaoSolver context
414 .  J - Matrix used for the jacobian
415 .  Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.  Only used if Jinv is NULL
416 .  Jinv - [optional] Matrix used to apply the inverse of the state jacobian. Use NULL to default to PETSc KSP solvers to apply the inverse.
417 .  jac - Jacobian evaluation routine
418 -  ctx - [optional] user-defined context for private data for the
419          Jacobian evaluation routine (may be NULL)
420 
421    Calling sequence of jac:
422 $    jac (TaoSolver tao,Vec x,Mat *J,Mat *Jpre,MatStructure *flag,void *ctx);
423 
424 +  tao - the TaoSolver  context
425 .  x - input vector
426 .  J - Jacobian matrix
427 .  Jpre - preconditioner matrix, usually the same as J
428 .  Jinv - inverse of J
429 .  flag - flag indicating information about the preconditioner matrix
430    structure (see below)
431 -  ctx - [optional] user-defined Jacobian context
432 
433 
434    Notes:
435    Because of the structure of the jacobian matrix,
436    It may be more efficient for a pde-constrained application to provide
437    its own Jinv matrix.
438 
439    The function jac() takes Mat * as the matrix arguments rather than Mat.
440    This allows the Jacobian evaluation routine to replace A and/or B with a
441    completely new new maitrix structure (not just different matrix elements)
442    when appropriate, for instance, if the nonzero structure is changing
443    throughout the global iterations.
444 
445    The flag can be used to eliminate unnecessary work in the preconditioner
446    during the repeated solution of linear systems of the same size.  The
447    available options are
448 $    SAME_PRECONDITIONER -
449 $      Jpre is identical during successive linear solves.
450 $      This option is intended for folks who are using
451 $      different Amat and Pmat matrices and want to reuse the
452 $      same preconditioner matrix.  For example, this option
453 $      saves work by not recomputing incomplete factorization
454 $      for ILU/ICC preconditioners.
455 $    SAME_NONZERO_PATTERN -
456 $      Jpre has the same nonzero structure during
457 $      successive linear solves.
458 $    DIFFERENT_NONZERO_PATTERN -
459 $      Jpre does not have the same nonzero structure.
460 
461    Caution:
462    If you specify SAME_NONZERO_PATTERN, the software believes your assertion
463    and does not check the structure of the matrix.  If you erroneously
464    claim that the structure is the same when it actually is not, the new
465    preconditioner will not function correctly.  Thus, use this optimization
466    feature carefully!
467 
468    If in doubt about whether your preconditioner matrix has changed
469    structure or not, use the flag DIFFERENT_NONZERO_PATTERN.
470 
471    Level: intermediate
472 .seealse: TaoComputeJacobianState(), TaoSetJacobianDesignRoutine(), TaoSetStateDesignIS()
473 @*/
474 PetscErrorCode TaoSetJacobianStateRoutine(TaoSolver tao, Mat J, Mat Jpre, Mat Jinv, PetscErrorCode (*func)(TaoSolver, Vec, Mat*, Mat *, Mat *, MatStructure *, void*), void *ctx)
475 {
476   PetscErrorCode ierr;
477   PetscFunctionBegin;
478   PetscValidHeaderSpecific(tao,TAOSOLVER_CLASSID,1);
479   if (J) {
480     PetscValidHeaderSpecific(J,MAT_CLASSID,2);
481     PetscCheckSameComm(tao,1,J,2);
482   }
483   if (Jpre) {
484     PetscValidHeaderSpecific(Jpre,MAT_CLASSID,3);
485     PetscCheckSameComm(tao,1,Jpre,3);
486   }
487   if (Jinv) {
488     PetscValidHeaderSpecific(Jinv,MAT_CLASSID,4);
489     PetscCheckSameComm(tao,1,Jinv,4);
490   }
491   if (ctx) {
492     tao->user_jac_stateP = ctx;
493   }
494   if (func) {
495     tao->ops->computejacobianstate = func;
496   }
497   if (J) {
498     ierr = PetscObjectReference((PetscObject)J);CHKERRQ(ierr);
499     ierr = MatDestroy(&tao->jacobian_state);CHKERRQ(ierr);
500     tao->jacobian_state = J;
501   }
502   if (Jpre) {
503     ierr = PetscObjectReference((PetscObject)Jpre);CHKERRQ(ierr);
504     ierr = MatDestroy(&tao->jacobian_state_pre);CHKERRQ(ierr);
505     tao->jacobian_state_pre=Jpre;
506   }
507   if (Jinv) {
508     ierr = PetscObjectReference((PetscObject)Jinv);CHKERRQ(ierr);
509     ierr = MatDestroy(&tao->jacobian_state_inv);CHKERRQ(ierr);
510     tao->jacobian_state_inv=Jinv;
511   }
512   PetscFunctionReturn(0);
513 }
514 
515 #undef __FUNCT__
516 #define __FUNCT__ "TaoSetJacobianDesignRoutine"
517 /*@C
518    TaoSetJacobianDesignRoutine - Sets the function to compute the Jacobian of
519    the constraint function with respect to the design variables.  Used only for
520    pde-constrained optimization.
521 
522    Logically collective on TaoSolver
523 
524    Input Parameters:
525 +  tao - the TaoSolver context
526 .  J - Matrix used for the jacobian
527 .  jac - Jacobian evaluation routine
528 -  ctx - [optional] user-defined context for private data for the
529          Jacobian evaluation routine (may be NULL)
530 
531    Calling sequence of jac:
532 $    jac (TaoSolver tao,Vec x,Mat *J,void *ctx);
533 
534 +  tao - the TaoSolver  context
535 .  x - input vector
536 .  J - Jacobian matrix
537 -  ctx - [optional] user-defined Jacobian context
538 
539 
540    Notes:
541 
542    The function jac() takes Mat * as the matrix arguments rather than Mat.
543    This allows the Jacobian evaluation routine to replace A and/or B with a
544    completely new new matrix structure (not just different matrix elements)
545    when appropriate, for instance, if the nonzero structure is changing
546    throughout the global iterations.
547 
548    Level: intermediate
549 .seealso: TaoComputeJacobianDesign(), TaoSetJacobianStateRoutine(), TaoSetStateDesignIS()
550 @*/
551 PetscErrorCode TaoSetJacobianDesignRoutine(TaoSolver tao, Mat J, PetscErrorCode (*func)(TaoSolver, Vec, Mat*, void*), void *ctx)
552 {
553   PetscErrorCode ierr;
554 
555   PetscFunctionBegin;
556   PetscValidHeaderSpecific(tao,TAOSOLVER_CLASSID,1);
557   if (J) {
558     PetscValidHeaderSpecific(J,MAT_CLASSID,2);
559     PetscCheckSameComm(tao,1,J,2);
560   }
561   if (ctx) {
562     tao->user_jac_designP = ctx;
563   }
564   if (func) {
565     tao->ops->computejacobiandesign = func;
566   }
567   if (J) {
568     ierr = PetscObjectReference((PetscObject)J);CHKERRQ(ierr);
569     ierr = MatDestroy(&tao->jacobian_design);CHKERRQ(ierr);
570     tao->jacobian_design = J;
571   }
572   PetscFunctionReturn(0);
573 }
574 
575 #undef __FUNCT__
576 #define __FUNCT__ "TaoSetStateDesignIS"
577 /*@
578    TaoSetStateDesignIS - Indicate to the TaoSolver which variables in the
579    solution vector are state variables and which are design.  Only applies to
580    pde-constrained optimization.
581 
582    Logically Collective on TaoSolver
583 
584    Input Parameters:
585 +  tao - The TaoSolver context
586 .  s_is - the index set corresponding to the state variables
587 -  d_is - the index set corresponding to the design variables
588 
589    Level: intermediate
590 
591 .seealso: TaoSetJacobianStateRoutine(), TaoSetJacobianDesignRoutine()
592 @*/
593 PetscErrorCode TaoSetStateDesignIS(TaoSolver tao, IS s_is, IS d_is)
594 {
595   PetscErrorCode ierr;
596 
597   PetscFunctionBegin;
598   ierr = PetscObjectReference((PetscObject)s_is);CHKERRQ(ierr);
599   ierr = ISDestroy(&tao->state_is);CHKERRQ(ierr);
600   tao->state_is = s_is;
601   ierr = PetscObjectReference((PetscObject)(d_is));CHKERRQ(ierr);
602   ierr = ISDestroy(&tao->design_is);CHKERRQ(ierr);
603   tao->design_is = d_is;
604   PetscFunctionReturn(0);
605 }
606 
607 #undef __FUNCT__
608 #define __FUNCT__ "TaoComputeJacobianEquality"
609 /*@C
610    TaoComputeJacobianEquality - Computes the Jacobian matrix that has been
611    set with TaoSetJacobianEqualityRoutine().
612 
613    Collective on TaoSolver
614 
615    Input Parameters:
616 +  solver - the TaoSolver solver context
617 -  xx - input vector
618 
619    Output Parameters:
620 +  H - Jacobian matrix
621 .  Hpre - Preconditioning matrix
622 -  flag - flag indicating matrix structure (SAME_NONZERO_PATTERN, DIFFERENT_NONZERO_PATTERN, or SAME_PRECONDITIONER)
623 
624    Notes:
625    Most users should not need to explicitly call this routine, as it
626    is used internally within the minimization solvers.
627 
628    Level: developer
629 
630 .seealso:  TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
631 
632 @*/
633 PetscErrorCode TaoComputeJacobianEquality(TaoSolver tao, Vec X, Mat *J, Mat *Jpre, MatStructure *flg)
634 {
635   PetscErrorCode ierr;
636 
637   PetscFunctionBegin;
638   PetscValidHeaderSpecific(tao,TAOSOLVER_CLASSID,1);
639   PetscValidHeaderSpecific(X, VEC_CLASSID,2);
640   PetscValidPointer(flg,5);
641   PetscCheckSameComm(tao,1,X,2);
642 
643   if (!tao->ops->computejacobianequality) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianEquality() first");
644   *flg = DIFFERENT_NONZERO_PATTERN;
645   ++tao->njac_equality;
646   ierr = PetscLogEventBegin(TaoSolver_JacobianEval,tao,X,*J,*Jpre);CHKERRQ(ierr);
647   PetscStackPush("TaoSolver user Jacobian(equality) function");
648   CHKMEMQ;
649   ierr = (*tao->ops->computejacobianequality)(tao,X,J,Jpre,flg,tao->user_jac_equalityP);CHKERRQ(ierr);
650   CHKMEMQ;
651   PetscStackPop;
652   ierr = PetscLogEventEnd(TaoSolver_JacobianEval,tao,X,*J,*Jpre);CHKERRQ(ierr);
653   PetscFunctionReturn(0);
654 }
655 
656 #undef __FUNCT__
657 #define __FUNCT__ "TaoComputeJacobianInequality"
658 /*@C
659    TaoComputeJacobianInequality - Computes the Jacobian matrix that has been
660    set with TaoSetJacobianInequalityRoutine().
661 
662    Collective on TaoSolver
663 
664    Input Parameters:
665 +  solver - the TaoSolver solver context
666 -  xx - input vector
667 
668    Output Parameters:
669 +  H - Jacobian matrix
670 .  Hpre - Preconditioning matrix
671 -  flag - flag indicating matrix structure (SAME_NONZERO_PATTERN, DIFFERENT_NONZERO_PATTERN, or SAME_PRECONDITIONER)
672 
673    Notes:
674    Most users should not need to explicitly call this routine, as it
675    is used internally within the minimization solvers.
676 
677    Level: developer
678 
679 .seealso:  TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
680 
681 @*/
682 PetscErrorCode TaoComputeJacobianInequality(TaoSolver tao, Vec X, Mat *J, Mat *Jpre, MatStructure *flg)
683 {
684   PetscErrorCode ierr;
685 
686   PetscFunctionBegin;
687   PetscValidHeaderSpecific(tao,TAOSOLVER_CLASSID,1);
688   PetscValidHeaderSpecific(X, VEC_CLASSID,2);
689   PetscValidPointer(flg,5);
690   PetscCheckSameComm(tao,1,X,2);
691 
692   if (!tao->ops->computejacobianinequality) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianInequality() first");
693   *flg = DIFFERENT_NONZERO_PATTERN;
694   ++tao->njac_inequality;
695   ierr = PetscLogEventBegin(TaoSolver_JacobianEval,tao,X,*J,*Jpre);CHKERRQ(ierr);
696   PetscStackPush("TaoSolver user Jacobian(inequality) function");
697   CHKMEMQ;
698   ierr = (*tao->ops->computejacobianinequality)(tao,X,J,Jpre,flg,tao->user_jac_inequalityP);CHKERRQ(ierr);
699   CHKMEMQ;
700   PetscStackPop;
701   ierr = PetscLogEventEnd(TaoSolver_JacobianEval,tao,X,*J,*Jpre);CHKERRQ(ierr);
702   PetscFunctionReturn(0);
703 }
704 
705 #undef __FUNCT__
706 #define __FUNCT__ "TaoSetJacobianEqualityRoutine"
707 /*@C
708    TaoSetJacobianEqualityRoutine - Sets the function to compute the Jacobian
709    (and its inverse) of the constraint function with respect to the equality variables.
710    Used only for pde-constrained optimization.
711 
712    Logically collective on TaoSolver
713 
714    Input Parameters:
715 +  tao - the TaoSolver context
716 .  J - Matrix used for the jacobian
717 .  Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.
718 .  jac - Jacobian evaluation routine
719 -  ctx - [optional] user-defined context for private data for the
720          Jacobian evaluation routine (may be NULL)
721 
722    Calling sequence of jac:
723 $    jac (TaoSolver tao,Vec x,Mat *J,Mat *Jpre,MatStructure *flag,void *ctx);
724 
725 +  tao - the TaoSolver  context
726 .  x - input vector
727 .  J - Jacobian matrix
728 .  Jpre - preconditioner matrix, usually the same as J
729 .  flag - flag indicating information about the preconditioner matrix
730    structure (see below)
731 -  ctx - [optional] user-defined Jacobian context
732 
733    Notes:
734    Because of the structure of the jacobian matrix,
735    It may be more efficient for a pde-constrained application to provide
736    its own Jinv matrix.
737 
738    The function jac() takes Mat * as the matrix arguments rather than Mat.
739    This allows the Jacobian evaluation routine to replace A and/or B with a
740    completely new new maitrix structure (not just different matrix elements)
741    when appropriate, for instance, if the nonzero structure is changing
742    throughout the global iterations.
743 
744    The flag can be used to eliminate unnecessary work in the preconditioner
745    during the repeated solution of linear systems of the same size.  The
746    available options are
747 $    SAME_PRECONDITIONER -
748 $      Jpre is identical during successive linear solves.
749 $      This option is intended for folks who are using
750 $      different Amat and Pmat matrices and want to reuse the
751 $      same preconditioner matrix.  For example, this option
752 $      saves work by not recomputing incomplete factorization
753 $      for ILU/ICC preconditioners.
754 $    SAME_NONZERO_PATTERN -
755 $      Jpre has the same nonzero structure during
756 $      successive linear solves.
757 $    DIFFERENT_NONZERO_PATTERN -
758 $      Jpre does not have the same nonzero structure.
759 
760    Caution:
761    If you specify SAME_NONZERO_PATTERN, the software believes your assertion
762    and does not check the structure of the matrix.  If you erroneously
763    claim that the structure is the same when it actually is not, the new
764    preconditioner will not function correctly.  Thus, use this optimization
765    feature carefully!
766 
767    If in doubt about whether your preconditioner matrix has changed
768    structure or not, use the flag DIFFERENT_NONZERO_PATTERN.
769 
770    Level: intermediate
771 .seealse: TaoComputeJacobianEquality(), TaoSetJacobianDesignRoutine(), TaoSetEqualityDesignIS()
772 @*/
773 PetscErrorCode TaoSetJacobianEqualityRoutine(TaoSolver tao, Mat J, Mat Jpre, PetscErrorCode (*func)(TaoSolver, Vec, Mat*, Mat *, MatStructure *, void*), void *ctx)
774 {
775   PetscErrorCode ierr;
776 
777   PetscFunctionBegin;
778   PetscValidHeaderSpecific(tao,TAOSOLVER_CLASSID,1);
779   if (J) {
780     PetscValidHeaderSpecific(J,MAT_CLASSID,2);
781     PetscCheckSameComm(tao,1,J,2);
782   }
783   if (Jpre) {
784     PetscValidHeaderSpecific(Jpre,MAT_CLASSID,3);
785     PetscCheckSameComm(tao,1,Jpre,3);
786   }
787   if (ctx) {
788     tao->user_jac_equalityP = ctx;
789   }
790   if (func) {
791     tao->ops->computejacobianequality = func;
792   }
793   if (J) {
794     ierr = PetscObjectReference((PetscObject)J);CHKERRQ(ierr);
795     ierr = MatDestroy(&tao->jacobian_equality);CHKERRQ(ierr);
796     tao->jacobian_equality = J;
797   }
798   if (Jpre) {
799     ierr = PetscObjectReference((PetscObject)Jpre);CHKERRQ(ierr);
800     ierr = MatDestroy(&tao->jacobian_equality_pre);CHKERRQ(ierr);
801     tao->jacobian_equality_pre=Jpre;
802   }
803   PetscFunctionReturn(0);
804 }
805 
806 #undef __FUNCT__
807 #define __FUNCT__ "TaoSetJacobianInequalityRoutine"
808 /*@C
809    TaoSetJacobianInequalityRoutine - Sets the function to compute the Jacobian
810    (and its inverse) of the constraint function with respect to the inequality variables.
811    Used only for pde-constrained optimization.
812 
813    Logically collective on TaoSolver
814 
815    Input Parameters:
816 +  tao - the TaoSolver context
817 .  J - Matrix used for the jacobian
818 .  Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.
819 .  jac - Jacobian evaluation routine
820 -  ctx - [optional] user-defined context for private data for the
821          Jacobian evaluation routine (may be NULL)
822 
823    Calling sequence of jac:
824 $    jac (TaoSolver tao,Vec x,Mat *J,Mat *Jpre,MatStructure *flag,void *ctx);
825 
826 +  tao - the TaoSolver  context
827 .  x - input vector
828 .  J - Jacobian matrix
829 .  Jpre - preconditioner matrix, usually the same as J
830 .  flag - flag indicating information about the preconditioner matrix
831    structure (see below)
832 -  ctx - [optional] user-defined Jacobian context
833 
834 
835    Notes:
836    Because of the structure of the jacobian matrix,
837    It may be more efficient for a pde-constrained application to provide
838    its own Jinv matrix.
839 
840    The function jac() takes Mat * as the matrix arguments rather than Mat.
841    This allows the Jacobian evaluation routine to replace A and/or B with a
842    completely new new maitrix structure (not just different matrix elements)
843    when appropriate, for instance, if the nonzero structure is changing
844    throughout the global iterations.
845 
846    The flag can be used to eliminate unnecessary work in the preconditioner
847    during the repeated solution of linear systems of the same size.  The
848    available options are
849 $    SAME_PRECONDITIONER -
850 $      Jpre is identical during successive linear solves.
851 $      This option is intended for folks who are using
852 $      different Amat and Pmat matrices and want to reuse the
853 $      same preconditioner matrix.  For example, this option
854 $      saves work by not recomputing incomplete factorization
855 $      for ILU/ICC preconditioners.
856 $    SAME_NONZERO_PATTERN -
857 $      Jpre has the same nonzero structure during
858 $      successive linear solves.
859 $    DIFFERENT_NONZERO_PATTERN -
860 $      Jpre does not have the same nonzero structure.
861 
862    Caution:
863    If you specify SAME_NONZERO_PATTERN, the software believes your assertion
864    and does not check the structure of the matrix.  If you erroneously
865    claim that the structure is the same when it actually is not, the new
866    preconditioner will not function correctly.  Thus, use this optimization
867    feature carefully!
868 
869    If in doubt about whether your preconditioner matrix has changed
870    structure or not, use the flag DIFFERENT_NONZERO_PATTERN.
871 
872    Level: intermediate
873 .seealse: TaoComputeJacobianInequality(), TaoSetJacobianDesignRoutine(), TaoSetInequalityDesignIS()
874 @*/
875 PetscErrorCode TaoSetJacobianInequalityRoutine(TaoSolver tao, Mat J, Mat Jpre, PetscErrorCode (*func)(TaoSolver, Vec, Mat*, Mat *, MatStructure *, void*), void *ctx)
876 {
877   PetscErrorCode ierr;
878   PetscFunctionBegin;
879   PetscValidHeaderSpecific(tao,TAOSOLVER_CLASSID,1);
880   if (J) {
881     PetscValidHeaderSpecific(J,MAT_CLASSID,2);
882     PetscCheckSameComm(tao,1,J,2);
883   }
884   if (Jpre) {
885     PetscValidHeaderSpecific(Jpre,MAT_CLASSID,3);
886     PetscCheckSameComm(tao,1,Jpre,3);
887   }
888   if (ctx) {
889     tao->user_jac_inequalityP = ctx;
890   }
891   if (func) {
892     tao->ops->computejacobianinequality = func;
893   }
894   if (J) {
895     ierr = PetscObjectReference((PetscObject)J);CHKERRQ(ierr);
896     ierr = MatDestroy(&tao->jacobian_inequality);CHKERRQ(ierr);
897     tao->jacobian_inequality = J;
898   }
899   if (Jpre) {
900     ierr = PetscObjectReference((PetscObject)Jpre);CHKERRQ(ierr);
901     ierr = MatDestroy(&tao->jacobian_inequality_pre);CHKERRQ(ierr);
902     tao->jacobian_inequality_pre=Jpre;
903   }
904   PetscFunctionReturn(0);
905 }
906