Actual source code: f90_cwrap.c

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

  3: /*@C

  5:    PetscMPIFortranDatatypeToC - Converts a `MPI_Fint` that contains a Fortran `MPI_Datatype` to its C `MPI_Datatype` equivalent

  7:    Not Collective, No Fortran Support

  9:    Input Parameter:
 10: .  unit - The Fortran `MPI_Datatype`

 12:    Output Parameter:
 13: .  dtype - the corresponding C `MPI_Datatype`

 15:    Level: developer

 17:    Developer Note:
 18:    The MPI documentation in multiple places says that one can never us
 19:    Fortran `MPI_Datatype`s in C (or vice-versa) but this is problematic since users could never
 20:    call C routines from Fortran that have `MPI_Datatype` arguments. Jed states that the Fortran
 21:    `MPI_Datatype`s will always be available in C if the MPI was built to support Fortran. This function
 22:    relies on this.
 23: @*/
 24: PetscErrorCode PetscMPIFortranDatatypeToC(MPI_Fint unit, MPI_Datatype *dtype)
 25: {
 26:   MPI_Datatype ftype;

 28:   PetscFunctionBegin;
 29:   ftype = MPI_Type_f2c(unit);
 30:   if (ftype == MPI_INTEGER || ftype == MPI_INT) *dtype = MPI_INT;
 31:   else if (ftype == MPI_INTEGER8 || ftype == MPIU_INT64) *dtype = MPIU_INT64;
 32:   else if (ftype == MPI_DOUBLE_PRECISION || ftype == MPI_DOUBLE) *dtype = MPI_DOUBLE;
 33: #if defined(PETSC_HAVE_COMPLEX)
 34:   else if (ftype == MPI_COMPLEX16 || ftype == MPI_C_DOUBLE_COMPLEX) *dtype = MPI_C_DOUBLE_COMPLEX;
 35: #endif
 36:   else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unknown Fortran MPI_Datatype");
 37:   PetscFunctionReturn(PETSC_SUCCESS);
 38: }

 40: /*************************************************************************/

 42: #if defined(PETSC_HAVE_FORTRAN_CAPS)
 43:   #define f90array1dcreatescalar_       F90ARRAY1DCREATESCALAR
 44:   #define f90array1daccessscalar_       F90ARRAY1DACCESSSCALAR
 45:   #define f90array1ddestroyscalar_      F90ARRAY1DDESTROYSCALAR
 46:   #define f90array1dcreatereal_         F90ARRAY1DCREATEREAL
 47:   #define f90array1daccessreal_         F90ARRAY1DACCESSREAL
 48:   #define f90array1ddestroyreal_        F90ARRAY1DDESTROYREAL
 49:   #define f90array1dcreateint_          F90ARRAY1DCREATEINT
 50:   #define f90array1daccessint_          F90ARRAY1DACCESSINT
 51:   #define f90array1ddestroyint_         F90ARRAY1DDESTROYINT
 52:   #define f90array1dcreatempiint_       F90ARRAY1DCREATEMPIINT
 53:   #define f90array1daccessmpiint_       F90ARRAY1DACCESSMPIINT
 54:   #define f90array1ddestroympiint_      F90ARRAY1DDESTROYMPIINT
 55:   #define f90array1dcreatefortranaddr_  F90ARRAY1DCREATEFORTRANADDR
 56:   #define f90array1daccessfortranaddr_  F90ARRAY1DACCESSFORTRANADDR
 57:   #define f90array1ddestroyfortranaddr_ F90ARRAY1DDESTROYFORTRANADDR
 58: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 59:   #define f90array1dcreatescalar_       f90array1dcreatescalar
 60:   #define f90array1daccessscalar_       f90array1daccessscalar
 61:   #define f90array1ddestroyscalar_      f90array1ddestroyscalar
 62:   #define f90array1dcreatereal_         f90array1dcreatereal
 63:   #define f90array1daccessreal_         f90array1daccessreal
 64:   #define f90array1ddestroyreal_        f90array1ddestroyreal
 65:   #define f90array1dcreateint_          f90array1dcreateint
 66:   #define f90array1daccessint_          f90array1daccessint
 67:   #define f90array1ddestroyint_         f90array1ddestroyint
 68:   #define f90array1dcreatempiint_       f90array1dcreatempiint
 69:   #define f90array1daccessmpiint_       f90array1daccessmpiint
 70:   #define f90array1ddestroympiint_      f90array1ddestroympiint
 71:   #define f90array1dcreatefortranaddr_  f90array1dcreatefortranaddr
 72:   #define f90array1daccessfortranaddr_  f90array1daccessfortranaddr
 73:   #define f90array1ddestroyfortranaddr_ f90array1ddestroyfortranaddr
 74: #endif

 76: PETSC_EXTERN void f90array1dcreatescalar_(void *, PetscInt *, PetscInt *, F90Array1d *PETSC_F90_2PTR_PROTO_NOVAR);
 77: PETSC_EXTERN void f90array1daccessscalar_(F90Array1d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
 78: PETSC_EXTERN void f90array1ddestroyscalar_(F90Array1d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
 79: PETSC_EXTERN void f90array1dcreatereal_(void *, PetscInt *, PetscInt *, F90Array1d *PETSC_F90_2PTR_PROTO_NOVAR);
 80: PETSC_EXTERN void f90array1daccessreal_(F90Array1d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
 81: PETSC_EXTERN void f90array1ddestroyreal_(F90Array1d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
 82: PETSC_EXTERN void f90array1dcreateint_(void *, PetscInt *, PetscInt *, F90Array1d *PETSC_F90_2PTR_PROTO_NOVAR);
 83: PETSC_EXTERN void f90array1daccessint_(F90Array1d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
 84: PETSC_EXTERN void f90array1ddestroyint_(F90Array1d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
 85: PETSC_EXTERN void f90array1dcreatempiint_(void *, PetscInt *, PetscInt *, F90Array1d *PETSC_F90_2PTR_PROTO_NOVAR);
 86: PETSC_EXTERN void f90array1daccessmpiint_(F90Array1d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
 87: PETSC_EXTERN void f90array1ddestroympiint_(F90Array1d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
 88: PETSC_EXTERN void f90array1dcreatefortranaddr_(void *, PetscInt *, PetscInt *, F90Array1d *PETSC_F90_2PTR_PROTO_NOVAR);
 89: PETSC_EXTERN void f90array1daccessfortranaddr_(F90Array1d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
 90: PETSC_EXTERN void f90array1ddestroyfortranaddr_(F90Array1d *ptr PETSC_F90_2PTR_PROTO_NOVAR);

 92: PetscErrorCode F90Array1dCreate(void *array, MPI_Datatype type, PetscInt start, PetscInt len, F90Array1d *ptr PETSC_F90_2PTR_PROTO(ptrd))
 93: {
 94:   PetscFunctionBegin;
 95:   if (type == MPIU_SCALAR) {
 96:     if (!len) array = PETSC_NULL_SCALAR_Fortran;
 97:     f90array1dcreatescalar_(array, &start, &len, ptr PETSC_F90_2PTR_PARAM(ptrd));
 98:   } else if (type == MPIU_REAL) {
 99:     if (!len) array = PETSC_NULL_REAL_Fortran;
100:     f90array1dcreatereal_(array, &start, &len, ptr PETSC_F90_2PTR_PARAM(ptrd));
101:   } else if (type == MPIU_INT) {
102:     if (!len) array = PETSC_NULL_INTEGER_Fortran;
103:     f90array1dcreateint_(array, &start, &len, ptr PETSC_F90_2PTR_PARAM(ptrd));
104:   } else if (type == MPI_INT) {
105:     /* PETSC_NULL_MPIINT_Fortran is not needed since there is no petsc APIs allowing NULL in place of 'PetscMPIInt *' arguments.
106:        At this line, we only need to assign 'array' a valid address when len is 0, thus PETSC_NULL_INTEGER_Fortran is enough.
107:     */
108:     if (!len) array = PETSC_NULL_INTEGER_Fortran;
109:     f90array1dcreatempiint_(array, &start, &len, ptr PETSC_F90_2PTR_PARAM(ptrd));
110:   } else if (type == MPIU_FORTRANADDR) {
111:     f90array1dcreatefortranaddr_(array, &start, &len, ptr PETSC_F90_2PTR_PARAM(ptrd));
112:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
113:   PetscFunctionReturn(PETSC_SUCCESS);
114: }

116: PetscErrorCode F90Array1dAccess(F90Array1d *ptr, MPI_Datatype type, void **array PETSC_F90_2PTR_PROTO(ptrd))
117: {
118:   PetscFunctionBegin;
119:   if (type == MPIU_SCALAR) {
120:     f90array1daccessscalar_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
121:     if (*array == PETSC_NULL_SCALAR_Fortran) *array = NULL;
122:   } else if (type == MPIU_REAL) {
123:     f90array1daccessreal_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
124:     if (*array == PETSC_NULL_REAL_Fortran) *array = NULL;
125:   } else if (type == MPIU_INT) {
126:     f90array1daccessint_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
127:     if (*array == PETSC_NULL_INTEGER_Fortran) *array = NULL;
128:   } else if (type == MPI_INT) {
129:     f90array1daccessmpiint_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
130:     if (*array == PETSC_NULL_INTEGER_Fortran) *array = NULL;
131:   } else if (type == MPIU_FORTRANADDR) {
132:     f90array1daccessfortranaddr_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
133:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
134:   PetscFunctionReturn(PETSC_SUCCESS);
135: }

137: PetscErrorCode F90Array1dDestroy(F90Array1d *ptr, MPI_Datatype type PETSC_F90_2PTR_PROTO(ptrd))
138: {
139:   PetscFunctionBegin;
140:   if (type == MPIU_SCALAR) {
141:     f90array1ddestroyscalar_(ptr PETSC_F90_2PTR_PARAM(ptrd));
142:   } else if (type == MPIU_REAL) {
143:     f90array1ddestroyreal_(ptr PETSC_F90_2PTR_PARAM(ptrd));
144:   } else if (type == MPIU_INT) {
145:     f90array1ddestroyint_(ptr PETSC_F90_2PTR_PARAM(ptrd));
146:   } else if (type == MPI_INT) {
147:     f90array1ddestroympiint_(ptr PETSC_F90_2PTR_PARAM(ptrd));
148:   } else if (type == MPIU_FORTRANADDR) {
149:     f90array1ddestroyfortranaddr_(ptr PETSC_F90_2PTR_PARAM(ptrd));
150:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
151:   PetscFunctionReturn(PETSC_SUCCESS);
152: }

154: /*************************************************************************/

156: #if defined(PETSC_HAVE_FORTRAN_CAPS)
157:   #define f90array2dcreatescalar_       F90ARRAY2DCREATESCALAR
158:   #define f90array2daccessscalar_       F90ARRAY2DACCESSSCALAR
159:   #define f90array2ddestroyscalar_      F90ARRAY2DDESTROYSCALAR
160:   #define f90array2dcreatereal_         F90ARRAY2DCREATEREAL
161:   #define f90array2daccessreal_         F90ARRAY2DACCESSREAL
162:   #define f90array2ddestroyreal_        F90ARRAY2DDESTROYREAL
163:   #define f90array2dcreateint_          F90ARRAY2DCREATEINT
164:   #define f90array2daccessint_          F90ARRAY2DACCESSINT
165:   #define f90array2ddestroyint_         F90ARRAY2DDESTROYINT
166:   #define f90array2dcreatefortranaddr_  F90ARRAY2DCREATEFORTRANADDR
167:   #define f90array2daccessfortranaddr_  F90ARRAY2DACCESSFORTRANADDR
168:   #define f90array2ddestroyfortranaddr_ F90ARRAY2DDESTROYFORTRANADDR
169: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
170:   #define f90array2dcreatescalar_       f90array2dcreatescalar
171:   #define f90array2daccessscalar_       f90array2daccessscalar
172:   #define f90array2ddestroyscalar_      f90array2ddestroyscalar
173:   #define f90array2dcreatereal_         f90array2dcreatereal
174:   #define f90array2daccessreal_         f90array2daccessreal
175:   #define f90array2ddestroyreal_        f90array2ddestroyreal
176:   #define f90array2dcreateint_          f90array2dcreateint
177:   #define f90array2daccessint_          f90array2daccessint
178:   #define f90array2ddestroyint_         f90array2ddestroyint
179:   #define f90array2dcreatefortranaddr_  f90array2dcreatefortranaddr
180:   #define f90array2daccessfortranaddr_  f90array2daccessfortranaddr
181:   #define f90array2ddestroyfortranaddr_ f90array2ddestroyfortranaddr
182: #endif

184: PETSC_EXTERN void f90array2dcreatescalar_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array2d *PETSC_F90_2PTR_PROTO_NOVAR);
185: PETSC_EXTERN void f90array2daccessscalar_(F90Array2d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
186: PETSC_EXTERN void f90array2ddestroyscalar_(F90Array2d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
187: PETSC_EXTERN void f90array2dcreatereal_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array2d *PETSC_F90_2PTR_PROTO_NOVAR);
188: PETSC_EXTERN void f90array2daccessreal_(F90Array2d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
189: PETSC_EXTERN void f90array2ddestroyreal_(F90Array2d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
190: PETSC_EXTERN void f90array2dcreateint_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array2d *PETSC_F90_2PTR_PROTO_NOVAR);
191: PETSC_EXTERN void f90array2daccessint_(F90Array2d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
192: PETSC_EXTERN void f90array2ddestroyint_(F90Array2d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
193: PETSC_EXTERN void f90array2dcreatefortranaddr_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array2d *PETSC_F90_2PTR_PROTO_NOVAR);
194: PETSC_EXTERN void f90array2daccessfortranaddr_(F90Array2d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
195: PETSC_EXTERN void f90array2ddestroyfortranaddr_(F90Array2d *ptr PETSC_F90_2PTR_PROTO_NOVAR);

197: PetscErrorCode F90Array2dCreate(void *array, MPI_Datatype type, PetscInt start1, PetscInt len1, PetscInt start2, PetscInt len2, F90Array2d *ptr PETSC_F90_2PTR_PROTO(ptrd))
198: {
199:   PetscFunctionBegin;
200:   if (type == MPIU_SCALAR) {
201:     f90array2dcreatescalar_(array, &start1, &len1, &start2, &len2, ptr PETSC_F90_2PTR_PARAM(ptrd));
202:   } else if (type == MPIU_REAL) {
203:     f90array2dcreatereal_(array, &start1, &len1, &start2, &len2, ptr PETSC_F90_2PTR_PARAM(ptrd));
204:   } else if (type == MPIU_INT) {
205:     f90array2dcreateint_(array, &start1, &len1, &start2, &len2, ptr PETSC_F90_2PTR_PARAM(ptrd));
206:   } else if (type == MPIU_FORTRANADDR) {
207:     f90array2dcreatefortranaddr_(array, &start1, &len1, &start2, &len2, ptr PETSC_F90_2PTR_PARAM(ptrd));
208:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
209:   PetscFunctionReturn(PETSC_SUCCESS);
210: }

212: PetscErrorCode F90Array2dAccess(F90Array2d *ptr, MPI_Datatype type, void **array PETSC_F90_2PTR_PROTO(ptrd))
213: {
214:   PetscFunctionBegin;
215:   if (type == MPIU_SCALAR) {
216:     f90array2daccessscalar_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
217:   } else if (type == MPIU_REAL) {
218:     f90array2daccessreal_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
219:   } else if (type == MPIU_INT) {
220:     f90array2daccessint_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
221:   } else if (type == MPIU_FORTRANADDR) {
222:     f90array2daccessfortranaddr_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
223:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
224:   PetscFunctionReturn(PETSC_SUCCESS);
225: }

227: PetscErrorCode F90Array2dDestroy(F90Array2d *ptr, MPI_Datatype type PETSC_F90_2PTR_PROTO(ptrd))
228: {
229:   PetscFunctionBegin;
230:   if (type == MPIU_SCALAR) {
231:     f90array2ddestroyscalar_(ptr PETSC_F90_2PTR_PARAM(ptrd));
232:   } else if (type == MPIU_REAL) {
233:     f90array2ddestroyreal_(ptr PETSC_F90_2PTR_PARAM(ptrd));
234:   } else if (type == MPIU_INT) {
235:     f90array2ddestroyint_(ptr PETSC_F90_2PTR_PARAM(ptrd));
236:   } else if (type == MPIU_FORTRANADDR) {
237:     f90array2ddestroyfortranaddr_(ptr PETSC_F90_2PTR_PARAM(ptrd));
238:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
239:   PetscFunctionReturn(PETSC_SUCCESS);
240: }

242: /*************************************************************************/

244: #if defined(PETSC_HAVE_FORTRAN_CAPS)
245:   #define f90array3dcreatescalar_       F90ARRAY3DCREATESCALAR
246:   #define f90array3daccessscalar_       F90ARRAY3DACCESSSCALAR
247:   #define f90array3ddestroyscalar_      F90ARRAY3DDESTROYSCALAR
248:   #define f90array3dcreatereal_         F90ARRAY3DCREATEREAL
249:   #define f90array3daccessreal_         F90ARRAY3DACCESSREAL
250:   #define f90array3ddestroyreal_        F90ARRAY3DDESTROYREAL
251:   #define f90array3dcreateint_          F90ARRAY3DCREATEINT
252:   #define f90array3daccessint_          F90ARRAY3DACCESSINT
253:   #define f90array3ddestroyint_         F90ARRAY3DDESTROYINT
254:   #define f90array3dcreatefortranaddr_  F90ARRAY3DCREATEFORTRANADDR
255:   #define f90array3daccessfortranaddr_  F90ARRAY3DACCESSFORTRANADDR
256:   #define f90array3ddestroyfortranaddr_ F90ARRAY3DDESTROYFORTRANADDR
257: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
258:   #define f90array3dcreatescalar_       f90array3dcreatescalar
259:   #define f90array3daccessscalar_       f90array3daccessscalar
260:   #define f90array3ddestroyscalar_      f90array3ddestroyscalar
261:   #define f90array3dcreatereal_         f90array3dcreatereal
262:   #define f90array3daccessreal_         f90array3daccessreal
263:   #define f90array3ddestroyreal_        f90array3ddestroyreal
264:   #define f90array3dcreateint_          f90array3dcreateint
265:   #define f90array3daccessint_          f90array3daccessint
266:   #define f90array3ddestroyint_         f90array3ddestroyint
267:   #define f90array3dcreatefortranaddr_  f90array3dcreatefortranaddr
268:   #define f90array3daccessfortranaddr_  f90array3daccessfortranaddr
269:   #define f90array3ddestroyfortranaddr_ f90array3ddestroyfortranaddr
270: #endif

272: PETSC_EXTERN void f90array3dcreatescalar_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array3d *PETSC_F90_2PTR_PROTO_NOVAR);
273: PETSC_EXTERN void f90array3daccessscalar_(F90Array3d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
274: PETSC_EXTERN void f90array3ddestroyscalar_(F90Array3d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
275: PETSC_EXTERN void f90array3dcreatereal_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array3d *PETSC_F90_2PTR_PROTO_NOVAR);
276: PETSC_EXTERN void f90array3daccessreal_(F90Array3d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
277: PETSC_EXTERN void f90array3ddestroyreal_(F90Array3d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
278: PETSC_EXTERN void f90array3dcreateint_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array3d *PETSC_F90_2PTR_PROTO_NOVAR);
279: PETSC_EXTERN void f90array3daccessint_(F90Array3d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
280: PETSC_EXTERN void f90array3ddestroyint_(F90Array3d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
281: PETSC_EXTERN void f90array3dcreatefortranaddr_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array3d *PETSC_F90_2PTR_PROTO_NOVAR);
282: PETSC_EXTERN void f90array3daccessfortranaddr_(F90Array3d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
283: PETSC_EXTERN void f90array3ddestroyfortranaddr_(F90Array3d *ptr PETSC_F90_2PTR_PROTO_NOVAR);

285: PetscErrorCode F90Array3dCreate(void *array, MPI_Datatype type, PetscInt start1, PetscInt len1, PetscInt start2, PetscInt len2, PetscInt start3, PetscInt len3, F90Array3d *ptr PETSC_F90_2PTR_PROTO(ptrd))
286: {
287:   PetscFunctionBegin;
288:   if (type == MPIU_SCALAR) {
289:     f90array3dcreatescalar_(array, &start1, &len1, &start2, &len2, &start3, &len3, ptr PETSC_F90_2PTR_PARAM(ptrd));
290:   } else if (type == MPIU_REAL) {
291:     f90array3dcreatereal_(array, &start1, &len1, &start2, &len2, &start3, &len3, ptr PETSC_F90_2PTR_PARAM(ptrd));
292:   } else if (type == MPIU_INT) {
293:     f90array3dcreateint_(array, &start1, &len1, &start2, &len2, &start3, &len3, ptr PETSC_F90_2PTR_PARAM(ptrd));
294:   } else if (type == MPIU_FORTRANADDR) {
295:     f90array3dcreatefortranaddr_(array, &start1, &len1, &start2, &len2, &start3, &len3, ptr PETSC_F90_2PTR_PARAM(ptrd));
296:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
297:   PetscFunctionReturn(PETSC_SUCCESS);
298: }

300: PetscErrorCode F90Array3dAccess(F90Array3d *ptr, MPI_Datatype type, void **array PETSC_F90_2PTR_PROTO(ptrd))
301: {
302:   PetscFunctionBegin;
303:   if (type == MPIU_SCALAR) {
304:     f90array3daccessscalar_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
305:   } else if (type == MPIU_REAL) {
306:     f90array3daccessreal_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
307:   } else if (type == MPIU_INT) {
308:     f90array3daccessint_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
309:   } else if (type == MPIU_FORTRANADDR) {
310:     f90array3daccessfortranaddr_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
311:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
312:   PetscFunctionReturn(PETSC_SUCCESS);
313: }

315: PetscErrorCode F90Array3dDestroy(F90Array3d *ptr, MPI_Datatype type PETSC_F90_2PTR_PROTO(ptrd))
316: {
317:   PetscFunctionBegin;
318:   if (type == MPIU_SCALAR) {
319:     f90array3ddestroyscalar_(ptr PETSC_F90_2PTR_PARAM(ptrd));
320:   } else if (type == MPIU_REAL) {
321:     f90array3ddestroyreal_(ptr PETSC_F90_2PTR_PARAM(ptrd));
322:   } else if (type == MPIU_INT) {
323:     f90array3ddestroyint_(ptr PETSC_F90_2PTR_PARAM(ptrd));
324:   } else if (type == MPIU_FORTRANADDR) {
325:     f90array3ddestroyfortranaddr_(ptr PETSC_F90_2PTR_PARAM(ptrd));
326:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
327:   PetscFunctionReturn(PETSC_SUCCESS);
328: }

330: /*************************************************************************/
331: #if defined(PETSC_HAVE_FORTRAN_CAPS)
332:   #define f90array4dcreatescalar_       F90ARRAY4DCREATESCALAR
333:   #define f90array4daccessscalar_       F90ARRAY4DACCESSSCALAR
334:   #define f90array4ddestroyscalar_      F90ARRAY4DDESTROYSCALAR
335:   #define f90array4dcreatereal_         F90ARRAY4DCREATEREAL
336:   #define f90array4daccessreal_         F90ARRAY4DACCESSREAL
337:   #define f90array4ddestroyreal_        F90ARRAY4DDESTROYREAL
338:   #define f90array4dcreateint_          F90ARRAY4DCREATEINT
339:   #define f90array4daccessint_          F90ARRAY4DACCESSINT
340:   #define f90array4ddestroyint_         F90ARRAY4DDESTROYINT
341:   #define f90array4dcreatefortranaddr_  F90ARRAY4DCREATEFORTRANADDR
342:   #define f90array4daccessfortranaddr_  F90ARRAY4DACCESSFORTRANADDR
343:   #define f90array4ddestroyfortranaddr_ F90ARRAY4DDESTROYFORTRANADDR
344: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
345:   #define f90array4dcreatescalar_       f90array4dcreatescalar
346:   #define f90array4daccessscalar_       f90array4daccessscalar
347:   #define f90array4ddestroyscalar_      f90array4ddestroyscalar
348:   #define f90array4dcreatereal_         f90array4dcreatereal
349:   #define f90array4daccessreal_         f90array4daccessreal
350:   #define f90array4ddestroyreal_        f90array4ddestroyreal
351:   #define f90array4dcreateint_          f90array4dcreateint
352:   #define f90array4daccessint_          f90array4daccessint
353:   #define f90array4ddestroyint_         f90array4ddestroyint
354:   #define f90array4dcreatefortranaddr_  f90array4dcreatefortranaddr
355:   #define f90array4daccessfortranaddr_  f90array4daccessfortranaddr
356:   #define f90array4ddestroyfortranaddr_ f90array4ddestroyfortranaddr
357: #endif

359: PETSC_EXTERN void f90array4dcreatescalar_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array4d *PETSC_F90_2PTR_PROTO_NOVAR);
360: PETSC_EXTERN void f90array4daccessscalar_(F90Array4d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
361: PETSC_EXTERN void f90array4ddestroyscalar_(F90Array4d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
362: PETSC_EXTERN void f90array4dcreatereal_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array4d *PETSC_F90_2PTR_PROTO_NOVAR);
363: PETSC_EXTERN void f90array4daccessreal_(F90Array4d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
364: PETSC_EXTERN void f90array4ddestroyreal_(F90Array4d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
365: PETSC_EXTERN void f90array4dcreateint_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array4d *PETSC_F90_2PTR_PROTO_NOVAR);
366: PETSC_EXTERN void f90array4daccessint_(F90Array4d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
367: PETSC_EXTERN void f90array4ddestroyint_(F90Array4d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
368: PETSC_EXTERN void f90array4dcreatefortranaddr_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array4d *PETSC_F90_2PTR_PROTO_NOVAR);
369: PETSC_EXTERN void f90array4daccessfortranaddr_(F90Array4d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
370: PETSC_EXTERN void f90array4ddestroyfortranaddr_(F90Array4d *ptr PETSC_F90_2PTR_PROTO_NOVAR);

372: PetscErrorCode F90Array4dCreate(void *array, MPI_Datatype type, PetscInt start1, PetscInt len1, PetscInt start2, PetscInt len2, PetscInt start3, PetscInt len3, PetscInt start4, PetscInt len4, F90Array4d *ptr PETSC_F90_2PTR_PROTO(ptrd))
373: {
374:   PetscFunctionBegin;
375:   if (type == MPIU_SCALAR) {
376:     f90array4dcreatescalar_(array, &start1, &len1, &start2, &len2, &start3, &len3, &start4, &len4, ptr PETSC_F90_2PTR_PARAM(ptrd));
377:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
378:   PetscFunctionReturn(PETSC_SUCCESS);
379: }

381: PetscErrorCode F90Array4dAccess(F90Array4d *ptr, MPI_Datatype type, void **array PETSC_F90_2PTR_PROTO(ptrd))
382: {
383:   PetscFunctionBegin;
384:   if (type == MPIU_SCALAR) {
385:     f90array4daccessscalar_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
386:   } else if (type == MPIU_REAL) {
387:     f90array4daccessreal_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
388:   } else if (type == MPIU_INT) {
389:     f90array4daccessint_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
390:   } else if (type == MPIU_FORTRANADDR) {
391:     f90array4daccessfortranaddr_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
392:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
393:   PetscFunctionReturn(PETSC_SUCCESS);
394: }

396: PetscErrorCode F90Array4dDestroy(F90Array4d *ptr, MPI_Datatype type PETSC_F90_2PTR_PROTO(ptrd))
397: {
398:   PetscFunctionBegin;
399:   if (type == MPIU_SCALAR) {
400:     f90array4ddestroyscalar_(ptr PETSC_F90_2PTR_PARAM(ptrd));
401:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
402:   PetscFunctionReturn(PETSC_SUCCESS);
403: }

405: /*************************************************************************/
406: #if defined(PETSC_HAVE_FORTRAN_CAPS)
407:   #define f90array1dgetaddrscalar_      F90ARRAY1DGETADDRSCALAR
408:   #define f90array1dgetaddrreal_        F90ARRAY1DGETADDRREAL
409:   #define f90array1dgetaddrint_         F90ARRAY1DGETADDRINT
410:   #define f90array1dgetaddrmpiint_      F90ARRAY1DGETADDRMPIINT
411:   #define f90array1dgetaddrfortranaddr_ F90ARRAY1DGETADDRFORTRANADDR
412: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
413:   #define f90array1dgetaddrscalar_      f90array1dgetaddrscalar
414:   #define f90array1dgetaddrreal_        f90array1dgetaddrreal
415:   #define f90array1dgetaddrint_         f90array1dgetaddrint
416:   #define f90array1dgetaddrmpiint_      f90array1dgetaddrmpiint
417:   #define f90array1dgetaddrfortranaddr_ f90array1dgetaddrfortranaddr
418: #endif

420: PETSC_EXTERN void f90array1dgetaddrscalar_(void *array, PetscFortranAddr *address)
421: {
422:   *address = (PetscFortranAddr)array;
423: }
424: PETSC_EXTERN void f90array1dgetaddrreal_(void *array, PetscFortranAddr *address)
425: {
426:   *address = (PetscFortranAddr)array;
427: }
428: PETSC_EXTERN void f90array1dgetaddrint_(void *array, PetscFortranAddr *address)
429: {
430:   *address = (PetscFortranAddr)array;
431: }
432: PETSC_EXTERN void f90array1dgetaddrmpiint_(void *array, PetscFortranAddr *address)
433: {
434:   *address = (PetscFortranAddr)array;
435: }
436: PETSC_EXTERN void f90array1dgetaddrfortranaddr_(void *array, PetscFortranAddr *address)
437: {
438:   *address = (PetscFortranAddr)array;
439: }

441: /*************************************************************************/
442: #if defined(PETSC_HAVE_FORTRAN_CAPS)
443:   #define f90array2dgetaddrscalar_      F90ARRAY2DGETADDRSCALAR
444:   #define f90array2dgetaddrreal_        F90ARRAY2DGETADDRREAL
445:   #define f90array2dgetaddrint_         F90ARRAY2DGETADDRINT
446:   #define f90array2dgetaddrfortranaddr_ F90ARRAY2DGETADDRFORTRANADDR
447: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
448:   #define f90array2dgetaddrscalar_      f90array2dgetaddrscalar
449:   #define f90array2dgetaddrreal_        f90array2dgetaddrreal
450:   #define f90array2dgetaddrint_         f90array2dgetaddrint
451:   #define f90array2dgetaddrfortranaddr_ f90array2dgetaddrfortranaddr
452: #endif

454: PETSC_EXTERN void f90array2dgetaddrscalar_(void *array, PetscFortranAddr *address)
455: {
456:   *address = (PetscFortranAddr)array;
457: }
458: PETSC_EXTERN void f90array2dgetaddrreal_(void *array, PetscFortranAddr *address)
459: {
460:   *address = (PetscFortranAddr)array;
461: }
462: PETSC_EXTERN void f90array2dgetaddrint_(void *array, PetscFortranAddr *address)
463: {
464:   *address = (PetscFortranAddr)array;
465: }
466: PETSC_EXTERN void f90array2dgetaddrfortranaddr_(void *array, PetscFortranAddr *address)
467: {
468:   *address = (PetscFortranAddr)array;
469: }

471: /*************************************************************************/
472: #if defined(PETSC_HAVE_FORTRAN_CAPS)
473:   #define f90array3dgetaddrscalar_      F90ARRAY3DGETADDRSCALAR
474:   #define f90array3dgetaddrreal_        F90ARRAY3DGETADDRREAL
475:   #define f90array3dgetaddrint_         F90ARRAY3DGETADDRINT
476:   #define f90array3dgetaddrfortranaddr_ F90ARRAY3DGETADDRFORTRANADDR
477: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
478:   #define f90array3dgetaddrscalar_      f90array3dgetaddrscalar
479:   #define f90array3dgetaddrreal_        f90array3dgetaddrreal
480:   #define f90array3dgetaddrint_         f90array3dgetaddrint
481:   #define f90array3dgetaddrfortranaddr_ f90array3dgetaddrfortranaddr
482: #endif

484: PETSC_EXTERN void f90array3dgetaddrscalar_(void *array, PetscFortranAddr *address)
485: {
486:   *address = (PetscFortranAddr)array;
487: }
488: PETSC_EXTERN void f90array3dgetaddrreal_(void *array, PetscFortranAddr *address)
489: {
490:   *address = (PetscFortranAddr)array;
491: }
492: PETSC_EXTERN void f90array3dgetaddrint_(void *array, PetscFortranAddr *address)
493: {
494:   *address = (PetscFortranAddr)array;
495: }
496: PETSC_EXTERN void f90array3dgetaddrfortranaddr_(void *array, PetscFortranAddr *address)
497: {
498:   *address = (PetscFortranAddr)array;
499: }

501: /*************************************************************************/
502: #if defined(PETSC_HAVE_FORTRAN_CAPS)
503:   #define f90array4dgetaddrscalar_      F90ARRAY4DGETADDRSCALAR
504:   #define f90array4dgetaddrreal_        F90ARRAY4DGETADDRREAL
505:   #define f90array4dgetaddrint_         F90ARRAY4DGETADDRINT
506:   #define f90array4dgetaddrfortranaddr_ F90ARRAY4DGETADDRFORTRANADDR
507: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
508:   #define f90array4dgetaddrscalar_      f90array4dgetaddrscalar
509:   #define f90array4dgetaddrreal_        f90array4dgetaddrreal
510:   #define f90array4dgetaddrint_         f90array4dgetaddrint
511:   #define f90array4dgetaddrfortranaddr_ f90array4dgetaddrfortranaddr
512: #endif

514: PETSC_EXTERN void f90array4dgetaddrscalar_(void *array, PetscFortranAddr *address)
515: {
516:   *address = (PetscFortranAddr)array;
517: }
518: PETSC_EXTERN void f90array4dgetaddrreal_(void *array, PetscFortranAddr *address)
519: {
520:   *address = (PetscFortranAddr)array;
521: }
522: PETSC_EXTERN void f90array4dgetaddrint_(void *array, PetscFortranAddr *address)
523: {
524:   *address = (PetscFortranAddr)array;
525: }
526: PETSC_EXTERN void f90array4dgetaddrfortranaddr_(void *array, PetscFortranAddr *address)
527: {
528:   *address = (PetscFortranAddr)array;
529: }