Actual source code: zsf.c

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

  4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
  5:   #define petscsfrestoregraph_     PETSCSFRESTOREGRAPH
  6:   #define petscsfgetgraph_         PETSCSFGETGRAPH
  7:   #define petscsfbcastbegin_       PETSCSFBCASTBEGIN
  8:   #define petscsfbcastend_         PETSCSFBCASTEND
  9:   #define petscsfreducebegin_      PETSCSFREDUCEBEGIN
 10:   #define petscsfreduceend_        PETSCSFREDUCEEND
 11:   #define petscsfgetleafranks_     PETSCSFGETLEAFRANKS
 12:   #define petscsfgetrootranks_     PETSCSFGETROOTRANKS
 13:   #define f90array1dcreatesfnode_  F90ARRAY1DCREATESFNODE
 14:   #define f90array1ddestroysfnode_ F90ARRAY1DDESTROYSFNODE
 15: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 16:   #define petscsfrestoregraph_     petscsfrestoregraph
 17:   #define petscsfgetgraph_         petscsfgetgraph
 18:   #define petscsfbcastbegin_       petscsfbcastbegin
 19:   #define petscsfbcastend_         petscsfbcastend
 20:   #define petscsfreducebegin_      petscsfreducebegin
 21:   #define petscsfreduceend_        petscsfreduceend
 22:   #define petscsfgetleafranks_     petscsfgetleafranks
 23:   #define petscsfgetrootranks_     petscsfgetrootranks
 24:   #define f90array1dcreatesfnode_  f90array1dcreatesfnode
 25:   #define f90array1ddestroysfnode_ f90array1ddestroysfnode
 26: #endif

 28: PETSC_EXTERN void f90array1dcreatesfnode_(const PetscSFNode *, PetscInt *, PetscInt *, void *PETSC_F90_2PTR_PROTO_NOVAR);
 29: PETSC_EXTERN void f90array1ddestroysfnode_(void *PETSC_F90_2PTR_PROTO_NOVAR);

 31: 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))
 32: {
 33:   const PetscInt    *ilocal;
 34:   const PetscSFNode *iremote;
 35:   PetscInt           nl, one = 1;

 37:   *ierr = PetscSFGetGraph(*sf, nroots, nleaves, &ilocal, &iremote);
 38:   if (*ierr) return;
 39:   nl = *nleaves;
 40:   if (!ilocal) nl = 0;
 41:   *ierr = F90Array1dCreate((void *)ilocal, MPIU_INT, 1, nl, ailocal PETSC_F90_2PTR_PARAM(pilocal));
 42:   /* this creates a memory leak */
 43:   f90array1dcreatesfnode_(iremote, &one, nleaves, airemote PETSC_F90_2PTR_PARAM(piremote));
 44: }

 46: PETSC_EXTERN void petscsfrestoregraph_(PetscSF *sf, PetscInt *nroots, PetscInt *nleaves, F90Array1d *ailocal, F90Array1d *airemote, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(pilocal) PETSC_F90_2PTR_PROTO(piremote))
 47: {
 48:   *ierr = F90Array1dDestroy(ailocal, MPIU_INT PETSC_F90_2PTR_PARAM(pilocal));
 49:   if (*ierr) return;
 50:   f90array1ddestroysfnode_(airemote PETSC_F90_2PTR_PARAM(piremote));
 51: }

 53: PETSC_EXTERN void petscsfgetleafranks_(PetscSF *sf, PetscMPIInt *niranks, F90Array1d *airanks, F90Array1d *aioffset, F90Array1d *airootloc, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(piranks) PETSC_F90_2PTR_PROTO(pioffset) PETSC_F90_2PTR_PROTO(pirootloc))
 54: {
 55:   const PetscMPIInt *iranks   = NULL;
 56:   const PetscInt    *ioffset  = NULL;
 57:   const PetscInt    *irootloc = NULL;

 59:   *ierr = PetscSFGetLeafRanks(*sf, niranks, &iranks, &ioffset, &irootloc);
 60:   if (*ierr) return;
 61:   *ierr = F90Array1dCreate((void *)irootloc, MPIU_INT, 1, ioffset[*niranks], airootloc PETSC_F90_2PTR_PARAM(pirootloc));
 62:   if (*ierr) return;
 63:   *ierr = F90Array1dCreate((void *)iranks, MPI_INT, 1, *niranks, airanks PETSC_F90_2PTR_PARAM(piranks));
 64:   if (*ierr) return;
 65:   *ierr = F90Array1dCreate((void *)ioffset, MPIU_INT, 1, *niranks + 1, aioffset PETSC_F90_2PTR_PARAM(pioffset));
 66:   if (*ierr) return;
 67: }

 69: PETSC_EXTERN void petscsfgetrootranks_(PetscSF *sf, PetscMPIInt *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))
 70: {
 71:   const PetscMPIInt *ranks   = NULL;
 72:   const PetscInt    *roffset = NULL;
 73:   const PetscInt    *rmine   = NULL;
 74:   const PetscInt    *rremote = NULL;

 76:   *ierr = PetscSFGetRootRanks(*sf, nranks, &ranks, &roffset, &rmine, &rremote);
 77:   if (*ierr) return;
 78:   *ierr = F90Array1dCreate((void *)ranks, MPI_INT, 1, *nranks, aranks PETSC_F90_2PTR_PARAM(pranks));
 79:   if (*ierr) return;
 80:   *ierr = F90Array1dCreate((void *)roffset, MPIU_INT, 1, *nranks + 1, aroffset PETSC_F90_2PTR_PARAM(proffset));
 81:   if (*ierr) return;
 82:   *ierr = F90Array1dCreate((void *)rmine, MPIU_INT, 1, roffset[*nranks], armine PETSC_F90_2PTR_PARAM(prmine));
 83:   if (*ierr) return;
 84:   *ierr = F90Array1dCreate((void *)rremote, MPIU_INT, 1, roffset[*nranks], arremote PETSC_F90_2PTR_PARAM(prremote));
 85:   if (*ierr) return;
 86: }

 88: #if defined(PETSC_HAVE_F90_ASSUMED_TYPE_NOT_PTR)
 89: PETSC_EXTERN void petscsfbcastbegin_(PetscSF *sf, MPI_Fint *unit, const void *rptr, void *lptr, MPI_Fint *op, PetscErrorCode *ierr)
 90: {
 91:   MPI_Datatype dtype;
 92:   MPI_Op       cop = MPI_Op_f2c(*op);

 94:   *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype);
 95:   if (*ierr) return;
 96:   *ierr = PetscSFBcastBegin(*sf, dtype, rptr, lptr, cop);
 97: }

 99: PETSC_EXTERN void petscsfbcastend_(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 = PetscSFBcastEnd(*sf, dtype, rptr, lptr, cop);
107: }

