Actual source code: zts.c

  1: /*$Id: zts.c,v 1.39 2001/09/25 14:32:39 balay Exp $*/

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

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

 52: EXTERN_C_BEGIN

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

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

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

 72: void PETSC_STDCALL tssetproblemtype_(TS *ts,TSProblemType *t,int *ierr)
 73: {
 74:   *TSSetProblemType(*ts,*t);
 75: }

 77: void PETSC_STDCALL tsgetrhsmatrix_(TS *ts,Mat *J,Mat *M,void **ctx,int *ierr)
 78: {
 79:   *TSGetRHSMatrix(*ts,J,M,ctx);
 80: }

 82: void PETSC_STDCALL tsview_(TS *ts,PetscViewer *viewer, int *ierr)
 83: {
 84:   PetscViewer v;
 85:   PetscPatchDefaultViewers_Fortran(viewer,v);
 86:   *TSView(*ts,v);
 87: }

 89: /* function */
 90: void tsdefaultcomputejacobian_(TS *ts,PetscReal *t,Vec *xx1,Mat *J,Mat *B,MatStructure *flag,void *ctx,int *ierr)
 91: {
 92:   *TSDefaultComputeJacobian(*ts,*t,*xx1,J,B,flag,ctx);
 93: }

 95: /* function */
 96: void tsdefaultcomputejacobiancolor_(TS *ts,PetscReal *t,Vec *xx1,Mat *J,Mat *B,MatStructure *flag,void *ctx,int *ierr)
 97: {
 98:   *TSDefaultComputeJacobianColor(*ts,*t,*xx1,J,B,flag,*(MatFDColoring*)ctx);
 99: }

101: void PETSC_STDCALL tssettype_(TS *ts,CHAR type PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
102: {
103:   char *t;

105:   FIXCHAR(type,len,t);
106:   *TSSetType(*ts,t);
107:   FREECHAR(type,t);
108: }

110: static int ourtsfunction(TS ts,PetscReal d,Vec x,Vec f,void *ctx)
111: {
112:   int 0;
113:   (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Vec*,Vec*,void*,int*))(((PetscObject)ts)->fortran_func_pointers[1]))(&ts,&d,&x,&f,ctx,&ierr);
114:   return 0;
115: }

117: void PETSC_STDCALL tssetrhsfunction_(TS *ts,int (PETSC_STDCALL *f)(TS*,PetscReal*,Vec*,Vec*,void*,int*),void*fP,int *ierr)
118: {
119:   ((PetscObject)*ts)->fortran_func_pointers[1] = (void(*)(void))f;
120:   *TSSetRHSFunction(*ts,ourtsfunction,fP);
121: }


124: /* ---------------------------------------------------------*/
125: static int ourtsmatrix(TS ts,PetscReal d,Mat* m,Mat* p,MatStructure* type,void*ctx)
126: {
127:   int 0;
128:   (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Mat*,Mat*,MatStructure*,void*,int*))(((PetscObject)ts)->fortran_func_pointers[2]))(&ts,&d,m,p,type,ctx,&ierr);
129:   return 0;
130: }

132: void PETSC_STDCALL tssetrhsmatrix_(TS *ts,Mat *A,Mat *B,int (PETSC_STDCALL *f)(TS*,PetscReal*,Mat*,Mat*,MatStructure*,
133:                                                    void*,int *),void*fP,int *ierr)
134: {
135:   if (FORTRANNULLFUNCTION(f)) {
136:     *TSSetRHSMatrix(*ts,*A,*B,PETSC_NULL,fP);
137:   } else {
138:     ((PetscObject)*ts)->fortran_func_pointers[2] = (void(*)(void))f;
139:     *TSSetRHSMatrix(*ts,*A,*B,ourtsmatrix,fP);
140:   }
141: }

143: /* ---------------------------------------------------------*/
144: static int ourtsjacobian(TS ts,PetscReal d,Vec x,Mat* m,Mat* p,MatStructure* type,void*ctx)
145: {
146:   int 0;
147:   (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Vec*,Mat*,Mat*,MatStructure*,void*,int*))(((PetscObject)ts)->fortran_func_pointers[3]))(&ts,&d,&x,m,p,type,ctx,&ierr);
148:   return 0;
149: }

151: void PETSC_STDCALL tssetrhsjacobian_(TS *ts,Mat *A,Mat *B,void (PETSC_STDCALL *f)(TS*,PetscReal*,Vec*,Mat*,Mat*,MatStructure*,
152:                void*,int*),void*fP,int *ierr)
153: {
154:   if (FORTRANNULLFUNCTION(f)) {
155:     *TSSetRHSJacobian(*ts,*A,*B,PETSC_NULL,fP);
156:   } else if ((void(*)(void))f == (void(*)(void))tsdefaultcomputejacobian_) {
157:     *TSSetRHSJacobian(*ts,*A,*B,TSDefaultComputeJacobian,fP);
158:   } else if ((void(*)(void))f == (void(*)(void))tsdefaultcomputejacobiancolor_) {
159:     *TSSetRHSJacobian(*ts,*A,*B,TSDefaultComputeJacobianColor,*(MatFDColoring*)fP);
160:   } else {
161:   ((PetscObject)*ts)->fortran_func_pointers[3] = (void(*)(void))f;
162:     *TSSetRHSJacobian(*ts,*A,*B,ourtsjacobian,fP);
163:   }
164: }

