Actual source code: zgmres2f.c

  1: #include <petsc/private/ftnimpl.h>
  2: #include <petscksp.h>

  4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
  5:   #define kspgmressetorthogonalization_                  KSPGMRESSETORTHOGONALIZATION
  6:   #define kspgmresmodifiedgramschmidtorthogonalization_  KSPGMRESMODIFIEDGRAMSCHMIDTORTHOGONALIZATION
  7:   #define kspgmresclassicalgramschmidtorthogonalization_ KSPGMRESCLASSICALGRAMSCHMIDTORTHOGONALIZATION
  8: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
  9:   #define kspgmressetorthogonalization_                  kspgmressetorthogonalization
 10:   #define kspgmresmodifiedgramschmidtorthogonalization_  kspgmresmodifiedgramschmidtorthogonalization
 11:   #define kspgmresclassicalgramschmidtorthogonalization_ kspgmresclassicalgramschmidtorthogonalization
 12: #endif

 14: static struct {
 15:   PetscFortranCallbackId orthog;
 16: } _cb;

 18: PETSC_EXTERN void kspgmresmodifiedgramschmidtorthogonalization_(KSP *, PetscInt *, PetscErrorCode *);
 19: PETSC_EXTERN void kspgmresclassicalgramschmidtorthogonalization_(KSP *, PetscInt *, PetscErrorCode *);

 21: static PetscErrorCode ourorthog(KSP ksp, PetscInt n)
 22: {
 23:   PetscObjectUseFortranCallback(ksp, _cb.orthog, (KSP *, PetscInt *, PetscErrorCode *), (&ksp, &n, &ierr));
 24: }

 26: PETSC_EXTERN void kspgmressetorthogonalization_(KSP *ksp, void (*orthog)(KSP *, PetscInt *, PetscErrorCode *), PetscErrorCode *ierr)
 27: {
 28:   if ((PetscVoidFn *)orthog == (PetscVoidFn *)kspgmresmodifiedgramschmidtorthogonalization_) {
 29:     *ierr = KSPGMRESSetOrthogonalization(*ksp, KSPGMRESModifiedGramSchmidtOrthogonalization);
 30:   } else if ((PetscVoidFn *)orthog == (PetscVoidFn *)kspgmresclassicalgramschmidtorthogonalization_) {
 31:     *ierr = KSPGMRESSetOrthogonalization(*ksp, KSPGMRESClassicalGramSchmidtOrthogonalization);
 32:   } else {
 33:     *ierr = PetscObjectSetFortranCallback((PetscObject)*ksp, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.orthog, (PetscVoidFn *)orthog, NULL);
 34:     if (*ierr) return;
 35:     *ierr = KSPGMRESSetOrthogonalization(*ksp, ourorthog);
 36:   }
 37: }