Actual source code: zsfutilsf90.c
1: #include <petsc/private/ftnimpl.h>
2: #include <petscsf.h>
3: #include <petscsection.h>
5: #if defined(PETSC_HAVE_FORTRAN_CAPS)
6: #define petscsfdestroyremoteoffsets_ PETSCSFDESTROYREMOTEOFFSETS
7: #define petscsfdistributesection_ PETSCSFDISTRIBUTESECTION
8: #define petscsfcreateremoteoffsets_ PETSCSFCREATEREMOTEOFFSETS
9: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
10: #define petscsfdestroyremoteoffsets_ petscsfdestroyremoteoffsets
11: #define petscsfdistributesection_ petscsfdistributesection
12: #define petscsfcreateremoteoffsets_ petscsfcreateremoteoffsets
13: #endif
15: PETSC_EXTERN void petscsfdestroyremoteoffsets_(F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
16: {
17: PetscInt *fa;
18: *ierr = F90Array1dAccess(ptr, MPIU_INT, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
19: if (*ierr) return;
20: *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
21: if (*ierr) return;
22: *ierr = PetscFree(fa);
23: }
25: PETSC_EXTERN void petscsfdistributesection_(PetscSF *sf, PetscSection *rootSection, F90Array1d *ptr, PetscSection *leafSection, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
26: {
27: if (ptr == PETSC_NULL_INTEGER_POINTER_Fortran) {
28: *ierr = PetscSFDistributeSection(*sf, *rootSection, NULL, *leafSection);
29: } else {
30: PetscInt *fa;
31: PetscInt lpStart, lpEnd;
33: *ierr = PetscSFDistributeSection(*sf, *rootSection, &fa, *leafSection);
34: if (*ierr) return;
35: *ierr = PetscSectionGetChart(*leafSection, &lpStart, &lpEnd);
36: if (*ierr) return;
37: *ierr = F90Array1dCreate((void *)fa, MPIU_INT, 1, lpEnd - lpStart, ptr PETSC_F90_2PTR_PARAM(ptrd));
38: }
39: }
41: PETSC_EXTERN void petscsfcreateremoteoffsets_(PetscSF *pointSF, PetscSection *rootSection, PetscSection *leafSection, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
42: {
43: PetscInt *fa;
44: PetscInt lpStart, lpEnd;
46: *ierr = PetscSFCreateRemoteOffsets(*pointSF, *rootSection, *leafSection, &fa);
47: if (*ierr) return;
48: *ierr = PetscSectionGetChart(*leafSection, &lpStart, &lpEnd);
49: if (*ierr) return;
50: *ierr = F90Array1dCreate((void *)fa, MPIU_INT, 1, lpEnd - lpStart, ptr PETSC_F90_2PTR_PARAM(ptrd));
51: }