Actual source code: zts.c

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

  5: #ifdef PETSC_HAVE_FORTRAN_CAPS
  6: #define tssetproblemtype_                    TSSETPROBLEMTYPE
  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 tsgetksp_                           TSGETKSP
 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 tssetproblemtype_                    tssetproblemtype
 29: #define tsdefaultcomputejacobian_            tsdefaultcomputejacobian
 30: #define tsdefaultcomputejacobiancolor_       tsdefaultcomputejacobiancolor
 31: #define tspvodegetiterations_                tspvodegetiterations
 32: #define tssetrhsfunction_                    tssetrhsfunction
 33: #define tssetrhsmatrix_                      tssetrhsmatrix
 34: #define tssetrhsjacobian_                    tssetrhsjacobian
 35: #define tscreate_                            tscreate
 36: #define tsgetsolution_                       tsgetsolution
 37: #define tsgetsnes_                           tsgetsnes
 38: #define tsgetksp_                           tsgetksp
 39: #define tsgettype_                           tsgettype
 40: #define tsdestroy_                           tsdestroy
 41: #define tssetmonitor_                        tssetmonitor
 42: #define tssettype_                           tssettype
 43: #define tsgetoptionsprefix_                  tsgetoptionsprefix
 44: #define tsdefaultmonitor_                    tsdefaultmonitor
 45: #define tsview_                              tsview
 46: #define tsgetrhsjacobian_                    tsgetrhsjacobian
 47: #define tsgetrhsmatrix_                      tsgetrhsmatrix
 48: #define tssetrhsboundaryconditions_          tssetrhsboundaryconditions
 49: #endif


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

 59: static PetscErrorCode ourtsfunction(TS ts,PetscReal d,Vec x,Vec f,void *ctx)
 60: {
 61:   PetscErrorCode 0;
 62:   (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Vec*,Vec*,void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[1]))(&ts,&d,&x,&f,ctx,&ierr);
 63:   return 0;
 64: }

 66: static PetscErrorCode ourtsmatrix(TS ts,PetscReal d,Mat* m,Mat* p,MatStructure* type,void*ctx)
 67: {
 68:   PetscErrorCode 0;
 69:   (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Mat*,Mat*,MatStructure*,void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[2]))(&ts,&d,m,p,type,ctx,&ierr);
 70:   return 0;
 71: }

 73: static PetscErrorCode ourtsjacobian(TS ts,PetscReal d,Vec x,Mat* m,Mat* p,MatStructure* type,void*ctx)
 74: {
 75:   PetscErrorCode 0;
 76:   (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Vec*,Mat*,Mat*,MatStructure*,void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[3]))(&ts,&d,&x,m,p,type,ctx,&ierr);
 77:   return 0;
 78: }

 80: /*
 81:    Note ctx is the same as ts so we need to get the Fortran context out of the TS
 82: */
 83: static PetscErrorCode ourtsmonitor(TS ts,PetscInt i,PetscReal d,Vec v,void*ctx)
 84: {
 85:   PetscErrorCode 0;
 86:   void       (*mctx)(void) = ((PetscObject)ts)->fortran_func_pointers[6];
 87:   (*(void (PETSC_STDCALL *)(TS*,PetscInt*,PetscReal*,Vec*,FCNVOID,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[4]))(&ts,&i,&d,&v,mctx,&ierr);
 88:   return 0;
 89: }

 91: static PetscErrorCode ourtsdestroy(void *ctx)
 92: {
 93:   PetscErrorCode 0;
 94:   TS          ts = (TS)ctx;
 95:   void        (*mctx)(void) = ((PetscObject)ts)->fortran_func_pointers[6];
 96:   (*(void (PETSC_STDCALL *)(FCNVOID,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[5]))(mctx,&ierr);
 97:   return 0;
 98: }


102: void PETSC_STDCALL tssetrhsboundaryconditions_(TS *ts,PetscErrorCode (PETSC_STDCALL *f)(TS*,PetscReal*,Vec*,void*,PetscErrorCode*),void *ctx,PetscErrorCode *ierr)
103: {
104:   ((PetscObject)*ts)->fortran_func_pointers[0] = (FCNVOID)f;
105:   *TSSetRHSBoundaryConditions(*ts,ourtsbcfunction,ctx);
106: }

108: void PETSC_STDCALL tsgetrhsjacobian_(TS *ts,Mat *J,Mat *M,void **ctx,PetscErrorCode *ierr)
109: {
110:   *TSGetRHSJacobian(*ts,J,M,ctx);
111: }

113: void PETSC_STDCALL tssetproblemtype_(TS *ts,TSProblemType *t,PetscErrorCode *ierr)
114: {
115:   *TSSetProblemType(*ts,*t);
116: }

118: void PETSC_STDCALL tsgetrhsmatrix_(TS *ts,Mat *J,Mat *M,void **ctx,PetscErrorCode *ierr)
119: {
120:   *TSGetRHSMatrix(*ts,J,M,ctx);
121: }

123: void PETSC_STDCALL tsview_(TS *ts,PetscViewer *viewer, PetscErrorCode *ierr)
124: {
125:   PetscViewer v;
126:   PetscPatchDefaultViewers_Fortran(viewer,v);
127:   *TSView(*ts,v);
128: }

130: /* function */
131: void tsdefaultcomputejacobian_(TS *ts,PetscReal *t,Vec *xx1,Mat *J,Mat *B,MatStructure *flag,void *ctx,PetscErrorCode *ierr)
132: {
133:   *TSDefaultComputeJacobian(*ts,*t,*xx1,J,B,flag,ctx);
134: }

136: /* function */
137: void tsdefaultcomputejacobiancolor_(TS *ts,PetscReal *t,Vec *xx1,Mat *J,Mat *B,MatStructure *flag,void *ctx,PetscErrorCode *ierr)
138: {
139:   *TSDefaultComputeJacobianColor(*ts,*t,*xx1,J,B,flag,*(MatFDColoring*)ctx);
140: }

142: void PETSC_STDCALL tssettype_(TS *ts,CHAR type PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len))
143: {
144:   char *t;

146:   FIXCHAR(type,len,t);
147:   *TSSetType(*ts,t);
148:   FREECHAR(type,t);
149: }


152: void PETSC_STDCALL tssetrhsfunction_(TS *ts,PetscErrorCode (PETSC_STDCALL *f)(TS*,PetscReal*,Vec*,Vec*,void*,PetscErrorCode*),void*fP,PetscErrorCode *ierr)
153: {
154:   ((PetscObject)*ts)->fortran_func_pointers[1] = (FCNVOID)f;
155:   *TSSetRHSFunction(*ts,ourtsfunction,fP);
156: }


159: /* ---------------------------------------------------------*/

161: void PETSC_STDCALL tssetrhsmatrix_(TS *ts,Mat *A,Mat *B,PetscErrorCode (PETSC_STDCALL *f)(TS*,PetscReal*,Mat*,Mat*,MatStructure*,
162:                                                    void*,PetscInt *),void*fP,PetscErrorCode *ierr)
163: {
164:   if (FORTRANNULLFUNCTION(f)) {
165:     *TSSetRHSMatrix(*ts,*A,*B,PETSC_NULL,fP);
166:   } else {
167:     ((PetscObject)*ts)->fortran_func_pointers[2] = (FCNVOID)f;
168:     *TSSetRHSMatrix(*ts,*A,*B,ourtsmatrix,fP);
169:   }
170: }

172: /* ---------------------------------------------------------*/

174: void PETSC_STDCALL tssetrhsjacobian_(TS *ts,Mat *A,Mat *B,void (PETSC_STDCALL *f)(TS*,PetscReal*,Vec*,Mat*,Mat*,MatStructure*,
175:                void*,PetscErrorCode*),void*fP,PetscErrorCode *ierr)
176: {
177:   if (FORTRANNULLFUNCTION(f)) {
178:     *TSSetRHSJacobian(*ts,*A,*B,PETSC_NULL,fP);
179:   } else if ((FCNVOID)f == (FCNVOID)tsdefaultcomputejacobian_) {
180:     *TSSetRHSJacobian(*ts,*A,*B,TSDefaultComputeJacobian,fP);
181:   } else if ((FCNVOID)f == (FCNVOID)tsdefaultcomputejacobiancolor_) {
182:     *TSSetRHSJacobian(*ts,*A,*B,TSDefaultComputeJacobianColor,*(MatFDColoring*)fP);
183:   } else {
184:   ((PetscObject)*ts)->fortran_func_pointers[3] = (FCNVOID)f;
185:     *TSSetRHSJacobian(*ts,*A,*B,ourtsjacobian,fP);
186:   }
187: }

189: void PETSC_STDCALL tsgetsolution_(TS *ts,Vec *v,PetscErrorCode *ierr)
190: {
191:   *TSGetSolution(*ts,v);
192: }

194: void PETSC_STDCALL tscreate_(MPI_Comm *comm,TS *outts,PetscErrorCode *ierr)
195: {
196:   *TSCreate((MPI_Comm)PetscToPointerComm(*comm),outts);
197:   *PetscMalloc(7*sizeof(void*),&((PetscObject)*outts)->fortran_func_pointers);
198: }

200: void PETSC_STDCALL tsgetsnes_(TS *ts,SNES *snes,PetscErrorCode *ierr)
201: {
202:   *TSGetSNES(*ts,snes);
203: }

205: void PETSC_STDCALL tsgetksp_(TS *ts,KSP *ksp,PetscErrorCode *ierr)
206: {
207:   *TSGetKSP(*ts,ksp);
208: }

210: void PETSC_STDCALL tsgettype_(TS *ts,CHAR name PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len))
211: {
212:   char *tname;

214:   *TSGetType(*ts,(TSType *)&tname);
215: #if defined(PETSC_USES_CPTOFCD)
216:   {
217:     char *t = _fcdtocp(name); int len1 = _fcdlen(name);
218:     *PetscStrncpy(t,tname,len1);
219:   }
220: #else
221:   *PetscStrncpy(name,tname,len);
222: #endif
223:   FIXRETURNCHAR(name,len);
224: }

