Actual source code: zsvdf.c

slepc-3.22.2 2024-12-02
Report Typos and Errors
  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/fortranimpl.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 svdmonitorconvergeddestroy_       SVDMONITORCONVERGEDDESTROY
 22: #define svdconvergedabsolute_             SVDCONVERGEDABSOLUTE
 23: #define svdconvergedrelative_             SVDCONVERGEDRELATIVE
 24: #define svdconvergednorm_                 SVDCONVERGEDNORM
 25: #define svdconvergedmaxit_                SVDCONVERGEDMAXIT
 26: #define svdsetconvergencetestfunction_    SVDSETCONVERGENCETESTFUNCTION
 27: #define svdsetstoppingtestfunction_       SVDSETSTOPPINGTESTFUNCTION
 28: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 29: #define svdmonitorset_                    svdmonitorset
 30: #define svdmonitorall_                    svdmonitorall
 31: #define svdmonitorfirst_                  svdmonitorfirst
 32: #define svdmonitorconditioning_           svdmonitorconditioning
 33: #define svdmonitorconverged_              svdmonitorconverged
 34: #define svdmonitorconvergedcreate_        svdmonitorconvergedcreate
 35: #define svdmonitorconvergeddestroy_       svdmonitorconvergeddestroy
 36: #define svdconvergedabsolute_             svdconvergedabsolute
 37: #define svdconvergedrelative_             svdconvergedrelative
 38: #define svdconvergednorm_                 svdconvergednorm
 39: #define svdconvergedmaxit_                svdconvergedmaxit
 40: #define svdsetconvergencetestfunction_    svdsetconvergencetestfunction
 41: #define svdsetstoppingtestfunction_       svdsetstoppingtestfunction
 42: #endif

 44: /*
 45:    These are not usually called from Fortran but allow Fortran users
 46:    to transparently set these monitors from .F code
 47: */
 48: SLEPC_EXTERN void svdmonitorall_(SVD *svd,PetscInt *it,PetscInt *nconv,PetscReal *sigma,PetscReal *errest,PetscInt *nest,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
 49: {
 50:   *ierr = SVDMonitorAll(*svd,*it,*nconv,sigma,errest,*nest,*vf);
 51: }

 53: SLEPC_EXTERN void svdmonitorfirst_(SVD *svd,PetscInt *it,PetscInt *nconv,PetscReal *sigma,PetscReal *errest,PetscInt *nest,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
 54: {
 55:   *ierr = SVDMonitorFirst(*svd,*it,*nconv,sigma,errest,*nest,*vf);
 56: }

 58: SLEPC_EXTERN void svdmonitorconditioning_(SVD *svd,PetscInt *it,PetscInt *nconv,PetscReal *sigma,PetscReal *errest,PetscInt *nest,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
 59: {
 60:   *ierr = SVDMonitorConditioning(*svd,*it,*nconv,sigma,errest,*nest,*vf);
 61: }

 63: SLEPC_EXTERN void svdmonitorconverged_(SVD *svd,PetscInt *it,PetscInt *nconv,PetscReal *sigma,PetscReal *errest,PetscInt *nest,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
 64: {
 65:   *ierr = SVDMonitorConverged(*svd,*it,*nconv,sigma,errest,*nest,*vf);
 66: }

 68: SLEPC_EXTERN void svdmonitorconvergedcreate_(PetscViewer *vin,PetscViewerFormat *format,void *ctx,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
 69: {
 70:   PetscViewer v;
 71:   PetscPatchDefaultViewers_Fortran(vin,v);
 72:   CHKFORTRANNULLOBJECT(ctx);
 73:   *ierr = SVDMonitorConvergedCreate(v,*format,ctx,vf);
 74: }

 76: SLEPC_EXTERN void svdmonitorconvergeddestroy_(PetscViewerAndFormat **vf,PetscErrorCode *ierr)
 77: {
 78:   *ierr = SVDMonitorConvergedDestroy(vf);
 79: }

 81: static struct {
 82:   PetscFortranCallbackId monitor;
 83:   PetscFortranCallbackId monitordestroy;
 84:   PetscFortranCallbackId convergence;
 85:   PetscFortranCallbackId convdestroy;
 86:   PetscFortranCallbackId stopping;
 87:   PetscFortranCallbackId stopdestroy;
 88: } _cb;

 90: /* These are not extern C because they are passed into non-extern C user level functions */
 91: static PetscErrorCode ourmonitor(SVD svd,PetscInt i,PetscInt nc,PetscReal *sigma,PetscReal *d,PetscInt l,void* ctx)
 92: {
 93:   PetscObjectUseFortranCallback(svd,_cb.monitor,(SVD*,PetscInt*,PetscInt*,PetscReal*,PetscReal*,PetscInt*,void*,PetscErrorCode*),(&svd,&i,&nc,sigma,d,&l,_ctx,&ierr));
 94: }

 96: static PetscErrorCode ourdestroy(void** ctx)
 97: {
 98:   SVD svd = (SVD)*ctx;
 99:   PetscObjectUseFortranCallback(svd,_cb.monitordestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
100: }

102: static PetscErrorCode ourconvergence(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
103: {
104:   PetscObjectUseFortranCallback(svd,_cb.convergence,(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*),(&svd,&sigma,&res,errest,_ctx,&ierr));
105: }

107: static PetscErrorCode ourconvdestroy(void *ctx)
108: {
109:   SVD svd = (SVD)ctx;
110:   PetscObjectUseFortranCallback(svd,_cb.convdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
111: }

113: static PetscErrorCode ourstopping(SVD svd,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nsv,SVDConvergedReason *reason,void *ctx)
114: {
115:   PetscObjectUseFortranCallback(svd,_cb.stopping,(SVD*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,SVDConvergedReason*,void*,PetscErrorCode*),(&svd,&its,&max_it,&nconv,&nsv,reason,_ctx,&ierr));
116: }

118: static PetscErrorCode ourstopdestroy(void *ctx)
119: {
120:   SVD svd = (SVD)ctx;
121:   PetscObjectUseFortranCallback(svd,_cb.stopdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
122: }

124: SLEPC_EXTERN void svdmonitorset_(SVD *svd,void (*monitor)(SVD*,PetscInt*,PetscInt*,PetscReal*,PetscReal*,PetscInt*,void*,PetscErrorCode*),void *mctx,void (*monitordestroy)(void *,PetscErrorCode*),PetscErrorCode *ierr)
125: {
126:   CHKFORTRANNULLOBJECT(mctx);
127:   CHKFORTRANNULLFUNCTION(monitordestroy);
128:   if ((PetscVoidFunction)monitor == (PetscVoidFunction)svdmonitorall_) {
129:     *ierr = SVDMonitorSet(*svd,(PetscErrorCode (*)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*))SVDMonitorAll,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
130:   } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)svdmonitorconverged_) {
131:     *ierr = SVDMonitorSet(*svd,(PetscErrorCode (*)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*))SVDMonitorConverged,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void**))SVDMonitorConvergedDestroy);
132:   } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)svdmonitorfirst_) {
133:     *ierr = SVDMonitorSet(*svd,(PetscErrorCode (*)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*))SVDMonitorFirst,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
134:   } else {
135:     *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)monitor,mctx); if (*ierr) return;
136:     *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitordestroy,(PetscVoidFunction)monitordestroy,mctx); if (*ierr) return;
137:     *ierr = SVDMonitorSet(*svd,ourmonitor,*svd,ourdestroy);
138:   }
139: }

