Actual source code: zlmef.c

slepc-3.21.1 2024-04-26
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 <slepclme.h>

 14: #if defined(PETSC_HAVE_FORTRAN_CAPS)
 15: #define lmedestroy_                       LMEDESTROY
 16: #define lmeview_                          LMEVIEW
 17: #define lmeviewfromoptions_               LMEVIEWFROMOPTIONS
 18: #define lmeconvergedreasonview_           LMECONVERGEDREASONVIEW
 19: #define lmesetoptionsprefix_              LMESETOPTIONSPREFIX
 20: #define lmeappendoptionsprefix_           LMEAPPENDOPTIONSPREFIX
 21: #define lmegetoptionsprefix_              LMEGETOPTIONSPREFIX
 22: #define lmesettype_                       LMESETTYPE
 23: #define lmegettype_                       LMEGETTYPE
 24: #define lmemonitordefault_                LMEMONITORDEFAULT
 25: #define lmemonitorset_                    LMEMONITORSET
 26: #define lmegettolerances00_               LMEGETTOLERANCES00
 27: #define lmegettolerances10_               LMEGETTOLERANCES10
 28: #define lmegettolerances01_               LMEGETTOLERANCES01
 29: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 30: #define lmedestroy_                       lmedestroy
 31: #define lmeview_                          lmeview
 32: #define lmeviewfromoptions_               lmeviewfromoptions
 33: #define lmeconvergedreasonview_           lmeconvergedreasonview
 34: #define lmesetoptionsprefix_              lmesetoptionsprefix
 35: #define lmeappendoptionsprefix_           lmeappendoptionsprefix
 36: #define lmegetoptionsprefix_              lmegetoptionsprefix
 37: #define lmesettype_                       lmesettype
 38: #define lmegettype_                       lmegettype
 39: #define lmemonitordefault_                lmemonitordefault
 40: #define lmemonitorset_                    lmemonitorset
 41: #define lmegettolerances00_               lmegettolerances00
 42: #define lmegettolerances10_               lmegettolerances10
 43: #define lmegettolerances01_               lmegettolerances01
 44: #endif

 46: /*
 47:    These are not usually called from Fortran but allow Fortran users
 48:    to transparently set these monitors from .F code
 49: */
 50: SLEPC_EXTERN void lmemonitordefault_(LME *lme,PetscInt *it,PetscReal *errest,PetscViewerAndFormat **ctx,PetscErrorCode *ierr)
 51: {
 52:   *ierr = LMEMonitorDefault(*lme,*it,*errest,*ctx);
 53: }

 55: static struct {
 56:   PetscFortranCallbackId monitor;
 57:   PetscFortranCallbackId monitordestroy;
 58: } _cb;

 60: /* These are not extern C because they are passed into non-extern C user level functions */
 61: static PetscErrorCode ourmonitor(LME lme,PetscInt i,PetscReal d,void* ctx)
 62: {
 63:   PetscObjectUseFortranCallback(lme,_cb.monitor,(LME*,PetscInt*,PetscReal*,void*,PetscErrorCode*),(&lme,&i,&d,_ctx,&ierr));
 64: }

 66: static PetscErrorCode ourdestroy(void** ctx)
 67: {
 68:   LME lme = (LME)*ctx;
 69:   PetscObjectUseFortranCallback(lme,_cb.monitordestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
 70: }

 72: SLEPC_EXTERN void lmedestroy_(LME *lme,PetscErrorCode *ierr)
 73: {
 74:   PETSC_FORTRAN_OBJECT_F_DESTROYED_TO_C_NULL(lme);
 75:   *ierr = LMEDestroy(lme); if (*ierr) return;
 76:   PETSC_FORTRAN_OBJECT_C_NULL_TO_F_DESTROYED(lme);
 77: }

 79: SLEPC_EXTERN void lmeview_(LME *lme,PetscViewer *viewer,PetscErrorCode *ierr)
 80: {
 81:   PetscViewer v;
 82:   PetscPatchDefaultViewers_Fortran(viewer,v);
 83:   *ierr = LMEView(*lme,v);
 84: }

 86: SLEPC_EXTERN void lmeviewfromoptions_(LME *lme,PetscObject obj,char* type,PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len)
 87: {
 88:   char *t;

 90:   FIXCHAR(type,len,t);
 91:   CHKFORTRANNULLOBJECT(obj);
 92:   *ierr = LMEViewFromOptions(*lme,obj,t);if (*ierr) return;
 93:   FREECHAR(type,t);
 94: }

 96: SLEPC_EXTERN void lmeconvergedreasonview_(LME *lme,PetscViewer *viewer,PetscErrorCode *ierr)
 97: {
 98:   PetscViewer v;
 99:   PetscPatchDefaultViewers_Fortran(viewer,v);
100:   *ierr = LMEConvergedReasonView(*lme,v);
101: }

103: SLEPC_EXTERN void lmesettype_(LME *lme,char *type,PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len)
104: {
105:   char *t;

107:   FIXCHAR(type,len,t);
108:   *ierr = LMESetType(*lme,t);if (*ierr) return;
109:   FREECHAR(type,t);
110: }

112: SLEPC_EXTERN void lmegettype_(LME *lme,char *name,PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len)
113: {
114:   LMEType tname;

116:   *ierr = LMEGetType(*lme,&tname);if (*ierr) return;
117:   *ierr = PetscStrncpy(name,tname,len);if (*ierr) return;
118:   FIXRETURNCHAR(PETSC_TRUE,name,len);
119: }

121: SLEPC_EXTERN void lmesetoptionsprefix_(LME *lme,char *prefix,PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len)
122: {
123:   char *t;

125:   FIXCHAR(prefix,len,t);
126:   *ierr = LMESetOptionsPrefix(*lme,t);if (*ierr) return;
127:   FREECHAR(prefix,t);
128: }

130: SLEPC_EXTERN void lmeappendoptionsprefix_(LME *lme,char *prefix,PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len)
131: {
132:   char *t;

134:   FIXCHAR(prefix,len,t);
135:   *ierr = LMEAppendOptionsPrefix(*lme,t);if (*ierr) return;
136:   FREECHAR(prefix,t);
137: }

139: SLEPC_EXTERN void lmegetoptionsprefix_(LME *lme,char *prefix,PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len)
140: {
141:   const char *tname;

143:   *ierr = LMEGetOptionsPrefix(*lme,&tname); if (*ierr) return;
144:   *ierr = PetscStrncpy(prefix,tname,len);if (*ierr) return;
145:   FIXRETURNCHAR(PETSC_TRUE,prefix,len);
146: }

148: SLEPC_EXTERN void lmemonitorset_(LME *lme,void (*monitor)(LME*,PetscInt*,PetscReal*,void*,PetscErrorCode*),void *mctx,void (*monitordestroy)(void *,PetscErrorCode*),PetscErrorCode *ierr)
149: {
150:   CHKFORTRANNULLOBJECT(mctx);
151:   CHKFORTRANNULLFUNCTION(monitordestroy);
152:   if ((PetscVoidFunction)monitor == (PetscVoidFunction)lmemonitordefault_) {
153:     *ierr = LMEMonitorSet(*lme,(PetscErrorCode (*)(LME,PetscInt,PetscReal,void*))LMEMonitorDefault,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
154:   } else {
155:     *ierr = PetscObjectSetFortranCallback((PetscObject)*lme,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)monitor,mctx); if (*ierr) return;
156:     *ierr = PetscObjectSetFortranCallback((PetscObject)*lme,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitordestroy,(PetscVoidFunction)monitordestroy,mctx); if (*ierr) return;
157:     *ierr = LMEMonitorSet(*lme,ourmonitor,*lme,ourdestroy);
158:   }
159: }

161: SLEPC_EXTERN void lmegettolerances_(LME *lme,PetscReal *tol,PetscInt *maxits,PetscErrorCode *ierr)
162: {
163:   CHKFORTRANNULLREAL(tol);
164:   CHKFORTRANNULLINTEGER(maxits);
165:   *ierr = LMEGetTolerances(*lme,tol,maxits);
166: }

168: SLEPC_EXTERN void lmegettolerances00_(LME *lme,PetscReal *tol,PetscInt *maxits,PetscErrorCode *ierr)
169: {
170:   lmegettolerances_(lme,tol,maxits,ierr);
171: }

173: SLEPC_EXTERN void lmegettolerances10_(LME *lme,PetscReal *tol,PetscInt *maxits,PetscErrorCode *ierr)
174: {
175:   lmegettolerances_(lme,tol,maxits,ierr);
176: }

178: SLEPC_EXTERN void lmegettolerances01_(LME *lme,PetscReal *tol,PetscInt *maxits,PetscErrorCode *ierr)
179: {
180:   lmegettolerances_(lme,tol,maxits,ierr);
181: }