Actual source code: taosolver_hj.c
petsc-3.8.4 2018-03-24
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: . hess - 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 hess:
17: $ hess (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
27: @*/
28: PetscErrorCode TaoSetHessianRoutine(Tao tao, Mat H, Mat Hpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, void*), void *ctx)
29: {
34: if (H) {
37: }
38: if (Hpre) {
41: }
42: if (ctx) {
43: tao->user_hessP = ctx;
44: }
45: if (func) {
46: tao->ops->computehessian = func;
47: }
48: if (H) {
49: PetscObjectReference((PetscObject)H);
50: MatDestroy(&tao->hessian);
51: tao->hessian = H;
52: }
53: if (Hpre) {
54: PetscObjectReference((PetscObject)Hpre);
55: MatDestroy(&tao->hessian_pre);
56: tao->hessian_pre = Hpre;
57: }
58: return(0);
59: }
61: /*@C
62: TaoComputeHessian - Computes the Hessian matrix that has been
63: set with TaoSetHessianRoutine().
65: Collective on Tao
67: Input Parameters:
68: + solver - the Tao solver context
69: - xx - input vector
71: Output Parameters:
72: + H - Hessian matrix
73: - Hpre - Preconditioning matrix
75: Notes:
76: Most users should not need to explicitly call this routine, as it
77: is used internally within the minimization solvers.
79: TaoComputeHessian() is typically used within minimization
80: implementations, so most users would not generally call this routine
81: themselves.
83: Level: developer
85: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetHessian()
87: @*/
88: PetscErrorCode TaoComputeHessian(Tao tao, Vec X, Mat H, Mat Hpre)
89: {
97: if (!tao->ops->computehessian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetHessian() first");
98: ++tao->nhess;
99: PetscLogEventBegin(Tao_HessianEval,tao,X,H,Hpre);
100: PetscStackPush("Tao user Hessian function");
101: (*tao->ops->computehessian)(tao,X,H,Hpre,tao->user_hessP);
102: PetscStackPop;
103: PetscLogEventEnd(Tao_HessianEval,tao,X,H,Hpre);
104: return(0);
105: }
107: /*@C
108: TaoComputeJacobian - Computes the Jacobian matrix that has been
109: set with TaoSetJacobianRoutine().
111: Collective on Tao
113: Input Parameters:
114: + solver - the Tao solver context
115: - xx - input vector
117: Output Parameters:
118: + H - Jacobian matrix
119: - Hpre - Preconditioning matrix
121: Notes:
122: Most users should not need to explicitly call this routine, as it
123: is used internally within the minimization solvers.
125: TaoComputeJacobian() is typically used within minimization
126: implementations, so most users would not generally call this routine
127: themselves.
129: Level: developer
131: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobian()
133: @*/
134: PetscErrorCode TaoComputeJacobian(Tao tao, Vec X, Mat J, Mat Jpre)
135: {
143: if (!tao->ops->computejacobian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobian() first");
144: ++tao->njac;
145: PetscLogEventBegin(Tao_JacobianEval,tao,X,J,Jpre);
146: PetscStackPush("Tao user Jacobian function");
147: (*tao->ops->computejacobian)(tao,X,J,Jpre,tao->user_jacP);
148: PetscStackPop;
149: PetscLogEventEnd(Tao_JacobianEval,tao,X,J,Jpre);
150: return(0);
151: }
153: /*@C
154: TaoComputeJacobianState - Computes the Jacobian matrix that has been
155: set with TaoSetJacobianStateRoutine().
157: Collective on Tao
159: Input Parameters:
160: + solver - the Tao solver context
161: - xx - input vector
163: Output Parameters:
164: + H - Jacobian matrix
165: - Hpre - Preconditioning matrix
167: Notes:
168: Most users should not need to explicitly call this routine, as it
169: is used internally within the minimization solvers.
171: TaoComputeJacobianState() is typically used within minimization
172: implementations, so most users would not generally call this routine
173: themselves.
175: Level: developer
177: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
179: @*/
180: PetscErrorCode TaoComputeJacobianState(Tao tao, Vec X, Mat J, Mat Jpre, Mat Jinv)
181: {
189: if (!tao->ops->computejacobianstate) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianState() first");
190: ++tao->njac_state;
191: PetscLogEventBegin(Tao_JacobianEval,tao,X,J,Jpre);
192: PetscStackPush("Tao user Jacobian(state) function");
193: (*tao->ops->computejacobianstate)(tao,X,J,Jpre,Jinv,tao->user_jac_stateP);
194: PetscStackPop;
195: PetscLogEventEnd(Tao_JacobianEval,tao,X,J,Jpre);
196: return(0);
197: }
199: /*@C
200: TaoComputeJacobianDesign - Computes the Jacobian matrix that has been
201: set with TaoSetJacobianDesignRoutine().
203: Collective on Tao
205: Input Parameters:
206: + solver - the Tao solver context
207: - xx - input vector
209: Output Parameters:
210: . H - Jacobian matrix
212: Notes:
213: Most users should not need to explicitly call this routine, as it
214: is used internally within the minimization solvers.
216: TaoComputeJacobianDesign() is typically used within minimization
217: implementations, so most users would not generally call this routine
218: themselves.
220: Level: developer
222: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianDesignRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
224: @*/
225: PetscErrorCode TaoComputeJacobianDesign(Tao tao, Vec X, Mat J)
226: {
234: if (!tao->ops->computejacobiandesign) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianDesign() first");
235: ++tao->njac_design;
236: PetscLogEventBegin(Tao_JacobianEval,tao,X,J,NULL);
237: PetscStackPush("Tao user Jacobian(design) function");
238: (*tao->ops->computejacobiandesign)(tao,X,J,tao->user_jac_designP);
239: PetscStackPop;
240: PetscLogEventEnd(Tao_JacobianEval,tao,X,J,NULL);
241: return(0);
242: }
244: /*@C
245: TaoSetJacobianRoutine - Sets the function to compute the Jacobian as well as the location to store the matrix.
247: Logically collective on Tao
249: Input Parameters:
250: + tao - the Tao context
251: . J - Matrix used for the jacobian
252: . Jpre - Matrix that will be used operated on by preconditioner, can be same as J
253: . jac - Jacobian evaluation routine
254: - ctx - [optional] user-defined context for private data for the
255: Jacobian evaluation routine (may be NULL)
257: Calling sequence of jac:
258: $ jac (Tao tao,Vec x,Mat *J,Mat *Jpre,void *ctx);
260: + tao - the Tao context
261: . x - input vector
262: . J - Jacobian matrix
263: . Jpre - preconditioner matrix, usually the same as J
264: - ctx - [optional] user-defined Jacobian context
266: Level: intermediate
268: @*/
269: PetscErrorCode TaoSetJacobianRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, void*), void *ctx)
270: {
274: if (J) {
277: }
278: if (Jpre) {
281: }
282: if (ctx) {
283: tao->user_jacP = ctx;
284: }
285: if (func) {
286: tao->ops->computejacobian = func;
287: }
288: if (J) {
289: PetscObjectReference((PetscObject)J);
290: MatDestroy(&tao->jacobian);
291: tao->jacobian = J;
292: }
293: if (Jpre) {
294: PetscObjectReference((PetscObject)Jpre);
295: MatDestroy(&tao->jacobian_pre);
296: tao->jacobian_pre=Jpre;
297: }
298: return(0);
299: }
301: /*@C
302: TaoSetJacobianStateRoutine - Sets the function to compute the Jacobian
303: (and its inverse) of the constraint function with respect to the state variables.
304: Used only for pde-constrained optimization.
306: Logically collective on Tao
308: Input Parameters:
309: + tao - the Tao context
310: . J - Matrix used for the jacobian
311: . Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J. Only used if Jinv is NULL
312: . Jinv - [optional] Matrix used to apply the inverse of the state jacobian. Use NULL to default to PETSc KSP solvers to apply the inverse.
313: . jac - Jacobian evaluation routine
314: - ctx - [optional] user-defined context for private data for the
315: Jacobian evaluation routine (may be NULL)
317: Calling sequence of jac:
318: $ jac (Tao tao,Vec x,Mat *J,Mat *Jpre,void *ctx);
320: + tao - the Tao context
321: . x - input vector
322: . J - Jacobian matrix
323: . Jpre - preconditioner matrix, usually the same as J
324: . Jinv - inverse of J
325: - ctx - [optional] user-defined Jacobian context
327: Level: intermediate
328: .seealse: TaoComputeJacobianState(), TaoSetJacobianDesignRoutine(), TaoSetStateDesignIS()
329: @*/
330: PetscErrorCode TaoSetJacobianStateRoutine(Tao tao, Mat J, Mat Jpre, Mat Jinv, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, Mat,void*), void *ctx)
331: {
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: . jac - 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 jac:
386: $ jac (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
394: Notes:
396: The function jac() takes Mat * as the matrix arguments rather than Mat.
397: This allows the Jacobian evaluation routine to replace A and/or B with a
398: completely new new matrix structure (not just different matrix elements)
399: when appropriate, for instance, if the nonzero structure is changing
400: throughout the global iterations.
402: Level: intermediate
403: .seealso: TaoComputeJacobianDesign(), TaoSetJacobianStateRoutine(), TaoSetStateDesignIS()
404: @*/
405: PetscErrorCode TaoSetJacobianDesignRoutine(Tao tao, Mat J, PetscErrorCode (*func)(Tao, Vec, Mat, void*), void *ctx)
406: {
411: if (J) {
414: }
415: if (ctx) {
416: tao->user_jac_designP = ctx;
417: }
418: if (func) {
419: tao->ops->computejacobiandesign = func;
420: }
421: if (J) {
422: PetscObjectReference((PetscObject)J);
423: MatDestroy(&tao->jacobian_design);
424: tao->jacobian_design = J;
425: }
426: return(0);
427: }
429: /*@
430: TaoSetStateDesignIS - Indicate to the Tao which variables in the
431: solution vector are state variables and which are design. Only applies to
432: pde-constrained optimization.
434: Logically Collective on Tao
436: Input Parameters:
437: + tao - The Tao context
438: . s_is - the index set corresponding to the state variables
439: - d_is - the index set corresponding to the design variables
441: Level: intermediate
443: .seealso: TaoSetJacobianStateRoutine(), TaoSetJacobianDesignRoutine()
444: @*/
445: PetscErrorCode TaoSetStateDesignIS(Tao tao, IS s_is, IS d_is)
446: {
450: PetscObjectReference((PetscObject)s_is);
451: ISDestroy(&tao->state_is);
452: tao->state_is = s_is;
453: PetscObjectReference((PetscObject)(d_is));
454: ISDestroy(&tao->design_is);
455: tao->design_is = d_is;
456: return(0);
457: }
459: /*@C
460: TaoComputeJacobianEquality - Computes the Jacobian matrix that has been
461: set with TaoSetJacobianEqualityRoutine().
463: Collective on Tao
465: Input Parameters:
466: + solver - the Tao solver context
467: - xx - input vector
469: Output Parameters:
470: + H - Jacobian matrix
471: - Hpre - Preconditioning matrix
473: Notes:
474: Most users should not need to explicitly call this routine, as it
475: is used internally within the minimization solvers.
477: Level: developer
479: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
481: @*/
482: PetscErrorCode TaoComputeJacobianEquality(Tao tao, Vec X, Mat J, Mat Jpre)
483: {
491: if (!tao->ops->computejacobianequality) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianEquality() first");
492: ++tao->njac_equality;
493: PetscLogEventBegin(Tao_JacobianEval,tao,X,J,Jpre);
494: PetscStackPush("Tao user Jacobian(equality) function");
495: (*tao->ops->computejacobianequality)(tao,X,J,Jpre,tao->user_jac_equalityP);
496: PetscStackPop;
497: PetscLogEventEnd(Tao_JacobianEval,tao,X,J,Jpre);
498: return(0);
499: }
501: /*@C
502: TaoComputeJacobianInequality - Computes the Jacobian matrix that has been
503: set with TaoSetJacobianInequalityRoutine().
505: Collective on Tao
507: Input Parameters:
508: + solver - the Tao solver context
509: - xx - input vector
511: Output Parameters:
512: + H - Jacobian matrix
513: - Hpre - Preconditioning matrix
515: Notes:
516: Most users should not need to explicitly call this routine, as it
517: is used internally within the minimization solvers.
519: Level: developer
521: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
523: @*/
524: PetscErrorCode TaoComputeJacobianInequality(Tao tao, Vec X, Mat J, Mat Jpre)
525: {
533: if (!tao->ops->computejacobianinequality) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianInequality() first");
534: ++tao->njac_inequality;
535: PetscLogEventBegin(Tao_JacobianEval,tao,X,J,Jpre);
536: PetscStackPush("Tao user Jacobian(inequality) function");
537: (*tao->ops->computejacobianinequality)(tao,X,J,Jpre,tao->user_jac_inequalityP);
538: PetscStackPop;
539: PetscLogEventEnd(Tao_JacobianEval,tao,X,J,Jpre);
540: return(0);
541: }
543: /*@C
544: TaoSetJacobianEqualityRoutine - Sets the function to compute the Jacobian
545: (and its inverse) of the constraint function with respect to the equality variables.
546: Used only for pde-constrained optimization.
548: Logically collective on Tao
550: Input Parameters:
551: + tao - the Tao context
552: . J - Matrix used for the jacobian
553: . Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.
554: . jac - Jacobian evaluation routine
555: - ctx - [optional] user-defined context for private data for the
556: Jacobian evaluation routine (may be NULL)
558: Calling sequence of jac:
559: $ jac (Tao tao,Vec x,Mat *J,Mat *Jpre,void *ctx);
561: + tao - the Tao context
562: . x - input vector
563: . J - Jacobian matrix
564: . Jpre - preconditioner matrix, usually the same as J
565: - ctx - [optional] user-defined Jacobian context
567: Level: intermediate
568: .seealse: TaoComputeJacobianEquality(), TaoSetJacobianDesignRoutine(), TaoSetEqualityDesignIS()
569: @*/
570: PetscErrorCode TaoSetJacobianEqualityRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat,void*), void *ctx)
571: {
576: if (J) {
579: }
580: if (Jpre) {
583: }
584: if (ctx) {
585: tao->user_jac_equalityP = ctx;
586: }
587: if (func) {
588: tao->ops->computejacobianequality = func;
589: }
590: if (J) {
591: PetscObjectReference((PetscObject)J);
592: MatDestroy(&tao->jacobian_equality);
593: tao->jacobian_equality = J;
594: }
595: if (Jpre) {
596: PetscObjectReference((PetscObject)Jpre);
597: MatDestroy(&tao->jacobian_equality_pre);
598: tao->jacobian_equality_pre=Jpre;
599: }
600: return(0);
601: }
603: /*@C
604: TaoSetJacobianInequalityRoutine - Sets the function to compute the Jacobian
605: (and its inverse) of the constraint function with respect to the inequality variables.
606: Used only for pde-constrained optimization.
608: Logically collective on Tao
610: Input Parameters:
611: + tao - the Tao context
612: . J - Matrix used for the jacobian
613: . Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.
614: . jac - Jacobian evaluation routine
615: - ctx - [optional] user-defined context for private data for the
616: Jacobian evaluation routine (may be NULL)
618: Calling sequence of jac:
619: $ jac (Tao tao,Vec x,Mat *J,Mat *Jpre,void *ctx);
621: + tao - the Tao context
622: . x - input vector
623: . J - Jacobian matrix
624: . Jpre - preconditioner matrix, usually the same as J
625: - ctx - [optional] user-defined Jacobian context
627: Level: intermediate
628: .seealse: TaoComputeJacobianInequality(), TaoSetJacobianDesignRoutine(), TaoSetInequalityDesignIS()
629: @*/
630: PetscErrorCode TaoSetJacobianInequalityRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat,void*), void *ctx)
631: {
635: if (J) {
638: }
639: if (Jpre) {
642: }
643: if (ctx) {
644: tao->user_jac_inequalityP = ctx;
645: }
646: if (func) {
647: tao->ops->computejacobianinequality = func;
648: }
649: if (J) {
650: PetscObjectReference((PetscObject)J);
651: MatDestroy(&tao->jacobian_inequality);
652: tao->jacobian_inequality = J;
653: }
654: if (Jpre) {
655: PetscObjectReference((PetscObject)Jpre);
656: MatDestroy(&tao->jacobian_inequality_pre);
657: tao->jacobian_inequality_pre=Jpre;
658: }
659: return(0);
660: }