Actual source code: zsvdf.c
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
11: #include <petsc/private/ftnimpl.h>
12: #include <slepcsvd.h>
14: #if defined(PETSC_HAVE_FORTRAN_CAPS)
15: #define svdmonitorset_ SVDMONITORSET
16: #define svdmonitorall_ SVDMONITORALL
17: #define svdmonitorfirst_ SVDMONITORFIRST
18: #define svdmonitorconditioning_ SVDMONITORCONDITIONING
19: #define svdmonitorconverged_ SVDMONITORCONVERGED
20: #define svdmonitorconvergedcreate_ SVDMONITORCONVERGEDCREATE
21: #define svdconvergedabsolute_ SVDCONVERGEDABSOLUTE
22: #define svdconvergedrelative_ SVDCONVERGEDRELATIVE
23: #define svdconvergednorm_ SVDCONVERGEDNORM
24: #define svdconvergedmaxit_ SVDCONVERGEDMAXIT
25: #define svdsetconvergencetestfunction_ SVDSETCONVERGENCETESTFUNCTION
26: #define svdstoppingbasic_ SVDSTOPPINGBASIC
27: #define svdstoppingthreshold_ SVDSTOPPINGTHRESHOLD
28: #define svdsetstoppingtestfunction_ SVDSETSTOPPINGTESTFUNCTION
29: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
30: #define svdmonitorset_ svdmonitorset
31: #define svdmonitorall_ svdmonitorall
32: #define svdmonitorfirst_ svdmonitorfirst
33: #define svdmonitorconditioning_ svdmonitorconditioning
34: #define svdmonitorconverged_ svdmonitorconverged
35: #define svdmonitorconvergedcreate_ svdmonitorconvergedcreate
36: #define svdconvergedabsolute_ svdconvergedabsolute
37: #define svdconvergedrelative_ svdconvergedrelative
38: #define svdconvergednorm_ svdconvergednorm
39: #define svdconvergedmaxit_ svdconvergedmaxit
40: #define svdsetconvergencetestfunction_ svdsetconvergencetestfunction
41: #define svdstoppingbasic_ svdstoppingbasic
42: #define svdstoppingthreshold_ svdstoppingthreshold
43: #define svdsetstoppingtestfunction_ svdsetstoppingtestfunction
44: #endif
46: /*
47: These cannot be called from Fortran but allow Fortran users
48: to transparently set these monitors from .F code
49: */
50: SLEPC_EXTERN void svdmonitorall_(SVD*,PetscInt*,PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
51: SLEPC_EXTERN void svdmonitorfirst_(SVD*,PetscInt*,PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
52: SLEPC_EXTERN void svdmonitorconditioning_(SVD*,PetscInt*,PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
53: SLEPC_EXTERN void svdmonitorconverged_(SVD*,PetscInt*,PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
55: SLEPC_EXTERN void svdmonitorconvergedcreate_(PetscViewer *vin,PetscViewerFormat *format,void *ctx,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
56: {
57: PetscViewer v;
58: PetscPatchDefaultViewers_Fortran(vin,v);
59: CHKFORTRANNULLOBJECT(ctx);
60: *ierr = SVDMonitorConvergedCreate(v,*format,ctx,vf);
61: }
63: static struct {
64: PetscFortranCallbackId monitor;
65: PetscFortranCallbackId monitordestroy;
66: PetscFortranCallbackId convergence;
67: PetscFortranCallbackId convdestroy;
68: PetscFortranCallbackId stopping;
69: PetscFortranCallbackId stopdestroy;
70: } _cb;
72: /* These are not extern C because they are passed into non-extern C user level functions */
73: static PetscErrorCode ourmonitor(SVD svd,PetscInt i,PetscInt nc,PetscReal *sigma,PetscReal *d,PetscInt l,void *ctx)
74: {
75: PetscObjectUseFortranCallback(svd,_cb.monitor,(SVD*,PetscInt*,PetscInt*,PetscReal*,PetscReal*,PetscInt*,void*,PetscErrorCode*),(&svd,&i,&nc,sigma,d,&l,_ctx,&ierr));
76: }
78: static PetscErrorCode ourdestroy(PetscCtxRt ctx)
79: {
80: SVD svd = *(SVD*)ctx;
81: PetscObjectUseFortranCallback(svd,_cb.monitordestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
82: }
84: static PetscErrorCode ourconvergence(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
85: {
86: PetscObjectUseFortranCallback(svd,_cb.convergence,(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*),(&svd,&sigma,&res,errest,_ctx,&ierr));
87: }
89: static PetscErrorCode ourconvdestroy(PetscCtxRt ctx)
90: {
91: SVD svd = *(SVD*)ctx;
92: PetscObjectUseFortranCallback(svd,_cb.convdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
93: }
95: static PetscErrorCode ourstopping(SVD svd,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nsv,SVDConvergedReason *reason,void *ctx)
96: {
97: PetscObjectUseFortranCallback(svd,_cb.stopping,(SVD*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,SVDConvergedReason*,void*,PetscErrorCode*),(&svd,&its,&max_it,&nconv,&nsv,reason,_ctx,&ierr));
98: }
100: static PetscErrorCode ourstopdestroy(PetscCtxRt ctx)
101: {
102: SVD svd = *(SVD*)ctx;
103: PetscObjectUseFortranCallback(svd,_cb.stopdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
104: }
106: SLEPC_EXTERN void svdmonitorset_(SVD *svd,void (*monitor)(SVD*,PetscInt*,PetscInt*,PetscReal*,PetscReal*,PetscInt*,void*,PetscErrorCode*),void *mctx,void (*monitordestroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
107: {
108: CHKFORTRANNULLOBJECT(mctx);
109: CHKFORTRANNULLFUNCTION(monitordestroy);
110: if ((PetscFortranCallbackFn*)monitor == (PetscFortranCallbackFn*)svdmonitorall_) {
111: *ierr = SVDMonitorSet(*svd,(SVDMonitorFn*)SVDMonitorAll,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)PetscViewerAndFormatDestroy);
112: } else if ((PetscFortranCallbackFn*)monitor == (PetscFortranCallbackFn*)svdmonitorconverged_) {
113: *ierr = SVDMonitorSet(*svd,(SVDMonitorFn*)SVDMonitorConverged,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)PetscViewerAndFormatDestroy);
114: } else if ((PetscFortranCallbackFn*)monitor == (PetscFortranCallbackFn*)svdmonitorfirst_) {
115: *ierr = SVDMonitorSet(*svd,(SVDMonitorFn*)SVDMonitorFirst,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)PetscViewerAndFormatDestroy);
116: } else {
117: *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscFortranCallbackFn*)monitor,mctx); if (*ierr) return;
118: *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitordestroy,(PetscFortranCallbackFn*)monitordestroy,mctx); if (*ierr) return;
119: *ierr = SVDMonitorSet(*svd,ourmonitor,*svd,ourdestroy);
120: }
121: }
123: SLEPC_EXTERN void svdconvergedabsolute_(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*);
124: SLEPC_EXTERN void svdconvergedrelative_(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*);
125: SLEPC_EXTERN void svdconvergednorm_(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*);
126: SLEPC_EXTERN void svdconvergedmaxit_(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*);
128: SLEPC_EXTERN void svdsetconvergencetestfunction_(SVD *svd,void (*func)(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*),void *ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
129: {
130: CHKFORTRANNULLOBJECT(ctx);
131: CHKFORTRANNULLFUNCTION(destroy);
132: if (func == svdconvergedabsolute_) {
133: *ierr = SVDSetConvergenceTest(*svd,SVD_CONV_ABS);
134: } else if (func == svdconvergedrelative_) {
135: *ierr = SVDSetConvergenceTest(*svd,SVD_CONV_REL);
136: } else if (func == svdconvergednorm_) {
137: *ierr = SVDSetConvergenceTest(*svd,SVD_CONV_NORM);
138: } else if (func == svdconvergedmaxit_) {
139: *ierr = SVDSetConvergenceTest(*svd,SVD_CONV_MAXIT);
140: } else {
141: *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convergence,(PetscFortranCallbackFn*)func,ctx); if (*ierr) return;
142: *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convdestroy,(PetscFortranCallbackFn*)destroy,ctx); if (*ierr) return;
143: *ierr = SVDSetConvergenceTestFunction(*svd,ourconvergence,*svd,ourconvdestroy);
144: }
145: }
147: SLEPC_EXTERN void svdstoppingbasic_(SVD*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,SVDConvergedReason*,void*,PetscErrorCode*);
148: SLEPC_EXTERN void svdstoppingthreshold_(SVD*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,SVDConvergedReason*,void*,PetscErrorCode*);
150: SLEPC_EXTERN void svdsetstoppingtestfunction_(SVD *svd,void (*func)(SVD*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,SVDConvergedReason*,void*,PetscErrorCode*),void *ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
151: {
152: CHKFORTRANNULLOBJECT(ctx);
153: CHKFORTRANNULLFUNCTION(destroy);
154: if (func == svdstoppingbasic_) {
155: *ierr = SVDSetStoppingTest(*svd,SVD_STOP_BASIC);
156: } else if (func == svdstoppingthreshold_) {
157: *ierr = SVDSetStoppingTest(*svd,SVD_STOP_THRESHOLD);
158: } else {
159: *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopping,(PetscFortranCallbackFn*)func,ctx); if (*ierr) return;
160: *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopdestroy,(PetscFortranCallbackFn*)destroy,ctx); if (*ierr) return;
161: *ierr = SVDSetStoppingTestFunction(*svd,ourstopping,*svd,ourstopdestroy);
162: }
163: }