Actual source code: zgasmf.c
1: #include <petsc/private/ftnimpl.h>
2: #include <petscksp.h>
4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
5: #define pcgasmdestroysubdomains_ PCGASMDESTROYSUBDOMAINS
6: #define pcgasmgetsubksp_ PCGASMGETSUBKSP
7: #define pcgasmrestoresubksp_ PCGASMRESTORESUBKSP
8: #define pcgasmcreatesubdomains2d_ PCGASMCREATESUBDOMAINS2D
9: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
10: #define pcgasmdestroysubdomains_ pcgasmdestroysubdomains
11: #define pcgasmgetsubksp_ pcgasmgetsubksp
12: #define pcgasmrestoresubksp_ pcgasmrestoresubksp
13: #define pcgasmcreatesubdomains2d_ pcgasmcreatesubdomains2d
14: #endif
16: PETSC_EXTERN void pcgasmdestroysubdomains_(PetscInt *n, F90Array1d *is1, F90Array1d *is2, int *ierr PETSC_F90_2PTR_PROTO(ptrd1) PETSC_F90_2PTR_PROTO(ptrd2))
17: {
18: IS *isa, *isb;
20: *ierr = F90Array1dAccess(is1, MPIU_FORTRANADDR, (void **)&isa PETSC_F90_2PTR_PARAM(ptrd1));
21: if (*ierr) return;
22: *ierr = F90Array1dAccess(is2, MPIU_FORTRANADDR, (void **)&isb PETSC_F90_2PTR_PARAM(ptrd2));
23: if (*ierr) return;
24: *ierr = F90Array1dDestroy(is1, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd1));
25: if (*ierr) return;
26: *ierr = F90Array1dDestroy(is2, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd2));
27: if (*ierr) return;
28: *ierr = PCGASMDestroySubdomains(*n, &isa, &isb);
29: }
31: 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))
32: {
33: IS *iis, *iisl;
34: *ierr = PCGASMCreateSubdomains2D(*pc, *m, *n, *M, *N, *dof, *overlap, Nsub, &iis, &iisl);
35: if (*ierr) return;
36: *ierr = F90Array1dCreate(iis, MPIU_FORTRANADDR, 1, *Nsub, is1 PETSC_F90_2PTR_PARAM(ptrd1));
37: if (*ierr) return;
38: *ierr = F90Array1dCreate(iisl, MPIU_FORTRANADDR, 1, *Nsub, is2 PETSC_F90_2PTR_PARAM(ptrd2));
39: if (*ierr) return;
40: }
42: PETSC_EXTERN void pcgasmgetsubksp_(PC *pc, PetscInt *n_local, PetscInt *first_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
43: {
44: KSP *tksp;
45: PetscInt nloc, flocal;
47: CHKFORTRANNULLINTEGER(n_local);
48: CHKFORTRANNULLINTEGER(first_local);
49: *ierr = PCGASMGetSubKSP(*pc, &nloc, &flocal, &tksp);
50: if (n_local) *n_local = nloc;
51: if (first_local) *first_local = flocal;
52: *ierr = F90Array1dCreate(tksp, MPIU_FORTRANADDR, 1, nloc, ksp PETSC_F90_2PTR_PARAM(ptrd));
53: }
55: PETSC_EXTERN void pcgasmrestoresubksp_(PC *pc, PetscInt *n_local, PetscInt *first_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
56: {
57: *ierr = F90Array1dDestroy(ksp, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd));
58: }