Actual source code: f90_cwrap.c
1: #include <petsc/private/ftnimpl.h>
3: /*@C
4: PetscMPIFortranDatatypeToC - Converts a `MPI_Fint` that contains a Fortran `MPI_Datatype` to its C `MPI_Datatype` equivalent
6: Not Collective, No Fortran Support
8: Input Parameter:
9: . unit - The Fortran `MPI_Datatype`
11: Output Parameter:
12: . dtype - the corresponding C `MPI_Datatype`
14: Level: developer
16: Developer Note:
17: The MPI documentation in multiple places says that one can never us
18: Fortran `MPI_Datatype`s in C (or vice-versa) but this is problematic since users could never
19: call C routines from Fortran that have `MPI_Datatype` arguments. Jed states that the Fortran
20: `MPI_Datatype`s will always be available in C if the MPI was built to support Fortran. This function
21: relies on this.
23: .seealso: `MPI_Fint`, `MPI_Datatype`
24: @*/
25: PetscErrorCode PetscMPIFortranDatatypeToC(MPI_Fint unit, MPI_Datatype *dtype)
26: {
27: MPI_Datatype ftype;
29: PetscFunctionBegin;
30: ftype = MPI_Type_f2c(unit);
31: if (ftype == MPI_INTEGER || ftype == MPI_INT) *dtype = MPI_INT;
32: else if (ftype == MPI_INTEGER8 || ftype == MPIU_INT64) *dtype = MPIU_INT64;
33: else if (ftype == MPI_DOUBLE_PRECISION || ftype == MPI_DOUBLE) *dtype = MPI_DOUBLE;
34: else if (ftype == MPI_FLOAT) *dtype = MPI_FLOAT;
35: #if defined(PETSC_HAVE_COMPLEX)
36: else if (ftype == MPI_COMPLEX16 || ftype == MPI_C_DOUBLE_COMPLEX) *dtype = MPI_C_DOUBLE_COMPLEX;
37: #endif
38: #if defined(PETSC_HAVE_REAL___FLOAT128)
39: else if (ftype == MPIU___FLOAT128) *dtype = MPIU___FLOAT128;
40: #endif
41: else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unknown Fortran MPI_Datatype");
42: PetscFunctionReturn(PETSC_SUCCESS);
43: }
45: /*************************************************************************/
47: #if defined(PETSC_HAVE_FORTRAN_CAPS)
48: #define f90array1dcreatescalar_ F90ARRAY1DCREATESCALAR
49: #define f90array1daccessscalar_ F90ARRAY1DACCESSSCALAR
50: #define f90array1ddestroyscalar_ F90ARRAY1DDESTROYSCALAR
51: #define f90array1dcreatereal_ F90ARRAY1DCREATEREAL
52: #define f90array1daccessreal_ F90ARRAY1DACCESSREAL
53: #define f90array1ddestroyreal_ F90ARRAY1DDESTROYREAL
54: #define f90array1dcreateint_ F90ARRAY1DCREATEINT
55: #define f90array1daccessint_ F90ARRAY1DACCESSINT
56: #define f90array1ddestroyint_ F90ARRAY1DDESTROYINT
57: #define f90array1dcreatempiint_ F90ARRAY1DCREATEMPIINT
58: #define f90array1daccessmpiint_ F90ARRAY1DACCESSMPIINT
59: #define f90array1ddestroympiint_ F90ARRAY1DDESTROYMPIINT
60: #define f90array1dcreatefortranaddr_ F90ARRAY1DCREATEFORTRANADDR
61: #define f90array1daccessfortranaddr_ F90ARRAY1DACCESSFORTRANADDR
62: #define f90array1ddestroyfortranaddr_ F90ARRAY1DDESTROYFORTRANADDR
63: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
64: #define f90array1dcreatescalar_ f90array1dcreatescalar
65: #define f90array1daccessscalar_ f90array1daccessscalar
66: #define f90array1ddestroyscalar_ f90array1ddestroyscalar
67: #define f90array1dcreatereal_ f90array1dcreatereal
68: #define f90array1daccessreal_ f90array1daccessreal
69: #define f90array1ddestroyreal_ f90array1ddestroyreal
70: #define f90array1dcreateint_ f90array1dcreateint
71: #define f90array1daccessint_ f90array1daccessint
72: #define f90array1ddestroyint_ f90array1ddestroyint
73: #define f90array1dcreatempiint_ f90array1dcreatempiint
74: #define f90array1daccessmpiint_ f90array1daccessmpiint
75: #define f90array1ddestroympiint_ f90array1ddestroympiint
76: #define f90array1dcreatefortranaddr_ f90array1dcreatefortranaddr
77: #define f90array1daccessfortranaddr_ f90array1daccessfortranaddr
78: #define f90array1ddestroyfortranaddr_ f90array1ddestroyfortranaddr
79: #endif
81: PETSC_EXTERN void f90array1dcreatescalar_(void *, PetscInt *, PetscInt *, F90Array1d *PETSC_F90_2PTR_PROTO_NOVAR);
82: PETSC_EXTERN void f90array1daccessscalar_(F90Array1d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
83: PETSC_EXTERN void f90array1ddestroyscalar_(F90Array1d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
84: PETSC_EXTERN void f90array1dcreatereal_(void *, PetscInt *, PetscInt *, F90Array1d *PETSC_F90_2PTR_PROTO_NOVAR);
85: PETSC_EXTERN void f90array1daccessreal_(F90Array1d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
86: PETSC_EXTERN void f90array1ddestroyreal_(F90Array1d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
87: PETSC_EXTERN void f90array1dcreateint_(void *, PetscInt *, PetscInt *, F90Array1d *PETSC_F90_2PTR_PROTO_NOVAR);
88: PETSC_EXTERN void f90array1daccessint_(F90Array1d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
89: PETSC_EXTERN void f90array1ddestroyint_(F90Array1d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
90: PETSC_EXTERN void f90array1dcreatempiint_(void *, PetscInt *, PetscInt *, F90Array1d *PETSC_F90_2PTR_PROTO_NOVAR);
91: PETSC_EXTERN void f90array1daccessmpiint_(F90Array1d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
92: PETSC_EXTERN void f90array1ddestroympiint_(F90Array1d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
93: PETSC_EXTERN void f90array1dcreatefortranaddr_(void *, PetscInt *, PetscInt *, F90Array1d *PETSC_F90_2PTR_PROTO_NOVAR);
94: PETSC_EXTERN void f90array1daccessfortranaddr_(F90Array1d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
95: PETSC_EXTERN void f90array1ddestroyfortranaddr_(F90Array1d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
97: /*@C
98: F90Array1dCreate - given a `F90Array1d` passed from Fortran associate with it a C array, its starting index and length
100: Not Collective, No Fortran Support
102: Input Parameters:
103: + array - the C address pointer
104: . type - the MPI datatype of the array
105: . start - the first index of the array
106: . len - the length of the array
107: . ptr - the `F90Array1d` passed from Fortran
108: - ptrd - an extra pointer passed by some Fortran compilers
110: Level: developer
112: Developer Notes:
113: This is used in PETSc Fortran stubs that are used to pass C arrays to Fortran, for example `VecGetArray()`
115: This doesn't actually create the `F90Array1d()`, it just associates a C pointer with it.
117: There are equivalent routines for 2, 3, and 4 dimensional Fortran arrays.
119: .seealso: `F90Array1d`, `F90Array1dAccess()`, `F90Array1dDestroy()`, `F90Array2dCreate()`, `F90Array2dAccess()`, `F90Array2dDestroy()`
120: @*/
121: PetscErrorCode F90Array1dCreate(void *array, MPI_Datatype type, PetscInt start, PetscInt len, F90Array1d *ptr PETSC_F90_2PTR_PROTO(ptrd))
122: {
123: PetscFunctionBegin;
124: if (type == MPIU_SCALAR) {
125: if (!len) array = PETSC_NULL_SCALAR_Fortran;
126: f90array1dcreatescalar_(array, &start, &len, ptr PETSC_F90_2PTR_PARAM(ptrd));
127: } else if (type == MPIU_REAL) {
128: if (!len) array = PETSC_NULL_REAL_Fortran;
129: f90array1dcreatereal_(array, &start, &len, ptr PETSC_F90_2PTR_PARAM(ptrd));
130: } else if (type == MPIU_INT) {
131: if (!len) array = PETSC_NULL_INTEGER_Fortran;
132: f90array1dcreateint_(array, &start, &len, ptr PETSC_F90_2PTR_PARAM(ptrd));
133: } else if (type == MPI_INT) {
134: /* PETSC_NULL_MPIINT_Fortran is not needed since there is no PETSc APIs allowing NULL in place of 'PetscMPIInt *' arguments.
135: At this line, we only need to assign 'array' a valid address when len is 0, thus PETSC_NULL_INTEGER_Fortran is enough.
136: */
137: if (!len) array = PETSC_NULL_INTEGER_Fortran;
138: f90array1dcreatempiint_(array, &start, &len, ptr PETSC_F90_2PTR_PARAM(ptrd));
139: } else if (type == MPIU_FORTRANADDR) {
140: f90array1dcreatefortranaddr_(array, &start, &len, ptr PETSC_F90_2PTR_PARAM(ptrd));
141: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
142: PetscFunctionReturn(PETSC_SUCCESS);
143: }
145: /*@C
146: F90Array1dAccess - given a `F90Array1d` passed from Fortran, accesses from it the associate C array that was provided with `F90Array1dCreate()`
148: Not Collective, No Fortran Support
150: Input Parameters:
151: + ptr - the `F90Array1d` passed from Fortran
152: . type - the MPI datatype of the array
153: - ptrd - an extra pointer passed by some Fortran compilers
155: Output Parameter:
156: . array - the C address pointer
158: Level: developer
160: Developer Note:
161: This is used in PETSc Fortran stubs that access C arrays inside Fortran pointer arrays to Fortran. It is usually used in `XXXRestore()`` Fortran stubs.
163: There are equivalent routines for 2, 3, and 4 dimensional Fortran arrays.
165: .seealso: `F90Array1d`, `F90Array1dCreate()`, `F90Array1dDestroy()`, `F90Array2dCreate()`, `F90Array2dAccess()`, `F90Array2dDestroy()`
166: @*/
167: PetscErrorCode F90Array1dAccess(F90Array1d *ptr, MPI_Datatype type, void **array PETSC_F90_2PTR_PROTO(ptrd))
168: {
169: PetscFunctionBegin;
170: if (type == MPIU_SCALAR) {
171: f90array1daccessscalar_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
172: if (*array == PETSC_NULL_SCALAR_Fortran) *array = NULL;
173: } else if (type == MPIU_REAL) {
174: f90array1daccessreal_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
175: if (*array == PETSC_NULL_REAL_Fortran) *array = NULL;
176: } else if (type == MPIU_INT) {
177: f90array1daccessint_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
178: if (*array == PETSC_NULL_INTEGER_Fortran) *array = NULL;
179: } else if (type == MPI_INT) {
180: f90array1daccessmpiint_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
181: if (*array == PETSC_NULL_INTEGER_Fortran) *array = NULL;
182: } else if (type == MPIU_FORTRANADDR) {
183: f90array1daccessfortranaddr_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
184: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
185: PetscFunctionReturn(PETSC_SUCCESS);
186: }
188: /*@C
189: F90Array1dDestroy - given a `F90Array1d` passed from Fortran removes the C array associate with it with `F90Array1dCreate()`
191: Not Collective, No Fortran Support
193: Input Parameters:
194: + ptr - the `F90Array1d` passed from Fortran
195: . type - the MPI datatype of the array
196: - ptrd - an extra pointer passed by some Fortran compilers
198: Level: developer
200: Developer Notes:
201: This is used in PETSc Fortran stubs that are used to end access to C arrays from Fortran, for example `VecRestoreArray()`
203: This doesn't actually destroy the `F90Array1d()`, it just removes the associated C pointer from it.
205: There are equivalent routines for 2, 3, and 4 dimensional Fortran arrays.
207: .seealso: `F90Array1d`, `F90Array1dAccess()`, `F90Array1dCreate()`, `F90Array2dCreate()`, `F90Array2dAccess()`, `F90Array2dDestroy()`
208: @*/
209: PetscErrorCode F90Array1dDestroy(F90Array1d *ptr, MPI_Datatype type PETSC_F90_2PTR_PROTO(ptrd))
210: {
211: PetscFunctionBegin;
212: if (type == MPIU_SCALAR) {
213: f90array1ddestroyscalar_(ptr PETSC_F90_2PTR_PARAM(ptrd));
214: } else if (type == MPIU_REAL) {
215: f90array1ddestroyreal_(ptr PETSC_F90_2PTR_PARAM(ptrd));
216: } else if (type == MPIU_INT) {
217: f90array1ddestroyint_(ptr PETSC_F90_2PTR_PARAM(ptrd));
218: } else if (type == MPI_INT) {
219: f90array1ddestroympiint_(ptr PETSC_F90_2PTR_PARAM(ptrd));
220: } else if (type == MPIU_FORTRANADDR) {
221: f90array1ddestroyfortranaddr_(ptr PETSC_F90_2PTR_PARAM(ptrd));
222: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
223: *(void **)ptr = (void *)1;
224: PetscFunctionReturn(PETSC_SUCCESS);
225: }
227: /*MC
228: F90Array1d - a PETSc C representation of a Fortran `XXX, pointer :: array(:)` object
230: Not Collective, No Fortran Support
232: Level: developer
234: Developer Notes:
235: This is used in PETSc Fortran stubs that are used to control access to C arrays from Fortran, for example `VecGetArray()`
237: PETSc does not require any information about the format of this object, all operations on the object are performed by calling Fortran routines.
239: There are equivalent objects for 2, 3, and 4 dimensional Fortran arrays.
241: .seealso: `F90Array1dAccess()`, `F90Array1dCreate()`, , `F90Array1dDestroy()`, `F90Array2dCreate()`, `F90Array2dAccess()`, `F90Array2dDestroy()`
242: M*/
244: #if defined(PETSC_HAVE_FORTRAN_CAPS)
245: #define f90array2dcreatescalar_ F90ARRAY2DCREATESCALAR
246: #define f90array2daccessscalar_ F90ARRAY2DACCESSSCALAR
247: #define f90array2ddestroyscalar_ F90ARRAY2DDESTROYSCALAR
248: #define f90array2dcreatereal_ F90ARRAY2DCREATEREAL
249: #define f90array2daccessreal_ F90ARRAY2DACCESSREAL
250: #define f90array2ddestroyreal_ F90ARRAY2DDESTROYREAL
251: #define f90array2dcreateint_ F90ARRAY2DCREATEINT
252: #define f90array2daccessint_ F90ARRAY2DACCESSINT
253: #define f90array2ddestroyint_ F90ARRAY2DDESTROYINT
254: #define f90array2dcreatefortranaddr_ F90ARRAY2DCREATEFORTRANADDR
255: #define f90array2daccessfortranaddr_ F90ARRAY2DACCESSFORTRANADDR
256: #define f90array2ddestroyfortranaddr_ F90ARRAY2DDESTROYFORTRANADDR
257: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
258: #define f90array2dcreatescalar_ f90array2dcreatescalar
259: #define f90array2daccessscalar_ f90array2daccessscalar
260: #define f90array2ddestroyscalar_ f90array2ddestroyscalar
261: #define f90array2dcreatereal_ f90array2dcreatereal
262: #define f90array2daccessreal_ f90array2daccessreal
263: #define f90array2ddestroyreal_ f90array2ddestroyreal
264: #define f90array2dcreateint_ f90array2dcreateint
265: #define f90array2daccessint_ f90array2daccessint
266: #define f90array2ddestroyint_ f90array2ddestroyint
267: #define f90array2dcreatefortranaddr_ f90array2dcreatefortranaddr
268: #define f90array2daccessfortranaddr_ f90array2daccessfortranaddr
269: #define f90array2ddestroyfortranaddr_ f90array2ddestroyfortranaddr
270: #endif
272: PETSC_EXTERN void f90array2dcreatescalar_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array2d *PETSC_F90_2PTR_PROTO_NOVAR);
273: PETSC_EXTERN void f90array2daccessscalar_(F90Array2d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
274: PETSC_EXTERN void f90array2ddestroyscalar_(F90Array2d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
275: PETSC_EXTERN void f90array2dcreatereal_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array2d *PETSC_F90_2PTR_PROTO_NOVAR);
276: PETSC_EXTERN void f90array2daccessreal_(F90Array2d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
277: PETSC_EXTERN void f90array2ddestroyreal_(F90Array2d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
278: PETSC_EXTERN void f90array2dcreateint_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array2d *PETSC_F90_2PTR_PROTO_NOVAR);
279: PETSC_EXTERN void f90array2daccessint_(F90Array2d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
280: PETSC_EXTERN void f90array2ddestroyint_(F90Array2d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
281: PETSC_EXTERN void f90array2dcreatefortranaddr_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array2d *PETSC_F90_2PTR_PROTO_NOVAR);
282: PETSC_EXTERN void f90array2daccessfortranaddr_(F90Array2d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
283: PETSC_EXTERN void f90array2ddestroyfortranaddr_(F90Array2d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
285: PetscErrorCode F90Array2dCreate(void *array, MPI_Datatype type, PetscInt start1, PetscInt len1, PetscInt start2, PetscInt len2, F90Array2d *ptr PETSC_F90_2PTR_PROTO(ptrd))
286: {
287: PetscFunctionBegin;
288: if (type == MPIU_SCALAR) {
289: f90array2dcreatescalar_(array, &start1, &len1, &start2, &len2, ptr PETSC_F90_2PTR_PARAM(ptrd));
290: } else if (type == MPIU_REAL) {
291: f90array2dcreatereal_(array, &start1, &len1, &start2, &len2, ptr PETSC_F90_2PTR_PARAM(ptrd));
292: } else if (type == MPIU_INT) {
293: f90array2dcreateint_(array, &start1, &len1, &start2, &len2, ptr PETSC_F90_2PTR_PARAM(ptrd));
294: } else if (type == MPIU_FORTRANADDR) {
295: f90array2dcreatefortranaddr_(array, &start1, &len1, &start2, &len2, 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 F90Array2dAccess(F90Array2d *ptr, MPI_Datatype type, void **array PETSC_F90_2PTR_PROTO(ptrd))
301: {
302: PetscFunctionBegin;
303: if (type == MPIU_SCALAR) {
304: f90array2daccessscalar_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
305: } else if (type == MPIU_REAL) {
306: f90array2daccessreal_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
307: } else if (type == MPIU_INT) {
308: f90array2daccessint_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
309: } else if (type == MPIU_FORTRANADDR) {
310: f90array2daccessfortranaddr_(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 F90Array2dDestroy(F90Array2d *ptr, MPI_Datatype type PETSC_F90_2PTR_PROTO(ptrd))
316: {
317: PetscFunctionBegin;
318: if (type == MPIU_SCALAR) {
319: f90array2ddestroyscalar_(ptr PETSC_F90_2PTR_PARAM(ptrd));
320: } else if (type == MPIU_REAL) {
321: f90array2ddestroyreal_(ptr PETSC_F90_2PTR_PARAM(ptrd));
322: } else if (type == MPIU_INT) {
323: f90array2ddestroyint_(ptr PETSC_F90_2PTR_PARAM(ptrd));
324: } else if (type == MPIU_FORTRANADDR) {
325: f90array2ddestroyfortranaddr_(ptr PETSC_F90_2PTR_PARAM(ptrd));
326: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
327: PetscFunctionReturn(PETSC_SUCCESS);
328: }
330: #if defined(PETSC_HAVE_FORTRAN_CAPS)
331: #define f90array3dcreatescalar_ F90ARRAY3DCREATESCALAR
332: #define f90array3daccessscalar_ F90ARRAY3DACCESSSCALAR
333: #define f90array3ddestroyscalar_ F90ARRAY3DDESTROYSCALAR
334: #define f90array3dcreatereal_ F90ARRAY3DCREATEREAL
335: #define f90array3daccessreal_ F90ARRAY3DACCESSREAL
336: #define f90array3ddestroyreal_ F90ARRAY3DDESTROYREAL
337: #define f90array3dcreateint_ F90ARRAY3DCREATEINT
338: #define f90array3daccessint_ F90ARRAY3DACCESSINT
339: #define f90array3ddestroyint_ F90ARRAY3DDESTROYINT
340: #define f90array3dcreatefortranaddr_ F90ARRAY3DCREATEFORTRANADDR
341: #define f90array3daccessfortranaddr_ F90ARRAY3DACCESSFORTRANADDR
342: #define f90array3ddestroyfortranaddr_ F90ARRAY3DDESTROYFORTRANADDR
343: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
344: #define f90array3dcreatescalar_ f90array3dcreatescalar
345: #define f90array3daccessscalar_ f90array3daccessscalar
346: #define f90array3ddestroyscalar_ f90array3ddestroyscalar
347: #define f90array3dcreatereal_ f90array3dcreatereal
348: #define f90array3daccessreal_ f90array3daccessreal
349: #define f90array3ddestroyreal_ f90array3ddestroyreal
350: #define f90array3dcreateint_ f90array3dcreateint
351: #define f90array3daccessint_ f90array3daccessint
352: #define f90array3ddestroyint_ f90array3ddestroyint
353: #define f90array3dcreatefortranaddr_ f90array3dcreatefortranaddr
354: #define f90array3daccessfortranaddr_ f90array3daccessfortranaddr
355: #define f90array3ddestroyfortranaddr_ f90array3ddestroyfortranaddr
356: #endif
358: PETSC_EXTERN void f90array3dcreatescalar_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array3d *PETSC_F90_2PTR_PROTO_NOVAR);
359: PETSC_EXTERN void f90array3daccessscalar_(F90Array3d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
360: PETSC_EXTERN void f90array3ddestroyscalar_(F90Array3d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
361: PETSC_EXTERN void f90array3dcreatereal_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array3d *PETSC_F90_2PTR_PROTO_NOVAR);
362: PETSC_EXTERN void f90array3daccessreal_(F90Array3d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
363: PETSC_EXTERN void f90array3ddestroyreal_(F90Array3d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
364: PETSC_EXTERN void f90array3dcreateint_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array3d *PETSC_F90_2PTR_PROTO_NOVAR);
365: PETSC_EXTERN void f90array3daccessint_(F90Array3d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
366: PETSC_EXTERN void f90array3ddestroyint_(F90Array3d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
367: PETSC_EXTERN void f90array3dcreatefortranaddr_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array3d *PETSC_F90_2PTR_PROTO_NOVAR);
368: PETSC_EXTERN void f90array3daccessfortranaddr_(F90Array3d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
369: PETSC_EXTERN void f90array3ddestroyfortranaddr_(F90Array3d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
371: 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))
372: {
373: PetscFunctionBegin;
374: if (type == MPIU_SCALAR) {
375: f90array3dcreatescalar_(array, &start1, &len1, &start2, &len2, &start3, &len3, ptr PETSC_F90_2PTR_PARAM(ptrd));
376: } else if (type == MPIU_REAL) {
377: f90array3dcreatereal_(array, &start1, &len1, &start2, &len2, &start3, &len3, ptr PETSC_F90_2PTR_PARAM(ptrd));
378: } else if (type == MPIU_INT) {
379: f90array3dcreateint_(array, &start1, &len1, &start2, &len2, &start3, &len3, ptr PETSC_F90_2PTR_PARAM(ptrd));
380: } else if (type == MPIU_FORTRANADDR) {
381: f90array3dcreatefortranaddr_(array, &start1, &len1, &start2, &len2, &start3, &len3, ptr PETSC_F90_2PTR_PARAM(ptrd));
382: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
383: PetscFunctionReturn(PETSC_SUCCESS);
384: }
386: PetscErrorCode F90Array3dAccess(F90Array3d *ptr, MPI_Datatype type, void **array PETSC_F90_2PTR_PROTO(ptrd))
387: {
388: PetscFunctionBegin;
389: if (type == MPIU_SCALAR) {
390: f90array3daccessscalar_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
391: } else if (type == MPIU_REAL) {
392: f90array3daccessreal_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
393: } else if (type == MPIU_INT) {
394: f90array3daccessint_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
395: } else if (type == MPIU_FORTRANADDR) {
396: f90array3daccessfortranaddr_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
397: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
398: PetscFunctionReturn(PETSC_SUCCESS);
399: }
401: PetscErrorCode F90Array3dDestroy(F90Array3d *ptr, MPI_Datatype type PETSC_F90_2PTR_PROTO(ptrd))
402: {
403: PetscFunctionBegin;
404: if (type == MPIU_SCALAR) {
405: f90array3ddestroyscalar_(ptr PETSC_F90_2PTR_PARAM(ptrd));
406: } else if (type == MPIU_REAL) {
407: f90array3ddestroyreal_(ptr PETSC_F90_2PTR_PARAM(ptrd));
408: } else if (type == MPIU_INT) {
409: f90array3ddestroyint_(ptr PETSC_F90_2PTR_PARAM(ptrd));
410: } else if (type == MPIU_FORTRANADDR) {
411: f90array3ddestroyfortranaddr_(ptr PETSC_F90_2PTR_PARAM(ptrd));
412: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
413: PetscFunctionReturn(PETSC_SUCCESS);
414: }
416: #if defined(PETSC_HAVE_FORTRAN_CAPS)
417: #define f90array4dcreatescalar_ F90ARRAY4DCREATESCALAR
418: #define f90array4daccessscalar_ F90ARRAY4DACCESSSCALAR
419: #define f90array4ddestroyscalar_ F90ARRAY4DDESTROYSCALAR
420: #define f90array4dcreatereal_ F90ARRAY4DCREATEREAL
421: #define f90array4daccessreal_ F90ARRAY4DACCESSREAL
422: #define f90array4ddestroyreal_ F90ARRAY4DDESTROYREAL
423: #define f90array4dcreateint_ F90ARRAY4DCREATEINT
424: #define f90array4daccessint_ F90ARRAY4DACCESSINT
425: #define f90array4ddestroyint_ F90ARRAY4DDESTROYINT
426: #define f90array4dcreatefortranaddr_ F90ARRAY4DCREATEFORTRANADDR
427: #define f90array4daccessfortranaddr_ F90ARRAY4DACCESSFORTRANADDR
428: #define f90array4ddestroyfortranaddr_ F90ARRAY4DDESTROYFORTRANADDR
429: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
430: #define f90array4dcreatescalar_ f90array4dcreatescalar
431: #define f90array4daccessscalar_ f90array4daccessscalar
432: #define f90array4ddestroyscalar_ f90array4ddestroyscalar
433: #define f90array4dcreatereal_ f90array4dcreatereal
434: #define f90array4daccessreal_ f90array4daccessreal
435: #define f90array4ddestroyreal_ f90array4ddestroyreal
436: #define f90array4dcreateint_ f90array4dcreateint
437: #define f90array4daccessint_ f90array4daccessint
438: #define f90array4ddestroyint_ f90array4ddestroyint
439: #define f90array4dcreatefortranaddr_ f90array4dcreatefortranaddr
440: #define f90array4daccessfortranaddr_ f90array4daccessfortranaddr
441: #define f90array4ddestroyfortranaddr_ f90array4ddestroyfortranaddr
442: #endif
444: PETSC_EXTERN void f90array4dcreatescalar_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array4d *PETSC_F90_2PTR_PROTO_NOVAR);
445: PETSC_EXTERN void f90array4daccessscalar_(F90Array4d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
446: PETSC_EXTERN void f90array4ddestroyscalar_(F90Array4d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
447: PETSC_EXTERN void f90array4dcreatereal_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array4d *PETSC_F90_2PTR_PROTO_NOVAR);
448: PETSC_EXTERN void f90array4daccessreal_(F90Array4d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
449: PETSC_EXTERN void f90array4ddestroyreal_(F90Array4d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
450: PETSC_EXTERN void f90array4dcreateint_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array4d *PETSC_F90_2PTR_PROTO_NOVAR);
451: PETSC_EXTERN void f90array4daccessint_(F90Array4d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
452: PETSC_EXTERN void f90array4ddestroyint_(F90Array4d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
453: PETSC_EXTERN void f90array4dcreatefortranaddr_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array4d *PETSC_F90_2PTR_PROTO_NOVAR);
454: PETSC_EXTERN void f90array4daccessfortranaddr_(F90Array4d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
455: PETSC_EXTERN void f90array4ddestroyfortranaddr_(F90Array4d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
457: 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))
458: {
459: PetscFunctionBegin;
460: PetscCheck(type == MPIU_SCALAR, PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
461: f90array4dcreatescalar_(array, &start1, &len1, &start2, &len2, &start3, &len3, &start4, &len4, ptr PETSC_F90_2PTR_PARAM(ptrd));
462: PetscFunctionReturn(PETSC_SUCCESS);
463: }
465: PetscErrorCode F90Array4dAccess(F90Array4d *ptr, MPI_Datatype type, void **array PETSC_F90_2PTR_PROTO(ptrd))
466: {
467: PetscFunctionBegin;
468: if (type == MPIU_SCALAR) {
469: f90array4daccessscalar_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
470: } else if (type == MPIU_REAL) {
471: f90array4daccessreal_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
472: } else if (type == MPIU_INT) {
473: f90array4daccessint_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
474: } else if (type == MPIU_FORTRANADDR) {
475: f90array4daccessfortranaddr_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
476: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
477: PetscFunctionReturn(PETSC_SUCCESS);
478: }
480: PetscErrorCode F90Array4dDestroy(F90Array4d *ptr, MPI_Datatype type PETSC_F90_2PTR_PROTO(ptrd))
481: {
482: PetscFunctionBegin;
483: PetscCheck(type == MPIU_SCALAR, PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
484: f90array4ddestroyscalar_(ptr PETSC_F90_2PTR_PARAM(ptrd));
485: PetscFunctionReturn(PETSC_SUCCESS);
486: }
488: #if defined(PETSC_HAVE_FORTRAN_CAPS)
489: #define f90array1dgetaddrscalar_ F90ARRAY1DGETADDRSCALAR
490: #define f90array1dgetaddrreal_ F90ARRAY1DGETADDRREAL
491: #define f90array1dgetaddrint_ F90ARRAY1DGETADDRINT
492: #define f90array1dgetaddrmpiint_ F90ARRAY1DGETADDRMPIINT
493: #define f90array1dgetaddrfortranaddr_ F90ARRAY1DGETADDRFORTRANADDR
494: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
495: #define f90array1dgetaddrscalar_ f90array1dgetaddrscalar
496: #define f90array1dgetaddrreal_ f90array1dgetaddrreal
497: #define f90array1dgetaddrint_ f90array1dgetaddrint
498: #define f90array1dgetaddrmpiint_ f90array1dgetaddrmpiint
499: #define f90array1dgetaddrfortranaddr_ f90array1dgetaddrfortranaddr
500: #endif
502: PETSC_EXTERN void f90array1dgetaddrscalar_(void *array, PetscFortranAddr *address)
503: {
504: *address = (PetscFortranAddr)array;
505: }
506: PETSC_EXTERN void f90array1dgetaddrreal_(void *array, PetscFortranAddr *address)
507: {
508: *address = (PetscFortranAddr)array;
509: }
510: PETSC_EXTERN void f90array1dgetaddrint_(void *array, PetscFortranAddr *address)
511: {
512: *address = (PetscFortranAddr)array;
513: }
514: PETSC_EXTERN void f90array1dgetaddrmpiint_(void *array, PetscFortranAddr *address)
515: {
516: *address = (PetscFortranAddr)array;
517: }
518: PETSC_EXTERN void f90array1dgetaddrfortranaddr_(void *array, PetscFortranAddr *address)
519: {
520: *address = (PetscFortranAddr)array;
521: }
523: #if defined(PETSC_HAVE_FORTRAN_CAPS)
524: #define f90array2dgetaddrscalar_ F90ARRAY2DGETADDRSCALAR
525: #define f90array2dgetaddrreal_ F90ARRAY2DGETADDRREAL
526: #define f90array2dgetaddrint_ F90ARRAY2DGETADDRINT
527: #define f90array2dgetaddrfortranaddr_ F90ARRAY2DGETADDRFORTRANADDR
528: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
529: #define f90array2dgetaddrscalar_ f90array2dgetaddrscalar
530: #define f90array2dgetaddrreal_ f90array2dgetaddrreal
531: #define f90array2dgetaddrint_ f90array2dgetaddrint
532: #define f90array2dgetaddrfortranaddr_ f90array2dgetaddrfortranaddr
533: #endif
535: PETSC_EXTERN void f90array2dgetaddrscalar_(void *array, PetscFortranAddr *address)
536: {
537: *address = (PetscFortranAddr)array;
538: }
539: PETSC_EXTERN void f90array2dgetaddrreal_(void *array, PetscFortranAddr *address)
540: {
541: *address = (PetscFortranAddr)array;
542: }
543: PETSC_EXTERN void f90array2dgetaddrint_(void *array, PetscFortranAddr *address)
544: {
545: *address = (PetscFortranAddr)array;
546: }
547: PETSC_EXTERN void f90array2dgetaddrfortranaddr_(void *array, PetscFortranAddr *address)
548: {
549: *address = (PetscFortranAddr)array;
550: }
552: #if defined(PETSC_HAVE_FORTRAN_CAPS)
553: #define f90array3dgetaddrscalar_ F90ARRAY3DGETADDRSCALAR
554: #define f90array3dgetaddrreal_ F90ARRAY3DGETADDRREAL
555: #define f90array3dgetaddrint_ F90ARRAY3DGETADDRINT
556: #define f90array3dgetaddrfortranaddr_ F90ARRAY3DGETADDRFORTRANADDR
557: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
558: #define f90array3dgetaddrscalar_ f90array3dgetaddrscalar
559: #define f90array3dgetaddrreal_ f90array3dgetaddrreal
560: #define f90array3dgetaddrint_ f90array3dgetaddrint
561: #define f90array3dgetaddrfortranaddr_ f90array3dgetaddrfortranaddr
562: #endif
564: PETSC_EXTERN void f90array3dgetaddrscalar_(void *array, PetscFortranAddr *address)
565: {
566: *address = (PetscFortranAddr)array;
567: }
568: PETSC_EXTERN void f90array3dgetaddrreal_(void *array, PetscFortranAddr *address)
569: {
570: *address = (PetscFortranAddr)array;
571: }
572: PETSC_EXTERN void f90array3dgetaddrint_(void *array, PetscFortranAddr *address)
573: {
574: *address = (PetscFortranAddr)array;
575: }
576: PETSC_EXTERN void f90array3dgetaddrfortranaddr_(void *array, PetscFortranAddr *address)
577: {
578: *address = (PetscFortranAddr)array;
579: }
581: #if defined(PETSC_HAVE_FORTRAN_CAPS)
582: #define f90array4dgetaddrscalar_ F90ARRAY4DGETADDRSCALAR
583: #define f90array4dgetaddrreal_ F90ARRAY4DGETADDRREAL
584: #define f90array4dgetaddrint_ F90ARRAY4DGETADDRINT
585: #define f90array4dgetaddrfortranaddr_ F90ARRAY4DGETADDRFORTRANADDR
586: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
587: #define f90array4dgetaddrscalar_ f90array4dgetaddrscalar
588: #define f90array4dgetaddrreal_ f90array4dgetaddrreal
589: #define f90array4dgetaddrint_ f90array4dgetaddrint
590: #define f90array4dgetaddrfortranaddr_ f90array4dgetaddrfortranaddr
591: #endif
593: PETSC_EXTERN void f90array4dgetaddrscalar_(void *array, PetscFortranAddr *address)
594: {
595: *address = (PetscFortranAddr)array;
596: }
597: PETSC_EXTERN void f90array4dgetaddrreal_(void *array, PetscFortranAddr *address)
598: {
599: *address = (PetscFortranAddr)array;
600: }
601: PETSC_EXTERN void f90array4dgetaddrint_(void *array, PetscFortranAddr *address)
602: {
603: *address = (PetscFortranAddr)array;
604: }
605: PETSC_EXTERN void f90array4dgetaddrfortranaddr_(void *array, PetscFortranAddr *address)
606: {
607: *address = (PetscFortranAddr)array;
608: }