Actual source code: taosolver_hj.c
1: #include <petsc/private/taoimpl.h>
3: /*@C
4: TaoSetHessianRoutine - Sets the function to compute the Hessian as well as the location to store the matrix.
6: Logically collective on Tao
8: Input Parameters:
9: + tao - the Tao context
10: . H - Matrix used for the hessian
11: . Hpre - Matrix that will be used operated on by preconditioner, can be same as H
12: . func - Hessian evaluation routine
13: - ctx - [optional] user-defined context for private data for the
14: Hessian evaluation routine (may be NULL)
16: Calling sequence of func:
17: $ func(Tao tao,Vec x,Mat H,Mat Hpre,void *ctx);
19: + tao - the Tao context
20: . x - input vector
21: . H - Hessian matrix
22: . Hpre - preconditioner matrix, usually the same as H
23: - ctx - [optional] user-defined Hessian context
25: Level: beginner
26: @*/
27: PetscErrorCode TaoSetHessianRoutine(Tao tao, Mat H, Mat Hpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, void*), void *ctx)
28: {
33: if (H) {
36: }
37: if (Hpre) {
40: }
41: if (ctx) {
42: tao->user_hessP = ctx;
43: }
44: if (func) {
45: tao->ops->computehessian = func;
46: }
47: if (H) {
48: PetscObjectReference((PetscObject)H);
49: MatDestroy(&tao->hessian);
50: tao->hessian = H;
51: }
52: if (Hpre) {
53: PetscObjectReference((PetscObject)Hpre);
54: MatDestroy(&tao->hessian_pre);
55: tao->hessian_pre = Hpre;
56: }
57: return(0);
58: }
60: PetscErrorCode TaoTestHessian(Tao tao)
61: {
62: Mat A,B,C,D,hessian;
63: Vec x = tao->solution;
64: PetscErrorCode ierr;
65: PetscReal nrm,gnorm;
66: PetscReal threshold = 1.e-5;
67: PetscInt m,n,M,N;
68: PetscBool complete_print = PETSC_FALSE,test = PETSC_FALSE,flg;
69: PetscViewer viewer,mviewer;
70: MPI_Comm comm;
71: PetscInt tabs;
72: static PetscBool directionsprinted = PETSC_FALSE;
73: PetscViewerFormat format;
76: PetscObjectOptionsBegin((PetscObject)tao);
77: PetscOptionsName("-tao_test_hessian","Compare hand-coded and finite difference Hessians","None",&test);
78: PetscOptionsReal("-tao_test_hessian", "Threshold for element difference between hand-coded and finite difference being meaningful","None",threshold,&threshold,NULL);
79: PetscOptionsViewer("-tao_test_hessian_view","View difference between hand-coded and finite difference Hessians element entries","None",&mviewer,&format,&complete_print);
80: PetscOptionsEnd();
81: if (!test) return(0);
83: PetscObjectGetComm((PetscObject)tao,&comm);
84: PetscViewerASCIIGetStdout(comm,&viewer);
85: PetscViewerASCIIGetTab(viewer, &tabs);
86: PetscViewerASCIISetTab(viewer, ((PetscObject)tao)->tablevel);
87: PetscViewerASCIIPrintf(viewer," ---------- Testing Hessian -------------\n");
88: if (!complete_print && !directionsprinted) {
89: PetscViewerASCIIPrintf(viewer," Run with -tao_test_hessian_view and optionally -tao_test_hessian <threshold> to show difference\n");
90: PetscViewerASCIIPrintf(viewer," of hand-coded and finite difference Hessian entries greater than <threshold>.\n");
91: }
92: if (!directionsprinted) {
93: PetscViewerASCIIPrintf(viewer," Testing hand-coded Hessian, if (for double precision runs) ||J - Jfd||_F/||J||_F is\n");
94: PetscViewerASCIIPrintf(viewer," O(1.e-8), the hand-coded Hessian is probably correct.\n");
95: directionsprinted = PETSC_TRUE;
96: }
97: if (complete_print) {
98: PetscViewerPushFormat(mviewer,format);
99: }
101: PetscObjectTypeCompare((PetscObject)tao->hessian,MATMFFD,&flg);
102: if (!flg) hessian = tao->hessian;
103: else hessian = tao->hessian_pre;
105: while (hessian) {
106: PetscObjectBaseTypeCompareAny((PetscObject)hessian,&flg,MATSEQAIJ,MATMPIAIJ,MATSEQDENSE,MATMPIDENSE,MATSEQBAIJ,MATMPIBAIJ,MATSEQSBAIJ,MATMPIBAIJ,"");
107: if (flg) {
108: A = hessian;
109: PetscObjectReference((PetscObject)A);
110: } else {
111: MatComputeOperator(hessian,MATAIJ,&A);
112: }
114: MatCreate(PetscObjectComm((PetscObject)A),&B);
115: MatGetSize(A,&M,&N);
116: MatGetLocalSize(A,&m,&n);
117: MatSetSizes(B,m,n,M,N);
118: MatSetType(B,((PetscObject)A)->type_name);
119: MatSetUp(B);
120: MatSetOption(B,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
122: TaoDefaultComputeHessian(tao,x,B,B,NULL);
124: MatDuplicate(B,MAT_COPY_VALUES,&D);
125: MatAYPX(D,-1.0,A,DIFFERENT_NONZERO_PATTERN);
126: MatNorm(D,NORM_FROBENIUS,&nrm);
127: MatNorm(A,NORM_FROBENIUS,&gnorm);
128: MatDestroy(&D);
129: if (!gnorm) gnorm = 1; /* just in case */
130: PetscViewerASCIIPrintf(viewer," ||H - Hfd||_F/||H||_F = %g, ||H - Hfd||_F = %g\n",(double)(nrm/gnorm),(double)nrm);
132: if (complete_print) {
133: PetscViewerASCIIPrintf(viewer," Hand-coded Hessian ----------\n");
134: MatView(A,mviewer);
135: PetscViewerASCIIPrintf(viewer," Finite difference Hessian ----------\n");
136: MatView(B,mviewer);
137: }
139: if (complete_print) {
140: PetscInt Istart, Iend, *ccols, bncols, cncols, j, row;
141: PetscScalar *cvals;
142: const PetscInt *bcols;
143: const PetscScalar *bvals;
145: MatAYPX(B,-1.0,A,DIFFERENT_NONZERO_PATTERN);
146: MatCreate(PetscObjectComm((PetscObject)A),&C);
147: MatSetSizes(C,m,n,M,N);
148: MatSetType(C,((PetscObject)A)->type_name);
149: MatSetUp(C);
150: MatSetOption(C,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
151: MatGetOwnershipRange(B,&Istart,&Iend);
153: for (row = Istart; row < Iend; row++) {
154: MatGetRow(B,row,&bncols,&bcols,&bvals);
155: PetscMalloc2(bncols,&ccols,bncols,&cvals);
156: for (j = 0, cncols = 0; j < bncols; j++) {
157: if (PetscAbsScalar(bvals[j]) > threshold) {
158: ccols[cncols] = bcols[j];
159: cvals[cncols] = bvals[j];
160: cncols += 1;
161: }
162: }
163: if (cncols) {
164: MatSetValues(C,1,&row,cncols,ccols,cvals,INSERT_VALUES);
165: }
166: MatRestoreRow(B,row,&bncols,&bcols,&bvals);
167: PetscFree2(ccols,cvals);
168: }
169: MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);
170: MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);
171: PetscViewerASCIIPrintf(viewer," Finite-difference minus hand-coded Hessian with tolerance %g ----------\n",(double)threshold);
172: MatView(C,mviewer);
173: MatDestroy(&C);
174: }
175: MatDestroy(&A);
176: MatDestroy(&B);
178: if (hessian != tao->hessian_pre) {
179: hessian = tao->hessian_pre;
180: PetscViewerASCIIPrintf(viewer," ---------- Testing Hessian for preconditioner -------------\n");
181: } else hessian = NULL;
182: }
183: if (complete_print) {
184: PetscViewerPopFormat(mviewer);
185: PetscViewerDestroy(&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: + J - Jacobian matrix
351: . Jpre - Preconditioning matrix
352: - Jinv -
354: Notes:
355: Most users should not need to explicitly call this routine, as it
356: is used internally within the minimization solvers.
358: TaoComputeJacobianState() is typically used within minimization
359: implementations, so most users would not generally call this routine
360: themselves.
362: Level: developer
364: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
365: @*/
366: PetscErrorCode TaoComputeJacobianState(Tao tao, Vec X, Mat J, Mat Jpre, Mat Jinv)
367: {
374: if (!tao->ops->computejacobianstate) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianState() first");
375: ++tao->njac_state;
376: VecLockReadPush(X);
377: PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
378: PetscStackPush("Tao user Jacobian(state) function");
379: (*tao->ops->computejacobianstate)(tao,X,J,Jpre,Jinv,tao->user_jac_stateP);
380: PetscStackPop;
381: PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
382: VecLockReadPop(X);
383: return(0);
384: }
386: /*@C
387: TaoComputeJacobianDesign - Computes the Jacobian matrix that has been
388: set with TaoSetJacobianDesignRoutine().
390: Collective on Tao
392: Input Parameters:
393: + tao - the Tao solver context
394: - X - input vector
396: Output Parameters:
397: . J - Jacobian matrix
399: Notes:
400: Most users should not need to explicitly call this routine, as it
401: is used internally within the minimization solvers.
403: TaoComputeJacobianDesign() is typically used within minimization
404: implementations, so most users would not generally call this routine
405: themselves.
407: Level: developer
409: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianDesignRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
410: @*/
411: PetscErrorCode TaoComputeJacobianDesign(Tao tao, Vec X, Mat J)
412: {
419: if (!tao->ops->computejacobiandesign) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianDesign() first");
420: ++tao->njac_design;
421: VecLockReadPush(X);
422: PetscLogEventBegin(TAO_JacobianEval,tao,X,J,NULL);
423: PetscStackPush("Tao user Jacobian(design) function");
424: (*tao->ops->computejacobiandesign)(tao,X,J,tao->user_jac_designP);
425: PetscStackPop;
426: PetscLogEventEnd(TAO_JacobianEval,tao,X,J,NULL);
427: VecLockReadPop(X);
428: return(0);
429: }
431: /*@C
432: TaoSetJacobianRoutine - Sets the function to compute the Jacobian as well as the location to store the matrix.
434: Logically collective on Tao
436: Input Parameters:
437: + tao - the Tao context
438: . J - Matrix used for the jacobian
439: . Jpre - Matrix that will be used operated on by preconditioner, can be same as J
440: . func - Jacobian evaluation routine
441: - ctx - [optional] user-defined context for private data for the
442: Jacobian evaluation routine (may be NULL)
444: Calling sequence of func:
445: $ func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);
447: + tao - the Tao context
448: . x - input vector
449: . J - Jacobian matrix
450: . Jpre - preconditioning matrix, usually the same as J
451: - ctx - [optional] user-defined Jacobian context
453: Level: intermediate
454: @*/
455: PetscErrorCode TaoSetJacobianRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, void*), void *ctx)
456: {
461: if (J) {
464: }
465: if (Jpre) {
468: }
469: if (ctx) {
470: tao->user_jacP = ctx;
471: }
472: if (func) {
473: tao->ops->computejacobian = func;
474: }
475: if (J) {
476: PetscObjectReference((PetscObject)J);
477: MatDestroy(&tao->jacobian);
478: tao->jacobian = J;
479: }
480: if (Jpre) {
481: PetscObjectReference((PetscObject)Jpre);
482: MatDestroy(&tao->jacobian_pre);
483: tao->jacobian_pre=Jpre;
484: }
485: return(0);
486: }
488: /*@C
489: TaoSetJacobianResidualRoutine - Sets the function to compute the least-squares residual Jacobian as well as the
490: location to store the matrix.
492: Logically collective on Tao
494: Input Parameters:
495: + tao - the Tao context
496: . J - Matrix used for the jacobian
497: . Jpre - Matrix that will be used operated on by preconditioner, can be same as J
498: . func - Jacobian evaluation routine
499: - ctx - [optional] user-defined context for private data for the
500: Jacobian evaluation routine (may be NULL)
502: Calling sequence of func:
503: $ func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);
505: + tao - the Tao context
506: . x - input vector
507: . J - Jacobian matrix
508: . Jpre - preconditioning matrix, usually the same as J
509: - ctx - [optional] user-defined Jacobian context
511: Level: intermediate
512: @*/
513: PetscErrorCode TaoSetJacobianResidualRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, void*), void *ctx)
514: {
519: if (J) {
522: }
523: if (Jpre) {
526: }
527: if (ctx) {
528: tao->user_lsjacP = ctx;
529: }
530: if (func) {
531: tao->ops->computeresidualjacobian = func;
532: }
533: if (J) {
534: PetscObjectReference((PetscObject)J);
535: MatDestroy(&tao->ls_jac);
536: tao->ls_jac = J;
537: }
538: if (Jpre) {
539: PetscObjectReference((PetscObject)Jpre);
540: MatDestroy(&tao->ls_jac_pre);
541: tao->ls_jac_pre=Jpre;
542: }
543: return(0);
544: }
546: /*@C
547: TaoSetJacobianStateRoutine - Sets the function to compute the Jacobian
548: (and its inverse) of the constraint function with respect to the state variables.
549: Used only for pde-constrained optimization.
551: Logically collective on Tao
553: Input Parameters:
554: + tao - the Tao context
555: . J - Matrix used for the jacobian
556: . Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J. Only used if Jinv is NULL
557: . Jinv - [optional] Matrix used to apply the inverse of the state jacobian. Use NULL to default to PETSc KSP solvers to apply the inverse.
558: . func - Jacobian evaluation routine
559: - ctx - [optional] user-defined context for private data for the
560: Jacobian evaluation routine (may be NULL)
562: Calling sequence of func:
563: $ func(Tao tao,Vec x,Mat J,Mat Jpre,Mat Jinv,void *ctx);
565: + tao - the Tao context
566: . x - input vector
567: . J - Jacobian matrix
568: . Jpre - preconditioner matrix, usually the same as J
569: . Jinv - inverse of J
570: - ctx - [optional] user-defined Jacobian context
572: Level: intermediate
573: .seealso: TaoComputeJacobianState(), TaoSetJacobianDesignRoutine(), TaoSetStateDesignIS()
574: @*/
575: PetscErrorCode TaoSetJacobianStateRoutine(Tao tao, Mat J, Mat Jpre, Mat Jinv, PetscErrorCode (*func)(Tao, Vec, Mat, Mat, Mat, void*), void *ctx)
576: {
581: if (J) {
584: }
585: if (Jpre) {
588: }
589: if (Jinv) {
592: }
593: if (ctx) {
594: tao->user_jac_stateP = ctx;
595: }
596: if (func) {
597: tao->ops->computejacobianstate = func;
598: }
599: if (J) {
600: PetscObjectReference((PetscObject)J);
601: MatDestroy(&tao->jacobian_state);
602: tao->jacobian_state = J;
603: }
604: if (Jpre) {
605: PetscObjectReference((PetscObject)Jpre);
606: MatDestroy(&tao->jacobian_state_pre);
607: tao->jacobian_state_pre=Jpre;
608: }
609: if (Jinv) {
610: PetscObjectReference((PetscObject)Jinv);
611: MatDestroy(&tao->jacobian_state_inv);
612: tao->jacobian_state_inv=Jinv;
613: }
614: return(0);
615: }
617: /*@C
618: TaoSetJacobianDesignRoutine - Sets the function to compute the Jacobian of
619: the constraint function with respect to the design variables. Used only for
620: pde-constrained optimization.
622: Logically collective on Tao
624: Input Parameters:
625: + tao - the Tao context
626: . J - Matrix used for the jacobian
627: . func - Jacobian evaluation routine
628: - ctx - [optional] user-defined context for private data for the
629: Jacobian evaluation routine (may be NULL)
631: Calling sequence of func:
632: $ func(Tao tao,Vec x,Mat J,void *ctx);
634: + tao - the Tao context
635: . x - input vector
636: . J - Jacobian matrix
637: - ctx - [optional] user-defined Jacobian context
639: Level: intermediate
641: .seealso: TaoComputeJacobianDesign(), TaoSetJacobianStateRoutine(), TaoSetStateDesignIS()
642: @*/
643: PetscErrorCode TaoSetJacobianDesignRoutine(Tao tao, Mat J, PetscErrorCode (*func)(Tao, Vec, Mat, void*), void *ctx)
644: {
649: if (J) {
652: }
653: if (ctx) {
654: tao->user_jac_designP = ctx;
655: }
656: if (func) {
657: tao->ops->computejacobiandesign = func;
658: }
659: if (J) {
660: PetscObjectReference((PetscObject)J);
661: MatDestroy(&tao->jacobian_design);
662: tao->jacobian_design = J;
663: }
664: return(0);
665: }
667: /*@
668: TaoSetStateDesignIS - Indicate to the Tao which variables in the
669: solution vector are state variables and which are design. Only applies to
670: pde-constrained optimization.
672: Logically Collective on Tao
674: Input Parameters:
675: + tao - The Tao context
676: . s_is - the index set corresponding to the state variables
677: - d_is - the index set corresponding to the design variables
679: Level: intermediate
681: .seealso: TaoSetJacobianStateRoutine(), TaoSetJacobianDesignRoutine()
682: @*/
683: PetscErrorCode TaoSetStateDesignIS(Tao tao, IS s_is, IS d_is)
684: {
688: PetscObjectReference((PetscObject)s_is);
689: ISDestroy(&tao->state_is);
690: tao->state_is = s_is;
691: PetscObjectReference((PetscObject)(d_is));
692: ISDestroy(&tao->design_is);
693: tao->design_is = d_is;
694: return(0);
695: }
697: /*@C
698: TaoComputeJacobianEquality - Computes the Jacobian matrix that has been
699: set with TaoSetJacobianEqualityRoutine().
701: Collective on Tao
703: Input Parameters:
704: + tao - the Tao solver context
705: - X - input vector
707: Output Parameters:
708: + J - Jacobian matrix
709: - Jpre - Preconditioning matrix
711: Notes:
712: Most users should not need to explicitly call this routine, as it
713: is used internally within the minimization solvers.
715: Level: developer
717: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
718: @*/
719: PetscErrorCode TaoComputeJacobianEquality(Tao tao, Vec X, Mat J, Mat Jpre)
720: {
727: if (!tao->ops->computejacobianequality) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianEquality() first");
728: ++tao->njac_equality;
729: VecLockReadPush(X);
730: PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
731: PetscStackPush("Tao user Jacobian(equality) function");
732: (*tao->ops->computejacobianequality)(tao,X,J,Jpre,tao->user_jac_equalityP);
733: PetscStackPop;
734: PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
735: VecLockReadPop(X);
736: return(0);
737: }
739: /*@C
740: TaoComputeJacobianInequality - Computes the Jacobian matrix that has been
741: set with TaoSetJacobianInequalityRoutine().
743: Collective on Tao
745: Input Parameters:
746: + tao - the Tao solver context
747: - X - input vector
749: Output Parameters:
750: + J - Jacobian matrix
751: - Jpre - Preconditioning matrix
753: Notes:
754: Most users should not need to explicitly call this routine, as it
755: is used internally within the minimization solvers.
757: Level: developer
759: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianStateRoutine(), TaoComputeJacobianDesign(), TaoSetStateDesignIS()
760: @*/
761: PetscErrorCode TaoComputeJacobianInequality(Tao tao, Vec X, Mat J, Mat Jpre)
762: {
769: if (!tao->ops->computejacobianinequality) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobianInequality() first");
770: ++tao->njac_inequality;
771: VecLockReadPush(X);
772: PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
773: PetscStackPush("Tao user Jacobian(inequality) function");
774: (*tao->ops->computejacobianinequality)(tao,X,J,Jpre,tao->user_jac_inequalityP);
775: PetscStackPop;
776: PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
777: VecLockReadPop(X);
778: return(0);
779: }
781: /*@C
782: TaoSetJacobianEqualityRoutine - Sets the function to compute the Jacobian
783: (and its inverse) of the constraint function with respect to the equality variables.
784: Used only for pde-constrained optimization.
786: Logically collective on Tao
788: Input Parameters:
789: + tao - the Tao context
790: . J - Matrix used for the jacobian
791: . Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.
792: . func - Jacobian evaluation routine
793: - ctx - [optional] user-defined context for private data for the
794: Jacobian evaluation routine (may be NULL)
796: Calling sequence of func:
797: $ func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);
799: + tao - the Tao context
800: . x - input vector
801: . J - Jacobian matrix
802: . Jpre - preconditioner matrix, usually the same as J
803: - ctx - [optional] user-defined Jacobian context
805: Level: intermediate
807: .seealso: TaoComputeJacobianEquality(), TaoSetJacobianDesignRoutine(), TaoSetEqualityDesignIS()
808: @*/
809: PetscErrorCode TaoSetJacobianEqualityRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat,void*), void *ctx)
810: {
815: if (J) {
818: }
819: if (Jpre) {
822: }
823: if (ctx) {
824: tao->user_jac_equalityP = ctx;
825: }
826: if (func) {
827: tao->ops->computejacobianequality = func;
828: }
829: if (J) {
830: PetscObjectReference((PetscObject)J);
831: MatDestroy(&tao->jacobian_equality);
832: tao->jacobian_equality = J;
833: }
834: if (Jpre) {
835: PetscObjectReference((PetscObject)Jpre);
836: MatDestroy(&tao->jacobian_equality_pre);
837: tao->jacobian_equality_pre=Jpre;
838: }
839: return(0);
840: }
842: /*@C
843: TaoSetJacobianInequalityRoutine - Sets the function to compute the Jacobian
844: (and its inverse) of the constraint function with respect to the inequality variables.
845: Used only for pde-constrained optimization.
847: Logically collective on Tao
849: Input Parameters:
850: + tao - the Tao context
851: . J - Matrix used for the jacobian
852: . Jpre - Matrix that will be used operated on by PETSc preconditioner, can be same as J.
853: . func - Jacobian evaluation routine
854: - ctx - [optional] user-defined context for private data for the
855: Jacobian evaluation routine (may be NULL)
857: Calling sequence of func:
858: $ func(Tao tao,Vec x,Mat J,Mat Jpre,void *ctx);
860: + tao - the Tao context
861: . x - input vector
862: . J - Jacobian matrix
863: . Jpre - preconditioner matrix, usually the same as J
864: - ctx - [optional] user-defined Jacobian context
866: Level: intermediate
868: .seealso: TaoComputeJacobianInequality(), TaoSetJacobianDesignRoutine(), TaoSetInequalityDesignIS()
869: @*/
870: PetscErrorCode TaoSetJacobianInequalityRoutine(Tao tao, Mat J, Mat Jpre, PetscErrorCode (*func)(Tao, Vec, Mat, Mat,void*), void *ctx)
871: {
876: if (J) {
879: }
880: if (Jpre) {
883: }
884: if (ctx) {
885: tao->user_jac_inequalityP = ctx;
886: }
887: if (func) {
888: tao->ops->computejacobianinequality = func;
889: }
890: if (J) {
891: PetscObjectReference((PetscObject)J);
892: MatDestroy(&tao->jacobian_inequality);
893: tao->jacobian_inequality = J;
894: }
895: if (Jpre) {
896: PetscObjectReference((PetscObject)Jpre);
897: MatDestroy(&tao->jacobian_inequality_pre);
898: tao->jacobian_inequality_pre=Jpre;
899: }
900: return(0);
901: }