Actual source code: taosolver_hj.c
petsc-3.11.4 2019-09-28
1: #include <petsc/private/taoimpl.h>
4: /*@C
5: TaoSetHessianRoutine - Sets the function to compute the Hessian as well as the location to store the matrix.
7: Logically collective on Tao
9: Input Parameters:
10: + tao - the Tao context
11: . H - Matrix used for the hessian
12: . Hpre - Matrix that will be used operated on by preconditioner, can be same as H
13: . func - Hessian evaluation routine
14: - ctx - [optional] user-defined context for private data for the
15: Hessian evaluation routine (may be NULL)
17: Calling sequence of func:
18: $ func(Tao tao,Vec x,Mat H,Mat Hpre,void *ctx);
20: + tao - the Tao context
21: . x - input vector
22: . H - Hessian matrix
23: . Hpre - preconditioner matrix, usually the same as H
24: - ctx - [optional] user-defined Hessian context
26: 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: PetscErrorCode TaoTestHessian(Tao tao)
62: {
63: Mat A,B,C,D,hessian;
64: Vec x = tao->solution;
65: PetscErrorCode ierr;
66: PetscReal nrm,gnorm;
67: PetscReal threshold = 1.e-5;
68: PetscInt m,n,M,N;
69: PetscBool complete_print = PETSC_FALSE,test = PETSC_FALSE,flg;
70: PetscViewer viewer,mviewer;
71: MPI_Comm comm;
72: PetscInt tabs;
73: static PetscBool directionsprinted = PETSC_FALSE;
74: PetscViewerFormat format;
77: PetscObjectOptionsBegin((PetscObject)tao);
78: PetscOptionsName("-tao_test_hessian","Compare hand-coded and finite difference Hessians","None",&test);
79: PetscOptionsReal("-tao_test_hessian", "Threshold for element difference between hand-coded and finite difference being meaningful","None",threshold,&threshold,NULL);
80: PetscOptionsViewer("-tao_test_hessian_view","View difference between hand-coded and finite difference Hessians element entries","None",&mviewer,&format,&complete_print);
81: PetscOptionsEnd();
82: if (!test) return(0);
84: PetscObjectGetComm((PetscObject)tao,&comm);
85: PetscViewerASCIIGetStdout(comm,&viewer);
86: PetscViewerASCIIGetTab(viewer, &tabs);
87: PetscViewerASCIISetTab(viewer, ((PetscObject)tao)->tablevel);
88: PetscViewerASCIIPrintf(viewer," ---------- Testing Hessian -------------\n");
89: if (!complete_print && !directionsprinted) {
90: PetscViewerASCIIPrintf(viewer," Run with -tao_test_hessian_view and optionally -tao_test_hessian <threshold> to show difference\n");
91: PetscViewerASCIIPrintf(viewer," of hand-coded and finite difference Hessian entries greater than <threshold>.\n");
92: }
93: if (!directionsprinted) {
94: PetscViewerASCIIPrintf(viewer," Testing hand-coded Hessian, if (for double precision runs) ||J - Jfd||_F/||J||_F is\n");
95: PetscViewerASCIIPrintf(viewer," O(1.e-8), the hand-coded Hessian is probably correct.\n");
96: directionsprinted = PETSC_TRUE;
97: }
98: if (complete_print) {
99: PetscViewerPushFormat(mviewer,format);
100: }
102: PetscObjectTypeCompare((PetscObject)tao->hessian,MATMFFD,&flg);
103: if (!flg) hessian = tao->hessian;
104: else hessian = tao->hessian_pre;
106: while (hessian) {
107: PetscObjectTypeCompareAny((PetscObject)hessian,&flg,MATSEQAIJ,MATMPIAIJ,MATSEQDENSE,MATMPIDENSE,MATSEQBAIJ,MATMPIBAIJ,MATSEQSBAIJ,MATMPIBAIJ,"");
108: if (flg) {
109: A = hessian;
110: PetscObjectReference((PetscObject)A);
111: } else {
112: MatComputeExplicitOperator(hessian,&A);
113: }
115: MatCreate(PetscObjectComm((PetscObject)A),&B);
116: MatGetSize(A,&M,&N);
117: MatGetLocalSize(A,&m,&n);
118: MatSetSizes(B,m,n,M,N);
119: MatSetType(B,((PetscObject)A)->type_name);
120: MatSetUp(B);
121: MatSetOption(B,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
123: TaoDefaultComputeHessian(tao,x,B,B,NULL);
125: MatDuplicate(B,MAT_COPY_VALUES,&D);
126: MatAYPX(D,-1.0,A,DIFFERENT_NONZERO_PATTERN);
127: MatNorm(D,NORM_FROBENIUS,&nrm);
128: MatNorm(A,NORM_FROBENIUS,&gnorm);
129: MatDestroy(&D);
130: if (!gnorm) gnorm = 1; /* just in case */
131: PetscViewerASCIIPrintf(viewer," ||H - Hfd||_F/||H||_F = %g, ||H - Hfd||_F = %g\n",(double)(nrm/gnorm),(double)nrm);
133: if (complete_print) {
134: PetscViewerASCIIPrintf(viewer," Hand-coded Hessian ----------\n");
135: MatView(A,mviewer);
136: PetscViewerASCIIPrintf(viewer," Finite difference Hessian ----------\n");
137: MatView(B,mviewer);
138: }
140: if (complete_print) {
141: PetscInt Istart, Iend, *ccols, bncols, cncols, j, row;
142: PetscScalar *cvals;
143: const PetscInt *bcols;
144: const PetscScalar *bvals;
146: MatAYPX(B,-1.0,A,DIFFERENT_NONZERO_PATTERN);
147: MatCreate(PetscObjectComm((PetscObject)A),&C);
148: MatSetSizes(C,m,n,M,N);
149: MatSetType(C,((PetscObject)A)->type_name);
150: MatSetUp(C);
151: MatSetOption(C,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
152: MatGetOwnershipRange(B,&Istart,&Iend);
154: for (row = Istart; row < Iend; row++) {
155: MatGetRow(B,row,&bncols,&bcols,&bvals);
156: PetscMalloc2(bncols,&ccols,bncols,&cvals);
157: for (j = 0, cncols = 0; j < bncols; j++) {
158: if (PetscAbsScalar(bvals[j]) > threshold) {
159: ccols[cncols] = bcols[j];
160: cvals[cncols] = bvals[j];
161: cncols += 1;
162: }
163: }
164: if (cncols) {
165: MatSetValues(C,1,&row,cncols,ccols,cvals,INSERT_VALUES);
166: }
167: MatRestoreRow(B,row,&bncols,&bcols,&bvals);
168: PetscFree2(ccols,cvals);
169: }
170: MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);
171: MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);
172: PetscViewerASCIIPrintf(viewer," Finite-difference minus hand-coded Hessian with tolerance %g ----------\n",(double)threshold);
173: MatView(C,mviewer);
174: MatDestroy(&C);
175: }
176: MatDestroy(&A);
177: MatDestroy(&B);
179: if (hessian != tao->hessian_pre) {
180: hessian = tao->hessian_pre;
181: PetscViewerASCIIPrintf(viewer," ---------- Testing Hessian for preconditioner -------------\n");
182: } else hessian = NULL;
183: }
184: if (complete_print) {
185: PetscViewerPopFormat(mviewer);
186: }
187: PetscViewerASCIISetTab(viewer,tabs);
188: return(0);
189: }
191: /*@C
192: TaoComputeHessian - Computes the Hessian matrix that has been
193: set with TaoSetHessianRoutine().
195: Collective on Tao
197: Input Parameters:
198: + tao - the Tao solver context
199: - X - input vector
201: Output Parameters:
202: + H - Hessian matrix
203: - Hpre - Preconditioning matrix
205: Options Database Keys:
206: + -tao_test_hessian - compare the user provided Hessian with one compute via finite differences to check for errors
207: . -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
208: - -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
210: Notes:
211: Most users should not need to explicitly call this routine, as it
212: is used internally within the minimization solvers.
214: TaoComputeHessian() is typically used within minimization
215: implementations, so most users would not generally call this routine
216: themselves.
218: Developer Notes:
219: The Hessian test mechanism follows SNESTestJacobian().
221: Level: developer
223: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetHessianRoutine()
224: @*/
225: PetscErrorCode TaoComputeHessian(Tao tao, Vec X, Mat H, Mat Hpre)
226: {
233: if (!tao->ops->computehessian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetHessianRoutine() first");
234: ++tao->nhess;
235: VecLockReadPush(X);
236: PetscLogEventBegin(TAO_HessianEval,tao,X,H,Hpre);
237: PetscStackPush("Tao user Hessian function");
238: (*tao->ops->computehessian)(tao,X,H,Hpre,tao->user_hessP);
239: PetscStackPop;
240: PetscLogEventEnd(TAO_HessianEval,tao,X,H,Hpre);
241: VecLockReadPop(X);
243: TaoTestHessian(tao);
244: return(0);
245: }
247: /*@C
248: TaoComputeJacobian - Computes the Jacobian matrix that has been
249: set with TaoSetJacobianRoutine().
251: Collective on Tao
253: Input Parameters:
254: + tao - the Tao solver context
255: - X - input vector
257: Output Parameters:
258: + J - Jacobian matrix
259: - Jpre - Preconditioning matrix
261: Notes:
262: Most users should not need to explicitly call this routine, as it
263: is used internally within the minimization solvers.
265: TaoComputeJacobian() is typically used within minimization
266: implementations, so most users would not generally call this routine
267: themselves.
269: Level: developer
271: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianRoutine()
272: @*/
273: PetscErrorCode TaoComputeJacobian(Tao tao, Vec X, Mat J, Mat Jpre)
274: {
281: if (!tao->ops->computejacobian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobian() first");
282: ++tao->njac;
283: VecLockReadPush(X);
284: PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
285: PetscStackPush("Tao user Jacobian function");
286: (*tao->ops->computejacobian)(tao,X,J,Jpre,tao->user_jacP);
287: PetscStackPop;
288: PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
289: VecLockReadPop(X);
290: return(0);
291: }
293: /*@C
294: TaoComputeResidualJacobian - Computes the least-squares residual Jacobian matrix that has been
295: set with TaoSetJacobianResidual().
297: Collective on Tao
299: Input Parameters:
300: + tao - the Tao solver context
301: - X - input vector
303: Output Parameters:
304: + J - Jacobian matrix
305: - Jpre - Preconditioning matrix
307: Notes:
308: Most users should not need to explicitly call this routine, as it
309: is used internally within the minimization solvers.
311: TaoComputeResidualJacobian() is typically used within least-squares
312: implementations, so most users would not generally call this routine
313: themselves.
315: Level: developer
317: .seealso: TaoComputeResidual(), TaoSetJacobianResidual()
318: @*/
319: PetscErrorCode TaoComputeResidualJacobian(Tao tao, Vec X, Mat J, Mat Jpre)
320: {
327: if (!tao->ops->computeresidualjacobian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetResidualJacobian() first");
328: ++tao->njac;
329: VecLockReadPush(X);
330: PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
331: PetscStackPush("Tao user least-squares residual Jacobian function");
332: (*tao->ops->computeresidualjacobian)(tao,X,J,Jpre,tao->user_lsjacP);
333: PetscStackPop;
334: PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
335: VecLockReadPop(X);
336: return(0);
337: }
339: /*@C
340: TaoComputeJacobianState - Computes the Jacobian matrix that has been
341: set with TaoSetJacobianStateRoutine().
343: Collective on Tao
345: Input Parameters:
346: + tao - the Tao solver context
347: - X - input vector
349: Output Parameters:
350: + Jpre - Jacobian matrix
351: - Jinv - Preconditioning matrix
353: Notes:
354: Most users should not need to explicitly call this routine, as it
355: is used internally within the minimization solvers.
357: TaoComputeJacobianState() is typically used within minimization
358: implementations, so most users would not generally call this routine
359: themselves.
361: Level: developer
363: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
364: @*/
365: PetscErrorCode TaoComputeJacobianState(Tao tao, Vec X, Mat J, Mat Jpre, Mat Jinv)
366: {
373: if (!tao->ops->computejacobianstate) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianState() first");
374: ++tao->njac_state;
375: VecLockReadPush(X);
376: PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
377: PetscStackPush("Tao user Jacobian(state) function");
378: (*tao->ops->computejacobianstate)(tao,X,J,Jpre,Jinv,tao->user_jac_stateP);
379: PetscStackPop;
380: PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
381: VecLockReadPop(X);
382: return(0);
383: }
385: /*@C
386: TaoComputeJacobianDesign - Computes the Jacobian matrix that has been
387: set with TaoSetJacobianDesignRoutine().
389: Collective on Tao
391: Input Parameters:
392: + tao - the Tao solver context
393: - X - input vector
395: Output Parameters:
396: . J - Jacobian matrix
398: Notes:
399: Most users should not need to explicitly call this routine, as it
400: is used internally within the minimization solvers.
402: TaoComputeJacobianDesign() is typically used within minimization
403: implementations, so most users would not generally call this routine
404: themselves.
406: Level: developer
408: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianDesignRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
409: @*/
410: PetscErrorCode TaoComputeJacobianDesign(Tao tao, Vec X, Mat J)
411: {
418: if (!tao->ops->computejacobiandesign) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianDesign() first");
419: ++tao->njac_design;
420: VecLockReadPush(X);
421: PetscLogEventBegin(TAO_JacobianEval,tao,X,J,NULL);
422: PetscStackPush("Tao user Jacobian(design) function");
423: (*tao->ops->computejacobiandesign)(tao,X,J,tao->user_jac_designP);
424: PetscStackPop;
425: PetscLogEventEnd(TAO_JacobianEval,tao,X,J,NULL);
426: VecLockReadPop(X);
427: return(0);
428: }
430: /*@C
431: TaoSetJacobianRoutine - Sets the function to compute the Jacobian as well as the location to store the matrix.
433: Logically collective on Tao
435: Input Parameters:
436: + tao - the Tao context
437: . J - Matrix used for the jacobian
438: . Jpre - Matrix that will be used operated on by preconditioner, can be same as J
439: . func - Jacobian evaluation routine
440: - ctx - [optional] user-defined context for private data for the
441: Jacobian evaluation routine (may be NULL)
443: Calling sequence of func:
444: $ func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);
446: + tao - the Tao context
447: . x - input vector
448: . J - Jacobian matrix
449: . Jpre - preconditioning matrix, usually the same as J
450: - ctx - [optional] user-defined Jacobian context
452: Level: intermediate
453: @*/
454: PetscErrorCode TaoSetJacobianRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, void*), void *ctx)
455: {
460: if (J) {
463: }
464: if (Jpre) {
467: }
468: if (ctx) {
469: tao->user_jacP = ctx;
470: }
471: if (func) {
472: tao->ops->computejacobian = func;
473: }
474: if (J) {
475: PetscObjectReference((PetscObject)J);
476: MatDestroy(&tao->jacobian);
477: tao->jacobian = J;
478: }
479: if (Jpre) {
480: PetscObjectReference((PetscObject)Jpre);
481: MatDestroy(&tao->jacobian_pre);
482: tao->jacobian_pre=Jpre;
483: }
484: return(0);
485: }
487: /*@C
488: TaoSetJacobianResidualRoutine - Sets the function to compute the least-squares residual Jacobian as well as the
489: location to store the matrix.
491: Logically collective on Tao
493: Input Parameters:
494: + tao - the Tao context
495: . J - Matrix used for the jacobian
496: . Jpre - Matrix that will be used operated on by preconditioner, can be same as J
497: . func - Jacobian evaluation routine
498: - ctx - [optional] user-defined context for private data for the
499: Jacobian evaluation routine (may be NULL)
501: Calling sequence of func:
502: $ func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);
504: + tao - the Tao context
505: . x - input vector
506: . J - Jacobian matrix
507: . Jpre - preconditioning matrix, usually the same as J
508: - ctx - [optional] user-defined Jacobian context
510: Level: intermediate
511: @*/
512: PetscErrorCode TaoSetJacobianResidualRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, void*), void *ctx)
513: {
518: if (J) {
521: }
522: if (Jpre) {
525: }
526: if (ctx) {
527: tao->user_lsjacP = ctx;
528: }
529: if (func) {
530: tao->ops->computeresidualjacobian = func;
531: }
532: if (J) {
533: PetscObjectReference((PetscObject)J);
534: MatDestroy(&tao->ls_jac);
535: tao->ls_jac = J;
536: }
537: if (Jpre) {
538: PetscObjectReference((PetscObject)Jpre);
539: MatDestroy(&tao->ls_jac_pre);
540: tao->ls_jac_pre=Jpre;
541: }
542: return(0);
543: }
545: /*@C
546: TaoSetJacobianStateRoutine - Sets the function to compute the Jacobian
547: (and its inverse) of the constraint function with respect to the state variables.
548: Used only for pde-constrained optimization.
550: Logically collective on Tao
552: Input Parameters:
553: + tao - the Tao context
554: . J - Matrix used for the jacobian
555: . Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J. Only used if Jinv is NULL
556: . Jinv - [optional] Matrix used to apply the inverse of the state jacobian. Use NULL to default to PETSc KSP solvers to apply the inverse.
557: . func - Jacobian evaluation routine
558: - ctx - [optional] user-defined context for private data for the
559: Jacobian evaluation routine (may be NULL)
561: Calling sequence of func:
562: $ func(Tao tao,Vec x,Mat J,Mat Jpre,Mat Jinv,void *ctx);
564: + tao - the Tao context
565: . x - input vector
566: . J - Jacobian matrix
567: . Jpre - preconditioner matrix, usually the same as J
568: . Jinv - inverse of J
569: - ctx - [optional] user-defined Jacobian context
571: Level: intermediate
572: .seealso: TaoComputeJacobianState(), TaoSetJacobianDesignRoutine(), TaoSetStateDesignIS()
573: @*/
574: PetscErrorCode TaoSetJacobianStateRoutine(Tao tao, Mat J, Mat Jpre, Mat Jinv, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, Mat, void*), void *ctx)
575: {
580: if (J) {
583: }
584: if (Jpre) {
587: }
588: if (Jinv) {
591: }
592: if (ctx) {
593: tao->user_jac_stateP = ctx;
594: }
595: if (func) {
596: tao->ops->computejacobianstate = func;
597: }
598: if (J) {
599: PetscObjectReference((PetscObject)J);
600: MatDestroy(&tao->jacobian_state);
601: tao->jacobian_state = J;
602: }
603: if (Jpre) {
604: PetscObjectReference((PetscObject)Jpre);
605: MatDestroy(&tao->jacobian_state_pre);
606: tao->jacobian_state_pre=Jpre;
607: }
608: if (Jinv) {
609: PetscObjectReference((PetscObject)Jinv);
610: MatDestroy(&tao->jacobian_state_inv);
611: tao->jacobian_state_inv=Jinv;
612: }
613: return(0);
614: }
616: /*@C
617: TaoSetJacobianDesignRoutine - Sets the function to compute the Jacobian of
618: the constraint function with respect to the design variables. Used only for
619: pde-constrained optimization.
621: Logically collective on Tao
623: Input Parameters:
624: + tao - the Tao context
625: . J - Matrix used for the jacobian
626: . func - Jacobian evaluation routine
627: - ctx - [optional] user-defined context for private data for the
628: Jacobian evaluation routine (may be NULL)
630: Calling sequence of func:
631: $ func(Tao tao,Vec x,Mat J,void *ctx);
633: + tao - the Tao context
634: . x - input vector
635: . J - Jacobian matrix
636: - ctx - [optional] user-defined Jacobian context
638: Level: intermediate
640: .seealso: TaoComputeJacobianDesign(), TaoSetJacobianStateRoutine(), TaoSetStateDesignIS()
641: @*/
642: PetscErrorCode TaoSetJacobianDesignRoutine(Tao tao, Mat J, PetscErrorCode (*func)(Tao, Vec, Mat, void*), void *ctx)
643: {
648: if (J) {
651: }
652: if (ctx) {
653: tao->user_jac_designP = ctx;
654: }
655: if (func) {
656: tao->ops->computejacobiandesign = func;
657: }
658: if (J) {
659: PetscObjectReference((PetscObject)J);
660: MatDestroy(&tao->jacobian_design);
661: tao->jacobian_design = J;
662: }
663: return(0);
664: }
666: /*@
667: TaoSetStateDesignIS - Indicate to the Tao which variables in the
668: solution vector are state variables and which are design. Only applies to
669: pde-constrained optimization.
671: Logically Collective on Tao
673: Input Parameters:
674: + tao - The Tao context
675: . s_is - the index set corresponding to the state variables
676: - d_is - the index set corresponding to the design variables
678: Level: intermediate
680: .seealso: TaoSetJacobianStateRoutine(), TaoSetJacobianDesignRoutine()
681: @*/
682: PetscErrorCode TaoSetStateDesignIS(Tao tao, IS s_is, IS d_is)
683: {
687: PetscObjectReference((PetscObject)s_is);
688: ISDestroy(&tao->state_is);
689: tao->state_is = s_is;
690: PetscObjectReference((PetscObject)(d_is));
691: ISDestroy(&tao->design_is);
692: tao->design_is = d_is;
693: return(0);
694: }
696: /*@C
697: TaoComputeJacobianEquality - Computes the Jacobian matrix that has been
698: set with TaoSetJacobianEqualityRoutine().
700: Collective on Tao
702: Input Parameters:
703: + tao - the Tao solver context
704: - X - input vector
706: Output Parameters:
707: + J - Jacobian matrix
708: - Jpre - Preconditioning matrix
710: Notes:
711: Most users should not need to explicitly call this routine, as it
712: is used internally within the minimization solvers.
714: Level: developer
716: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
717: @*/
718: PetscErrorCode TaoComputeJacobianEquality(Tao tao, Vec X, Mat J, Mat Jpre)
719: {
726: if (!tao->ops->computejacobianequality) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianEquality() first");
727: ++tao->njac_equality;
728: VecLockReadPush(X);
729: PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
730: PetscStackPush("Tao user Jacobian(equality) function");
731: (*tao->ops->computejacobianequality)(tao,X,J,Jpre,tao->user_jac_equalityP);
732: PetscStackPop;
733: PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
734: VecLockReadPop(X);
735: return(0);
736: }
738: /*@C
739: TaoComputeJacobianInequality - Computes the Jacobian matrix that has been
740: set with TaoSetJacobianInequalityRoutine().
742: Collective on Tao
744: Input Parameters:
745: + tao - the Tao solver context
746: - X - input vector
748: Output Parameters:
749: + J - Jacobian matrix
750: - Jpre - Preconditioning matrix
752: Notes:
753: Most users should not need to explicitly call this routine, as it
754: is used internally within the minimization solvers.
756: Level: developer
758: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
759: @*/
760: PetscErrorCode TaoComputeJacobianInequality(Tao tao, Vec X, Mat J, Mat Jpre)
761: {
768: if (!tao->ops->computejacobianinequality) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianInequality() first");
769: ++tao->njac_inequality;
770: VecLockReadPush(X);
771: PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
772: PetscStackPush("Tao user Jacobian(inequality) function");
773: (*tao->ops->computejacobianinequality)(tao,X,J,Jpre,tao->user_jac_inequalityP);
774: PetscStackPop;
775: PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
776: VecLockReadPop(X);
777: return(0);
778: }
780: /*@C
781: TaoSetJacobianEqualityRoutine - Sets the function to compute the Jacobian
782: (and its inverse) of the constraint function with respect to the equality variables.
783: Used only for pde-constrained optimization.
785: Logically collective on Tao
787: Input Parameters:
788: + tao - the Tao context
789: . J - Matrix used for the jacobian
790: . Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.
791: . func - Jacobian evaluation routine
792: - ctx - [optional] user-defined context for private data for the
793: Jacobian evaluation routine (may be NULL)
795: Calling sequence of func:
796: $ func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);
798: + tao - the Tao context
799: . x - input vector
800: . J - Jacobian matrix
801: . Jpre - preconditioner matrix, usually the same as J
802: - ctx - [optional] user-defined Jacobian context
804: Level: intermediate
806: .seealso: TaoComputeJacobianEquality(), TaoSetJacobianDesignRoutine(), TaoSetEqualityDesignIS()
807: @*/
808: PetscErrorCode TaoSetJacobianEqualityRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat,void*), void *ctx)
809: {
814: if (J) {
817: }
818: if (Jpre) {
821: }
822: if (ctx) {
823: tao->user_jac_equalityP = ctx;
824: }
825: if (func) {
826: tao->ops->computejacobianequality = func;
827: }
828: if (J) {
829: PetscObjectReference((PetscObject)J);
830: MatDestroy(&tao->jacobian_equality);
831: tao->jacobian_equality = J;
832: }
833: if (Jpre) {
834: PetscObjectReference((PetscObject)Jpre);
835: MatDestroy(&tao->jacobian_equality_pre);
836: tao->jacobian_equality_pre=Jpre;
837: }
838: return(0);
839: }
841: /*@C
842: TaoSetJacobianInequalityRoutine - Sets the function to compute the Jacobian
843: (and its inverse) of the constraint function with respect to the inequality variables.
844: Used only for pde-constrained optimization.
846: Logically collective on Tao
848: Input Parameters:
849: + tao - the Tao context
850: . J - Matrix used for the jacobian
851: . Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.
852: . func - Jacobian evaluation routine
853: - ctx - [optional] user-defined context for private data for the
854: Jacobian evaluation routine (may be NULL)
856: Calling sequence of func:
857: $ func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);
859: + tao - the Tao context
860: . x - input vector
861: . J - Jacobian matrix
862: . Jpre - preconditioner matrix, usually the same as J
863: - ctx - [optional] user-defined Jacobian context
865: Level: intermediate
867: .seealso: TaoComputeJacobianInequality(), TaoSetJacobianDesignRoutine(), TaoSetInequalityDesignIS()
868: @*/
869: PetscErrorCode TaoSetJacobianInequalityRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat,void*), void *ctx)
870: {
875: if (J) {
878: }
879: if (Jpre) {
882: }
883: if (ctx) {
884: tao->user_jac_inequalityP = ctx;
885: }
886: if (func) {
887: tao->ops->computejacobianinequality = func;
888: }
889: if (J) {
890: PetscObjectReference((PetscObject)J);
891: MatDestroy(&tao->jacobian_inequality);
892: tao->jacobian_inequality = J;
893: }
894: if (Jpre) {
895: PetscObjectReference((PetscObject)Jpre);
896: MatDestroy(&tao->jacobian_inequality_pre);
897: tao->jacobian_inequality_pre=Jpre;
898: }
899: return(0);
900: }