Actual source code: zsf.c

  1: #include <petsc/private/f90impl.h>
  2: #include <petsc/private/sfimpl.h>

  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 petscsfreducebegin_     PETSCSFREDUCEBEGIN
 10:   #define petscsfreduceend_       PETSCSFREDUCEEND
 11:   #define f90arraysfnodecreate_   F90ARRAYSFNODECREATE
 12:   #define petscsfviewfromoptions_ PETSCSFVIEWFROMOPTIONS
 13:   #define petscsfdestroy_         PETSCSFDESTROY
 14:   #define petscsfsetgraph_        PETSCSFSETGRAPH
 15:   #define petscsfgetleafranks_    PETSCSFGETLEAFRANKS
 16:   #define petscsfgetrootranks_    PETSCSFGETROOTRANKS
 17: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 18:   #define petscsfgetgraph_        petscsfgetgraph
 19:   #define petscsfview_            petscsfview
 20:   #define petscsfbcastbegin_      petscsfbcastbegin
 21:   #define petscsfbcastend_        petscsfbcastend
 22:   #define petscsfreducebegin_     petscsfreducebegin
 23:   #define petscsfreduceend_       petscsfreduceend
 24:   #define f90arraysfnodecreate_   f90arraysfnodecreate
 25:   #define petscsfviewfromoptions_ petscsfviewfromoptions
 26:   #define petscsfdestroy_         petscsfdestroy
 27:   #define petscsfsetgraph_        petscsfsetgraph
 28:   #define petscsfgetleafranks_    petscsfgetleafranks
 29:   #define petscsfgetrootranks_    petscsfgetrootranks
 30: #endif

 32: PETSC_EXTERN void f90arraysfnodecreate_(const PetscInt *, PetscInt *, void *PETSC_F90_2PTR_PROTO_NOVAR);

 34: PETSC_EXTERN void petscsfsetgraph_(PetscSF *sf, PetscInt *nroots, PetscInt *nleaves, PetscInt *ilocal, PetscCopyMode *localmode, PetscSFNode *iremote, PetscCopyMode *remotemode, int *ierr)
 35: {
 36:   if (ilocal == PETSC_NULL_INTEGER_Fortran) ilocal = NULL;
 37:   *ierr = PetscSFSetGraph(*sf, *nroots, *nleaves, ilocal, *localmode, iremote, *remotemode);
 38: }

 40: PETSC_EXTERN void petscsfview_(PetscSF *sf, PetscViewer *vin, PetscErrorCode *ierr)
 41: {
 42:   PetscViewer v;

 44:   PetscPatchDefaultViewers_Fortran(vin, v);
 45:   *ierr = PetscSFView(*sf, v);
 46: }

 48: 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))
 49: {
 50:   const PetscInt    *ilocal;
 51:   const PetscSFNode *iremote;
 52:   PetscInt           nl;

 54:   *ierr = PetscSFGetGraph(*sf, nroots, nleaves, &ilocal, &iremote);
 55:   if (*ierr) return;
 56:   nl = *nleaves;
 57:   if (!ilocal) nl = 0;
 58:   *ierr = F90Array1dCreate((void *)ilocal, MPIU_INT, 1, nl, ailocal PETSC_F90_2PTR_PARAM(pilocal));
 59:   /* this creates a memory leak */
 60:   f90arraysfnodecreate_((PetscInt *)iremote, nleaves, airemote PETSC_F90_2PTR_PARAM(piremote));
 61: }

 63: 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))
 64: {
 65:   const PetscMPIInt *iranks   = NULL;
 66:   const PetscInt    *ioffset  = NULL;
 67:   const PetscInt    *irootloc = NULL;

 69:   *ierr = PetscSFGetLeafRanks(*sf, niranks, &iranks, &ioffset, &irootloc);
 70:   if (*ierr) return;
 71:   *ierr = F90Array1dCreate((void *)irootloc, MPIU_INT, 1, ioffset[*niranks], airootloc PETSC_F90_2PTR_PARAM(pirootloc));
 72:   if (*ierr) return;
 73:   *ierr = F90Array1dCreate((void *)iranks, MPI_INT, 1, *niranks, airanks PETSC_F90_2PTR_PARAM(piranks));
 74:   if (*ierr) return;
 75:   *ierr = F90Array1dCreate((void *)ioffset, MPIU_INT, 1, *niranks + 1, aioffset PETSC_F90_2PTR_PARAM(pioffset));
 76:   if (*ierr) return;
 77: }

 79: 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))
 80: {
 81:   const PetscMPIInt *ranks   = NULL;
 82:   const PetscInt    *roffset = NULL;
 83:   const PetscInt    *rmine   = NULL;
 84:   const PetscInt    *rremote = NULL;

 86:   *ierr = PetscSFGetRootRanks(*sf, nranks, &ranks, &roffset, &rmine, &rremote);
 87:   if (*ierr) return;
 88:   *ierr = F90Array1dCreate((void *)ranks, MPI_INT, 1, *nranks, aranks PETSC_F90_2PTR_PARAM(pranks));
 89:   if (*ierr) return;
 90:   *ierr = F90Array1dCreate((void *)roffset, MPIU_INT, 1, *nranks + 1, aroffset PETSC_F90_2PTR_PARAM(proffset));
 91:   if (*ierr) return;
 92:   *ierr = F90Array1dCreate((void *)rmine, MPIU_INT, 1, roffset[*nranks], armine PETSC_F90_2PTR_PARAM(prmine));
 93:   if (*ierr) return;
 94:   *ierr = F90Array1dCreate((void *)rremote, MPIU_INT, 1, roffset[*nranks], arremote PETSC_F90_2PTR_PARAM(prremote));
 95:   if (*ierr) return;
 96: }

 98: #if defined(PETSC_HAVE_F90_ASSUMED_TYPE_NOT_PTR)
 99: PETSC_EXTERN void petscsfbcastbegin_(PetscSF *sf, MPI_Fint *unit, const void *rptr, void *lptr, MPI_Fint *op, PetscErrorCode *ierr)
