Actual source code: taosolver_hj.c
petsc-3.9.4 2018-09-11
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: /*@C
61: TaoComputeHessian - Computes the Hessian matrix that has been
62: set with TaoSetHessianRoutine().
64: Collective on Tao
66: Input Parameters:
67: + tao - the Tao solver context
68: - X - input vector
70: Output Parameters:
71: + H - Hessian matrix
72: - Hpre - Preconditioning matrix
74: Notes:
75: Most users should not need to explicitly call this routine, as it
76: is used internally within the minimization solvers.
78: TaoComputeHessian() is typically used within minimization
79: implementations, so most users would not generally call this routine
80: themselves.
82: Level: developer
84: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetHessianRoutine()
85: @*/
86: PetscErrorCode TaoComputeHessian(Tao tao, Vec X, Mat H, Mat Hpre)
87: {
94: if (!tao->ops->computehessian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetHessianRoutine() first");
95: ++tao->nhess;
96: VecLockPush(X);
97: PetscLogEventBegin(Tao_HessianEval,tao,X,H,Hpre);
98: PetscStackPush("Tao user Hessian function");
99: (*tao->ops->computehessian)(tao,X,H,Hpre,tao->user_hessP);
100: PetscStackPop;
101: PetscLogEventEnd(Tao_HessianEval,tao,X,H,Hpre);
102: VecLockPop(X);
103: return(0);
104: }
106: /*@C
107: TaoComputeJacobian - Computes the Jacobian matrix that has been
108: set with TaoSetJacobianRoutine().
110: Collective on Tao
112: Input Parameters:
113: + tao - the Tao solver context
114: - X - input vector
116: Output Parameters:
117: + J - Jacobian matrix
118: - Jpre - Preconditioning matrix
120: Notes:
121: Most users should not need to explicitly call this routine, as it
122: is used internally within the minimization solvers.
124: TaoComputeJacobian() is typically used within minimization
125: implementations, so most users would not generally call this routine
126: themselves.
128: Level: developer
130: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianRoutine()
131: @*/
132: PetscErrorCode TaoComputeJacobian(Tao tao, Vec X, Mat J, Mat Jpre)
133: {
140: if (!tao->ops->computejacobian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobian() first");
141: ++tao->njac;
142: VecLockPush(X);
143: PetscLogEventBegin(Tao_JacobianEval,tao,X,J,Jpre);
144: PetscStackPush("Tao user Jacobian function");
145: (*tao->ops->computejacobian)(tao,X,J,Jpre,tao->user_jacP);
146: PetscStackPop;
147: PetscLogEventEnd(Tao_JacobianEval,tao,X,J,Jpre);
148: VecLockPop(X);
149: return(0);
150: }
152: /*@C
153: TaoComputeJacobianState - Computes the Jacobian matrix that has been
154: set with TaoSetJacobianStateRoutine().
156: Collective on Tao
158: Input Parameters:
159: + tao - the Tao solver context
160: - X - input vector
162: Output Parameters:
163: + Jpre - Jacobian matrix
164: - Jinv - Preconditioning matrix
166: Notes:
167: Most users should not need to explicitly call this routine, as it
168: is used internally within the minimization solvers.
170: TaoComputeJacobianState() is typically used within minimization
171: implementations, so most users would not generally call this routine
172: themselves.
174: Level: developer
176: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
177: @*/
178: PetscErrorCode TaoComputeJacobianState(Tao tao, Vec X, Mat J, Mat Jpre, Mat Jinv)
179: {
186: if (!tao->ops->computejacobianstate) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianState() first");
187: ++tao->njac_state;
188: VecLockPush(X);
189: PetscLogEventBegin(Tao_JacobianEval,tao,X,J,Jpre);
190: PetscStackPush("Tao user Jacobian(state) function");
191: (*tao->ops->computejacobianstate)(tao,X,J,Jpre,Jinv,tao->user_jac_stateP);
192: PetscStackPop;
193: PetscLogEventEnd(Tao_JacobianEval,tao,X,J,Jpre);
194: VecLockPop(X);
195: return(0);
196: }
198: /*@C
199: TaoComputeJacobianDesign - Computes the Jacobian matrix that has been
200: set with TaoSetJacobianDesignRoutine().
202: Collective on Tao
204: Input Parameters:
205: + tao - the Tao solver context
206: - X - input vector
208: Output Parameters:
209: . J - Jacobian matrix
211: Notes:
212: Most users should not need to explicitly call this routine, as it
213: is used internally within the minimization solvers.
215: TaoComputeJacobianDesign() is typically used within minimization
216: implementations, so most users would not generally call this routine
217: themselves.
219: Level: developer
221: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianDesignRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
222: @*/
223: PetscErrorCode TaoComputeJacobianDesign(Tao tao, Vec X, Mat J)
224: {
231: if (!tao->ops->computejacobiandesign) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianDesign() first");
232: ++tao->njac_design;
233: VecLockPush(X);
234: PetscLogEventBegin(Tao_JacobianEval,tao,X,J,NULL);
235: PetscStackPush("Tao user Jacobian(design) function");
236: (*tao->ops->computejacobiandesign)(tao,X,J,tao->user_jac_designP);
237: PetscStackPop;
238: PetscLogEventEnd(Tao_JacobianEval,tao,X,J,NULL);
239: VecLockPop(X);
240: return(0);
241: }
243: /*@C
244: TaoSetJacobianRoutine - Sets the function to compute the Jacobian as well as the location to store the matrix.
246: Logically collective on Tao
248: Input Parameters:
249: + tao - the Tao context
250: . J - Matrix used for the jacobian
251: . Jpre - Matrix that will be used operated on by preconditioner, can be same as J
252: . func - Jacobian evaluation routine
253: - ctx - [optional] user-defined context for private data for the
254: Jacobian evaluation routine (may be NULL)
256: Calling sequence of func:
257: $ func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);
259: + tao - the Tao context
260: . x - input vector
261: . J - Jacobian matrix
262: . Jpre - preconditioning matrix, usually the same as J
263: - ctx - [optional] user-defined Jacobian context
265: Level: intermediate
266: @*/
267: PetscErrorCode TaoSetJacobianRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, void*), void *ctx)
268: {
273: if (J) {
276: }
277: if (Jpre) {
280: }
281: if (ctx) {
282: tao->user_jacP = ctx;
283: }
284: if (func) {
285: tao->ops->computejacobian = func;
286: }
287: if (J) {
288: PetscObjectReference((PetscObject)J);
289: MatDestroy(&tao->jacobian);
290: tao->jacobian = J;
291: }
292: if (Jpre) {
293: PetscObjectReference((PetscObject)Jpre);
294: MatDestroy(&tao->jacobian_pre);
295: tao->jacobian_pre=Jpre;
296: }
297: return(0);
298: }
300: /*@C
301: TaoSetJacobianStateRoutine - Sets the function to compute the Jacobian
302: (and its inverse) of the constraint function with respect to the state variables.
303: Used only for pde-constrained optimization.
305: Logically collective on Tao
307: Input Parameters:
308: + tao - the Tao context
309: . J - Matrix used for the jacobian
310: . Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J. Only used if Jinv is NULL
311: . Jinv - [optional] Matrix used to apply the inverse of the state jacobian. Use NULL to default to PETSc KSP solvers to apply the inverse.
312: . func - Jacobian evaluation routine
313: - ctx - [optional] user-defined context for private data for the
314: Jacobian evaluation routine (may be NULL)
316: Calling sequence of func:
317: $ func(Tao tao,Vec x,Mat J,Mat Jpre,Mat Jinv,void *ctx);
319: + tao - the Tao context
320: . x - input vector
321: . J - Jacobian matrix
322: . Jpre - preconditioner matrix, usually the same as J
323: . Jinv - inverse of J
324: - ctx - [optional] user-defined Jacobian context
326: Level: intermediate
327: .seealso: TaoComputeJacobianState(), TaoSetJacobianDesignRoutine(), TaoSetStateDesignIS()
328: @*/
329: PetscErrorCode TaoSetJacobianStateRoutine(Tao tao, Mat J, Mat Jpre, Mat Jinv, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, Mat, void*), void *ctx)
330: {
335: if (J) {
338: }
339: if (Jpre) {
342: }
343: if (Jinv) {
346: }
347: if (ctx) {
348: tao->user_jac_stateP = ctx;
349: }
350: if (func) {
351: tao->ops->computejacobianstate = func;
352: }
353: if (J) {
354: PetscObjectReference((PetscObject)J);
355: MatDestroy(&tao->jacobian_state);
356: tao->jacobian_state = J;
357: }
358: if (Jpre) {
359: PetscObjectReference((PetscObject)Jpre);
360: MatDestroy(&tao->jacobian_state_pre);
361: tao->jacobian_state_pre=Jpre;
362: }
363: if (Jinv) {
364: PetscObjectReference((PetscObject)Jinv);
365: MatDestroy(&tao->jacobian_state_inv);
366: tao->jacobian_state_inv=Jinv;
367: }
368: return(0);
369: }
371: /*@C
372: TaoSetJacobianDesignRoutine - Sets the function to compute the Jacobian of
373: the constraint function with respect to the design variables. Used only for
374: pde-constrained optimization.
376: Logically collective on Tao
378: Input Parameters:
379: + tao - the Tao context
380: . J - Matrix used for the jacobian
381: . func - Jacobian evaluation routine
382: - ctx - [optional] user-defined context for private data for the
383: Jacobian evaluation routine (may be NULL)
385: Calling sequence of func:
386: $ func(Tao tao,Vec x,Mat J,void *ctx);
388: + tao - the Tao context
389: . x - input vector
390: . J - Jacobian matrix
391: - ctx - [optional] user-defined Jacobian context
393: Level: intermediate
395: .seealso: TaoComputeJacobianDesign(), TaoSetJacobianStateRoutine(), TaoSetStateDesignIS()
396: @*/
397: PetscErrorCode TaoSetJacobianDesignRoutine(Tao tao, Mat J, PetscErrorCode (*func)(Tao, Vec, Mat, void*), void *ctx)
398: {
403: if (J) {
406: }
407: if (ctx) {
408: tao->user_jac_designP = ctx;
409: }
410: if (func) {
411: tao->ops->computejacobiandesign = func;
412: }
413: if (J) {
414: PetscObjectReference((PetscObject)J);
415: MatDestroy(&tao->jacobian_design);
416: tao->jacobian_design = J;
417: }
418: return(0);
419: }
421: /*@
422: TaoSetStateDesignIS - Indicate to the Tao which variables in the
423: solution vector are state variables and which are design. Only applies to
424: pde-constrained optimization.
426: Logically Collective on Tao
428: Input Parameters:
429: + tao - The Tao context
430: . s_is - the index set corresponding to the state variables
431: - d_is - the index set corresponding to the design variables
433: Level: intermediate
435: .seealso: TaoSetJacobianStateRoutine(), TaoSetJacobianDesignRoutine()
436: @*/
437: PetscErrorCode TaoSetStateDesignIS(Tao tao, IS s_is, IS d_is)
438: {
442: PetscObjectReference((PetscObject)s_is);
443: ISDestroy(&tao->state_is);
444: tao->state_is = s_is;
445: PetscObjectReference((PetscObject)(d_is));
446: ISDestroy(&tao->design_is);
447: tao->design_is = d_is;
448: return(0);
449: }
451: /*@C
452: TaoComputeJacobianEquality - Computes the Jacobian matrix that has been
453: set with TaoSetJacobianEqualityRoutine().
455: Collective on Tao
457: Input Parameters:
458: + tao - the Tao solver context
459: - X - input vector
461: Output Parameters:
462: + J - Jacobian matrix
463: - Jpre - Preconditioning matrix
465: Notes:
466: Most users should not need to explicitly call this routine, as it
467: is used internally within the minimization solvers.
469: Level: developer
471: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
472: @*/
473: PetscErrorCode TaoComputeJacobianEquality(Tao tao, Vec X, Mat J, Mat Jpre)
474: {
481: if (!tao->ops->computejacobianequality) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianEquality() first");
482: ++tao->njac_equality;
483: VecLockPush(X);
484: PetscLogEventBegin(Tao_JacobianEval,tao,X,J,Jpre);
485: PetscStackPush("Tao user Jacobian(equality) function");
486: (*tao->ops->computejacobianequality)(tao,X,J,Jpre,tao->user_jac_equalityP);
487: PetscStackPop;
488: PetscLogEventEnd(Tao_JacobianEval,tao,X,J,Jpre);
489: VecLockPop(X);
490: return(0);
491: }
493: /*@C
494: TaoComputeJacobianInequality - Computes the Jacobian matrix that has been
495: set with TaoSetJacobianInequalityRoutine().
497: Collective on Tao
499: Input Parameters:
500: + tao - the Tao solver context
501: - X - input vector
503: Output Parameters:
504: + J - Jacobian matrix
505: - Jpre - Preconditioning matrix
507: Notes:
508: Most users should not need to explicitly call this routine, as it
509: is used internally within the minimization solvers.
511: Level: developer
513: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
514: @*/
515: PetscErrorCode TaoComputeJacobianInequality(Tao tao, Vec X, Mat J, Mat Jpre)
516: {
523: if (!tao->ops->computejacobianinequality) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianInequality() first");
524: ++tao->njac_inequality;
525: VecLockPush(X);
526: PetscLogEventBegin(Tao_JacobianEval,tao,X,J,Jpre);
527: PetscStackPush("Tao user Jacobian(inequality) function");
528: (*tao->ops->computejacobianinequality)(tao,X,J,Jpre,tao->user_jac_inequalityP);
529: PetscStackPop;
530: PetscLogEventEnd(Tao_JacobianEval,tao,X,J,Jpre);
531: VecLockPop(X);
532: return(0);
533: }
535: /*@C
536: TaoSetJacobianEqualityRoutine - Sets the function to compute the Jacobian
537: (and its inverse) of the constraint function with respect to the equality variables.
538: Used only for pde-constrained optimization.
540: Logically collective on Tao
542: Input Parameters:
543: + tao - the Tao context
544: . J - Matrix used for the jacobian
545: . Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.
546: . func - Jacobian evaluation routine
547: - ctx - [optional] user-defined context for private data for the
548: Jacobian evaluation routine (may be NULL)
550: Calling sequence of func:
551: $ func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);
553: + tao - the Tao context
554: . x - input vector
555: . J - Jacobian matrix
556: . Jpre - preconditioner matrix, usually the same as J
557: - ctx - [optional] user-defined Jacobian context
559: Level: intermediate
561: .seealso: TaoComputeJacobianEquality(), TaoSetJacobianDesignRoutine(), TaoSetEqualityDesignIS()
562: @*/
563: PetscErrorCode TaoSetJacobianEqualityRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat,void*), void *ctx)
564: {
569: if (J) {
572: }
573: if (Jpre) {
576: }
577: if (ctx) {
578: tao->user_jac_equalityP = ctx;
579: }
580: if (func) {
581: tao->ops->computejacobianequality = func;
582: }
583: if (J) {
584: PetscObjectReference((PetscObject)J);
585: MatDestroy(&tao->jacobian_equality);
586: tao->jacobian_equality = J;
587: }
588: if (Jpre) {
589: PetscObjectReference((PetscObject)Jpre);
590: MatDestroy(&tao->jacobian_equality_pre);
591: tao->jacobian_equality_pre=Jpre;
592: }
593: return(0);
594: }
596: /*@C
597: TaoSetJacobianInequalityRoutine - Sets the function to compute the Jacobian
598: (and its inverse) of the constraint function with respect to the inequality variables.
599: Used only for pde-constrained optimization.
601: Logically collective on Tao
603: Input Parameters:
604: + tao - the Tao context
605: . J - Matrix used for the jacobian
606: . Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.
607: . func - Jacobian evaluation routine
608: - ctx - [optional] user-defined context for private data for the
609: Jacobian evaluation routine (may be NULL)
611: Calling sequence of func:
612: $ func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);
614: + tao - the Tao context
615: . x - input vector
616: . J - Jacobian matrix
617: . Jpre - preconditioner matrix, usually the same as J
618: - ctx - [optional] user-defined Jacobian context
620: Level: intermediate
622: .seealso: TaoComputeJacobianInequality(), TaoSetJacobianDesignRoutine(), TaoSetInequalityDesignIS()
623: @*/
624: PetscErrorCode TaoSetJacobianInequalityRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat,void*), void *ctx)
625: {
630: if (J) {
633: }
634: if (Jpre) {
637: }
638: if (ctx) {
639: tao->user_jac_inequalityP = ctx;
640: }
641: if (func) {
642: tao->ops->computejacobianinequality = func;
643: }
644: if (J) {
645: PetscObjectReference((PetscObject)J);
646: MatDestroy(&tao->jacobian_inequality);
647: tao->jacobian_inequality = J;
648: }
649: if (Jpre) {
650: PetscObjectReference((PetscObject)Jpre);
651: MatDestroy(&tao->jacobian_inequality_pre);
652: tao->jacobian_inequality_pre=Jpre;
653: }
654: return(0);
655: }