Actual source code: zts.c

  1: /*$Id: zts.c,v 1.37 2001/03/28 19:43:08 balay Exp $*/

  3: #include "src/fortran/custom/zpetsc.h"
  4: #include "petscts.h"

  6: #ifdef PETSC_HAVE_FORTRAN_CAPS
  7: #define tssetrhsfunction_                    TSSETRHSFUNCTION
  8: #define tssetrhsmatrix_                      TSSETRHSMATRIX
  9: #define tssetrhsjacobian_                    TSSETRHSJACOBIAN
 10: #define tscreate_                            TSCREATE
 11: #define tsgetsolution_                       TSGETSOLUTION
 12: #define tsgetsnes_                           TSGETSNES
 13: #define tsgetsles_                           TSGETSLES
 14: #define tsgettype_                           TSGETTYPE
 15: #define tsdestroy_                           TSDESTROY
 16: #define tssetmonitor_                        TSSETMONITOR
 17: #define tssettype_                           TSSETTYPE
 18: #define tspvodegetiterations_                TSPVODEGETITERATIONS
 19: #define tsdefaultcomputejacobian_            TSDEFAULTCOMPUTEJACOBIAN
 20: #define tsdefaultcomputejacobiancolor_       TSDEFAULTCOMPUTEJACOBIANCOLOR
 21: #define tsgetoptionsprefix_                  TSGETOPTIONSPREFIX
 22: #define tsdefaultmonitor_                    TSDEFAULTMONITOR
 23: #define tsview_                              TSVIEW
 24: #define tsgetrhsjacobian_                    TSGETRHSJACOBIAN
 25: #define tsgetrhsmatrix_                      TSGETRHSMATRIX
 26: #define tssetrhsboundaryconditions_          TSSETRHSBOUNDARYCONDITIONS
 27: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 28: #define tsdefaultcomputejacobian_            tsdefaultcomputejacobian
 29: #define tsdefaultcomputejacobiancolor_       tsdefaultcomputejacobiancolor
 30: #define tspvodegetiterations_                tspvodegetiterations
 31: #define tssetrhsfunction_                    tssetrhsfunction
 32: #define tssetrhsmatrix_                      tssetrhsmatrix
 33: #define tssetrhsjacobian_                    tssetrhsjacobian
 34: #define tscreate_                            tscreate
 35: #define tsgetsolution_                       tsgetsolution
 36: #define tsgetsnes_                           tsgetsnes
 37: #define tsgetsles_                           tsgetsles
 38: #define tsgettype_                           tsgettype
 39: #define tsdestroy_                           tsdestroy
 40: #define tssetmonitor_                        tssetmonitor
 41: #define tssettype_                           tssettype
 42: #define tsgetoptionsprefix_                  tsgetoptionsprefix
 43: #define tsdefaultmonitor_                    tsdefaultmonitor
 44: #define tsview_                              tsview
 45: #define tsgetrhsjacobian_                    tsgetrhsjacobian
 46: #define tsgetrhsmatrix_                      tsgetrhsmatrix
 47: #define tssetrhsboundaryconditions_          tssetrhsboundaryconditions
 48: #endif

 50: EXTERN_C_BEGIN

 52: static int ourtsbcfunction(TS ts,double d,Vec x,void *ctx)
 53: {
 54:   int 0;
 55:   (*(void (PETSC_STDCALL *)(TS*,double*,Vec*,void*,int*))(((PetscObject)ts)->fortran_func_pointers[0]))(&ts,&d,&x,ctx,&ierr);
 56:   return 0;
 57: }

 59: void PETSC_STDCALL tssetrhsboundaryconditions_(TS *ts,int (PETSC_STDCALL *f)(TS*,double*,Vec*,void*,int*),void *ctx,int *ierr)
 60: {
 61:   ((PetscObject)*ts)->fortran_func_pointers[0] = (void(*)())f;
 62:   *TSSetRHSBoundaryConditions(*ts,ourtsbcfunction,ctx);
 63: }

 65: void PETSC_STDCALL tsgetrhsjacobian_(TS *ts,Mat *J,Mat *M,void **ctx,int *ierr)
 66: {
 67:   *TSGetRHSJacobian(*ts,J,M,ctx);
 68: }

 70: void PETSC_STDCALL tsgetrhsmatrix_(TS *ts,Mat *J,Mat *M,void **ctx,int *ierr)
 71: {
 72:   *TSGetRHSMatrix(*ts,J,M,ctx);
 73: }

 75: void PETSC_STDCALL tsview_(TS *ts,PetscViewer *viewer, int *ierr)
 76: {
 77:   PetscViewer v;
 78:   PetscPatchDefaultViewers_Fortran(viewer,v);
 79:   *TSView(*ts,v);
 80: }

 82: /* function */
 83: void tsdefaultcomputejacobian_(TS *ts,double *t,Vec *xx1,Mat *J,Mat *B,MatStructure *flag,void *ctx,int *ierr)
 84: {
 85:   *TSDefaultComputeJacobian(*ts,*t,*xx1,J,B,flag,ctx);
 86: }

 88: /* function */
 89: void tsdefaultcomputejacobiancolor_(TS *ts,double *t,Vec *xx1,Mat *J,Mat *B,MatStructure *flag,void *ctx,int *ierr)
 90: {
 91:   *TSDefaultComputeJacobianColor(*ts,*t,*xx1,J,B,flag,*(MatFDColoring*)ctx);
 92: }

 94: void PETSC_STDCALL tssettype_(TS *ts,CHAR type PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
 95: {
 96:   char *t;

 98:   FIXCHAR(type,len,t);
 99:   *TSSetType(*ts,t);
100:   FREECHAR(type,t);
101: }

103: static int ourtsfunction(TS ts,double d,Vec x,Vec f,void *ctx)
104: {
105:   int 0;
106:   (*(void (PETSC_STDCALL *)(TS*,double*,Vec*,Vec*,void*,int*))(((PetscObject)ts)->fortran_func_pointers[1]))(&ts,&d,&x,&f,ctx,&ierr);
107:   return 0;
108: }

110: void PETSC_STDCALL tssetrhsfunction_(TS *ts,int (PETSC_STDCALL *f)(TS*,double*,Vec*,Vec*,void*,int*),void*fP,int *ierr)
111: {
112:   ((PetscObject)*ts)->fortran_func_pointers[1] = (void(*)())f;
113:   *TSSetRHSFunction(*ts,ourtsfunction,fP);
114: }


117: /* ---------------------------------------------------------*/
118: static int ourtsmatrix(TS ts,double d,Mat* m,Mat* p,MatStructure* type,void*ctx)
119: {
120:   int 0;
121:   (*(void (PETSC_STDCALL *)(TS*,double*,Mat*,Mat*,MatStructure*,void*,int*))(((PetscObject)ts)->fortran_func_pointers[2]))(&ts,&d,m,p,type,ctx,&ierr);
122:   return 0;
123: }

125: void PETSC_STDCALL tssetrhsmatrix_(TS *ts,Mat *A,Mat *B,int (PETSC_STDCALL *f)(TS*,double*,Mat*,Mat*,MatStructure*,
126:                                                    void*,int *),void*fP,int *ierr)
127: {
128:   if (FORTRANNULLFUNCTION(f)) {
129:     *TSSetRHSMatrix(*ts,*A,*B,PETSC_NULL,fP);
130:   } else {
131:     ((PetscObject)*ts)->fortran_func_pointers[2] = (void(*)())f;
132:     *TSSetRHSMatrix(*ts,*A,*B,ourtsmatrix,fP);
133:   }
134: }

136: /* ---------------------------------------------------------*/
137: static int ourtsjacobian(TS ts,double d,Vec x,Mat* m,Mat* p,MatStructure* type,void*ctx)
138: {
139:   int 0;
140:   (*(void (PETSC_STDCALL *)(TS*,double*,Vec*,Mat*,Mat*,MatStructure*,void*,int*))(((PetscObject)ts)->fortran_func_pointers[3]))(&ts,&d,&x,m,p,type,ctx,&ierr);
141:   return 0;
142: }

144: void PETSC_STDCALL tssetrhsjacobian_(TS *ts,Mat *A,Mat *B,void (PETSC_STDCALL *f)(TS*,double*,Vec*,Mat*,Mat*,MatStructure*,
145:                void*,int*),void*fP,int *ierr)
146: {
147:   if (FORTRANNULLFUNCTION(f)) {
148:     *TSSetRHSJacobian(*ts,*A,*B,PETSC_NULL,fP);
149:   } else if ((void(*)())f == (void(*)())tsdefaultcomputejacobian_) {
150:     *TSSetRHSJacobian(*ts,*A,*B,TSDefaultComputeJacobian,fP);
151:   } else if ((void(*)())f == (void(*)())tsdefaultcomputejacobiancolor_) {
152:     *TSSetRHSJacobian(*ts,*A,*B,TSDefaultComputeJacobianColor,*(MatFDColoring*)fP);
153:   } else {
154:   ((PetscObject)*ts)->fortran_func_pointers[3] = (void(*)())f;
155:     *TSSetRHSJacobian(*ts,*A,*B,ourtsjacobian,fP);
156:   }
157: }

159: void PETSC_STDCALL tsgetsolution_(TS *ts,Vec *v,int *ierr)
160: {
161:   *TSGetSolution(*ts,v);
162: }

164: void PETSC_STDCALL tscreate_(MPI_Comm *comm,TSProblemType *problemtype,TS *outts,int *ierr)
165: {
166:   *TSCreate((MPI_Comm)PetscToPointerComm(*comm),*problemtype,outts);
167:   *PetscMalloc(7*sizeof(void *),&((PetscObject)*outts)->fortran_func_pointers);
168: }

170: void PETSC_STDCALL tsgetsnes_(TS *ts,SNES *snes,int *ierr)
171: {
172:   *TSGetSNES(*ts,snes);
173: }

175: void PETSC_STDCALL tsgetsles_(TS *ts,SLES *sles,int *ierr)
176: {
177:   *TSGetSLES(*ts,sles);
178: }

180: void PETSC_STDCALL tsgettype_(TS *ts,CHAR name PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
181: {
182:   char *tname;

184:   *TSGetType(*ts,(TSType *)&tname);
185: #if defined(PETSC_USES_CPTOFCD)
186:   {
187:     char *t = _fcdtocp(name); int len1 = _fcdlen(name);
188:     *PetscStrncpy(t,tname,len1);
189:   }
190: #else
191:   *PetscStrncpy(name,tname,len);
192: #endif
193: }

195: #if defined(PETSC_HAVE_PVODE)  && !defined(__cplusplus)
196: void PETSC_STDCALL tspvodegetiterations_(TS *ts,int *nonlin,int *lin,int *ierr)
197: {
198:   if (FORTRANNULLINTEGER(nonlin)) nonlin = PETSC_NULL;
199:   if (FORTRANNULLINTEGER(lin))    lin    = PETSC_NULL;
200:   *TSPVodeGetIterations(*ts,nonlin,lin);
201: }
202: #endif

204: void PETSC_STDCALL tsdestroy_(TS *ts,int *ierr){
205:   *TSDestroy(*ts);
206: }

208: void PETSC_STDCALL tsdefaultmonitor_(TS *ts,int *step,double *dt,Vec *x,void *ctx,int *ierr)
209: {
210:   *TSDefaultMonitor(*ts,*step,*dt,*x,ctx);
211: }

213: /*
214:    Note ctx is the same as ts so we need to get the Fortran context out of the TS
215: */
216: static int ourtsmonitor(TS ts,int i,double d,Vec v,void*ctx)
217: {
218:   int        0;
219:   void       (*mctx)() = ((PetscObject)ts)->fortran_func_pointers[6];
220:   (*(void (PETSC_STDCALL *)(TS*,int*,double*,Vec*,void(*)(),int*))(((PetscObject)ts)->fortran_func_pointers[4]))(&ts,&i,&d,&v,mctx,&ierr);
221:   return 0;
222: }

224: static int ourtsdestroy(void *ctx)
225: {
226:   int         0;
227:   TS          ts = (TS)ctx;
228:   void        (*mctx)() = ((PetscObject)ts)->fortran_func_pointers[6];
229:   (*(void (PETSC_STDCALL *)(void(*)(),int*))(((PetscObject)ts)->fortran_func_pointers[5]))(mctx,&ierr);
230:   return 0;
231: }

233: void PETSC_STDCALL tssetmonitor_(TS *ts,void (PETSC_STDCALL *func)(TS*,int*,double*,Vec*,void*,int*),void (*mctx)(),void (PETSC_STDCALL *d)(void*,int*),int *ierr)
234: {
235:   if ((void(*)())func == (void(*)())tsdefaultmonitor_) {
236:     *TSSetMonitor(*ts,TSDefaultMonitor,0,0);
237:   } else {
238:     ((PetscObject)*ts)->fortran_func_pointers[4] = (void(*)())func;
239:     ((PetscObject)*ts)->fortran_func_pointers[5] = (void(*)())d;
240:     ((PetscObject)*ts)->fortran_func_pointers[6] = mctx;
241:     if (FORTRANNULLFUNCTION(d)) {
242:       *TSSetMonitor(*ts,ourtsmonitor,*ts,0);
243:     } else {
244:       *TSSetMonitor(*ts,ourtsmonitor,*ts,ourtsdestroy);
245:     }
246:   }
247: }

249: void PETSC_STDCALL tsgetoptionsprefix_(TS *ts,CHAR prefix PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
250: {
251:   char *tname;

253:   *TSGetOptionsPrefix(*ts,&tname);
254: #if defined(PETSC_USES_CPTOFCD)
255:   {
256:     char *t = _fcdtocp(prefix); int len1 = _fcdlen(prefix);
257:     *PetscStrncpy(t,tname,len1);
258:   }
259: #else
260:   *PetscStrncpy(prefix,tname,len);
261: #endif
262: }


265: EXTERN_C_END