Actual source code: zksp.c
1: /*$Id: zksp.c,v 1.52 2001/08/07 21:32:16 bsmith Exp $*/
3: #include src/fortran/custom/zpetsc.h
4: #include petscksp.h
6: #ifdef PETSC_HAVE_FORTRAN_CAPS
7: #define kspgetresidualnorm_ KSPGETRESIDUALNORM
8: #define kspgetconvergedreason_ KSPGETCONVERGEDREASON
9: #define kspfgmressetmodifypc_ KSPFGMRESSETMODIFYPC
10: #define kspfgmresmodifypcsles_ KSPFGMRESMODIFYPCSLES
11: #define kspfgmresmodifypcnochange_ KSPFGMRESMODIFYPCNOCHANGE
12: #define kspdefaultconverged_ KSPDEFAULTCONVERGED
13: #define kspskipconverged_ KSPSKIPCONVERGED
14: #define kspdefaultmonitor_ KSPDEFAULTMONITOR
15: #define ksptruemonitor_ KSPTRUEMONITOR
16: #define kspvecviewmonitor_ KSPVECVIEWMONITOR
17: #define ksplgmonitor_ KSPLGMONITOR
18: #define ksplgtruemonitor_ KSPLGTRUEMONITOR
19: #define kspsingularvaluemonitor_ KSPSINGULARVALUEMONITOR
20: #define kspregisterdestroy_ KSPREGISTERDESTROY
21: #define kspdestroy_ KSPDESTROY
22: #define ksplgmonitordestroy_ KSPLGMONITORDESTROY
23: #define ksplgmonitorcreate_ KSPLGMONITORCREATE
24: #define kspgetrhs_ KSPGETRHS
25: #define kspgetsolution_ KSPGETSOLUTION
26: #define kspgetpc_ KSPGETPC
27: #define kspsetmonitor_ KSPSETMONITOR
28: #define kspsetconvergencetest_ KSPSETCONVERGENCETEST
29: #define kspcreate_ KSPCREATE
30: #define kspsetoptionsprefix_ KSPSETOPTIONSPREFIX
31: #define kspappendoptionsprefix_ KSPAPPENDOPTIONSPREFIX
32: #define kspgettype_ KSPGETTYPE
33: #define kspgetpreconditionerside_ KSPGETPRECONDITIONERSIDE
34: #define kspbuildsolution_ KSPBUILDSOLUTION
35: #define kspsettype_ KSPSETTYPE
36: #define kspgetresidualhistory_ KSPGETRESIDUALHISTORY
37: #define kspgetoptionsprefix_ KSPGETOPTIONSPREFIX
38: #define kspview_ KSPVIEW
39: #define kspgmressetrestart_ KSPGMRESSETRESTART
40: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
41: #define kspgetconvergedreason_ kspgetconvergedreason
42: #define kspfgmressetmodifypc_ kspfgmressetmodifypc
43: #define kspfgmresmodifypcsles_ kspfgmresmodifypcsles
44: #define kspfgmresmodifypcnochange_ kspfgmresmodifypcnochange
45: #define kspdefaultconverged_ kspdefaultconverged
46: #define kspskipconverged_ kspskipconverged
47: #define kspsingularvaluemonitor_ kspsingularvaluemonitor
48: #define kspdefaultmonitor_ kspdefaultmonitor
49: #define ksptruemonitor_ ksptruemonitor
50: #define kspvecviewmonitor_ kspvecviewmonitor
51: #define ksplgmonitor_ ksplgmonitor
52: #define ksplgtruemonitor_ ksplgtruemonitor
53: #define kspgetresidualhistory_ kspgetresidualhistory
54: #define kspsettype_ kspsettype
55: #define kspregisterdestroy_ kspregisterdestroy
56: #define kspdestroy_ kspdestroy
57: #define ksplgmonitordestroy_ ksplgmonitordestroy
58: #define ksplgmonitorcreate_ ksplgmonitorcreate
59: #define kspgetrhs_ kspgetrhs
60: #define kspgetsolution_ kspgetsolution
61: #define kspgetpc_ kspgetpc
62: #define kspsetmonitor_ kspsetmonitor
63: #define kspsetconvergencetest_ kspsetconvergencetest
64: #define kspcreate_ kspcreate
65: #define kspsetoptionsprefix_ kspsetoptionsprefix
66: #define kspappendoptionsprefix_ kspappendoptionsprefix
67: #define kspgettype_ kspgettype
68: #define kspgetpreconditionerside_ kspgetpreconditionerside
69: #define kspbuildsolution_ kspbuildsolution
70: #define kspgetoptionsprefix_ kspgetoptionsprefix
71: #define kspview_ kspview
72: #define kspgetresidualnorm_ kspgetresidualnorm
73: #define kspgmressetrestart_ kspgmressetrestart
74: #endif
76: EXTERN_C_BEGIN
78: void PETSC_STDCALL kspgmressetrestart_(KSP *ksp,int *max_k, int *ierr )
79: {
80: *KSPGMRESSetRestart(*ksp,*max_k);
81: }
83: void PETSC_STDCALL kspgetresidualnorm_(KSP *ksp,PetscReal *rnorm,int *ierr)
84: {
85: *KSPGetResidualNorm(*ksp,rnorm);
86: }
88: void PETSC_STDCALL kspgetconvergedreason_(KSP *ksp,KSPConvergedReason *reason,int *ierr)
89: {
90: *KSPGetConvergedReason(*ksp,reason);
91: }
93: /* function */
94: void PETSC_STDCALL kspview_(KSP *ksp,PetscViewer *viewer, int *ierr)
95: {
96: PetscViewer v;
97: PetscPatchDefaultViewers_Fortran(viewer,v);
98: *KSPView(*ksp,v);
99: }
101: void kspdefaultconverged_(KSP *ksp,int *n,PetscReal *rnorm,KSPConvergedReason *flag,void *dummy,int *ierr)
102: {
103: CHKFORTRANNULLOBJECT(dummy);
104: *KSPDefaultConverged(*ksp,*n,*rnorm,flag,dummy);
105: }
107: void kspskipconverged_(KSP *ksp,int *n,PetscReal *rnorm,KSPConvergedReason *flag,void *dummy,int *ierr)
108: {
109: CHKFORTRANNULLOBJECT(dummy);
110: *KSPSkipConverged(*ksp,*n,*rnorm,flag,dummy);
111: }
113: void PETSC_STDCALL kspgetresidualhistory_(KSP *ksp,int *na,int *ierr)
114: {
115: *KSPGetResidualHistory(*ksp,PETSC_NULL,na);
116: }
118: void PETSC_STDCALL kspsettype_(KSP *ksp,CHAR type PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
119: {
120: char *t;
122: FIXCHAR(type,len,t);
123: *KSPSetType(*ksp,t);
124: FREECHAR(type,t);
125: }
127: void PETSC_STDCALL kspgettype_(KSP *ksp,CHAR name PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
128: {
129: char *tname;
131: *KSPGetType(*ksp,&tname);if (*ierr) return;
132: #if defined(PETSC_USES_CPTOFCD)
133: {
134: char *t = _fcdtocp(name); int len1 = _fcdlen(name);
135: *PetscStrncpy(t,tname,len1);
136: }
137: #else
138: *PetscStrncpy(name,tname,len);
139: #endif
140: }
142: void PETSC_STDCALL kspgetpreconditionerside_(KSP *ksp,PCSide *side,int *ierr){
143: *KSPGetPreconditionerSide(*ksp,side);
144: }
146: void PETSC_STDCALL kspsetoptionsprefix_(KSP *ksp,CHAR prefix PETSC_MIXED_LEN(len),
147: int *ierr PETSC_END_LEN(len))
148: {
149: char *t;
151: FIXCHAR(prefix,len,t);
152: *KSPSetOptionsPrefix(*ksp,t);
153: FREECHAR(prefix,t);
154: }
156: void PETSC_STDCALL kspappendoptionsprefix_(KSP *ksp,CHAR prefix PETSC_MIXED_LEN(len),
157: int *ierr PETSC_END_LEN(len))
158: {
159: char *t;
161: FIXCHAR(prefix,len,t);
162: *KSPAppendOptionsPrefix(*ksp,t);
163: FREECHAR(prefix,t);
164: }
166: void PETSC_STDCALL kspcreate_(MPI_Comm *comm,KSP *ksp,int *ierr){
167: *KSPCreate((MPI_Comm)PetscToPointerComm(*comm),ksp);
168: }
170: static void (PETSC_STDCALL *f2)(KSP*,int*,PetscReal*,KSPConvergedReason*,void*,int*);
171: static int ourtest(KSP ksp,int i,PetscReal d,KSPConvergedReason *reason,void* ctx)
172: {
174: (*f2)(&ksp,&i,&d,reason,ctx,&ierr);
175: return 0;
176: }
177: void PETSC_STDCALL kspsetconvergencetest_(KSP *ksp,
178: void (PETSC_STDCALL *converge)(KSP*,int*,PetscReal*,KSPConvergedReason*,void*,int*),void *cctx,int *ierr)
179: {
180: if ((void(*)(void))converge == (void(*)(void))kspdefaultconverged_) {
181: *KSPSetConvergenceTest(*ksp,KSPDefaultConverged,0);
182: } else if ((void(*)(void))converge == (void(*)(void))kspskipconverged_) {
183: *KSPSetConvergenceTest(*ksp,KSPSkipConverged,0);
184: } else {
185: f2 = converge;
186: *KSPSetConvergenceTest(*ksp,ourtest,cctx);
187: }
188: }
190: /*
191: These are not usually called from Fortran but allow Fortran users
192: to transparently set these monitors from .F code
193:
194: functions, hence no STDCALL
195: */
196: void kspdefaultmonitor_(KSP *ksp,int *it,PetscReal *norm,void *ctx,int *ierr)
197: {
198: *KSPDefaultMonitor(*ksp,*it,*norm,ctx);
199: }
200:
201: void kspsingularvaluemonitor_(KSP *ksp,int *it,PetscReal *norm,void *ctx,int *ierr)
202: {
203: *KSPSingularValueMonitor(*ksp,*it,*norm,ctx);
204: }
206: void ksplgmonitor_(KSP *ksp,int *it,PetscReal *norm,void *ctx,int *ierr)
207: {
208: *KSPLGMonitor(*ksp,*it,*norm,ctx);
209: }
211: void ksplgtruemonitor_(KSP *ksp,int *it,PetscReal *norm,void *ctx,int *ierr)
212: {
213: *KSPLGTrueMonitor(*ksp,*it,*norm,ctx);
214: }
216: void ksptruemonitor_(KSP *ksp,int *it,PetscReal *norm,void *ctx,int *ierr)
217: {
218: *KSPTrueMonitor(*ksp,*it,*norm,ctx);
219: }
221: void kspvecviewmonitor_(KSP *ksp,int *it,PetscReal *norm,void *ctx,int *ierr)
222: {
223: *KSPVecViewMonitor(*ksp,*it,*norm,ctx);
224: }
226: static void (PETSC_STDCALL *f1)(KSP*,int*,PetscReal*,void*,int*);
227: static int ourmonitor(KSP ksp,int i,PetscReal d,void* ctx)
228: {
229: int 0;
230: (*f1)(&ksp,&i,&d,ctx,&ierr);
231: return 0;
232: }
233: static void (PETSC_STDCALL *f21)(void*,int*);
234: static int ourdestroy(void* ctx)
235: {
236: int 0;
237: (*f21)(ctx,&ierr);
238: return 0;
239: }
241: void PETSC_STDCALL kspsetmonitor_(KSP *ksp,void (PETSC_STDCALL *monitor)(KSP*,int*,PetscReal*,void*,int*),
242: void *mctx,void (PETSC_STDCALL *monitordestroy)(void *,int *),int *ierr)
243: {
244: if ((void(*)(void))monitor == (void(*)(void))kspdefaultmonitor_) {
245: *KSPSetMonitor(*ksp,KSPDefaultMonitor,0,0);
246: } else if ((void(*)(void))monitor == (void(*)(void))ksplgmonitor_) {
247: *KSPSetMonitor(*ksp,KSPLGMonitor,0,0);
248: } else if ((void(*)(void))monitor == (void(*)(void))ksplgtruemonitor_) {
249: *KSPSetMonitor(*ksp,KSPLGTrueMonitor,0,0);
250: } else if ((void(*)(void))monitor == (void(*)(void))kspvecviewmonitor_) {
251: *KSPSetMonitor(*ksp,KSPVecViewMonitor,0,0);
252: } else if ((void(*)(void))monitor == (void(*)(void))ksptruemonitor_) {
253: *KSPSetMonitor(*ksp,KSPTrueMonitor,0,0);
254: } else if ((void(*)(void))monitor == (void(*)(void))kspsingularvaluemonitor_) {
255: *KSPSetMonitor(*ksp,KSPSingularValueMonitor,0,0);
256: } else {
257: f1 = monitor;
258: if (FORTRANNULLFUNCTION(monitordestroy)) {
259: *KSPSetMonitor(*ksp,ourmonitor,mctx,0);
260: } else {
261: f21 = monitordestroy;
262: *KSPSetMonitor(*ksp,ourmonitor,mctx,ourdestroy);
263: }
264: }
265: }
267: void PETSC_STDCALL kspgetpc_(KSP *ksp,PC *B,int *ierr)
268: {
269: *KSPGetPC(*ksp,B);
270: }
272: void PETSC_STDCALL kspgetsolution_(KSP *ksp,Vec *v,int *ierr)
273: {
274: *KSPGetSolution(*ksp,v);
275: }
277: void PETSC_STDCALL kspgetrhs_(KSP *ksp,Vec *r,int *ierr)
278: {
279: *KSPGetRhs(*ksp,r);
280: }
282: /*
283: Possible bleeds memory but cannot be helped.
284: */
285: void PETSC_STDCALL ksplgmonitorcreate_(CHAR host PETSC_MIXED_LEN(len1),
286: CHAR label PETSC_MIXED_LEN(len2),int *x,int *y,int *m,int *n,PetscDrawLG *ctx,
287: int *ierr PETSC_END_LEN(len1) PETSC_END_LEN(len2))
288: {
289: char *t1,*t2;
291: FIXCHAR(host,len1,t1);
292: FIXCHAR(label,len2,t2);
293: *KSPLGMonitorCreate(t1,t2,*x,*y,*m,*n,ctx);
294: }
296: void PETSC_STDCALL ksplgmonitordestroy_(PetscDrawLG *ctx,int *ierr)
297: {
298: *KSPLGMonitorDestroy(*ctx);
299: }
301: void PETSC_STDCALL kspdestroy_(KSP *ksp,int *ierr)
302: {
303: *KSPDestroy(*ksp);
304: }
306: void PETSC_STDCALL kspregisterdestroy_(int* ierr)
307: {
308: *KSPRegisterDestroy();
309: }
311: void PETSC_STDCALL kspbuildsolution_(KSP *ctx,Vec *v,Vec *V,int *ierr)
312: {
313: *KSPBuildSolution(*ctx,*v,V);
314: }
316: void PETSC_STDCALL kspbuildresidual_(KSP *ctx,Vec *t,Vec *v,Vec *V,int *ierr)
317: {
318: *KSPBuildResidual(*ctx,*t,*v,V);
319: }
321: void PETSC_STDCALL kspgetoptionsprefix_(KSP *ksp,CHAR prefix PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
322: {
323: char *tname;
325: *KSPGetOptionsPrefix(*ksp,&tname);
326: #if defined(PETSC_USES_CPTOFCD)
327: {
328: char *t = _fcdtocp(prefix); int len1 = _fcdlen(prefix);
329: *PetscStrncpy(t,tname,len1); if (*ierr) return;
330: }
331: #else
332: *PetscStrncpy(prefix,tname,len); if (*ierr) return;
333: #endif
334: }
336: static void (PETSC_STDCALL *f109)(KSP*,int*,int*,PetscReal*,void*,int*);
337: static int ourmodify(KSP ksp,int i,int i2,PetscReal d,void* ctx)
338: {
339: int 0;
340: (*f109)(&ksp,&i,&i2,&d,ctx,&ierr);
341: return 0;
342: }
344: static void (PETSC_STDCALL *f210)(void*,int*);
345: static int ourmoddestroy(void* ctx)
346: {
347: int 0;
348: (*f210)(ctx,&ierr);
349: return 0;
350: }
352: void PETSC_STDCALL kspfgmresmodifypcnochange_(KSP *ksp,int *total_its,int *loc_its,PetscReal *res_norm,void* dummy,int *ierr)
353: {
354: *KSPFGMRESModifyPCNoChange(*ksp,*total_its,*loc_its,*res_norm,dummy);
355: }
357: void PETSC_STDCALL kspfgmresmodifypcsles_(KSP *ksp,int *total_its,int *loc_its,PetscReal *res_norm,void*dummy,int *ierr)
358: {
359: *KSPFGMRESModifyPCSLES(*ksp,*total_its,*loc_its,*res_norm,dummy);
360: }
362: void PETSC_STDCALL kspfgmressetmodifypc_(KSP *ksp,void (PETSC_STDCALL *fcn)(KSP*,int*,int*,PetscReal*,void*,int*),void* ctx,void (PETSC_STDCALL *d)(void*,int*),int *ierr)
363: {
364: if ((void(*)(void))fcn == (void(*)(void))kspfgmresmodifypcsles_) {
365: *KSPFGMRESSetModifyPC(*ksp,KSPFGMRESModifyPCSLES,0,0);
366: } else if ((void(*)(void))fcn == (void(*)(void))kspfgmresmodifypcnochange_) {
367: *KSPFGMRESSetModifyPC(*ksp,KSPFGMRESModifyPCNoChange,0,0);
368: } else {
369: f109 = fcn;
370: if (FORTRANNULLFUNCTION(d)) {
371: *KSPFGMRESSetModifyPC(*ksp,ourmodify,ctx,0);
372: } else {
373: f210 = d;
374: *KSPFGMRESSetModifyPC(*ksp,ourmodify,ctx,ourmoddestroy);
375: }
376: }
377: }
379: EXTERN_C_END