100: {
101:   MPI_Datatype dtype;
102:   MPI_Op       cop = MPI_Op_f2c(*op);

104:   *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype);
105:   if (*ierr) return;
106:   *ierr = PetscSFBcastBegin(*sf, dtype, rptr, lptr, cop);
107: }

109: PETSC_EXTERN void petscsfbcastend_(PetscSF *sf, MPI_Fint *unit, const void *rptr, void *lptr, MPI_Fint *op, PetscErrorCode *ierr)
110: {
111:   MPI_Datatype dtype;
112:   MPI_Op       cop = MPI_Op_f2c(*op);

114:   *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype);
115:   if (*ierr) return;
116:   *ierr = PetscSFBcastEnd(*sf, dtype, rptr, lptr, cop);
117: }

119: PETSC_EXTERN void petscsfreducebegin_(PetscSF *sf, MPI_Fint *unit, const void *lptr, void *rptr, MPI_Fint *op, PetscErrorCode *ierr)
120: {
121:   MPI_Datatype dtype;
122:   MPI_Op       cop = MPI_Op_f2c(*op);

124:   *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype);
125:   if (*ierr) return;
126:   *ierr = PetscSFReduceBegin(*sf, dtype, lptr, rptr, cop);
127: }

129: PETSC_EXTERN void petscsfreduceend_(PetscSF *sf, MPI_Fint *unit, const void *lptr, void *rptr, MPI_Fint *op, PetscErrorCode *ierr)
130: {
131:   MPI_Datatype dtype;
132:   MPI_Op       cop = MPI_Op_f2c(*op);

134:   *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype);
135:   if (*ierr) return;
136:   *ierr = PetscSFReduceEnd(*sf, dtype, lptr, rptr, cop);
137: }

