1 #include <petsc/private/fortranimpl.h> 2 #include <petscksp.h> 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define pcfieldsplitgetsubksp_ PCFIELDSPLITGETSUBKSP 6 #define pcfieldsplitschurgetsubksp_ PCFIELDSPLITSCHURGETSUBKSP 7 #define pcfieldsplitsetis_ PCFIELDSPLITSETIS 8 #define pcfieldsplitgetis_ PCFIELDSPLITGETIS 9 #define pcfieldsplitsetfields_ PCFIELDSPLITSETFIELDS 10 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 11 #define pcfieldsplitgetsubksp_ pcfieldsplitgetsubksp 12 #define pcfieldsplitschurgetsubksp_ pcfieldsplitschurgetsubksp 13 #define pcfieldsplitsetis_ pcfieldsplitsetis 14 #define pcfieldsplitgetis_ pcfieldsplitgetis 15 #define pcfieldsplitsetfields_ pcfieldsplitsetfields 16 #endif 17 18 PETSC_EXTERN void pcfieldsplitschurgetsubksp_(PC *pc,PetscInt *n_local,KSP *ksp,PetscErrorCode *ierr) 19 { 20 KSP *tksp; 21 PetscInt i,nloc; 22 CHKFORTRANNULLINTEGER(n_local); 23 *ierr = PCFieldSplitSchurGetSubKSP(*pc,&nloc,&tksp); if (*ierr) return; 24 if (n_local) *n_local = nloc; 25 CHKFORTRANNULLOBJECT(ksp); 26 if (ksp) { 27 for (i=0; i<nloc; i++) ksp[i] = tksp[i]; 28 } 29 *ierr = PetscFree(tksp); 30 } 31 32 PETSC_EXTERN void pcfieldsplitgetsubksp_(PC *pc,PetscInt *n_local,KSP *ksp,PetscErrorCode *ierr) 33 { 34 KSP *tksp; 35 PetscInt i,nloc; 36 CHKFORTRANNULLINTEGER(n_local); 37 *ierr = PCFieldSplitGetSubKSP(*pc,&nloc,&tksp); if (*ierr) return; 38 if (n_local) *n_local = nloc; 39 CHKFORTRANNULLOBJECT(ksp); 40 if (ksp) { 41 for (i=0; i<nloc; i++) ksp[i] = tksp[i]; 42 } 43 *ierr = PetscFree(tksp); 44 } 45 46 PETSC_EXTERN void pcfieldsplitsetis_(PC *pc, char* splitname,IS *is, PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len) 47 { 48 char *t; 49 FIXCHAR(splitname,len,t); 50 *ierr = PCFieldSplitSetIS(*pc,t,*is);if (*ierr) return; 51 FREECHAR(splitname,t); 52 } 53 54 PETSC_EXTERN void pcfieldsplitgetis_(PC *pc, char* splitname,IS *is, PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len) 55 { 56 char *t; 57 FIXCHAR(splitname,len,t); 58 *ierr = PCFieldSplitGetIS(*pc,t,is);if (*ierr) return; 59 FREECHAR(splitname,t); 60 } 61 62 PETSC_EXTERN void pcfieldsplitsetfields_(PC *pc, char *splitname, PetscInt *n, const PetscInt *fields, const PetscInt *fields_col, PetscErrorCode *ierr,PETSC_FORTRAN_CHARLEN_T len) 63 { 64 char *t; 65 FIXCHAR(splitname,len,t); 66 *ierr = PCFieldSplitSetFields(*pc, t, *n, fields, fields_col);if (*ierr) return; 67 FREECHAR(splitname,t); 68 } 69