1 #include <petscis.h> 2 #include <petsc/private/ftnimpl.h> 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define islocaltoglobalmappinggetindices_ ISLOCALTOGLOBALMAPPINGGETINDICES 6 #define islocaltoglobalmappingrestoreindices_ ISLOCALTOGLOBALMAPPINGRESTOREINDICES 7 #define islocaltoglobalmappinggetblockindices_ ISLOCALTOGLOBALMAPPINGGETBLOCKINDICES 8 #define islocaltoglobalmappingrestoreblockindices_ ISLOCALTOGLOBALMAPPINGRESTOREBLOCKINDICES 9 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 10 #define islocaltoglobalmappinggetindices_ islocaltoglobalmappinggetindices 11 #define islocaltoglobalmappingrestoreindices_ islocaltoglobalmappingrestoreindices 12 #define islocaltoglobalmappinggetblockindices_ islocaltoglobalmappinggetblockindices 13 #define islocaltoglobalmappingrestoreblockindices_ islocaltoglobalmappingrestoreblockindices 14 #endif 15 16 PETSC_EXTERN void islocaltoglobalmappinggetindices_(ISLocalToGlobalMapping *da, F90Array1d *indices, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 17 { 18 const PetscInt *idx; 19 PetscInt n; 20 *ierr = ISLocalToGlobalMappingGetIndices(*da, &idx); 21 if (*ierr) return; 22 *ierr = ISLocalToGlobalMappingGetSize(*da, &n); 23 if (*ierr) return; 24 *ierr = F90Array1dCreate((void *)idx, MPIU_INT, 1, n, indices PETSC_F90_2PTR_PARAM(ptrd)); 25 } 26 27 PETSC_EXTERN void islocaltoglobalmappingrestoreindices_(ISLocalToGlobalMapping *da, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 28 { 29 const PetscInt *fa; 30 31 *ierr = F90Array1dAccess(ptr, MPIU_INT, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 32 if (*ierr) return; 33 *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd)); 34 if (*ierr) return; 35 *ierr = ISLocalToGlobalMappingRestoreIndices(*da, &fa); 36 if (*ierr) return; 37 } 38 39 PETSC_EXTERN void islocaltoglobalmappinggetblockindices_(ISLocalToGlobalMapping *da, F90Array1d *indices, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 40 { 41 const PetscInt *idx; 42 PetscInt n; 43 *ierr = ISLocalToGlobalMappingGetBlockIndices(*da, &idx); 44 if (*ierr) return; 45 *ierr = ISLocalToGlobalMappingGetSize(*da, &n); 46 if (*ierr) return; 47 *ierr = F90Array1dCreate((void *)idx, MPIU_INT, 1, n, indices PETSC_F90_2PTR_PARAM(ptrd)); 48 } 49 50 PETSC_EXTERN void islocaltoglobalmappingrestoreblockindices_(ISLocalToGlobalMapping *da, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 51 { 52 const PetscInt *fa; 53 54 *ierr = F90Array1dAccess(ptr, MPIU_INT, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 55 if (*ierr) return; 56 *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd)); 57 if (*ierr) return; 58 *ierr = ISLocalToGlobalMappingRestoreBlockIndices(*da, &fa); 59 if (*ierr) return; 60 } 61