Actual source code: zsvdf.c

slepc-3.23.1 2025-05-01
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/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 svdmonitorconvergeddestroy_       SVDMONITORCONVERGEDDESTROY
 22: #define svdconvergedabsolute_             SVDCONVERGEDABSOLUTE
 23: #define svdconvergedrelative_             SVDCONVERGEDRELATIVE
 24: #define svdconvergednorm_                 SVDCONVERGEDNORM
 25: #define svdconvergedmaxit_                SVDCONVERGEDMAXIT
 26: #define svdsetconvergencetestfunction_    SVDSETCONVERGENCETESTFUNCTION
 27: #define svdstoppingbasic_                 SVDSTOPPINGBASIC
 28: #define svdstoppingthreshold_             SVDSTOPPINGTHRESHOLD
 29: #define svdsetstoppingtestfunction_       SVDSETSTOPPINGTESTFUNCTION
 30: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 31: #define svdmonitorset_                    svdmonitorset
 32: #define svdmonitorall_                    svdmonitorall
 33: #define svdmonitorfirst_                  svdmonitorfirst
 34: #define svdmonitorconditioning_           svdmonitorconditioning
 35: #define svdmonitorconverged_              svdmonitorconverged
 36: #define svdmonitorconvergedcreate_        svdmonitorconvergedcreate
 37: #define svdmonitorconvergeddestroy_       svdmonitorconvergeddestroy
 38: #define svdconvergedabsolute_             svdconvergedabsolute
 39: #define svdconvergedrelative_             svdconvergedrelative
 40: #define svdconvergednorm_                 svdconvergednorm
 41: #define svdconvergedmaxit_                svdconvergedmaxit
 42: #define svdsetconvergencetestfunction_    svdsetconvergencetestfunction
 43: #define svdstoppingbasic_                 svdstoppingbasic
 44: #define svdstoppingthreshold_             svdstoppingthreshold
 45: #define svdsetstoppingtestfunction_       svdsetstoppingtestfunction
 46: #endif

 48: /*
 49:    These cannot be called from Fortran but allow Fortran users
 50:    to transparently set these monitors from .F code
 51: */
 52: SLEPC_EXTERN void svdmonitorall_(SVD*,PetscInt*,PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
 53: SLEPC_EXTERN void svdmonitorfirst_(SVD*,PetscInt*,PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
 54: SLEPC_EXTERN void svdmonitorconditioning_(SVD*,PetscInt*,PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
 55: SLEPC_EXTERN void svdmonitorconverged_(SVD*,PetscInt*,PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);

 57: SLEPC_EXTERN void svdmonitorconvergedcreate_(PetscViewer *vin,PetscViewerFormat *format,void *ctx,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
 58: {
 59:   PetscViewer v;
 60:   PetscPatchDefaultViewers_Fortran(vin,v);
 61:   CHKFORTRANNULLOBJECT(ctx);
 62:   *ierr = SVDMonitorConvergedCreate(v,*format,ctx,vf);
 63: }

 65: SLEPC_EXTERN void svdmonitorconvergeddestroy_(PetscViewerAndFormat **vf,PetscErrorCode *ierr)
 66: {
 67:   *ierr = SVDMonitorConvergedDestroy(vf);
 68: }

 70: static struct {
 71:   PetscFortranCallbackId monitor;
 72:   PetscFortranCallbackId monitordestroy;
 73:   PetscFortranCallbackId convergence;
 74:   PetscFortranCallbackId convdestroy;
 75:   PetscFortranCallbackId stopping;
 76:   PetscFortranCallbackId stopdestroy;
 77: } _cb;

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

 85: static PetscErrorCode ourdestroy(void **ctx)
 86: {
 87:   SVD svd = (SVD)*ctx;
 88:   PetscObjectUseFortranCallback(svd,_cb.monitordestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
 89: }

 91: static PetscErrorCode ourconvergence(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
 92: {
 93:   PetscObjectUseFortranCallback(svd,_cb.convergence,(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*),(&svd,&sigma,&res,errest,_ctx,&ierr));
 94: }

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

102: static PetscErrorCode ourstopping(SVD svd,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nsv,SVDConvergedReason *reason,void *ctx)
103: {
104:   PetscObjectUseFortranCallback(svd,_cb.stopping,(SVD*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,SVDConvergedReason*,void*,PetscErrorCode*),(&svd,&its,&max_it,&nconv,&nsv,reason,_ctx,&ierr));
105: }

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

113: SLEPC_EXTERN void svdmonitorset_(SVD *svd,void (*monitor)(SVD*,PetscInt*,PetscInt*,PetscReal*,PetscReal*,PetscInt*,void*,PetscErrorCode*),void *mctx,void (*monitordestroy)(void *,PetscErrorCode*),PetscErrorCode *ierr)
114: {
115:   CHKFORTRANNULLOBJECT(mctx);
116:   CHKFORTRANNULLFUNCTION(monitordestroy);
117:   if ((PetscVoidFunction)monitor == (PetscVoidFunction)svdmonitorall_) {
118:     *ierr = SVDMonitorSet(*svd,(PetscErrorCode (*)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*))SVDMonitorAll,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)PetscViewerAndFormatDestroy);
119:   } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)svdmonitorconverged_) {
120:     *ierr = SVDMonitorSet(*svd,(PetscErrorCode (*)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*))SVDMonitorConverged,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)SVDMonitorConvergedDestroy);
121:   } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)svdmonitorfirst_) {
122:     *ierr = SVDMonitorSet(*svd,(PetscErrorCode (*)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*))SVDMonitorFirst,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)PetscViewerAndFormatDestroy);
123:   } else {
124:     *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)monitor,mctx); if (*ierr) return;
125:     *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitordestroy,(PetscVoidFunction)monitordestroy,mctx); if (*ierr) return;
126:     *ierr = SVDMonitorSet(*svd,ourmonitor,*svd,ourdestroy);
127:   }
128: }