166: void PETSC_STDCALL tsgetsolution_(TS *ts,Vec *v,int *ierr)
167: {
168:   *TSGetSolution(*ts,v);
169: }

171: void PETSC_STDCALL tscreate_(MPI_Comm *comm,TS *outts,int *ierr)
172: {
173:   *TSCreate((MPI_Comm)PetscToPointerComm(*comm),outts);
174:   *PetscMalloc(7*sizeof(void *),&((PetscObject)*outts)->fortran_func_pointers);
175: }

177: void PETSC_STDCALL tsgetsnes_(TS *ts,SNES *snes,int *ierr)
178: {
179:   *TSGetSNES(*ts,snes);
180: }

182: void PETSC_STDCALL tsgetsles_(TS *ts,SLES *sles,int *ierr)
183: {
184:   *TSGetSLES(*ts,sles);
185: }

187: void PETSC_STDCALL tsgettype_(TS *ts,CHAR name PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
188: {
189:   char *tname;

191:   *TSGetType(*ts,(TSType *)&tname);
192: #if defined(PETSC_USES_CPTOFCD)
193:   {
194:     char *t = _fcdtocp(name); int len1 = _fcdlen(name);
195:     *PetscStrncpy(t,tname,len1);
196:   }
197: #else
198:   *PetscStrncpy(name,tname,len);
199: #endif
200: }

202: #if defined(PETSC_HAVE_PVODE)  && !defined(__cplusplus)
203: void PETSC_STDCALL tspvodegetiterations_(TS *ts,int *nonlin,int *lin,int *ierr)
204: {
205:   CHKFORTRANNULLINTEGER(nonlin);
206:   CHKFORTRANNULLINTEGER(lin);
207:   *TSPVodeGetIterations(*ts,nonlin,lin);
208: }
209: #endif

211: void PETSC_STDCALL tsdestroy_(TS *ts,int *ierr){
212:   *TSDestroy(*ts);
213: }

215: void PETSC_STDCALL tsdefaultmonitor_(TS *ts,int *step,PetscReal *dt,Vec *x,void *ctx,int *ierr)
216: {
217:   *TSDefaultMonitor(*ts,*step,*dt,*x,ctx);
218: }

220: /*
221:    Note ctx is the same as ts so we need to get the Fortran context out of the TS
222: */
223: static int ourtsmonitor(TS ts,int i,PetscReal d,Vec v,void*ctx)
224: {
225:   int        0;
226:   void       (*mctx)(void) = ((PetscObject)ts)->fortran_func_pointers[6];
227:   (*(void (PETSC_STDCALL *)(TS*,int*,PetscReal*,Vec*,void(*)(void),int*))(((PetscObject)ts)->fortran_func_pointers[4]))(&ts,&i,&d,&v,mctx,&ierr);
228:   return 0;
229: }

231: static int ourtsdestroy(void *ctx)
232: {
233:   int         0;
234:   TS          ts = (TS)ctx;
235:   void        (*mctx)(void) = ((PetscObject)ts)->fortran_func_pointers[6];
236:   (*(void (PETSC_STDCALL *)(void(*)(void),int*))(((PetscObject)ts)->fortran_func_pointers[5]))(mctx,&ierr);
237:   return 0;
238: }

240: void PETSC_STDCALL tssetmonitor_(TS *ts,void (PETSC_STDCALL *func)(TS*,int*,PetscReal*,Vec*,void*,int*),void (*mctx)(void),void (PETSC_STDCALL *d)(void*,int*),int *ierr)
241: {
242:   if ((void(*)(void))func == (void(*)(void))tsdefaultmonitor_) {
243:     *TSSetMonitor(*ts,TSDefaultMonitor,0,0);
244:   } else {
245:     ((PetscObject)*ts)->fortran_func_pointers[4] = (void(*)(void))func;
246:     ((PetscObject)*ts)->fortran_func_pointers[5] = (void(*)(void))d;
247:     ((PetscObject)*ts)->fortran_func_pointers[6] = mctx;
248:     if (FORTRANNULLFUNCTION(d)) {
249:       *TSSetMonitor(*ts,ourtsmonitor,*ts,0);
250:     } else {
251:       *TSSetMonitor(*ts,ourtsmonitor,*ts,ourtsdestroy);
252:     }
253:   }
254: }

256: void PETSC_STDCALL tsgetoptionsprefix_(TS *ts,CHAR prefix PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
257: {
258:   char *tname;

260:   *TSGetOptionsPrefix(*ts,&tname);
261: #if defined(PETSC_USES_CPTOFCD)
262:   {
263:     char *t = _fcdtocp(prefix); int len1 = _fcdlen(prefix);
264:     *PetscStrncpy(t,tname,len1);
265:   }
266: #else
267:   *PetscStrncpy(prefix,tname,len);
268: #endif
269: }


272: EXTERN_C_END