1 #include <petsc/private/ftnimpl.h>
2 #include <petscksp.h>
3
4 #if defined(PETSC_HAVE_FORTRAN_CAPS)
5 #define pcfieldsplitgetsubksp_ PCFIELDSPLITGETSUBKSP
6 #define pcfieldsplitschurgetsubksp_ PCFIELDSPLITSCHURGETSUBKSP
7 #define pcfieldsplitrestoresubksp_ PCFIELDSPLITRESTORESUBKSP
8 #define pcfieldsplitschurrestoresubksp_ PCFIELDSPLITSCHURRESTORESUBKSP
9 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
10 #define pcfieldsplitgetsubksp_ pcfieldsplitgetsubksp
11 #define pcfieldsplitschurgetsubksp_ pcfieldsplitschurgetsubksp
12 #define pcfieldsplitrestoresubksp_ pcfieldsplitrestoresubksp
13 #define pcfieldsplitschurrestoresubksp_ pcfieldsplitschurrestoresubksp
14 #endif
15
pcfieldsplitschurgetsubksp_(PC * pc,PetscInt * n_local,F90Array1d * ksp,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))16 PETSC_EXTERN void pcfieldsplitschurgetsubksp_(PC *pc, PetscInt *n_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
17 {
18 KSP *tksp;
19 PetscInt nloc;
20 CHKFORTRANNULLINTEGER(n_local);
21 *ierr = PCFieldSplitSchurGetSubKSP(*pc, &nloc, &tksp);
22 if (*ierr) return;
23 if (n_local) *n_local = nloc;
24 *ierr = F90Array1dCreate(tksp, MPIU_FORTRANADDR, 1, nloc, ksp PETSC_F90_2PTR_PARAM(ptrd));
25 }
26
pcfieldsplitgetsubksp_(PC * pc,PetscInt * n_local,F90Array1d * ksp,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))27 PETSC_EXTERN void pcfieldsplitgetsubksp_(PC *pc, PetscInt *n_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
28 {
29 KSP *tksp;
30 PetscInt nloc;
31 CHKFORTRANNULLINTEGER(n_local);
32 *ierr = PCFieldSplitGetSubKSP(*pc, &nloc, &tksp);
33 if (*ierr) return;
34 if (n_local) *n_local = nloc;
35 *ierr = F90Array1dCreate(tksp, MPIU_FORTRANADDR, 1, nloc, ksp PETSC_F90_2PTR_PARAM(ptrd));
36 }
37
pcfieldsplitrestoresubksp_(PC * pc,PetscInt * n_local,F90Array1d * ksp,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))38 PETSC_EXTERN void pcfieldsplitrestoresubksp_(PC *pc, PetscInt *n_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
39 {
40 void *array;
41 *ierr = F90Array1dAccess(ksp, MPIU_FORTRANADDR, (void **)&array PETSC_F90_2PTR_PARAM(ptrd));
42 if (*ierr) return;
43 *ierr = F90Array1dDestroy(ksp, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd));
44 if (*ierr) return;
45 *ierr = PetscFree(array);
46 }
47
pcfieldsplitschurerestoresubksp_(PC * pc,PetscInt * n_local,F90Array1d * ksp,PetscErrorCode * ierr PETSC_F90_2PTR_PROTO (ptrd))48 PETSC_EXTERN void pcfieldsplitschurerestoresubksp_(PC *pc, PetscInt *n_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
49 {
50 void *array;
51 *ierr = F90Array1dAccess(ksp, MPIU_FORTRANADDR, (void **)&array PETSC_F90_2PTR_PARAM(ptrd));
52 if (*ierr) return;
53 *ierr = F90Array1dDestroy(ksp, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd));
54 if (*ierr) return;
55 *ierr = PetscFree(array);
56 }
57