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 kspgmreskrylovmonitor_ KSPGMRESKRYLOVMONITOR
15: #define kspdefaultmonitor_ KSPDEFAULTMONITOR
16: #define ksptruemonitor_ KSPTRUEMONITOR
17: #define kspvecviewmonitor_ KSPVECVIEWMONITOR
18: #define ksplgmonitor_ KSPLGMONITOR
19: #define ksplgtruemonitor_ KSPLGTRUEMONITOR
20: #define kspsingularvaluemonitor_ KSPSINGULARVALUEMONITOR
21: #define kspregisterdestroy_ KSPREGISTERDESTROY
22: #define kspdestroy_ KSPDESTROY
23: #define ksplgmonitordestroy_ KSPLGMONITORDESTROY
24: #define ksplgmonitorcreate_ KSPLGMONITORCREATE
25: #define kspgetrhs_ KSPGETRHS
26: #define kspgetsolution_ KSPGETSOLUTION
27: #define kspgetpc_ KSPGETPC
28: #define kspsetmonitor_ KSPSETMONITOR
29: #define kspsetconvergencetest_ KSPSETCONVERGENCETEST
30: #define kspcreate_ KSPCREATE
31: #define kspsetoptionsprefix_ KSPSETOPTIONSPREFIX
32: #define kspappendoptionsprefix_ KSPAPPENDOPTIONSPREFIX
33: #define kspgettype_ KSPGETTYPE
34: #define kspgetpreconditionerside_ KSPGETPRECONDITIONERSIDE
35: #define kspbuildsolution_ KSPBUILDSOLUTION
36: #define kspsettype_ KSPSETTYPE
37: #define kspgetresidualhistory_ KSPGETRESIDUALHISTORY
38: #define kspgetoptionsprefix_ KSPGETOPTIONSPREFIX
39: #define kspview_ KSPVIEW
40: #define kspgmressetrestart_ KSPGMRESSETRESTART
41: #define kspsetnormtype_ KSPSETNORMTYPE
42: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
43: #define kspgetconvergedreason_ kspgetconvergedreason
44: #define kspfgmressetmodifypc_ kspfgmressetmodifypc
45: #define kspfgmresmodifypcsles_ kspfgmresmodifypcsles
46: #define kspfgmresmodifypcnochange_ kspfgmresmodifypcnochange
47: #define kspdefaultconverged_ kspdefaultconverged
48: #define kspskipconverged_ kspskipconverged
49: #define kspsingularvaluemonitor_ kspsingularvaluemonitor
50: #define kspgmreskrylovmonitor_ kspgmreskrylovmonitor
51: #define kspdefaultmonitor_ kspdefaultmonitor
52: #define ksptruemonitor_ ksptruemonitor
53: #define kspvecviewmonitor_ kspvecviewmonitor
54: #define ksplgmonitor_ ksplgmonitor
55: #define ksplgtruemonitor_ ksplgtruemonitor
56: #define kspgetresidualhistory_ kspgetresidualhistory
57: #define kspsettype_ kspsettype
58: #define kspregisterdestroy_ kspregisterdestroy
59: #define kspdestroy_ kspdestroy
60: #define ksplgmonitordestroy_ ksplgmonitordestroy
61: #define ksplgmonitorcreate_ ksplgmonitorcreate
62: #define kspgetrhs_ kspgetrhs
63: #define kspgetsolution_ kspgetsolution
64: #define kspgetpc_ kspgetpc
65: #define kspsetmonitor_ kspsetmonitor
66: #define kspsetconvergencetest_ kspsetconvergencetest
67: #define kspcreate_ kspcreate
68: #define kspsetoptionsprefix_ kspsetoptionsprefix
69: #define kspappendoptionsprefix_ kspappendoptionsprefix
70: #define kspgettype_ kspgettype
71: #define kspgetpreconditionerside_ kspgetpreconditionerside
72: #define kspbuildsolution_ kspbuildsolution
73: #define kspgetoptionsprefix_ kspgetoptionsprefix
74: #define kspview_ kspview
75: #define kspgetresidualnorm_ kspgetresidualnorm
76: #define kspgmressetrestart_ kspgmressetrestart
77: #define kspsetnormtype_ kspsetnormtype
78: #endif
80: EXTERN_C_BEGIN
82: void PETSC_STDCALL kspgmressetrestart_(KSP *ksp,int *max_k, int *ierr )
83: {
84: *KSPGMRESSetRestart(*ksp,*max_k);
85: }
87: void PETSC_STDCALL kspgetresidualnorm_(KSP *ksp,PetscReal *rnorm,int *ierr)
88: {
89: *KSPGetResidualNorm(*ksp,rnorm);
90: }
92: void PETSC_STDCALL kspgetconvergedreason_(KSP *ksp,KSPConvergedReason *reason,int *ierr)
93: {
94: *KSPGetConvergedReason(*ksp,reason);
95: }
97: /* function */
98: void PETSC_STDCALL kspview_(KSP *ksp,PetscViewer *viewer, int *ierr)
99: {
100: PetscViewer v;
101: PetscPatchDefaultViewers_Fortran(viewer,v);
102: *KSPView(*ksp,v);
103: }
105: void kspdefaultconverged_(KSP *ksp,int *n,PetscReal *rnorm,KSPConvergedReason *flag,void *dummy,int *ierr)
106: {
107: CHKFORTRANNULLOBJECT(dummy);
108: *KSPDefaultConverged(*ksp,*n,*rnorm,flag,dummy);
109: }
111: void kspskipconverged_(KSP *ksp,int *n,PetscReal *rnorm,KSPConvergedReason *flag,void *dummy,int *ierr)
112: {
113: CHKFORTRANNULLOBJECT(dummy);
114: *KSPSkipConverged(*ksp,*n,*rnorm,flag,dummy);
115: }
117: void PETSC_STDCALL kspgetresidualhistory_(KSP *ksp,int *na,int *ierr)
118: {
119: *KSPGetResidualHistory(*ksp,PETSC_NULL,na);
120: }
122: void PETSC_STDCALL kspsettype_(KSP *ksp,CHAR type PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
123: {
124: char *t;
126: FIXCHAR(type,len,t);
127: *KSPSetType(*ksp,t);
128: FREECHAR(type,t);
129: }
131: void PETSC_STDCALL kspgettype_(KSP *ksp,CHAR name PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
132: {
133: char *tname;
135: *KSPGetType(*ksp,&tname);if (*ierr) return;
136: #if defined(PETSC_USES_CPTOFCD)
137: {
138: char *t = _fcdtocp(name); int len1 = _fcdlen(name);
139: *PetscStrncpy(t,tname,len1);
140: }
141: #else
142: *PetscStrncpy(name,tname,len);
143: #endif
144: }
146: void PETSC_STDCALL kspgetpreconditionerside_(KSP *ksp,PCSide *side,int *ierr){
147: *KSPGetPreconditionerSide(*ksp,side);
148: }
150: void PETSC_STDCALL kspsetoptionsprefix_(KSP *ksp,CHAR prefix PETSC_MIXED_LEN(len),
151: int *ierr PETSC_END_LEN(len))
152: {
153: char *t;
155: FIXCHAR(prefix,len,t);
156: *KSPSetOptionsPrefix(*ksp,t);
157: FREECHAR(prefix,t);
158: }
160: void PETSC_STDCALL kspappendoptionsprefix_(KSP *ksp,CHAR prefix PETSC_MIXED_LEN(len),
161: int *ierr PETSC_END_LEN(len))
162: {
163: char *t;
165: FIXCHAR(prefix,len,t);
166: *KSPAppendOptionsPrefix(*ksp,t);
167: FREECHAR(prefix,t);
168: }
170: void PETSC_STDCALL kspcreate_(MPI_Comm *comm,KSP *ksp,int *ierr){
171: *KSPCreate((MPI_Comm)PetscToPointerComm(*comm),ksp);
172: }
174: static void (PETSC_STDCALL *f2)(KSP*,int*,PetscReal*,KSPConvergedReason*,void*,int*);
175: static int ourtest(KSP ksp,int i,PetscReal d,KSPConvergedReason *reason,void* ctx)
176: {
178: (*f2)(&ksp,&i,&d,reason,ctx,&ierr);
179: return 0;
180: }
181: void PETSC_STDCALL kspsetconvergencetest_(KSP *ksp,
182: void (PETSC_STDCALL *converge)(KSP*,int*,PetscReal*,KSPConvergedReason*,void*,int*),void *cctx,int *ierr)
183: {
184: if ((void(*)(void))converge == (void(*)(void))kspdefaultconverged_) {
185: *KSPSetConvergenceTest(*ksp,KSPDefaultConverged,0);
186: } else if ((void(*)(void))converge == (void(*)(void))kspskipconverged_) {
187: *KSPSetConvergenceTest(*ksp,KSPSkipConverged,0);
188: } else {
189: f2 = converge;
190: *KSPSetConvergenceTest(*ksp,ourtest,cctx);
191: }
192: }
194: /*
195: These are not usually called from Fortran but allow Fortran users
196: to transparently set these monitors from .F code
197:
198: functions, hence no STDCALL
199: */
200: void kspgmreskrylovmonitor_(KSP *ksp,int *it,PetscReal *norm,void *ctx,int *ierr)
201: {
202: *KSPGMRESKrylovMonitor(*ksp,*it,*norm,ctx);
203: }
205: void kspdefaultmonitor_(KSP *ksp,int *it,PetscReal *norm,void *ctx,int *ierr)
206: {
207: *KSPDefaultMonitor(*ksp,*it,*norm,ctx);
208: }
209:
210: void kspsingularvaluemonitor_(KSP *ksp,int *it,PetscReal *norm,void *ctx,int *ierr)
211: {
212: *KSPSingularValueMonitor(*ksp,*it,*norm,ctx);
213: }
215: void ksplgmonitor_(KSP *ksp,int *it,PetscReal *norm,void *ctx,int *ierr)
216: {
217: *KSPLGMonitor(*ksp,*it,*norm,ctx);
218: }
220: void ksplgtruemonitor_(KSP *ksp,int *it,PetscReal *norm,void *ctx,int *ierr)
221: {
222: *KSPLGTrueMonitor(*ksp,*it,*norm,ctx);
223: }
225: void ksptruemonitor_(KSP *ksp,int *it,PetscReal *norm,void *ctx,int *ierr)
226: {
227: *KSPTrueMonitor(*ksp,*it,*norm,ctx);
228: }
230: void kspvecviewmonitor_(KSP *ksp,int *it,PetscReal *norm,void *ctx,int *ierr)
231: {
232: *KSPVecViewMonitor(*ksp,*it,*norm,ctx);
233: }
235: static void (PETSC_STDCALL *f1)(KSP*,int*,PetscReal*,void*,int*);
236: static int ourmonitor(KSP ksp,int i,PetscReal d,void* ctx)
237: {
238: int 0;
239: (*f1)(&ksp,&i,&d,ctx,&ierr);
240: return 0;
241: }
242: static void (PETSC_STDCALL *f21)(void*,int*);
243: static int ourdestroy(void* ctx)
244: {
245: int 0;
246: (*f21)(ctx,&ierr);
247: return 0;
248: }
250: void PETSC_STDCALL kspsetmonitor_(KSP *ksp,void (PETSC_STDCALL *monitor)(KSP*,int*,PetscReal*,void*,int*),
251: void *mctx,void (PETSC_STDCALL *monitordestroy)(void *,int *),int *ierr)
252: {
253: if ((void(*)(void))monitor == (void(*)(void))kspdefaultmonitor_) {
254: *KSPSetMonitor(*ksp,KSPDefaultMonitor,0,0);
255: } else if ((void(*)(void))monitor == (void(*)(void))ksplgmonitor_) {
256: *KSPSetMonitor(*ksp,KSPLGMonitor,0,0);
257: } else if ((void(*)(void))monitor == (void(*)(void))ksplgtruemonitor_) {
258: *KSPSetMonitor(*ksp,KSPLGTrueMonitor,0,0);
259: } else if ((void(*)(void))monitor == (void(*)(void))kspvecviewmonitor_) {
260: *KSPSetMonitor(*ksp,KSPVecViewMonitor,0,0);
261: } else if ((void(*)(void))monitor == (void(*)(void))ksptruemonitor_) {
262: *KSPSetMonitor(*ksp,KSPTrueMonitor,0,0);
263: } else if ((void(*)(void))monitor == (void(*)(void))kspsingularvaluemonitor_) {
264: *KSPSetMonitor(*ksp,KSPSingularValueMonitor,0,0);
265: } else {
266: f1 = monitor;
267: if (FORTRANNULLFUNCTION(monitordestroy)) {
268: *KSPSetMonitor(*ksp,ourmonitor,mctx,0);
269: } else {
270: f21 = monitordestroy;
271: *KSPSetMonitor(*ksp,ourmonitor,mctx,ourdestroy);
272: }
273: }
274: }
276: void PETSC_STDCALL kspgetpc_(KSP *ksp,PC *B,int *ierr)
277: {
278: *KSPGetPC(*ksp,B);
279: }
281: void PETSC_STDCALL kspgetsolution_(KSP *ksp,Vec *v,int *ierr)
282: {
283: *KSPGetSolution(*ksp,v);
284: }
286: void PETSC_STDCALL kspgetrhs_(KSP *ksp,Vec *r,int *ierr)
287: {
288: *KSPGetRhs(*ksp,r);
289: }
291: /*
292: Possible bleeds memory but cannot be helped.
293: */
294: void PETSC_STDCALL ksplgmonitorcreate_(CHAR host PETSC_MIXED_LEN(len1),
295: CHAR label PETSC_MIXED_LEN(len2),int *x,int *y,int *m,int *n,PetscDrawLG *ctx,
296: int *ierr PETSC_END_LEN(len1) PETSC_END_LEN(len2))
297: {
298: char *t1,*t2;
300: FIXCHAR(host,len1,t1);
301: FIXCHAR(label,len2,t2);
302: *KSPLGMonitorCreate(t1,t2,*x,*y,*m,*n,ctx);
303: }
305: void PETSC_STDCALL ksplgmonitordestroy_(PetscDrawLG *ctx,int *ierr)
306: {
307: *KSPLGMonitorDestroy(*ctx);
308: }
310: void PETSC_STDCALL kspdestroy_(KSP *ksp,int *ierr)
311: {
312: *KSPDestroy(*ksp);
313: }
315: void PETSC_STDCALL kspregisterdestroy_(int* ierr)
316: {
317: *KSPRegisterDestroy();
318: }
320: void PETSC_STDCALL kspbuildsolution_(KSP *ctx,Vec *v,Vec *V,int *ierr)
321: {
322: *KSPBuildSolution(*ctx,*v,V);
323: }
325: void PETSC_STDCALL kspbuildresidual_(KSP *ctx,Vec *t,Vec *v,Vec *V,int *ierr)
326: {
327: *KSPBuildResidual(*ctx,*t,*v,V);
328: }
330: void PETSC_STDCALL kspgetoptionsprefix_(KSP *ksp,CHAR prefix PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
331: {
332: char *tname;
334: *KSPGetOptionsPrefix(*ksp,&tname);
335: #if defined(PETSC_USES_CPTOFCD)
336: {
337: char *t = _fcdtocp(prefix); int len1 = _fcdlen(prefix);
338: *PetscStrncpy(t,tname,len1); if (*ierr) return;
339: }
340: #else
341: *PetscStrncpy(prefix,tname,len); if (*ierr) return;
342: #endif
343: }
345: static void (PETSC_STDCALL *f109)(KSP*,int*,int*,PetscReal*,void*,int*);
346: static int ourmodify(KSP ksp,int i,int i2,PetscReal d,void* ctx)
347: {
348: int 0;
349: (*f109)(&ksp,&i,&i2,&d,ctx,&ierr);
350: return 0;
351: }
353: static void (PETSC_STDCALL *f210)(void*,int*);
354: static int ourmoddestroy(void* ctx)
355: {
356: int 0;
357: (*f210)(ctx,&ierr);
358: return 0;
359: }
361: void PETSC_STDCALL kspfgmresmodifypcnochange_(KSP *ksp,int *total_its,int *loc_its,PetscReal *res_norm,void* dummy,int *ierr)
362: {
363: *KSPFGMRESModifyPCNoChange(*ksp,*total_its,*loc_its,*res_norm,dummy);
364: }
366: void PETSC_STDCALL kspfgmresmodifypcsles_(KSP *ksp,int *total_its,int *loc_its,PetscReal *res_norm,void*dummy,int *ierr)
367: {
368: *KSPFGMRESModifyPCSLES(*ksp,*total_its,*loc_its,*res_norm,dummy);
369: }
371: 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)
372: {
373: if ((void(*)(void))fcn == (void(*)(void))kspfgmresmodifypcsles_) {
374: *KSPFGMRESSetModifyPC(*ksp,KSPFGMRESModifyPCSLES,0,0);
375: } else if ((void(*)(void))fcn == (void(*)(void))kspfgmresmodifypcnochange_) {
376: *KSPFGMRESSetModifyPC(*ksp,KSPFGMRESModifyPCNoChange,0,0);
377: } else {
378: f109 = fcn;
379: if (FORTRANNULLFUNCTION(d)) {
380: *KSPFGMRESSetModifyPC(*ksp,ourmodify,ctx,0);
381: } else {
382: f210 = d;
383: *KSPFGMRESSetModifyPC(*ksp,ourmodify,ctx,ourmoddestroy);
384: }
385: }
386: }
388: void PETSC_STDCALL kspsetnormtype_(KSP *ksp,KSPNormType *type,int *ierr)
389: {
390: *KSPSetNormType(*ksp,*type);
391: }
393: EXTERN_C_END