16dd63270SBarry Smith #include <petsc/private/ftnimpl.h> 2768de9afSRichard Tran Mills #include <petscksp.h> 3768de9afSRichard Tran Mills 4768de9afSRichard Tran Mills #if defined(PETSC_HAVE_FORTRAN_CAPS) 5b3598b19SBarry Smith #define pcgasmdestroysubdomains_ PCGASMDESTROYSUBDOMAINS 66f8503afSBarry Smith #define pcgasmgetsubksp_ PCGASMGETSUBKSP 7*36083efbSBarry Smith #define pcgasmrestoresubksp_ PCGASMRESTORESUBKSP 86141accfSBarry Smith #define pcgasmcreatesubdomains2d_ PCGASMCREATESUBDOMAINS2D 9768de9afSRichard Tran Mills #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 10b3598b19SBarry Smith #define pcgasmdestroysubdomains_ pcgasmdestroysubdomains 116f8503afSBarry Smith #define pcgasmgetsubksp_ pcgasmgetsubksp 12*36083efbSBarry Smith #define pcgasmrestoresubksp_ pcgasmrestoresubksp 136141accfSBarry Smith #define pcgasmcreatesubdomains2d_ pcgasmcreatesubdomains2d 14768de9afSRichard Tran Mills #endif 15768de9afSRichard Tran Mills 16ce78bad3SBarry Smith PETSC_EXTERN void pcgasmdestroysubdomains_(PetscInt *n, F90Array1d *is1, F90Array1d *is2, int *ierr PETSC_F90_2PTR_PROTO(ptrd1) PETSC_F90_2PTR_PROTO(ptrd2)) 17b3598b19SBarry Smith { 18ce78bad3SBarry Smith IS *isa, *isb; 19ce78bad3SBarry Smith 20ce78bad3SBarry Smith *ierr = F90Array1dAccess(is1, MPIU_FORTRANADDR, (void **)&isa PETSC_F90_2PTR_PARAM(ptrd1)); 21b3598b19SBarry Smith if (*ierr) return; 22ce78bad3SBarry Smith *ierr = F90Array1dAccess(is2, MPIU_FORTRANADDR, (void **)&isb PETSC_F90_2PTR_PARAM(ptrd2)); 23b3598b19SBarry Smith if (*ierr) return; 24ce78bad3SBarry Smith *ierr = F90Array1dDestroy(is1, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd1)); 25b3598b19SBarry Smith if (*ierr) return; 26ce78bad3SBarry Smith *ierr = F90Array1dDestroy(is2, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd2)); 27ce78bad3SBarry Smith if (*ierr) return; 28ce78bad3SBarry Smith *ierr = PCGASMDestroySubdomains(*n, &isa, &isb); 29b3598b19SBarry Smith } 30b3598b19SBarry Smith 31ce78bad3SBarry Smith PETSC_EXTERN void pcgasmcreatesubdomains2d_(PC *pc, PetscInt *m, PetscInt *n, PetscInt *M, PetscInt *N, PetscInt *dof, PetscInt *overlap, PetscInt *Nsub, F90Array1d *is1, F90Array1d *is2, int *ierr PETSC_F90_2PTR_PROTO(ptrd1) PETSC_F90_2PTR_PROTO(ptrd2)) 326141accfSBarry Smith { 336141accfSBarry Smith IS *iis, *iisl; 346141accfSBarry Smith *ierr = PCGASMCreateSubdomains2D(*pc, *m, *n, *M, *N, *dof, *overlap, Nsub, &iis, &iisl); 356141accfSBarry Smith if (*ierr) return; 36ce78bad3SBarry Smith *ierr = F90Array1dCreate(iis, MPIU_FORTRANADDR, 1, *Nsub, is1 PETSC_F90_2PTR_PARAM(ptrd1)); 376141accfSBarry Smith if (*ierr) return; 38ce78bad3SBarry Smith *ierr = F90Array1dCreate(iisl, MPIU_FORTRANADDR, 1, *Nsub, is2 PETSC_F90_2PTR_PARAM(ptrd2)); 396141accfSBarry Smith if (*ierr) return; 406141accfSBarry Smith } 416141accfSBarry Smith 42ce78bad3SBarry Smith PETSC_EXTERN void pcgasmgetsubksp_(PC *pc, PetscInt *n_local, PetscInt *first_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 43768de9afSRichard Tran Mills { 44768de9afSRichard Tran Mills KSP *tksp; 45ce78bad3SBarry Smith PetscInt nloc, flocal; 46ce78bad3SBarry Smith 47768de9afSRichard Tran Mills CHKFORTRANNULLINTEGER(n_local); 48768de9afSRichard Tran Mills CHKFORTRANNULLINTEGER(first_local); 49ce78bad3SBarry Smith *ierr = PCGASMGetSubKSP(*pc, &nloc, &flocal, &tksp); 50768de9afSRichard Tran Mills if (n_local) *n_local = nloc; 51ce78bad3SBarry Smith if (first_local) *first_local = flocal; 52ce78bad3SBarry Smith *ierr = F90Array1dCreate(tksp, MPIU_FORTRANADDR, 1, nloc, ksp PETSC_F90_2PTR_PARAM(ptrd)); 5359fc98c0SBarry Smith } 54*36083efbSBarry Smith 55*36083efbSBarry Smith PETSC_EXTERN void pcgasmrestoresubksp_(PC *pc, PetscInt *n_local, PetscInt *first_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 56*36083efbSBarry Smith { 57*36083efbSBarry Smith *ierr = F90Array1dDestroy(ksp, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd)); 58*36083efbSBarry Smith } 59