16dd63270SBarry Smith #include <petsc/private/ftnimpl.h> 21193e19dSBarry Smith #include <petscksp.h> 31193e19dSBarry Smith 41193e19dSBarry Smith #if defined(PETSC_HAVE_FORTRAN_CAPS) 51193e19dSBarry Smith #define pcfieldsplitgetsubksp_ PCFIELDSPLITGETSUBKSP 6285fb4e2SStefano Zampini #define pcfieldsplitschurgetsubksp_ PCFIELDSPLITSCHURGETSUBKSP 7*e41f517fSBarry Smith #define pcfieldsplitrestoresubksp_ PCFIELDSPLITRESTORESUBKSP 8*e41f517fSBarry Smith #define pcfieldsplitschurrestoresubksp_ PCFIELDSPLITSCHURRESTORESUBKSP 91193e19dSBarry Smith #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 101193e19dSBarry Smith #define pcfieldsplitgetsubksp_ pcfieldsplitgetsubksp 11285fb4e2SStefano Zampini #define pcfieldsplitschurgetsubksp_ pcfieldsplitschurgetsubksp 12*e41f517fSBarry Smith #define pcfieldsplitrestoresubksp_ pcfieldsplitrestoresubksp 13*e41f517fSBarry Smith #define pcfieldsplitschurrestoresubksp_ pcfieldsplitschurrestoresubksp 141193e19dSBarry Smith #endif 151193e19dSBarry Smith 16*e41f517fSBarry Smith PETSC_EXTERN void pcfieldsplitschurgetsubksp_(PC *pc, PetscInt *n_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 17285fb4e2SStefano Zampini { 18285fb4e2SStefano Zampini KSP *tksp; 19*e41f517fSBarry Smith PetscInt nloc; 20285fb4e2SStefano Zampini CHKFORTRANNULLINTEGER(n_local); 215975b3b6SBarry Smith *ierr = PCFieldSplitSchurGetSubKSP(*pc, &nloc, &tksp); 225975b3b6SBarry Smith if (*ierr) return; 23285fb4e2SStefano Zampini if (n_local) *n_local = nloc; 24*e41f517fSBarry Smith *ierr = F90Array1dCreate(tksp, MPIU_FORTRANADDR, 1, nloc, ksp PETSC_F90_2PTR_PARAM(ptrd)); 25285fb4e2SStefano Zampini } 26285fb4e2SStefano Zampini 27*e41f517fSBarry Smith PETSC_EXTERN void pcfieldsplitgetsubksp_(PC *pc, PetscInt *n_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 281193e19dSBarry Smith { 291193e19dSBarry Smith KSP *tksp; 30*e41f517fSBarry Smith PetscInt nloc; 311193e19dSBarry Smith CHKFORTRANNULLINTEGER(n_local); 325975b3b6SBarry Smith *ierr = PCFieldSplitGetSubKSP(*pc, &nloc, &tksp); 335975b3b6SBarry Smith if (*ierr) return; 341193e19dSBarry Smith if (n_local) *n_local = nloc; 35*e41f517fSBarry Smith *ierr = F90Array1dCreate(tksp, MPIU_FORTRANADDR, 1, nloc, ksp PETSC_F90_2PTR_PARAM(ptrd)); 361193e19dSBarry Smith } 37*e41f517fSBarry Smith 38*e41f517fSBarry Smith PETSC_EXTERN void pcfieldsplitrestoresubksp_(PC *pc, PetscInt *n_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 39*e41f517fSBarry Smith { 40*e41f517fSBarry Smith void *array; 41*e41f517fSBarry Smith *ierr = F90Array1dAccess(ksp, MPIU_FORTRANADDR, (void **)&array PETSC_F90_2PTR_PARAM(ptrd)); 42*e41f517fSBarry Smith if (*ierr) return; 43*e41f517fSBarry Smith *ierr = F90Array1dDestroy(ksp, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd)); 44*e41f517fSBarry Smith if (*ierr) return; 45*e41f517fSBarry Smith *ierr = PetscFree(array); 46*e41f517fSBarry Smith } 47*e41f517fSBarry Smith 48*e41f517fSBarry Smith PETSC_EXTERN void pcfieldsplitschurerestoresubksp_(PC *pc, PetscInt *n_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 49*e41f517fSBarry Smith { 50*e41f517fSBarry Smith void *array; 51*e41f517fSBarry Smith *ierr = F90Array1dAccess(ksp, MPIU_FORTRANADDR, (void **)&array PETSC_F90_2PTR_PARAM(ptrd)); 52*e41f517fSBarry Smith if (*ierr) return; 53*e41f517fSBarry Smith *ierr = F90Array1dDestroy(ksp, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd)); 54*e41f517fSBarry Smith if (*ierr) return; 55*e41f517fSBarry Smith *ierr = PetscFree(array); 561193e19dSBarry Smith } 57