Actual source code: taopetsc.c
1: #include "src/fortran/custom/zpetsc.h"
2: #include "petscksp.h"
3: #include tao.h
4: /*
5: #include "src/externaltao/petsctao/vector/taovec_petsc.h"
6: #include "src/externaltao/petsctao/matrix/taomat_petsc.h"
7: #include "src/externaltao/petsctao/indexset/taois_petsc.h"
8: #include "src/externaltao/petsctao/linearsolver/taolinearsolver_petsc.h"
9: #include "src/externaltao/petsctao/application/taoapp_petsc.h"
10: */
13: #ifdef PETSC_HAVE_FORTRAN_CAPS
15: #define taowrappetscvec_ TAOWRAPPETSCVEC
16: #define taowrappetscmat_ TAOWRAPPETSCMAT
17: #define taovecgetpetscvec_ TAOVECGETPETSCVEC
18: #define taomatgetpetscmat_ TAOMATGETPETSCMAT
20: #define taowrapis_ TAOWRAPIS
21: #define taowrapksp_ TAOWRAPKSP
23: #define taoindexsetgetis_ TAOINDEXSETGETIS
24: #define taolinearsolvergetpetscksp_ TAOLINEARSOLVERGETPETSCKSP
26: #define taopetscapplicationcreate_ TAOPETSCAPPLICATIONCREATE
27: #define taoapplicationcreate_ TAOAPPLICATIONCREATE
28: #define taoapplicationdestroy_ TAOAPPLICATIONDESTROY
29: #define taosetapplication_ TAOSETAPPLICATION
30: /* #define taosetfromoptions_ TAOSETFROMOPTIONS */
32: /*
33: #define taosetobjectivefunctionroutine_ TAOSETOBJECTIVEFUNCTIONROUTINE
34: #define taosetobjectivefunctiongradientroutine_ TAOSETOBJECTIVEFUNCTIONGRADIENTROUTINE
35: #define taosetgradientroutine_ TAOSETGRADIENTROUTINE
36: #define taosethessianroutine_ TAOSETHESSIANROUTINE
37: #define taosetvariableboundsroutine_ TAOSETVARIABLEBOUNDSROUTINE
38: */
40: #define taoappsetobjectiveandgradientroutine_ TAOAPPSETOBJECTIVEANDGRADIENTROUTINE
41: #define taoappsetobjectiveandgradientro_ TAOAPPSETOBJECTIVEANDGRADIENTRO
42: #define taoappsetobjectiveroutine_ TAOAPPSETOBJECTIVEROUTINE
43: #define taoappsetgradientroutine_ TAOAPPSETGRADIENTROUTINE
44: #define taoappsethessianroutine_ TAOAPPSETHESSIANROUTINE
45: #define taoappsetvariableboundsroutine_ TAOAPPSETVARIABLEBOUNDSROUTINE
46: #define taoappsetjacobianroutine_ TAOAPPSETJACOBIANROUTINE
47: #define taoappsetconstraintroutine_ TAOAPPSETCONSTRAINTROUTINE
49: /* Grid application */
50: #define daappsethessianroutine_ DAAPPSETHESSIANROUTINE
51: #define daappsetobjectiveandgradientroutine_ DAAPPSETOBJECTIVEANDGRADIENTROUTINE
52: #define daappsetobjectiveandgradientrou_ DAAPPSETOBJECTIVEANDGRADIENTROU
53: #define daappsetgradientroutine_ DAAPPSETGRADIENTROUTINE
54: #define daappsetobjectiveroutine_ DAAPPSETOBJECTIVEROUTINE
55: #define daappsetvariableboundsroutine_ DAAPPSETVARIABLEBOUNDSROUTINE
56: #define daappsetconstraintroutine_ DAAPPSETCONSTRAINTROUTINE
57: #define daappsetjacobianroutine_ DAAPPSETJACOBIANROUTINE
58: #define daappsetbeforemonitor_ DAAPPSETBEFOREMONITOR
59: #define daappsetaftermonitor_ DAAPPSETAFTERMONITOR
60: #define daappsetelementobjectiveandgradientroutine_ DAAPPSETELEMENTOBJECTIVEANDGRADIENTROUTINE
61: #define daappsetelementobjectiveandgrad_ DAAPPSETELEMENTOBJECTIVEANDGRAD
62: #define daappsetelementhessianroutine_ DAAPPSETELEMENTHESSIANROUTINE
65: /* deprecated */
66: #define taosetpetscfunction_ TAOSETPETSCFUNCTION
67: #define taosetpetscfunctiongradient_ TAOSETPETSCFUNCTIONGRADIENT
68: #define taosetpetscgradient_ TAOSETPETSCGRADIENT
69: #define taosetpetschessian_ TAOSETPETSCHESSIAN
72: #define taodefaultcomputehessian_ TAODEFAULTCOMPUTEHESSIAN
73: #define taodefaultcomputehessiancolor_ TAODEFAULTCOMPUTEHESSIANCOLOR
74: #define taosetpetscjacobian_ TAOSETPETSCJACOBIAN
75: #define taosetjacobianroutine_ TAOSETJACOBIANROUTINE
76: #define taosetconstraintroutine_ TAOSETCONSTRAINTROUTINE
77: #define taosetjacobianmat_ TAOSETJACOBIANMAT
78: #define taosetconstraintvec_ TAOSETCONSTRAINTVEC
80: #define taosetpetscvariablebounds_ TAOSETPETSCVARIABLEBOUNDS
81: #define taosetpetscinitialvector_ TAOSETPETSCINITIALVECTOR
82: #define taogetksp_ TAOGETKSP
84: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
86: #define taowrappetscvec_ taowrappetscvec
87: #define taowrappetscmat_ taowrappetscmat
88: #define taovecgetpetscvec_ taovecgetpetscvec
89: #define taomatgetpetscmat_ taomatgetpetscmat
91: #define taowrapis_ taowrapis
92: #define taowrapksp_ taowrapksp
94: #define taoindexsetgetpetscis_ taoindexsetgetpetscis
95: #define taolinearsolvergetksp_ taolinearsolvergetksp
97: #define taopetscapplicationcreate_ taopetscapplicationcreate
98: #define taoapplicationcreate_ taoapplicationcreate
99: #define taoapplicationdestroy_ taoapplicationdestroy
100: #define taosetapplication_ taosetapplication
101: /* #define taosetfromoptions_ taosetfromoptions */
103: /*
104: #define taosetobjectivefunctionroutine_ taosetobjectivefunctionroutine
105: #define taosetobjectivefunctiongradientroutine_ taosetobjectivefunctiongradientroutine
106: #define taosetgradientroutine_ taosetgradientroutine
107: #define taosethessianroutine_ taosethessianroutine
108: #define taosetvariableboundsroutine_ taosetvariableboundsroutine
109: */
111: #define taoappsetobjectiveandgradientroutine_ taoappsetobjectiveandgradientroutine
112: #define taoappsetobjectiveandgradientro_ taoappsetobjectiveandgradientro
113: #define taoappsetobjectiveroutine_ taoappsetobjectiveroutine
114: #define taoappsetgradientroutine_ taoappsetgradientroutine
115: #define taoappsethessianroutine_ taoappsethessianroutine
116: #define taoappsetvariableboundsroutine_ taoappsetvariableboundsroutine
117: #define taoappsetjacobianroutine_ taoappsetjacobianroutine
118: #define taoappsetconstraintroutine_ taoappsetconstraintroutine
121: /* Grid application */
122: #define daappsethessianroutine_ daappsethessianroutine
123: #define daappsetobjectiveandgradientroutine_ daappsetobjectiveandgradientroutine
124: #define daappsetobjectiveandgradientrou_ daappsetobjectiveandgradientrou
125: #define daappsetgradientroutine_ daappsetgradientroutine
126: #define daappsetobjectiveroutine_ daappsetobjectiveroutine
127: #define daappsetvariableboundsroutine_ daappsetvariableboundsroutine
128: #define daappsetconstraintroutine_ daappsetconstraintroutine
129: #define daappsetjacobianroutine_ daappsetjacobianroutine
130: #define daappsetbeforemonitor_ daappsetbeforemonitor
131: #define daappsetaftermonitor_ daappsetaftermonitor
132: #define daappsetelementobjectiveandgradientroutine_ daappsetelementobjectiveandgradientroutine
133: #define daappsetelementobjectiveandgrad_ daappsetelementobjectiveandgrad
134: #define daappsetelementhessianroutine_ daappsetelementhessianroutine
136: /* deprecated */
137: #define taosetpetscfunction_ taosetpetscfunction
138: #define taosetpetscfunctiongradient_ taosetpetscfunctiongradient
139: #define taosetpetscgradient_ taosetpetscgradient
140: #define taosetpetschessian_ taosetpetschessian
141: #define taodefaultcomputehessian_ taodefaultcomputehessian
142: #define taodefaultcomputehessiancolor_ taodefaultcomputehessiancolor
144: #define taosetpetscjacobian_ taosetpetscjacobian
145: #define taosetconstraintsfunction_ taosetconstraintsfunction
147: #define taosetjacobianroutine_ taosetjacobianroutine
148: #define taosetconstraintroutine_ taosetconstraintroutine
149: #define taosetjacobianmat_ taosetjacobianmat
150: #define taosetconstraintvec_ taosetconstraintvec
152: #define taosetpetscvariablebounds_ taosetpetscvariablebounds
153: #define taosetpetscinitialvector_ taosetpetscinitialvector
154: #define taogetksp_ taogetksp
156: #endif
159: void PETSC_STDCALL taowrappetscvec_(Vec *V, TaoVecPetsc **TV, int *ierr)
160: {
161: *TaoWrapPetscVec(*V,TV);
162: }
164: void PETSC_STDCALL taovecgetpetscvec_(TaoVec* *TV, Vec *V, int *ierr)
165: {
166: *ierr=TaoVecGetPetscVec(*TV,V);
167: }
170: void PETSC_STDCALL taowrappetscmat_(Mat *M, TaoMatPetsc **TM, int *ierr)
171: {
172: *TaoWrapPetscMat(*M,TM);
173: }
176: void PETSC_STDCALL taomatgetpetscmat_(TaoMat* *TM, Mat *M, int *ierr)
177: {
178: *ierr=TaoMatGetPetscMat(*TM,M);
179: }
182: void PETSC_STDCALL taowrappetscis_(IS *M, int *imax,TaoIndexSetPetsc **TM, int *ierr)
183: {
184: *TaoWrapPetscIS(*M,*imax,TM);
185: }
188: void PETSC_STDCALL taoindexsetgetpetscis_(TaoIndexSet* *TM, IS *M, int *ierr)
189: {
190: *ierr=TaoIndexSetGetPetscIS(*TM,M);
191: }
194: void PETSC_STDCALL taowrapksp_( KSP *M, TaoLinearSolverPetsc **TM, int *ierr)
195: {
196: *TaoWrapKSP(*M,TM);
197: }
199: void PETSC_STDCALL taolinearsolvergetksp_(TaoLinearSolver **TM, KSP *M, int *ierr)
200: {
201: *ierr=TaoLinearSolverGetKSP(*TM,M);
202: }
204: void PETSC_STDCALL taopetscapplicationcreate_(MPI_Comm *comm,TAO_APPLICATION *outtao,int *ierr){
206: *TaoPetscApplicationCreate((MPI_Comm)PetscToPointerComm(*comm),outtao);
208: }
210: void PETSC_STDCALL taoapplicationcreate_(MPI_Comm *comm,TAO_APPLICATION *outtao,int *ierr){
212: *TaoApplicationCreate((MPI_Comm)PetscToPointerComm(*comm),outtao);
214: }
216: void PETSC_STDCALL taoapplicationdestroy_(TAO_APPLICATION *outtao,int *ierr){
218: *TaoApplicationDestroy(*outtao);
220: }
222: /* Very Old routine must get rid of it . Here only for backward compatibility */
223: void PETSC_STDCALL taosetapplication_(TAO_SOLVER *tao, TAO_APPLICATION *outtao, int *ierr){
225: *TaoSetupApplicationSolver(*outtao, *tao);
227: }
229: /* Very Old routine must get rid of it . Here only for backward compatibility */
230: /*
231: void PETSC_STDCALL taosetfromoptions_(TAO_SOLVER *tao, int *ierr){
233: *TaoSetOptions(*tao);
235: }
236: */
237: /* -------------- Setting call-back routines ---------------- */
239: static void (PETSC_STDCALL *f1a)(TAO_SOLVER*,Vec*,double*,Vec*,void*,int*);
240: static void (PETSC_STDCALL *f1)(TAO_APPLICATION*,Vec*,double*,Vec*,void*,int*);
244: static int ourtaominfunctiongradientroutine(TAO_SOLVER tao,Vec x,double *f,Vec r,void *ctx)
245: {
246: int info = 0;
247: (*f1a)(&tao,&x,f,&r,ctx,&info);CHKERRQ(info);
248: return 0;
249: }
251: static int ourtaominfunctiongradientroutine2(TAO_APPLICATION taoapp,Vec x,double *f,Vec r,void *ctx)
252: {
253: int info = 0;
254: (*f1)(&taoapp,&x,f,&r,ctx,&info);CHKERRQ(info);
255: return 0;
256: }
261: void PETSC_STDCALL taosetpetscfunctiongradient_(TAO_APPLICATION *taoapp,Vec *x,Vec *r,void (PETSC_STDCALL *func)(TAO_SOLVER*,Vec*,double*,Vec*,void*,int*),void *ctx,int *info){
262: f1a = func;
263: *info = TaoSetPetscFunctionGradient(*taoapp,*x,*r,ourtaominfunctiongradientroutine,ctx);
264: }
266: void PETSC_STDCALL taoappsetobjectiveandgradientroutine_(TAO_APPLICATION *taoapp,void (PETSC_STDCALL *func)(TAO_APPLICATION*,Vec*,double*,Vec*,void*,int*),void *ctx,int *info){
267: f1 = func;
268: *info = TaoAppSetObjectiveAndGradientRoutine(*taoapp,ourtaominfunctiongradientroutine2,ctx);
269: }
270: void PETSC_STDCALL taoappsetobjectiveandgradientro_(TAO_APPLICATION *taoapp,void (PETSC_STDCALL *func)(TAO_APPLICATION*,Vec*,double*,Vec*,void*,int*),void *ctx,int *info){
271: f1 = func;
272: *info = TaoAppSetObjectiveAndGradientRoutine(*taoapp,ourtaominfunctiongradientroutine2,ctx);
273: }
276: static void (PETSC_STDCALL *f2a)(TAO_SOLVER*,Vec*,double*,void*,int*);
277: static void (PETSC_STDCALL *f2)(TAO_APPLICATION*,Vec*,double*,void*,int*);
279: static int ourtaominfunction(TAO_SOLVER tao,Vec x,double* d,void *ctx)
280: {
281: int info = 0;
282: (*f2a)(&tao,&x,d,ctx,&info);CHKERRQ(info);
283: return 0;
284: }
285: static int ourtaominfunction2(TAO_APPLICATION taoapp,Vec x,double* d,void *ctx)
286: {
287: int info = 0;
288: (*f2)(&taoapp,&x,d,ctx,&info);CHKERRQ(info);
289: return 0;
290: }
293: void PETSC_STDCALL taosetpetscfunction_(TAO_APPLICATION *taoapp,Vec *x,
294: void (PETSC_STDCALL *func)(TAO_SOLVER*,Vec*,double*,void*,int*),void *ctx,int *info){
295: f2a = func;
296: *info = TaoSetPetscFunction(*taoapp,*x,ourtaominfunction,ctx);
297: }
299: void PETSC_STDCALL taoappsetobjectiveroutine(TAO_APPLICATION *taoapp,
300: void (PETSC_STDCALL *func)(TAO_APPLICATION*,Vec*,double*,void*,int*),void *ctx,int *info){
301: f2 = func;
302: *info = TaoAppSetObjectiveRoutine(*taoapp,ourtaominfunction2,ctx);
303: }
307: static void (PETSC_STDCALL *f3a)(TAO_SOLVER*,Vec*,Vec *,void*,int*);
308: static void (PETSC_STDCALL *f3)(TAO_APPLICATION*,Vec*,Vec *,void*,int*);
311: static int ourtaogradientfunction(TAO_SOLVER tao,Vec x,Vec d,void *ctx)
312: {
313: int info = 0;
314: (*f3a)(&tao,&x,&d,ctx,&info);CHKERRQ(info);
315: return 0;
316: }
317: static int ourtaogradientfunction2(TAO_APPLICATION taoapp,Vec x,Vec d,void *ctx)
318: {
319: int info = 0;
320: (*f3)(&taoapp,&x,&d,ctx,&info);CHKERRQ(info);
321: return 0;
322: }
325: void PETSC_STDCALL taosetpetscgradient_(TAO_APPLICATION *taoapp,Vec *r,void (PETSC_STDCALL *func)(TAO_SOLVER*,Vec*,Vec*,void*,int*),void *ctx,int *info){
326: f3a = func;
327: *info = TaoSetPetscGradient(*taoapp,*r,ourtaogradientfunction,ctx);
328: }
330: void PETSC_STDCALL taoappsetgradientroutine_(TAO_APPLICATION *taoapp,void (PETSC_STDCALL *func)(TAO_APPLICATION*,Vec*,Vec*,void*,int*),void *ctx,int *info){
331: f3 = func;
332: *info = TaoAppSetGradientRoutine(*taoapp,ourtaogradientfunction2,ctx);
333: }
336: /* ---------------------------------------------------------*/
337: /*
338: taodefaultcomputehessian() and taodefaultcomputehessiancolor().
339: These routines can be used directly from Fortran, but the following is done
340: primarily so that a Fortran call to TaoSetPetscHessian() will properly handle the
341: defaults being passed in.
343: functions, hence no STDCALL
344: */
345: void taodefaultcomputehessian_(TAO_APPLICATION *tao,Vec *x,Mat *m,Mat *p,MatStructure* type,
346: void *ctx,int *info)
347: {
348: *info = TaoAppDefaultComputeHessian(*tao,*x,m,p,type,ctx);
349: }
351: void taodefaultcomputehessiancolor_(TAO_APPLICATION *tao,Vec *x,Mat *m,Mat *p,
352: MatStructure* type,void *ctx,int *info)
353: {
354: *info = TaoAppDefaultComputeHessianColor(*tao,*x,m,p,type,*(MatFDColoring*)ctx);
355: }
357: static void (PETSC_STDCALL *f4a)(TAO_SOLVER*,Vec*,Mat*,Mat*,MatStructure*,void*,int*);
358: static void (PETSC_STDCALL *f4)(TAO_APPLICATION*,Vec*,Mat*,Mat*,MatStructure*,void*,int*);
361: static int ourtaohessian(TAO_SOLVER tao,Vec x,Mat* m,Mat* p,MatStructure* type,void*ctx)
362: {
363: int info = 0;
364: (*f4a)(&tao,&x,m,p,type,ctx,&info);CHKERRQ(info);
365: return 0;
366: }
367: static int ourtaohessian2(TAO_APPLICATION taoapp,Vec x,Mat* m,Mat* p,MatStructure* type,void*ctx)
368: {
369: int info = 0;
370: (*f4)(&taoapp,&x,m,p,type,ctx,&info);CHKERRQ(info);
371: return 0;
372: }
376: void PETSC_STDCALL taosetpetschessian_(TAO_APPLICATION *taoapp,Mat *A,Mat *B,void (PETSC_STDCALL *func)(TAO_SOLVER*,Vec*,Mat*,Mat*,
377: MatStructure*,void*,int*),void *ctx,int *info)
378: {
379: f4a = func;
380: *info = TaoSetPetscHessian(*taoapp,*A,*B,ourtaohessian,ctx);
381: }
383: void PETSC_STDCALL taoappsethessianroutine_(TAO_APPLICATION *taoapp,void (PETSC_STDCALL *func)(TAO_APPLICATION*,Vec*,Mat*,Mat*,
384: MatStructure*,void*,int*),void *ctx,int *info)
385: {
386: if ((void(*)())func == (void(*)()) taodefaultcomputehessian_) {
387: *info = TaoAppSetHessianRoutine(*taoapp,TaoAppDefaultComputeHessian,ctx);
388: } else if ( (void(*)()) func == (void(*)()) taodefaultcomputehessiancolor_) {
389: *info = TaoAppSetHessianRoutine(*taoapp,TaoAppDefaultComputeHessianColor,*(MatFDColoring*)ctx);
390: } else {
391: f4 = func;
392: *info = TaoAppSetHessianRoutine(*taoapp,ourtaohessian2,ctx);
393: }
394: }
397: static void (PETSC_STDCALL *fb)(TAO_APPLICATION*, Vec*,Vec *,void*,int*);
401: static int ourtaovariableboundsroutine(TAO_APPLICATION tao, Vec x,Vec d,void *ctx)
402: {
403: int info = 0;
404: (*fb)(&tao,&x,&d,ctx,&info);CHKERRQ(info);
405: return 0;
406: }
411: void PETSC_STDCALL taoappsetvariableboundsroutine_(TAO_APPLICATION *taoapp,void (PETSC_STDCALL *func)(TAO_APPLICATION*,Vec*,Vec*,void*,int*),void *ctx,int *info){
412: fb = func;
413: *info = TaoAppSetVariableBoundsRoutine(*taoapp,ourtaovariableboundsroutine,ctx);
414: }
417: static void (PETSC_STDCALL *f6a)(TAO_SOLVER*,Vec*,Vec *,void*,int*);
418: static void (PETSC_STDCALL *f6)(TAO_APPLICATION*,Vec*,Vec *,void*,int*);
420: static int ourtaoconstraintsfunction(TAO_SOLVER tao,Vec x,Vec d,void *ctx)
421: {
422: int info = 0;
423: (*f6a)(&tao,&x,&d,ctx,&info);CHKERRQ(info);
424: return 0;
425: }
426: static int ourtaoconstraintsfunction2(TAO_APPLICATION taoapp,Vec x,Vec d,void *ctx)
427: {
428: int info = 0;
429: (*f6)(&taoapp,&x,&d,ctx,&info);CHKERRQ(info);
430: return 0;
431: }
434: void PETSC_STDCALL taosetpetscconstraintsfunction_(TAO_APPLICATION *taoapp,Vec *r,void (PETSC_STDCALL *func)(TAO_SOLVER*,Vec*,Vec*,void*,int*),void *ctx,int *info){
435: f6a = func;
436: *info = TaoSetPetscConstraintsFunction(*taoapp,*r,ourtaoconstraintsfunction,ctx);
437: }
439: void PETSC_STDCALL taoappsetconstraintroutine_(TAO_APPLICATION *taoapp,void (PETSC_STDCALL *func)(TAO_APPLICATION*,Vec*,Vec*,void*,int*),void *ctx,int *info){
440: f6 = func;
441: *info = TaoAppSetConstraintRoutine(*taoapp,ourtaoconstraintsfunction2,ctx);
442: }
444: void PETSC_STDCALL taosetconstraintvec_(TAO_APPLICATION *taoapp, Vec *r, int *info){
445: *info = TaoAppSetFunctionVec(*taoapp,*r);
446: }
448: static void (PETSC_STDCALL *f7a)(TAO_SOLVER*,Vec*,Mat *,void*,int*);
449: static void (PETSC_STDCALL *f7)(TAO_APPLICATION*,Vec*,Mat *,Mat*, MatStructure*,void*,int*);
451: static int ourtaojacobianfunction(TAO_SOLVER tao,Vec x,Mat *J,void *ctx)
452: {
453: int info = 0;
454: (*f7a)(&tao,&x,J,ctx,&info);CHKERRQ(info);
455: return 0;
456: }
457: static int ourtaojacobianfunction2(TAO_APPLICATION tao,Vec x,Mat *J,Mat *JP, MatStructure *flag,void *ctx)
458: {
459: int info = 0;
460: (*f7)(&tao,&x,J,JP,flag,ctx,&info);CHKERRQ(info);
461: return 0;
462: }
465: void PETSC_STDCALL taosetpetscjacobian_(TAO_APPLICATION *taoapp, Mat*J, void (PETSC_STDCALL *func)(TAO_SOLVER*,Vec*,Mat*,void*,int*),void *ctx,int *info){
466: f7a = func;
467: *info = TaoSetPetscJacobian(*taoapp,*J,ourtaojacobianfunction,ctx);
468: }
470: void PETSC_STDCALL taoappsetjacobianroutine_(TAO_APPLICATION *taoapp,void (PETSC_STDCALL *func)(TAO_APPLICATION*,Vec*,Mat*,Mat*,MatStructure*,void*,int*),void *ctx,int *info){
471: f7 = func;
472: *info = TaoAppSetJacobianRoutine(*taoapp,ourtaojacobianfunction2,ctx);
473: }
475: void PETSC_STDCALL taosetjacobianmat_(TAO_APPLICATION *taoapp,Mat *J, Mat *JP,int *info){
476: *info = TaoAppSetJacobianMat(*taoapp,*J,*JP);
477: }
480: void PETSC_STDCALL taogetksp_(TAO_SOLVER *tao,KSP *tksp, int *info ){
481: *info = TaoGetKSP(*tao,tksp);
482: }
486: #include taodaapplication.h
488: /* Grid Application routines */
490: static void (PETSC_STDCALL *grid1)(TAO_APPLICATION*, DA*, Vec*, Mat*, void*, int*);
491: static void (PETSC_STDCALL *grid2)(TAO_APPLICATION*, DA*, Vec*, double*, Vec*, void*, int*);
492: static void (PETSC_STDCALL *grid3)(TAO_APPLICATION*, DA*, Vec*, Vec*, void*, int*);
493: static void (PETSC_STDCALL *grid4)(TAO_APPLICATION*, DA*, Vec*, double*, void*, int*);
494: static void (PETSC_STDCALL *grid5)(TAO_APPLICATION*, DA*, Vec*, Vec*, void*, int*);
495: static void (PETSC_STDCALL *grid6)(TAO_APPLICATION*, DA*, Vec*, Vec*, void*, int*);
496: static void (PETSC_STDCALL *grid7)(TAO_APPLICATION*, DA*, Vec*, Mat*, void*, int*);
497: static void (PETSC_STDCALL *grid8)(TAO_APPLICATION*, DA*, int*, void*, int*);
498: static void (PETSC_STDCALL *grid9)(TAO_APPLICATION*, DA*, int*, void*, int*);
499: /* Grid Application Element routines */
500: static void (PETSC_STDCALL *grid10)(int[2], double[4], double*, double[4], void*, int*);
501: static void (PETSC_STDCALL *grid11)(int[2], double[4], double[4][4], void*, int*);
503: static int ourdaappsethessianroutine(TAO_APPLICATION daapp, DA da, Vec x, Mat H, void *ctx) {
504: int info = 0;
505: (*grid1)(&daapp, &da, &x, &H, ctx, &info);
506: return info;
507: }
508: static int ourdaappsetobjectiveandgradientroutine(TAO_APPLICATION daapp, DA da, Vec x, double *f, Vec G, void *ctx) {
509: int info = 0;
510: (*grid2)(&daapp, &da, &x, f, &G, ctx, &info);
511: return info;
512: }
513: static int ourdaappsetgradientroutine(TAO_APPLICATION daapp, DA da, Vec x, Vec g, void *ctx) {
514: int info = 0;
515: (*grid3)(&daapp, &da, &x, &g, ctx, &info);
516: return info;
517: }
518: static int ourdaappsetobjectiveroutine(TAO_APPLICATION daapp, DA da, Vec x, double *f, void *ctx) {
519: int info = 0;
520: (*grid4)(&daapp, &da, &x, f, ctx, &info);
521: return info;
522: }
523: static int ourdaappsetvariableboundsroutine(TAO_APPLICATION daapp, DA da, Vec L, Vec U, void *ctx) {
524: int info = 0;
525: (*grid5)(&daapp, &da, &L, &U, ctx, &info);
526: return info;
527: }
528: static int ourdaappsetconstraintroutine(TAO_APPLICATION daapp, DA da, Vec X, Vec R, void *ctx) {
529: int info = 0;
530: (*grid6)(&daapp, &da, &X, &R, ctx, &info);
531: return info;
532: }
533: static int ourdaappsetjacobianroutine(TAO_APPLICATION daapp, DA da, Vec X, Mat J, void *ctx) {
534: int info = 0;
535: (*grid7)(&daapp, &da, &X, &J, ctx, &info);
536: return info;
537: }
538: static int ourdaappsetbeforemonitor(TAO_APPLICATION daapp, DA da, int levels, void *ctx) {
539: int info = 0;
540: (*grid8)(&daapp, &da, &levels, ctx, &info);
541: return info;
542: }
543: static int ourdaappsetaftermonitor(TAO_APPLICATION daapp, DA da, int levels, void *ctx) {
544: int info = 0;
545: (*grid9)(&daapp, &da, &levels, ctx, &info);
546: return info;
547: }
548: static int ourdaappsetelementobjectiveandgradientroutine(int coord[2], double x[4], double *f, double g[4], void *ctx) {
549: int info = 0;
550: (*grid10)(coord, x, f, g, ctx, &info);
551: return info;
552: }
553: static int ourdaappsetelementhessianroutine(int coord[2], double x[4], double H[4][4], void *ctx) {
554: int info = 0;
555: (*grid11)(coord, x, H, ctx, &info);
556: return info;
557: }
559: void PETSC_STDCALL daappsethessianroutine_(TAO_APPLICATION *daapp, void (PETSC_STDCALL *func)(TAO_APPLICATION*,DA*,Vec*, Mat*, void*, int*),void *ctx, int *info) {
560: grid1 = func;
561: *info = DAAppSetHessianRoutine(*daapp, ourdaappsethessianroutine, ctx);
562: }
563: void PETSC_STDCALL daappsetobjectiveandgradientroutine_(TAO_APPLICATION *daapp, void (PETSC_STDCALL *func)(TAO_APPLICATION*, DA*, Vec*, double*, Vec*, void*, int*), void *ctx, int *info) {
564: grid2 = func;
565: *info = DAAppSetObjectiveAndGradientRoutine(*daapp, ourdaappsetobjectiveandgradientroutine, ctx);
566: }
567: void PETSC_STDCALL daappsetobjectiveandgradientrou_(TAO_APPLICATION *daapp, void (PETSC_STDCALL *func)(TAO_APPLICATION*, DA*, Vec*, double*, Vec*, void*, int*), void *ctx, int *info) {
568: grid2 = func;
569: *info = DAAppSetObjectiveAndGradientRoutine(*daapp, ourdaappsetobjectiveandgradientroutine, ctx);
570: }
572: void PETSC_STDCALL daappsetgradientroutine_(TAO_APPLICATION *daapp, void (PETSC_STDCALL *func)(TAO_APPLICATION*, DA*, Vec*, Vec*, void*, int*), void *ctx, int *info) {
573: grid3 = func;
574: *info = DAAppSetGradientRoutine(*daapp, ourdaappsetgradientroutine, ctx);
575: }
576: void PETSC_STDCALL daappsetobjectiveroutine_(TAO_APPLICATION *daapp, void (PETSC_STDCALL *func)(TAO_APPLICATION*, DA*, Vec*, double*, void*, int*), void *ctx, int *info) {
577: grid4 = func;
578: *info = DAAppSetObjectiveRoutine(*daapp, ourdaappsetobjectiveroutine, ctx);
579: }
580: void PETSC_STDCALL daappsetvariableboundsroutine_(TAO_APPLICATION *daapp, void (PETSC_STDCALL *func)(TAO_APPLICATION*, DA*, Vec*, Vec*, void*, int*), void *ctx, int *info) {
581: grid5 = func;
582: *info = DAAppSetVariableBoundsRoutine(*daapp, ourdaappsetvariableboundsroutine, ctx);
583: }
584: void PETSC_STDCALL daappsetconstraintroutine_(TAO_APPLICATION *daapp, void (PETSC_STDCALL *func)(TAO_APPLICATION*, DA*, Vec*, Vec*, void*, int*), void *ctx, int *info) {
585: grid6 = func;
586: *info = DAAppSetConstraintRoutine(*daapp, ourdaappsetconstraintroutine, ctx);
587: }
588: void PETSC_STDCALL daappsetjacobianroutine_(TAO_APPLICATION *daapp, void (PETSC_STDCALL *func)(TAO_APPLICATION*, DA*, Vec*, Mat*, void*, int*), void *ctx, int *info) {
589: grid7 = func;
590: *info = DAAppSetJacobianRoutine(*daapp, ourdaappsetjacobianroutine, ctx);
591: }
592: void PETSC_STDCALL daappsetbeforemonitor_(TAO_APPLICATION *daapp, void (PETSC_STDCALL *func)(TAO_APPLICATION*, DA*, int*, void*, int*), void *ctx, int *info) {
593: grid8 = func;
594: *info = DAAppSetBeforeMonitor(*daapp, ourdaappsetbeforemonitor, ctx);
595: }
596: void PETSC_STDCALL daappsetaftermonitor_(TAO_APPLICATION *daapp, void (PETSC_STDCALL *func)(TAO_APPLICATION*, DA*, int*, void*, int*), void *ctx, int *info) {
597: grid9 = func;
598: *info = DAAppSetAfterMonitor(*daapp, ourdaappsetaftermonitor, ctx);
599: }
600: void PETSC_STDCALL daappsetelementobjectiveandgradientroutine_(TAO_APPLICATION *daapp, void (PETSC_STDCALL *func)(int[2], double[4], double*, double[4], void*, int*), int *flops, void *ctx, int *info) {
601: grid10 = func;
602: *info = DAAppSetElementObjectiveAndGradientRoutine(*daapp, ourdaappsetelementobjectiveandgradientroutine, *flops, ctx);
603: }
604: void PETSC_STDCALL daappsetelementobjectiveandgrad_(TAO_APPLICATION *daapp, void (PETSC_STDCALL *func)(int[2], double[4], double*, double[4], void*, int*), int *flops, void *ctx, int *info) {
605: grid10 = func;
606: *info = DAAppSetElementObjectiveAndGradientRoutine(*daapp, ourdaappsetelementobjectiveandgradientroutine, *flops, ctx);
607: }
608: void PETSC_STDCALL daappsetelementhessianroutine_(TAO_APPLICATION *daapp, void (PETSC_STDCALL *func)(int[2], double[4], double[4][4], void*, int*), int *flops, void *ctx, int *info) {
609: grid11 = func;
610: *info = DAAppSetElementHessianRoutine(*daapp, ourdaappsetelementhessianroutine, *flops, ctx);
611: }