226: #if defined(PETSC_HAVE_PVODE)
227: void PETSC_STDCALL tspvodegetiterations_(TS *ts,PetscInt *nonlin,PetscInt *lin,PetscErrorCode *ierr)
228: {
229:   CHKFORTRANNULLINTEGER(nonlin);
230:   CHKFORTRANNULLINTEGER(lin);
231:   *TSPVodeGetIterations(*ts,nonlin,lin);
232: }
233: #endif

235: void PETSC_STDCALL tsdestroy_(TS *ts,PetscErrorCode *ierr){
236:   *TSDestroy(*ts);
237: }

239: void PETSC_STDCALL tsdefaultmonitor_(TS *ts,PetscInt *step,PetscReal *dt,Vec *x,void *ctx,PetscErrorCode *ierr)
240: {
241:   *TSDefaultMonitor(*ts,*step,*dt,*x,ctx);
242: }


245: void PETSC_STDCALL tssetmonitor_(TS *ts,void (PETSC_STDCALL *func)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*),void (*mctx)(void),void (PETSC_STDCALL *d)(void*,PetscErrorCode*),PetscErrorCode *ierr)
246: {
247:   if ((FCNVOID)func == (FCNVOID)tsdefaultmonitor_) {
248:     *TSSetMonitor(*ts,TSDefaultMonitor,0,0);
249:   } else {
250:     ((PetscObject)*ts)->fortran_func_pointers[4] = (FCNVOID)func;
251:     ((PetscObject)*ts)->fortran_func_pointers[5] = (FCNVOID)d;
252:     ((PetscObject)*ts)->fortran_func_pointers[6] = (FCNVOID)mctx;
253:     if (FORTRANNULLFUNCTION(d)) {
254:       *TSSetMonitor(*ts,ourtsmonitor,*ts,0);
255:     } else {
256:       *TSSetMonitor(*ts,ourtsmonitor,*ts,ourtsdestroy);
257:     }
258:   }
259: }

261: void PETSC_STDCALL tsgetoptionsprefix_(TS *ts,CHAR prefix PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len))
262: {
263:   const char *tname;

265:   *TSGetOptionsPrefix(*ts,&tname);
266: #if defined(PETSC_USES_CPTOFCD)
267:   {
268:     char *t = _fcdtocp(prefix); int len1 = _fcdlen(prefix);
269:     *PetscStrncpy(t,tname,len1);
270:   }
271: #else
272:   *PetscStrncpy(prefix,tname,len);
273: #endif
274: }