1 #include <petsc/private/f90impl.h> 2 #include <petsc/private/sfimpl.h> 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define petscsfview_ PETSCSFVIEW 6 #define petscsfgetgraph_ PETSCSFGETGRAPH 7 #define petscsfbcastbegin_ PETSCSFBCASTBEGIN 8 #define petscsfbcastend_ PETSCSFBCASTEND 9 #define f90arraysfnodecreate_ F90ARRAYSFNODECREATE 10 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 11 #define petscsfgetgraph_ petscsfgetgraph 12 #define petscsfview_ petscsfview 13 #define petscsfbcastbegin_ petscsfbcastbegin 14 #define petscsfbcastend_ petscsfbcastend 15 #define f90arraysfnodecreate_ f90arraysfnodecreate 16 #endif 17 18 PETSC_EXTERN void PETSC_STDCALL f90arraysfnodecreate_(const PetscInt *,PetscInt *,void * PETSC_F90_2PTR_PROTO_NOVAR); 19 20 PETSC_EXTERN void PETSC_STDCALL petscsfview_(PetscSF *sf, PetscViewer *vin, PetscErrorCode *ierr) 21 { 22 PetscViewer v; 23 24 PetscPatchDefaultViewers_Fortran(vin, v); 25 *ierr = PetscSFView(*sf, v); 26 } 27 28 29 PETSC_EXTERN void PETSC_STDCALL petscsfgetgraph_(PetscSF *sf,PetscInt *nroots,PetscInt *nleaves, F90Array1d *ailocal, F90Array1d *airemote, int *ierr PETSC_F90_2PTR_PROTO(pilocal) PETSC_F90_2PTR_PROTO(piremote)) 30 { 31 const PetscInt *ilocal; 32 const PetscSFNode *iremote; 33 34 *ierr = PetscSFGetGraph(*sf,nroots,nleaves,&ilocal,&iremote);if (*ierr) return; 35 *ierr = F90Array1dCreate((void*)ilocal,MPIU_INT,1,*nleaves, ailocal PETSC_F90_2PTR_PARAM(pilocal)); 36 /* this creates a memory leak */ 37 f90arraysfnodecreate_((PetscInt*)iremote,nleaves, airemote PETSC_F90_2PTR_PARAM(piremote)); 38 } 39 40 #if defined(PETSC_HAVE_F90_ASSUMED_TYPE_NOT_PTR) 41 PETSC_EXTERN void PETSC_STDCALL petscsfbcastbegin_(PetscSF *sf, MPI_Fint *unit, const void *rptr, void *lptr, int *ierr) 42 { 43 MPI_Datatype dtype; 44 45 *ierr = PetscMPIFortranDatatypeToC(*unit,&dtype);if (*ierr) return; 46 *ierr = PetscSFBcastBegin(*sf, dtype, rptr, lptr); 47 } 48 49 50 PETSC_EXTERN void PETSC_STDCALL petscsfbcastend_(PetscSF *sf, MPI_Fint *unit, const void *rptr, void *lptr, int *ierr) 51 { 52 MPI_Datatype dtype; 53 54 *ierr = PetscMPIFortranDatatypeToC(*unit,&dtype);if (*ierr) return; 55 *ierr = PetscSFBcastEnd(*sf, dtype, rptr, lptr); 56 } 57 58 #else 59 60 PETSC_EXTERN void PETSC_STDCALL petscsfbcastbegin_(PetscSF *sf, MPI_Fint *unit,F90Array1d *rptr, F90Array1d *lptr, int *ierr PETSC_F90_2PTR_PROTO(rptrd) PETSC_F90_2PTR_PROTO(lptrd)) 61 { 62 MPI_Datatype dtype; 63 const void *rootdata; 64 void *leafdata; 65 66 67 *ierr = PetscMPIFortranDatatypeToC(*unit,&dtype);if (*ierr) return; 68 *ierr = F90Array1dAccess(rptr, dtype, (void**) &rootdata PETSC_F90_2PTR_PARAM(rptrd));if (*ierr) return; 69 *ierr = F90Array1dAccess(lptr, dtype, (void**) &leafdata PETSC_F90_2PTR_PARAM(lptrd));if (*ierr) return; 70 *ierr = PetscSFBcastBegin(*sf, dtype, rootdata, leafdata); 71 } 72 73 PETSC_EXTERN void PETSC_STDCALL petscsfbcastend_(PetscSF *sf, MPI_Fint *unit,F90Array1d *rptr, F90Array1d *lptr, int *ierr PETSC_F90_2PTR_PROTO(rptrd) PETSC_F90_2PTR_PROTO(lptrd)) 74 { 75 MPI_Datatype dtype; 76 const void *rootdata; 77 void *leafdata; 78 79 *ierr = PetscMPIFortranDatatypeToC(*unit,&dtype);if (*ierr) return; 80 *ierr = F90Array1dAccess(rptr, dtype, (void**) &rootdata PETSC_F90_2PTR_PARAM(rptrd));if (*ierr) return; 81 *ierr = F90Array1dAccess(lptr, dtype, (void**) &leafdata PETSC_F90_2PTR_PARAM(lptrd));if (*ierr) return; 82 *ierr = PetscSFBcastEnd(*sf, dtype, rootdata, leafdata); 83 } 84 85 #endif 86