1 #include <petsc/private/f90impl.h> 2 #include <petsc/private/sfimpl.h> 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define petscsfgetgraph_ PETSCSFGETGRAPH 6 #define petscsfbcastbegin_ PETSCSFBCASTBEGIN 7 #define petscsfbcastend_ PETSCSFBCASTEND 8 #define petscsfreducebegin_ PETSCSFREDUCEBEGIN 9 #define petscsfreduceend_ PETSCSFREDUCEEND 10 #define f90arraysfnodecreate_ F90ARRAYSFNODECREATE 11 #define petscsfgetleafranks_ PETSCSFGETLEAFRANKS 12 #define petscsfgetrootranks_ PETSCSFGETROOTRANKS 13 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 14 #define petscsfgetgraph_ petscsfgetgraph 15 #define petscsfbcastbegin_ petscsfbcastbegin 16 #define petscsfbcastend_ petscsfbcastend 17 #define petscsfreducebegin_ petscsfreducebegin 18 #define petscsfreduceend_ petscsfreduceend 19 #define f90arraysfnodecreate_ f90arraysfnodecreate 20 #define petscsfgetleafranks_ petscsfgetleafranks 21 #define petscsfgetrootranks_ petscsfgetrootranks 22 #endif 23 24 PETSC_EXTERN void f90arraysfnodecreate_(const PetscInt *, PetscInt *, void *PETSC_F90_2PTR_PROTO_NOVAR); 25 26 PETSC_EXTERN void petscsfgetgraph_(PetscSF *sf, PetscInt *nroots, PetscInt *nleaves, F90Array1d *ailocal, F90Array1d *airemote, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(pilocal) PETSC_F90_2PTR_PROTO(piremote)) 27 { 28 const PetscInt *ilocal; 29 const PetscSFNode *iremote; 30 PetscInt nl; 31 32 *ierr = PetscSFGetGraph(*sf, nroots, nleaves, &ilocal, &iremote); 33 if (*ierr) return; 34 nl = *nleaves; 35 if (!ilocal) nl = 0; 36 *ierr = F90Array1dCreate((void *)ilocal, MPIU_INT, 1, nl, ailocal PETSC_F90_2PTR_PARAM(pilocal)); 37 /* this creates a memory leak */ 38 f90arraysfnodecreate_((PetscInt *)iremote, nleaves, airemote PETSC_F90_2PTR_PARAM(piremote)); 39 } 40 41 PETSC_EXTERN void petscsfgetleafranks_(PetscSF *sf, PetscInt *niranks, F90Array1d *airanks, F90Array1d *aioffset, F90Array1d *airootloc, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(piranks) PETSC_F90_2PTR_PROTO(pioffset) PETSC_F90_2PTR_PROTO(pirootloc)) 42 { 43 const PetscMPIInt *iranks = NULL; 44 const PetscInt *ioffset = NULL; 45 const PetscInt *irootloc = NULL; 46 47 *ierr = PetscSFGetLeafRanks(*sf, niranks, &iranks, &ioffset, &irootloc); 48 if (*ierr) return; 49 *ierr = F90Array1dCreate((void *)irootloc, MPIU_INT, 1, ioffset[*niranks], airootloc PETSC_F90_2PTR_PARAM(pirootloc)); 50 if (*ierr) return; 51 *ierr = F90Array1dCreate((void *)iranks, MPI_INT, 1, *niranks, airanks PETSC_F90_2PTR_PARAM(piranks)); 52 if (*ierr) return; 53 *ierr = F90Array1dCreate((void *)ioffset, MPIU_INT, 1, *niranks + 1, aioffset PETSC_F90_2PTR_PARAM(pioffset)); 54 if (*ierr) return; 55 } 56 57 PETSC_EXTERN void petscsfgetrootranks_(PetscSF *sf, PetscInt *nranks, F90Array1d *aranks, F90Array1d *aroffset, F90Array1d *armine, F90Array1d *arremote, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(pranks) PETSC_F90_2PTR_PROTO(proffset) PETSC_F90_2PTR_PROTO(prmine) PETSC_F90_2PTR_PROTO(prremote)) 58 { 59 const PetscMPIInt *ranks = NULL; 60 const PetscInt *roffset = NULL; 61 const PetscInt *rmine = NULL; 62 const PetscInt *rremote = NULL; 63 64 *ierr = PetscSFGetRootRanks(*sf, nranks, &ranks, &roffset, &rmine, &rremote); 65 if (*ierr) return; 66 *ierr = F90Array1dCreate((void *)ranks, MPI_INT, 1, *nranks, aranks PETSC_F90_2PTR_PARAM(pranks)); 67 if (*ierr) return; 68 *ierr = F90Array1dCreate((void *)roffset, MPIU_INT, 1, *nranks + 1, aroffset PETSC_F90_2PTR_PARAM(proffset)); 69 if (*ierr) return; 70 *ierr = F90Array1dCreate((void *)rmine, MPIU_INT, 1, roffset[*nranks], armine PETSC_F90_2PTR_PARAM(prmine)); 71 if (*ierr) return; 72 *ierr = F90Array1dCreate((void *)rremote, MPIU_INT, 1, roffset[*nranks], arremote PETSC_F90_2PTR_PARAM(prremote)); 73 if (*ierr) return; 74 } 75 76 #if defined(PETSC_HAVE_F90_ASSUMED_TYPE_NOT_PTR) 77 PETSC_EXTERN void petscsfbcastbegin_(PetscSF *sf, MPI_Fint *unit, const void *rptr, void *lptr, MPI_Fint *op, PetscErrorCode *ierr) 78 { 79 MPI_Datatype dtype; 80 MPI_Op cop = MPI_Op_f2c(*op); 81 82 *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype); 83 if (*ierr) return; 84 *ierr = PetscSFBcastBegin(*sf, dtype, rptr, lptr, cop); 85 } 86 87 PETSC_EXTERN void petscsfbcastend_(PetscSF *sf, MPI_Fint *unit, const void *rptr, void *lptr, MPI_Fint *op, PetscErrorCode *ierr) 88 { 89 MPI_Datatype dtype; 90 MPI_Op cop = MPI_Op_f2c(*op); 91 92 *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype); 93 if (*ierr) return; 94 *ierr = PetscSFBcastEnd(*sf, dtype, rptr, lptr, cop); 95 } 96 97 PETSC_EXTERN void petscsfreducebegin_(PetscSF *sf, MPI_Fint *unit, const void *lptr, void *rptr, MPI_Fint *op, PetscErrorCode *ierr) 98 { 99 MPI_Datatype dtype; 100 MPI_Op cop = MPI_Op_f2c(*op); 101 102 *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype); 103 if (*ierr) return; 104 *ierr = PetscSFReduceBegin(*sf, dtype, lptr, rptr, cop); 105 } 106 107 PETSC_EXTERN void petscsfreduceend_(PetscSF *sf, MPI_Fint *unit, const void *lptr, void *rptr, MPI_Fint *op, PetscErrorCode *ierr) 108 { 109 MPI_Datatype dtype; 110 MPI_Op cop = MPI_Op_f2c(*op); 111 112 *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype); 113 if (*ierr) return; 114 *ierr = PetscSFReduceEnd(*sf, dtype, lptr, rptr, cop); 115 } 116 117 #else 118 119 PETSC_EXTERN void petscsfbcastbegin_(PetscSF *sf, MPI_Fint *unit, F90Array1d *rptr, F90Array1d *lptr, MPI_Fint *op, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(rptrd) PETSC_F90_2PTR_PROTO(lptrd)) 120 { 121 MPI_Datatype dtype; 122 const void *rootdata; 123 void *leafdata; 124 MPI_Op cop = MPI_Op_f2c(*op); 125 126 *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype); 127 if (*ierr) return; 128 *ierr = F90Array1dAccess(rptr, dtype, (void **)&rootdata PETSC_F90_2PTR_PARAM(rptrd)); 129 if (*ierr) return; 130 *ierr = F90Array1dAccess(lptr, dtype, (void **)&leafdata PETSC_F90_2PTR_PARAM(lptrd)); 131 if (*ierr) return; 132 *ierr = PetscSFBcastBegin(*sf, dtype, rootdata, leafdata, cop); 133 } 134 135 PETSC_EXTERN void petscsfbcastend_(PetscSF *sf, MPI_Fint *unit, F90Array1d *rptr, F90Array1d *lptr, MPI_Fint *op, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(rptrd) PETSC_F90_2PTR_PROTO(lptrd)) 136 { 137 MPI_Datatype dtype; 138 const void *rootdata; 139 void *leafdata; 140 MPI_Op cop = MPI_Op_f2c(*op); 141 142 *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype); 143 if (*ierr) return; 144 *ierr = F90Array1dAccess(rptr, dtype, (void **)&rootdata PETSC_F90_2PTR_PARAM(rptrd)); 145 if (*ierr) return; 146 *ierr = F90Array1dAccess(lptr, dtype, (void **)&leafdata PETSC_F90_2PTR_PARAM(lptrd)); 147 if (*ierr) return; 148 *ierr = PetscSFBcastEnd(*sf, dtype, rootdata, leafdata, cop); 149 } 150 151 PETSC_EXTERN void petscsfreducebegin_(PetscSF *sf, MPI_Fint *unit, F90Array1d *lptr, F90Array1d *rptr, MPI_Fint *op, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(lptrd) PETSC_F90_2PTR_PROTO(rptrd)) 152 { 153 MPI_Datatype dtype; 154 const void *leafdata; 155 void *rootdata; 156 MPI_Op cop = MPI_Op_f2c(*op); 157 158 *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype); 159 if (*ierr) return; 160 *ierr = F90Array1dAccess(rptr, dtype, (void **)&rootdata PETSC_F90_2PTR_PARAM(rptrd)); 161 if (*ierr) return; 162 *ierr = F90Array1dAccess(lptr, dtype, (void **)&leafdata PETSC_F90_2PTR_PARAM(lptrd)); 163 if (*ierr) return; 164 *ierr = PetscSFReduceBegin(*sf, dtype, leafdata, rootdata, cop); 165 } 166 167 PETSC_EXTERN void petscsfreduceend_(PetscSF *sf, MPI_Fint *unit, F90Array1d *lptr, F90Array1d *rptr, MPI_Fint *op, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(lptrd) PETSC_F90_2PTR_PROTO(rptrd)) 168 { 169 MPI_Datatype dtype; 170 const void *leafdata; 171 void *rootdata; 172 MPI_Op cop = MPI_Op_f2c(*op); 173 174 *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype); 175 if (*ierr) return; 176 *ierr = F90Array1dAccess(rptr, dtype, (void **)&rootdata PETSC_F90_2PTR_PARAM(rptrd)); 177 if (*ierr) return; 178 *ierr = F90Array1dAccess(lptr, dtype, (void **)&leafdata PETSC_F90_2PTR_PARAM(lptrd)); 179 if (*ierr) return; 180 *ierr = PetscSFReduceEnd(*sf, dtype, leafdata, rootdata, cop); 181 } 182 #endif 183