141: SLEPC_EXTERN void svdconvergedabsolute_(SVD *svd,PetscReal *sigma,PetscReal *res,PetscReal *errest,void *ctx,PetscErrorCode *ierr)
142: {
143:   *ierr = SVDConvergedAbsolute(*svd,*sigma,*res,errest,ctx);
144: }

146: SLEPC_EXTERN void svdconvergedrelative_(SVD *svd,PetscReal *sigma,PetscReal *res,PetscReal *errest,void *ctx,PetscErrorCode *ierr)
147: {
148:   *ierr = SVDConvergedRelative(*svd,*sigma,*res,errest,ctx);
149: }

151: SLEPC_EXTERN void svdconvergednorm_(SVD *svd,PetscReal *sigma,PetscReal *res,PetscReal *errest,void *ctx,PetscErrorCode *ierr)
152: {
153:   *ierr = SVDConvergedNorm(*svd,*sigma,*res,errest,ctx);
154: }

156: SLEPC_EXTERN void svdconvergedmaxit_(SVD *svd,PetscReal *sigma,PetscReal *res,PetscReal *errest,void *ctx,PetscErrorCode *ierr)
157: {
158:   *ierr = SVDConvergedMaxIt(*svd,*sigma,*res,errest,ctx);
159: }

161: SLEPC_EXTERN void svdsetconvergencetestfunction_(SVD *svd,void (*func)(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
162: {
163:   CHKFORTRANNULLOBJECT(ctx);
164:   CHKFORTRANNULLFUNCTION(destroy);
165:   if ((PetscVoidFunction)func == (PetscVoidFunction)svdconvergedabsolute_) {
166:     *ierr = SVDSetConvergenceTest(*svd,SVD_CONV_ABS);
167:   } else if ((PetscVoidFunction)func == (PetscVoidFunction)svdconvergedrelative_) {
168:     *ierr = SVDSetConvergenceTest(*svd,SVD_CONV_REL);
169:   } else if ((PetscVoidFunction)func == (PetscVoidFunction)svdconvergednorm_) {
170:     *ierr = SVDSetConvergenceTest(*svd,SVD_CONV_NORM);
171:   } else if ((PetscVoidFunction)func == (PetscVoidFunction)svdconvergedmaxit_) {
172:     *ierr = SVDSetConvergenceTest(*svd,SVD_CONV_MAXIT);
173:   } else {
174:     *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convergence,(PetscVoidFunction)func,ctx); if (*ierr) return;
175:     *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
176:     *ierr = SVDSetConvergenceTestFunction(*svd,ourconvergence,*svd,ourconvdestroy);
177:   }
178: }

180: SLEPC_EXTERN void svdstoppingbasic_(SVD *svd,PetscInt *its,PetscInt *max_it,PetscInt *nconv,PetscInt *nsv,SVDConvergedReason *reason,void *ctx,PetscErrorCode *ierr)
181: {
182:   *ierr = SVDStoppingBasic(*svd,*its,*max_it,*nconv,*nsv,reason,ctx);
183: }

185: SLEPC_EXTERN void svdsetstoppingtestfunction_(SVD *svd,void (*func)(SVD*,PetscInt,PetscInt,PetscInt,PetscInt,SVDConvergedReason*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
186: {
187:   CHKFORTRANNULLOBJECT(ctx);
188:   CHKFORTRANNULLFUNCTION(destroy);
189:   if ((PetscVoidFunction)func == (PetscVoidFunction)svdstoppingbasic_) {
190:     *ierr = SVDSetStoppingTest(*svd,SVD_STOP_BASIC);
191:   } else {
192:     *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopping,(PetscVoidFunction)func,ctx); if (*ierr) return;
193:     *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
194:     *ierr = SVDSetStoppingTestFunction(*svd,ourstopping,*svd,ourstopdestroy);
195:   }
196: }