109: PETSC_EXTERN void petscsfreducebegin_(PetscSF *sf, MPI_Fint *unit, const void *lptr, void *rptr, 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 = PetscSFReduceBegin(*sf, dtype, lptr, rptr, cop);
117: }

119: PETSC_EXTERN void petscsfreduceend_(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 = PetscSFReduceEnd(*sf, dtype, lptr, rptr, cop);
127: }

129: #else

131: 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))
132: {
133:   MPI_Datatype dtype;
134:   const void  *rootdata;
135:   void        *leafdata;
136:   MPI_Op       cop = MPI_Op_f2c(*op);

138:   *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype);
139:   if (*ierr) return;
140:   *ierr = F90Array1dAccess(rptr, dtype, (void **)&rootdata PETSC_F90_2PTR_PARAM(rptrd));
141:   if (*ierr) return;
142:   *ierr = F90Array1dAccess(lptr, dtype, (void **)&leafdata PETSC_F90_2PTR_PARAM(lptrd));
143:   if (*ierr) return;
144:   *ierr = PetscSFBcastBegin(*sf, dtype, rootdata, leafdata, cop);
145: }

147: 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))
148: {
149:   MPI_Datatype dtype;
150:   const void  *rootdata;
151:   void        *leafdata;
152:   MPI_Op       cop = MPI_Op_f2c(*op);

154:   *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype);
155:   if (*ierr) return;
156:   *ierr = F90Array1dAccess(rptr, dtype, (void **)&rootdata PETSC_F90_2PTR_PARAM(rptrd));
157:   if (*ierr) return;
158:   *ierr = F90Array1dAccess(lptr, dtype, (void **)&leafdata PETSC_F90_2PTR_PARAM(lptrd));
159:   if (*ierr) return;
160:   *ierr = PetscSFBcastEnd(*sf, dtype, rootdata, leafdata, cop);
161: }

163: 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))
164: {
165:   MPI_Datatype dtype;
166:   const void  *leafdata;
167:   void        *rootdata;
168:   MPI_Op       cop = MPI_Op_f2c(*op);

170:   *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype);
171:   if (*ierr) return;
172:   *ierr = F90Array1dAccess(rptr, dtype, (void **)&rootdata PETSC_F90_2PTR_PARAM(rptrd));
173:   if (*ierr) return;
174:   *ierr = F90Array1dAccess(lptr, dtype, (void **)&leafdata PETSC_F90_2PTR_PARAM(lptrd));
175:   if (*ierr) return;
176:   *ierr = PetscSFReduceBegin(*sf, dtype, leafdata, rootdata, cop);
177: }

179: 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))
180: {
181:   MPI_Datatype dtype;
182:   const void  *leafdata;
183:   void        *rootdata;
184:   MPI_Op       cop = MPI_Op_f2c(*op);

186:   *ierr = PetscMPIFortranDatatypeToC(*unit, &dtype);
187:   if (*ierr) return;
188:   *ierr = F90Array1dAccess(rptr, dtype, (void **)&rootdata PETSC_F90_2PTR_PARAM(rptrd));
189:   if (*ierr) return;
190:   *ierr = F90Array1dAccess(lptr, dtype, (void **)&leafdata PETSC_F90_2PTR_PARAM(lptrd));
191:   if (*ierr) return;
192:   *ierr = PetscSFReduceEnd(*sf, dtype, leafdata, rootdata, cop);
193: }
194: #endif