Actual source code: taosolver_hj.c
petsc-3.14.6 2021-03-30
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: PetscObjectBaseTypeCompareAny((PetscObject)hessian,&flg,MATSEQAIJ,MATMPIAIJ,MATSEQDENSE,MATMPIDENSE,MATSEQBAIJ,MATMPIBAIJ,MATSEQSBAIJ,MATMPIBAIJ,"");
108: if (flg) {
109: A = hessian;
110: PetscObjectReference((PetscObject)A);
111: } else {
112: MatComputeOperator(hessian,MATAIJ,&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: PetscViewerDestroy(&mviewer);
187: }
188: PetscViewerASCIISetTab(viewer,tabs);
189: return(0);
190: }
192: /*@C
193: TaoComputeHessian - Computes the Hessian matrix that has been
194: set with TaoSetHessianRoutine().
196: Collective on Tao
198: Input Parameters:
199: + tao - the Tao solver context
200: - X - input vector
202: Output Parameters:
203: + H - Hessian matrix
204: - Hpre - Preconditioning matrix
206: Options Database Keys:
207: + -tao_test_hessian - compare the user provided Hessian with one compute via finite differences to check for errors
208: . -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
209: - -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
211: Notes:
212: Most users should not need to explicitly call this routine, as it
213: is used internally within the minimization solvers.
215: TaoComputeHessian() is typically used within minimization
216: implementations, so most users would not generally call this routine
217: themselves.
219: Developer Notes:
220: The Hessian test mechanism follows SNESTestJacobian().
222: Level: developer
224: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetHessianRoutine()
225: @*/
226: PetscErrorCode TaoComputeHessian(Tao tao, Vec X, Mat H, Mat Hpre)
227: {
234: if (!tao->ops->computehessian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetHessianRoutine() first");
235: ++tao->nhess;
236: VecLockReadPush(X);
237: PetscLogEventBegin(TAO_HessianEval,tao,X,H,Hpre);
238: PetscStackPush("Tao user Hessian function");
239: (*tao->ops->computehessian)(tao,X,H,Hpre,tao->user_hessP);
240: PetscStackPop;
241: PetscLogEventEnd(TAO_HessianEval,tao,X,H,Hpre);
242: VecLockReadPop(X);
244: TaoTestHessian(tao);
245: return(0);
246: }
248: /*@C
249: TaoComputeJacobian - Computes the Jacobian matrix that has been
250: set with TaoSetJacobianRoutine().
252: Collective on Tao
254: Input Parameters:
255: + tao - the Tao solver context
256: - X - input vector
258: Output Parameters:
259: + J - Jacobian matrix
260: - Jpre - Preconditioning matrix
262: Notes:
263: Most users should not need to explicitly call this routine, as it
264: is used internally within the minimization solvers.
266: TaoComputeJacobian() is typically used within minimization
267: implementations, so most users would not generally call this routine
268: themselves.
270: Level: developer
272: .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetJacobianRoutine()
273: @*/
274: PetscErrorCode TaoComputeJacobian(Tao tao, Vec X, Mat J, Mat Jpre)
275: {
282: if (!tao->ops->computejacobian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetJacobian() first");
283: ++tao->njac;
284: VecLockReadPush(X);
285: PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
286: PetscStackPush("Tao user Jacobian function");
287: (*tao->ops->computejacobian)(tao,X,J,Jpre,tao->user_jacP);
288: PetscStackPop;
289: PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
290: VecLockReadPop(X);
291: return(0);
292: }
294: /*@C
295: TaoComputeResidualJacobian - Computes the least-squares residual Jacobian matrix that has been
296: set with TaoSetJacobianResidual().
298: Collective on Tao
300: Input Parameters:
301: + tao - the Tao solver context
302: - X - input vector
304: Output Parameters:
305: + J - Jacobian matrix
306: - Jpre - Preconditioning matrix
308: Notes:
309: Most users should not need to explicitly call this routine, as it
310: is used internally within the minimization solvers.
312: TaoComputeResidualJacobian() is typically used within least-squares
313: implementations, so most users would not generally call this routine
314: themselves.
316: Level: developer
318: .seealso: TaoComputeResidual(), TaoSetJacobianResidual()
319: @*/
320: PetscErrorCode TaoComputeResidualJacobian(Tao tao, Vec X, Mat J, Mat Jpre)
321: {
328: if (!tao->ops->computeresidualjacobian) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call TaoSetResidualJacobian() first");
329: ++tao->njac;
330: VecLockReadPush(X);
331: PetscLogEventBegin(TAO_JacobianEval,tao,X,J,Jpre);
332: PetscStackPush("Tao user least-squares residual Jacobian function");
333: (*tao->ops->computeresidualjacobian)(tao,X,J,Jpre,tao->user_lsjacP);
334: PetscStackPop;
335: PetscLogEventEnd(TAO_JacobianEval,tao,X,J,Jpre);
336: VecLockReadPop(X);
337: return(0);
338: }
340: /*@C
341: TaoComputeJacobianState - Computes the Jacobian matrix that has been
342: set with TaoSetJacobianStateRoutine().
344: Collective on Tao
346: Input Parameters:
347: + tao - the Tao solver context
348: - X - input vector
350: Output Parameters:
351: + Jpre - Jacobian matrix
352: - Jinv - Preconditioning matrix
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: }