139: #else

141: 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))
142: {
143:   MPI_Datatype dtype;
144:   const void  *rootdata;
145:   void        *leafdata;
146:   MPI_Op       cop = MPI_Op_f2c(*op);

148:   *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype);
149:   if (*ierr) return;
150:   *ierr = F90Array1dAccess(rptr, dtype, (void **)&rootdata PETSC_F90_2PTR_PARAM(rptrd));
151:   if (*ierr) return;
152:   *ierr = F90Array1dAccess(lptr, dtype, (void **)&leafdata PETSC_F90_2PTR_PARAM(lptrd));
153:   if (*ierr) return;
154:   *ierr = PetscSFBcastBegin(*sf, dtype, rootdata, leafdata, cop);
155: }

157: 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))
158: {
159:   MPI_Datatype dtype;
160:   const void  *rootdata;
161:   void        *leafdata;
162:   MPI_Op       cop = MPI_Op_f2c(*op);

164:   *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype);
165:   if (*ierr) return;
166:   *ierr = F90Array1dAccess(rptr, dtype, (void **)&rootdata PETSC_F90_2PTR_PARAM(rptrd));
167:   if (*ierr) return;
168:   *ierr = F90Array1dAccess(lptr, dtype, (void **)&leafdata PETSC_F90_2PTR_PARAM(lptrd));
169:   if (*ierr) return;
170:   *ierr = PetscSFBcastEnd(*sf, dtype, rootdata, leafdata, cop);
171: }

173: 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))
174: {
175:   MPI_Datatype dtype;
176:   const void  *leafdata;
177:   void        *rootdata;
178:   MPI_Op       cop = MPI_Op_f2c(*op);

180:   *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype);
181:   if (*ierr) return;
182:   *ierr = F90Array1dAccess(rptr, dtype, (void **)&rootdata PETSC_F90_2PTR_PARAM(rptrd));
183:   if (*ierr) return;
184:   *ierr = F90Array1dAccess(lptr, dtype, (void **)&leafdata PETSC_F90_2PTR_PARAM(lptrd));
185:   if (*ierr) return;
186:   *ierr = PetscSFReduceBegin(*sf, dtype, leafdata, rootdata, cop);
187: }

189: 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))
190: {
191:   MPI_Datatype dtype;
192:   const void  *leafdata;
193:   void        *rootdata;
194:   MPI_Op       cop = MPI_Op_f2c(*op);

196:   *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype);
197:   if (*ierr) return;
198:   *ierr = F90Array1dAccess(rptr, dtype, (void **)&rootdata PETSC_F90_2PTR_PARAM(rptrd));
199:   if (*ierr) return;
200:   *ierr = F90Array1dAccess(lptr, dtype, (void **)&leafdata PETSC_F90_2PTR_PARAM(lptrd));
201:   if (*ierr) return;
202:   *ierr = PetscSFReduceEnd(*sf, dtype, leafdata, rootdata, cop);
203: }

205: PETSC_EXTERN void petscsfviewfromoptions_(PetscSF *ao, PetscObject obj, char *type, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T len)
206: {
207:   char *t;

209:   FIXCHAR(type, len, t);
210:   CHKFORTRANNULLOBJECT(obj);
211:   *ierr = PetscSFViewFromOptions(*ao, obj, t);
212:   if (*ierr) return;
213:   FREECHAR(type, t);
214: }

216: PETSC_EXTERN void petscsfdestroy_(PetscSF *x, int *ierr)
217: {
218:   PETSC_FORTRAN_OBJECT_F_DESTROYED_TO_C_NULL(x);
219:   *ierr = PetscSFDestroy(x);
220:   if (*ierr) return;
221:   PETSC_FORTRAN_OBJECT_C_NULL_TO_F_DESTROYED(x);
222: }

224: #endif