Actual source code: zksp.c
2: #include src/fortran/custom/zpetsc.h
3: #include petscksp.h
5: #ifdef PETSC_HAVE_FORTRAN_CAPS
6: #define kspsetnullspace_ KSPSETNULLSPACE
7: #define kspgetnullspace_ KSPGETNULLSPACE
8: #define kspgetresidualnorm_ KSPGETRESIDUALNORM
9: #define kspgetconvergedreason_ KSPGETCONVERGEDREASON
10: #define kspfgmressetmodifypc_ KSPFGMRESSETMODIFYPC
11: #define kspfgmresmodifypcksp_ KSPFGMRESMODIFYPCKSP
12: #define kspfgmresmodifypcnochange_ KSPFGMRESMODIFYPCNOCHANGE
13: #define kspdefaultconverged_ KSPDEFAULTCONVERGED
14: #define kspskipconverged_ KSPSKIPCONVERGED
15: #define kspgmreskrylovmonitor_ KSPGMRESKRYLOVMONITOR
16: #define kspdefaultmonitor_ KSPDEFAULTMONITOR
17: #define ksptruemonitor_ KSPTRUEMONITOR
18: #define kspvecviewmonitor_ KSPVECVIEWMONITOR
19: #define ksplgmonitor_ KSPLGMONITOR
20: #define ksplgtruemonitor_ KSPLGTRUEMONITOR
21: #define kspsingularvaluemonitor_ KSPSINGULARVALUEMONITOR
22: #define kspregisterdestroy_ KSPREGISTERDESTROY
23: #define kspdestroy_ KSPDESTROY
24: #define ksplgmonitordestroy_ KSPLGMONITORDESTROY
25: #define ksplgmonitorcreate_ KSPLGMONITORCREATE
26: #define kspgetrhs_ KSPGETRHS
27: #define kspgetsolution_ KSPGETSOLUTION
28: #define kspgetpc_ KSPGETPC
29: #define kspsetmonitor_ KSPSETMONITOR
30: #define kspsetconvergencetest_ KSPSETCONVERGENCETEST
31: #define kspcreate_ KSPCREATE
32: #define kspsetoptionsprefix_ KSPSETOPTIONSPREFIX
33: #define kspappendoptionsprefix_ KSPAPPENDOPTIONSPREFIX
34: #define kspgettype_ KSPGETTYPE
35: #define kspgetpreconditionerside_ KSPGETPRECONDITIONERSIDE
36: #define kspbuildsolution_ KSPBUILDSOLUTION
37: #define kspbuildresidual_ KSPBUILDRESIDUAL
38: #define kspsettype_ KSPSETTYPE
39: #define kspgetresidualhistory_ KSPGETRESIDUALHISTORY
40: #define kspgetoptionsprefix_ KSPGETOPTIONSPREFIX
41: #define kspview_ KSPVIEW
42: #define kspgmressetrestart_ KSPGMRESSETRESTART
43: #define kspsetnormtype_ KSPSETNORMTYPE
44: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
45: #define kspsetnullspace_ kspsetnullspace
46: #define kspgetnullspace_ kspgetnullspace
47: #define kspgetconvergedreason_ kspgetconvergedreason
48: #define kspfgmressetmodifypc_ kspfgmressetmodifypc
49: #define kspfgmresmodifypcksp_ kspfgmresmodifypcksp
50: #define kspfgmresmodifypcnochange_ kspfgmresmodifypcnochange
51: #define kspdefaultconverged_ kspdefaultconverged
52: #define kspskipconverged_ kspskipconverged
53: #define kspsingularvaluemonitor_ kspsingularvaluemonitor
54: #define kspgmreskrylovmonitor_ kspgmreskrylovmonitor
55: #define kspdefaultmonitor_ kspdefaultmonitor
56: #define ksptruemonitor_ ksptruemonitor
57: #define kspvecviewmonitor_ kspvecviewmonitor
58: #define ksplgmonitor_ ksplgmonitor
59: #define ksplgtruemonitor_ ksplgtruemonitor
60: #define kspgetresidualhistory_ kspgetresidualhistory
61: #define kspsettype_ kspsettype
62: #define kspregisterdestroy_ kspregisterdestroy
63: #define kspdestroy_ kspdestroy
64: #define ksplgmonitordestroy_ ksplgmonitordestroy
65: #define ksplgmonitorcreate_ ksplgmonitorcreate
66: #define kspgetrhs_ kspgetrhs
67: #define kspgetsolution_ kspgetsolution
68: #define kspgetpc_ kspgetpc
69: #define kspsetmonitor_ kspsetmonitor
70: #define kspsetconvergencetest_ kspsetconvergencetest
71: #define kspcreate_ kspcreate
72: #define kspsetoptionsprefix_ kspsetoptionsprefix
73: #define kspappendoptionsprefix_ kspappendoptionsprefix
74: #define kspgettype_ kspgettype
75: #define kspgetpreconditionerside_ kspgetpreconditionerside
76: #define kspbuildsolution_ kspbuildsolution
77: #define kspbuildresidual_ kspbuildresidual
78: #define kspgetoptionsprefix_ kspgetoptionsprefix
79: #define kspview_ kspview
80: #define kspgetresidualnorm_ kspgetresidualnorm
81: #define kspgmressetrestart_ kspgmressetrestart
82: #define kspsetnormtype_ kspsetnormtype
83: #endif
86: static void (PETSC_STDCALL *f2)(KSP*,PetscInt*,PetscReal*,KSPConvergedReason*,void*,PetscErrorCode*);
87: static void (PETSC_STDCALL *f1)(KSP*,PetscInt*,PetscReal*,void*,PetscErrorCode*);
88: static void (PETSC_STDCALL *f21)(void*,PetscErrorCode*);
89: static void (PETSC_STDCALL *f109)(KSP*,PetscInt*,PetscInt*,PetscReal*,void*,PetscErrorCode*);
90: static void (PETSC_STDCALL *f210)(void*,PetscErrorCode*);
94: static PetscErrorCode ourtest(KSP ksp,PetscInt i,PetscReal d,KSPConvergedReason *reason,void* ctx)
95: {
97: (*f2)(&ksp,&i,&d,reason,ctx,&ierr);
98: return 0;
99: }
101: static PetscErrorCode ourmonitor(KSP ksp,PetscInt i,PetscReal d,void* ctx)
102: {
103: PetscErrorCode 0;
104: (*f1)(&ksp,&i,&d,ctx,&ierr);
105: return 0;
106: }
108: static PetscErrorCode ourdestroy(void* ctx)
109: {
110: PetscErrorCode 0;
111: (*f21)(ctx,&ierr);
112: return 0;
113: }
115: static PetscErrorCode ourmodify(KSP ksp,PetscInt i,PetscInt i2,PetscReal d,void* ctx)
116: {
117: PetscErrorCode 0;
118: (*f109)(&ksp,&i,&i2,&d,ctx,&ierr);
119: return 0;
120: }
122: static PetscErrorCode ourmoddestroy(void* ctx)
123: {
124: PetscErrorCode 0;
125: (*f210)(ctx,&ierr);
126: return 0;
127: }
131: void PETSC_STDCALL kspgmressetrestart_(KSP *ksp,PetscInt *max_k, PetscErrorCode *ierr)
132: {
133: *KSPGMRESSetRestart(*ksp,*max_k);
134: }
136: void PETSC_STDCALL kspgetresidualnorm_(KSP *ksp,PetscReal *rnorm,PetscErrorCode *ierr)
137: {
138: *KSPGetResidualNorm(*ksp,rnorm);
139: }
141: void PETSC_STDCALL kspgetconvergedreason_(KSP *ksp,KSPConvergedReason *reason,PetscErrorCode *ierr)
142: {
143: *KSPGetConvergedReason(*ksp,reason);
144: }
146: /* function */
147: void PETSC_STDCALL kspsetnullspace_(KSP *ksp,MatNullSpace *sp, PetscErrorCode *ierr)
148: {
149: *KSPSetNullSpace(*ksp,*sp);
150: }
152: void PETSC_STDCALL kspgetnullspace_(KSP *ksp,MatNullSpace *sp, PetscErrorCode *ierr)
153: {
154: *KSPGetNullSpace(*ksp,sp);
155: }
157: void PETSC_STDCALL kspview_(KSP *ksp,PetscViewer *viewer, PetscErrorCode *ierr)
158: {
159: PetscViewer v;
160: PetscPatchDefaultViewers_Fortran(viewer,v);
161: *KSPView(*ksp,v);
162: }
164: void kspdefaultconverged_(KSP *ksp,PetscInt *n,PetscReal *rnorm,KSPConvergedReason *flag,void *dummy,PetscErrorCode *ierr)
165: {
166: CHKFORTRANNULLOBJECT(dummy);
167: *KSPDefaultConverged(*ksp,*n,*rnorm,flag,dummy);
168: }
170: void kspskipconverged_(KSP *ksp,PetscInt *n,PetscReal *rnorm,KSPConvergedReason *flag,void *dummy,PetscErrorCode *ierr)
171: {
172: CHKFORTRANNULLOBJECT(dummy);
173: *KSPSkipConverged(*ksp,*n,*rnorm,flag,dummy);
174: }
176: void PETSC_STDCALL kspgetresidualhistory_(KSP *ksp,PetscInt *na,PetscErrorCode *ierr)
177: {
178: *KSPGetResidualHistory(*ksp,PETSC_NULL,na);
179: }
181: void PETSC_STDCALL kspsettype_(KSP *ksp,CHAR type PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len))
182: {
183: char *t;
185: FIXCHAR(type,len,t);
186: *KSPSetType(*ksp,t);
187: FREECHAR(type,t);
188: }
190: void PETSC_STDCALL kspgettype_(KSP *ksp,CHAR name PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len))
191: {
192: const char *tname;
194: *KSPGetType(*ksp,&tname);if (*ierr) return;
195: #if defined(PETSC_USES_CPTOFCD)
196: {
197: char *t = _fcdtocp(name); int len1 = _fcdlen(name);
198: *PetscStrncpy(t,tname,len1);
199: }
200: #else
201: *PetscStrncpy(name,tname,len);
202: #endif
203: FIXRETURNCHAR(name,len);
205: }
207: void PETSC_STDCALL kspgetpreconditionerside_(KSP *ksp,PCSide *side,PetscErrorCode *ierr){
208: *KSPGetPreconditionerSide(*ksp,side);
209: }
211: void PETSC_STDCALL kspsetoptionsprefix_(KSP *ksp,CHAR prefix PETSC_MIXED_LEN(len),
212: PetscErrorCode *ierr PETSC_END_LEN(len))
213: {
214: char *t;
216: FIXCHAR(prefix,len,t);
217: *KSPSetOptionsPrefix(*ksp,t);
218: FREECHAR(prefix,t);
219: }
221: void PETSC_STDCALL kspappendoptionsprefix_(KSP *ksp,CHAR prefix PETSC_MIXED_LEN(len),
222: PetscErrorCode *ierr PETSC_END_LEN(len))
223: {
224: char *t;
226: FIXCHAR(prefix,len,t);
227: *KSPAppendOptionsPrefix(*ksp,t);
228: FREECHAR(prefix,t);
229: }
231: void PETSC_STDCALL kspcreate_(MPI_Comm *comm,KSP *ksp,PetscErrorCode *ierr){
232: *KSPCreate((MPI_Comm)PetscToPointerComm(*comm),ksp);
233: }
235: void PETSC_STDCALL kspsetconvergencetest_(KSP *ksp,
236: void (PETSC_STDCALL *converge)(KSP*,PetscInt*,PetscReal*,KSPConvergedReason*,void*,PetscErrorCode*),void *cctx,PetscErrorCode *ierr)
237: {
238: if ((FCNVOID)converge == (FCNVOID)kspdefaultconverged_) {
239: *KSPSetConvergenceTest(*ksp,KSPDefaultConverged,0);
240: } else if ((FCNVOID)converge == (FCNVOID)kspskipconverged_) {
241: *KSPSetConvergenceTest(*ksp,KSPSkipConverged,0);
242: } else {
243: f2 = converge;
244: *KSPSetConvergenceTest(*ksp,ourtest,cctx);
245: }
246: }
248: /*
249: These are not usually called from Fortran but allow Fortran users
250: to transparently set these monitors from .F code
251:
252: functions, hence no STDCALL
253: */
254: void kspgmreskrylovmonitor_(KSP *ksp,PetscInt *it,PetscReal *norm,void *ctx,PetscErrorCode *ierr)
255: {
256: *KSPGMRESKrylovMonitor(*ksp,*it,*norm,ctx);
257: }
259: void kspdefaultmonitor_(KSP *ksp,PetscInt *it,PetscReal *norm,void *ctx,PetscErrorCode *ierr)
260: {
261: *KSPDefaultMonitor(*ksp,*it,*norm,ctx);
262: }
263:
264: void kspsingularvaluemonitor_(KSP *ksp,PetscInt *it,PetscReal *norm,void *ctx,PetscErrorCode *ierr)
265: {
266: *KSPSingularValueMonitor(*ksp,*it,*norm,ctx);
267: }
269: void ksplgmonitor_(KSP *ksp,PetscInt *it,PetscReal *norm,void *ctx,PetscErrorCode *ierr)
270: {
271: *KSPLGMonitor(*ksp,*it,*norm,ctx);
272: }
274: void ksplgtruemonitor_(KSP *ksp,PetscInt *it,PetscReal *norm,void *ctx,PetscErrorCode *ierr)
275: {
276: *KSPLGTrueMonitor(*ksp,*it,*norm,ctx);
277: }
279: void ksptruemonitor_(KSP *ksp,PetscInt *it,PetscReal *norm,void *ctx,PetscErrorCode *ierr)
280: {
281: *KSPTrueMonitor(*ksp,*it,*norm,ctx);
282: }
284: void kspvecviewmonitor_(KSP *ksp,PetscInt *it,PetscReal *norm,void *ctx,PetscErrorCode *ierr)
285: {
286: *KSPVecViewMonitor(*ksp,*it,*norm,ctx);
287: }
290: void PETSC_STDCALL kspsetmonitor_(KSP *ksp,void (PETSC_STDCALL *monitor)(KSP*,PetscInt*,PetscReal*,void*,PetscErrorCode*),
291: void *mctx,void (PETSC_STDCALL *monitordestroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
292: {
293: if ((FCNVOID)monitor == (FCNVOID)kspdefaultmonitor_) {
294: *KSPSetMonitor(*ksp,KSPDefaultMonitor,0,0);
295: } else if ((FCNVOID)monitor == (FCNVOID)ksplgmonitor_) {
296: *KSPSetMonitor(*ksp,KSPLGMonitor,0,0);
297: } else if ((FCNVOID)monitor == (FCNVOID)ksplgtruemonitor_) {
298: *KSPSetMonitor(*ksp,KSPLGTrueMonitor,0,0);
299: } else if ((FCNVOID)monitor == (FCNVOID)kspvecviewmonitor_) {
300: *KSPSetMonitor(*ksp,KSPVecViewMonitor,0,0);
301: } else if ((FCNVOID)monitor == (FCNVOID)ksptruemonitor_) {
302: *KSPSetMonitor(*ksp,KSPTrueMonitor,0,0);
303: } else if ((FCNVOID)monitor == (FCNVOID)kspsingularvaluemonitor_) {
304: *KSPSetMonitor(*ksp,KSPSingularValueMonitor,0,0);
305: } else {
306: f1 = monitor;
307: if (FORTRANNULLFUNCTION(monitordestroy)) {
308: *KSPSetMonitor(*ksp,ourmonitor,mctx,0);
309: } else {
310: f21 = monitordestroy;
311: *KSPSetMonitor(*ksp,ourmonitor,mctx,ourdestroy);
312: }
313: }
314: }
316: void PETSC_STDCALL kspgetpc_(KSP *ksp,PC *B,PetscErrorCode *ierr)
317: {
318: *KSPGetPC(*ksp,B);
319: }
321: void PETSC_STDCALL kspgetsolution_(KSP *ksp,Vec *v,PetscErrorCode *ierr)
322: {
323: *KSPGetSolution(*ksp,v);
324: }
326: void PETSC_STDCALL kspgetrhs_(KSP *ksp,Vec *r,PetscErrorCode *ierr)
327: {
328: *KSPGetRhs(*ksp,r);
329: }
331: /*
332: Possible bleeds memory but cannot be helped.
333: */
334: void PETSC_STDCALL ksplgmonitorcreate_(CHAR host PETSC_MIXED_LEN(len1),
335: CHAR label PETSC_MIXED_LEN(len2),int *x,int *y,int *m,int *n,PetscDrawLG *ctx,
336: PetscErrorCode *ierr PETSC_END_LEN(len1) PETSC_END_LEN(len2))
337: {
338: char *t1,*t2;
340: FIXCHAR(host,len1,t1);
341: FIXCHAR(label,len2,t2);
342: *KSPLGMonitorCreate(t1,t2,*x,*y,*m,*n,ctx);
343: }
345: void PETSC_STDCALL ksplgmonitordestroy_(PetscDrawLG *ctx,PetscErrorCode *ierr)
346: {
347: *KSPLGMonitorDestroy(*ctx);
348: }
350: void PETSC_STDCALL kspdestroy_(KSP *ksp,PetscErrorCode *ierr)
351: {
352: *KSPDestroy(*ksp);
353: }
355: void PETSC_STDCALL kspregisterdestroy_(PetscErrorCode* ierr)
356: {
357: *KSPRegisterDestroy();
358: }
360: void PETSC_STDCALL kspbuildsolution_(KSP *ctx,Vec *v,Vec *V,PetscErrorCode *ierr)
361: {
362: *KSPBuildSolution(*ctx,*v,V);
363: }
365: void PETSC_STDCALL kspbuildresidual_(KSP *ctx,Vec *t,Vec *v,Vec *V,PetscErrorCode *ierr)
366: {
367: *KSPBuildResidual(*ctx,*t,*v,V);
368: }
370: void PETSC_STDCALL kspgetoptionsprefix_(KSP *ksp,CHAR prefix PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len))
371: {
372: const char *tname;
374: *KSPGetOptionsPrefix(*ksp,&tname);
375: #if defined(PETSC_USES_CPTOFCD)
376: {
377: char *t = _fcdtocp(prefix); int len1 = _fcdlen(prefix);
378: *PetscStrncpy(t,tname,len1); if (*ierr) return;
379: }
380: #else
381: *PetscStrncpy(prefix,tname,len); if (*ierr) return;
382: #endif
383: }
385: void PETSC_STDCALL kspfgmresmodifypcnochange_(KSP *ksp,PetscInt *total_its,PetscInt *loc_its,PetscReal *res_norm,void* dummy,PetscErrorCode *ierr)
386: {
387: *KSPFGMRESModifyPCNoChange(*ksp,*total_its,*loc_its,*res_norm,dummy);
388: }
390: void PETSC_STDCALL kspfgmresmodifypcksp_(KSP *ksp,PetscInt *total_its,PetscInt *loc_its,PetscReal *res_norm,void*dummy,PetscErrorCode *ierr)
391: {
392: *KSPFGMRESModifyPCKSP(*ksp,*total_its,*loc_its,*res_norm,dummy);
393: }
395: void PETSC_STDCALL kspfgmressetmodifypc_(KSP *ksp,void (PETSC_STDCALL *fcn)(KSP*,PetscInt*,PetscInt*,PetscReal*,void*,PetscErrorCode*),void* ctx,void (PETSC_STDCALL *d)(void*,PetscErrorCode*),PetscErrorCode *ierr)
396: {
397: if ((FCNVOID)fcn == (FCNVOID)kspfgmresmodifypcksp_) {
398: *KSPFGMRESSetModifyPC(*ksp,KSPFGMRESModifyPCKSP,0,0);
399: } else if ((FCNVOID)fcn == (FCNVOID)kspfgmresmodifypcnochange_) {
400: *KSPFGMRESSetModifyPC(*ksp,KSPFGMRESModifyPCNoChange,0,0);
401: } else {
402: f109 = fcn;
403: if (FORTRANNULLFUNCTION(d)) {
404: *KSPFGMRESSetModifyPC(*ksp,ourmodify,ctx,0);
405: } else {
406: f210 = d;
407: *KSPFGMRESSetModifyPC(*ksp,ourmodify,ctx,ourmoddestroy);
408: }
409: }
410: }
412: void PETSC_STDCALL kspsetnormtype_(KSP *ksp,KSPNormType *type,PetscErrorCode *ierr)
413: {
414: *KSPSetNormType(*ksp,*type);
415: }