130: SLEPC_EXTERN void svdconvergedabsolute_(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*);
131: SLEPC_EXTERN void svdconvergedrelative_(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*);
132: SLEPC_EXTERN void svdconvergednorm_(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*);
133: SLEPC_EXTERN void svdconvergedmaxit_(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*);

135: SLEPC_EXTERN void svdsetconvergencetestfunction_(SVD *svd,void (*func)(SVD*,PetscReal*,PetscReal*,PetscReal*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
136: {
137:   CHKFORTRANNULLOBJECT(ctx);
138:   CHKFORTRANNULLFUNCTION(destroy);
139:   if ((PetscVoidFunction)func == (PetscVoidFunction)svdconvergedabsolute_) {
140:     *ierr = SVDSetConvergenceTest(*svd,SVD_CONV_ABS);
141:   } else if ((PetscVoidFunction)func == (PetscVoidFunction)svdconvergedrelative_) {
142:     *ierr = SVDSetConvergenceTest(*svd,SVD_CONV_REL);
143:   } else if ((PetscVoidFunction)func == (PetscVoidFunction)svdconvergednorm_) {
144:     *ierr = SVDSetConvergenceTest(*svd,SVD_CONV_NORM);
145:   } else if ((PetscVoidFunction)func == (PetscVoidFunction)svdconvergedmaxit_) {
146:     *ierr = SVDSetConvergenceTest(*svd,SVD_CONV_MAXIT);
147:   } else {
148:     *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convergence,(PetscVoidFunction)func,ctx); if (*ierr) return;
149:     *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
150:     *ierr = SVDSetConvergenceTestFunction(*svd,ourconvergence,*svd,ourconvdestroy);
151:   }
152: }

154: SLEPC_EXTERN void svdstoppingbasic_(SVD*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,SVDConvergedReason*,void*,PetscErrorCode*);
155: SLEPC_EXTERN void svdstoppingthreshold_(SVD*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,SVDConvergedReason*,void*,PetscErrorCode*);

157: SLEPC_EXTERN void svdsetstoppingtestfunction_(SVD *svd,void (*func)(SVD*,PetscInt,PetscInt,PetscInt,PetscInt,SVDConvergedReason*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
158: {
159:   CHKFORTRANNULLOBJECT(ctx);
160:   CHKFORTRANNULLFUNCTION(destroy);
161:   if ((PetscVoidFunction)func == (PetscVoidFunction)svdstoppingbasic_) {
162:     *ierr = SVDSetStoppingTest(*svd,SVD_STOP_BASIC);
163:   } else if ((PetscVoidFunction)func == (PetscVoidFunction)svdstoppingthreshold_) {
164:     *ierr = SVDSetStoppingTest(*svd,SVD_STOP_THRESHOLD);
165:   } else {
166:     *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopping,(PetscVoidFunction)func,ctx); if (*ierr) return;
167:     *ierr = PetscObjectSetFortranCallback((PetscObject)*svd,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
168:     *ierr = SVDSetStoppingTestFunction(*svd,ourstopping,*svd,ourstopdestroy);
169:   }
170: }