1 #include <petsc/private/ftnimpl.h>
2
3 /*@C
4 PetscMPIFortranDatatypeToC - Converts a `MPI_Fint` that contains a Fortran `MPI_Datatype` to its C `MPI_Datatype` equivalent
5
6 Not Collective, No Fortran Support
7
8 Input Parameter:
9 . unit - The Fortran `MPI_Datatype`
10
11 Output Parameter:
12 . dtype - the corresponding C `MPI_Datatype`
13
14 Level: developer
15
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.
22
23 .seealso: `MPI_Fint`, `MPI_Datatype`
24 @*/
PetscMPIFortranDatatypeToC(MPI_Fint unit,MPI_Datatype * dtype)25 PetscErrorCode PetscMPIFortranDatatypeToC(MPI_Fint unit, MPI_Datatype *dtype)
26 {
27 MPI_Datatype ftype;
28
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 }
44
45 /*************************************************************************/
46
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
80
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);
96
97 /*@C
98 F90Array1dCreate - given a `F90Array1d` passed from Fortran associate with it a C array, its starting index and length
99
100 Not Collective, No Fortran Support
101
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
109
110 Level: developer
111
112 Developer Notes:
113 This is used in PETSc Fortran stubs that are used to pass C arrays to Fortran, for example `VecGetArray()`
114
115 This doesn't actually create the `F90Array1d()`, it just associates a C pointer with it.
116
117 There are equivalent routines for 2, 3, and 4 dimensional Fortran arrays.
118
119 .seealso: `F90Array1d`, `F90Array1dAccess()`, `F90Array1dDestroy()`, `F90Array2dCreate()`, `F90Array2dAccess()`, `F90Array2dDestroy()`
120 @*/
F90Array1dCreate(void * array,MPI_Datatype type,PetscInt start,PetscInt len,F90Array1d * ptr PETSC_F90_2PTR_PROTO (ptrd))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 }
144
145 /*@C
146 F90Array1dAccess - given a `F90Array1d` passed from Fortran, accesses from it the associate C array that was provided with `F90Array1dCreate()`
147
148 Not Collective, No Fortran Support
149
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
154
155 Output Parameter:
156 . array - the C address pointer
157
158 Level: developer
159
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.
162
163 There are equivalent routines for 2, 3, and 4 dimensional Fortran arrays.
164
165 .seealso: `F90Array1d`, `F90Array1dCreate()`, `F90Array1dDestroy()`, `F90Array2dCreate()`, `F90Array2dAccess()`, `F90Array2dDestroy()`
166 @*/
F90Array1dAccess(F90Array1d * ptr,MPI_Datatype type,void ** array PETSC_F90_2PTR_PROTO (ptrd))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 }
187
188 /*@C
189 F90Array1dDestroy - given a `F90Array1d` passed from Fortran removes the C array associate with it with `F90Array1dCreate()`
190
191 Not Collective, No Fortran Support
192
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
197
198 Level: developer
199
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()`
202
203 This doesn't actually destroy the `F90Array1d()`, it just removes the associated C pointer from it.
204
205 There are equivalent routines for 2, 3, and 4 dimensional Fortran arrays.
206
207 .seealso: `F90Array1d`, `F90Array1dAccess()`, `F90Array1dCreate()`, `F90Array2dCreate()`, `F90Array2dAccess()`, `F90Array2dDestroy()`
208 @*/
F90Array1dDestroy(F90Array1d * ptr,MPI_Datatype type PETSC_F90_2PTR_PROTO (ptrd))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 PetscFunctionReturn(PETSC_SUCCESS);
224 }
225
226 /*MC
227 F90Array1d - a PETSc C representation of a Fortran `XXX, pointer :: array(:)` object
228
229 Not Collective, No Fortran Support
230
231 Level: developer
232
233 Developer Notes:
234 This is used in PETSc Fortran stubs that are used to control access to C arrays from Fortran, for example `VecGetArray()`
235
236 PETSc does not require any information about the format of this object, all operations on the object are performed by calling Fortran routines.
237
238 There are equivalent objects for 2, 3, and 4 dimensional Fortran arrays.
239
240 .seealso: `F90Array1dAccess()`, `F90Array1dCreate()`, , `F90Array1dDestroy()`, `F90Array2dCreate()`, `F90Array2dAccess()`, `F90Array2dDestroy()`
241 M*/
242
243 #if defined(PETSC_HAVE_FORTRAN_CAPS)
244 #define f90array2dcreatescalar_ F90ARRAY2DCREATESCALAR
245 #define f90array2daccessscalar_ F90ARRAY2DACCESSSCALAR
246 #define f90array2ddestroyscalar_ F90ARRAY2DDESTROYSCALAR
247 #define f90array2dcreatereal_ F90ARRAY2DCREATEREAL
248 #define f90array2daccessreal_ F90ARRAY2DACCESSREAL
249 #define f90array2ddestroyreal_ F90ARRAY2DDESTROYREAL
250 #define f90array2dcreateint_ F90ARRAY2DCREATEINT
251 #define f90array2daccessint_ F90ARRAY2DACCESSINT
252 #define f90array2ddestroyint_ F90ARRAY2DDESTROYINT
253 #define f90array2dcreatefortranaddr_ F90ARRAY2DCREATEFORTRANADDR
254 #define f90array2daccessfortranaddr_ F90ARRAY2DACCESSFORTRANADDR
255 #define f90array2ddestroyfortranaddr_ F90ARRAY2DDESTROYFORTRANADDR
256 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
257 #define f90array2dcreatescalar_ f90array2dcreatescalar
258 #define f90array2daccessscalar_ f90array2daccessscalar
259 #define f90array2ddestroyscalar_ f90array2ddestroyscalar
260 #define f90array2dcreatereal_ f90array2dcreatereal
261 #define f90array2daccessreal_ f90array2daccessreal
262 #define f90array2ddestroyreal_ f90array2ddestroyreal
263 #define f90array2dcreateint_ f90array2dcreateint
264 #define f90array2daccessint_ f90array2daccessint
265 #define f90array2ddestroyint_ f90array2ddestroyint
266 #define f90array2dcreatefortranaddr_ f90array2dcreatefortranaddr
267 #define f90array2daccessfortranaddr_ f90array2daccessfortranaddr
268 #define f90array2ddestroyfortranaddr_ f90array2ddestroyfortranaddr
269 #endif
270
271 PETSC_EXTERN void f90array2dcreatescalar_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array2d *PETSC_F90_2PTR_PROTO_NOVAR);
272 PETSC_EXTERN void f90array2daccessscalar_(F90Array2d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
273 PETSC_EXTERN void f90array2ddestroyscalar_(F90Array2d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
274 PETSC_EXTERN void f90array2dcreatereal_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array2d *PETSC_F90_2PTR_PROTO_NOVAR);
275 PETSC_EXTERN void f90array2daccessreal_(F90Array2d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
276 PETSC_EXTERN void f90array2ddestroyreal_(F90Array2d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
277 PETSC_EXTERN void f90array2dcreateint_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array2d *PETSC_F90_2PTR_PROTO_NOVAR);
278 PETSC_EXTERN void f90array2daccessint_(F90Array2d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
279 PETSC_EXTERN void f90array2ddestroyint_(F90Array2d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
280 PETSC_EXTERN void f90array2dcreatefortranaddr_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array2d *PETSC_F90_2PTR_PROTO_NOVAR);
281 PETSC_EXTERN void f90array2daccessfortranaddr_(F90Array2d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
282 PETSC_EXTERN void f90array2ddestroyfortranaddr_(F90Array2d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
283
F90Array2dCreate(void * array,MPI_Datatype type,PetscInt start1,PetscInt len1,PetscInt start2,PetscInt len2,F90Array2d * ptr PETSC_F90_2PTR_PROTO (ptrd))284 PetscErrorCode F90Array2dCreate(void *array, MPI_Datatype type, PetscInt start1, PetscInt len1, PetscInt start2, PetscInt len2, F90Array2d *ptr PETSC_F90_2PTR_PROTO(ptrd))
285 {
286 PetscFunctionBegin;
287 if (type == MPIU_SCALAR) {
288 f90array2dcreatescalar_(array, &start1, &len1, &start2, &len2, ptr PETSC_F90_2PTR_PARAM(ptrd));
289 } else if (type == MPIU_REAL) {
290 f90array2dcreatereal_(array, &start1, &len1, &start2, &len2, ptr PETSC_F90_2PTR_PARAM(ptrd));
291 } else if (type == MPIU_INT) {
292 f90array2dcreateint_(array, &start1, &len1, &start2, &len2, ptr PETSC_F90_2PTR_PARAM(ptrd));
293 } else if (type == MPIU_FORTRANADDR) {
294 f90array2dcreatefortranaddr_(array, &start1, &len1, &start2, &len2, ptr PETSC_F90_2PTR_PARAM(ptrd));
295 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
296 PetscFunctionReturn(PETSC_SUCCESS);
297 }
298
F90Array2dAccess(F90Array2d * ptr,MPI_Datatype type,void ** array PETSC_F90_2PTR_PROTO (ptrd))299 PetscErrorCode F90Array2dAccess(F90Array2d *ptr, MPI_Datatype type, void **array PETSC_F90_2PTR_PROTO(ptrd))
300 {
301 PetscFunctionBegin;
302 if (type == MPIU_SCALAR) {
303 f90array2daccessscalar_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
304 } else if (type == MPIU_REAL) {
305 f90array2daccessreal_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
306 } else if (type == MPIU_INT) {
307 f90array2daccessint_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
308 } else if (type == MPIU_FORTRANADDR) {
309 f90array2daccessfortranaddr_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
310 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
311 PetscFunctionReturn(PETSC_SUCCESS);
312 }
313
F90Array2dDestroy(F90Array2d * ptr,MPI_Datatype type PETSC_F90_2PTR_PROTO (ptrd))314 PetscErrorCode F90Array2dDestroy(F90Array2d *ptr, MPI_Datatype type PETSC_F90_2PTR_PROTO(ptrd))
315 {
316 PetscFunctionBegin;
317 if (type == MPIU_SCALAR) {
318 f90array2ddestroyscalar_(ptr PETSC_F90_2PTR_PARAM(ptrd));
319 } else if (type == MPIU_REAL) {
320 f90array2ddestroyreal_(ptr PETSC_F90_2PTR_PARAM(ptrd));
321 } else if (type == MPIU_INT) {
322 f90array2ddestroyint_(ptr PETSC_F90_2PTR_PARAM(ptrd));
323 } else if (type == MPIU_FORTRANADDR) {
324 f90array2ddestroyfortranaddr_(ptr PETSC_F90_2PTR_PARAM(ptrd));
325 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
326 PetscFunctionReturn(PETSC_SUCCESS);
327 }
328
329 #if defined(PETSC_HAVE_FORTRAN_CAPS)
330 #define f90array3dcreatescalar_ F90ARRAY3DCREATESCALAR
331 #define f90array3daccessscalar_ F90ARRAY3DACCESSSCALAR
332 #define f90array3ddestroyscalar_ F90ARRAY3DDESTROYSCALAR
333 #define f90array3dcreatereal_ F90ARRAY3DCREATEREAL
334 #define f90array3daccessreal_ F90ARRAY3DACCESSREAL
335 #define f90array3ddestroyreal_ F90ARRAY3DDESTROYREAL
336 #define f90array3dcreateint_ F90ARRAY3DCREATEINT
337 #define f90array3daccessint_ F90ARRAY3DACCESSINT
338 #define f90array3ddestroyint_ F90ARRAY3DDESTROYINT
339 #define f90array3dcreatefortranaddr_ F90ARRAY3DCREATEFORTRANADDR
340 #define f90array3daccessfortranaddr_ F90ARRAY3DACCESSFORTRANADDR
341 #define f90array3ddestroyfortranaddr_ F90ARRAY3DDESTROYFORTRANADDR
342 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
343 #define f90array3dcreatescalar_ f90array3dcreatescalar
344 #define f90array3daccessscalar_ f90array3daccessscalar
345 #define f90array3ddestroyscalar_ f90array3ddestroyscalar
346 #define f90array3dcreatereal_ f90array3dcreatereal
347 #define f90array3daccessreal_ f90array3daccessreal
348 #define f90array3ddestroyreal_ f90array3ddestroyreal
349 #define f90array3dcreateint_ f90array3dcreateint
350 #define f90array3daccessint_ f90array3daccessint
351 #define f90array3ddestroyint_ f90array3ddestroyint
352 #define f90array3dcreatefortranaddr_ f90array3dcreatefortranaddr
353 #define f90array3daccessfortranaddr_ f90array3daccessfortranaddr
354 #define f90array3ddestroyfortranaddr_ f90array3ddestroyfortranaddr
355 #endif
356
357 PETSC_EXTERN void f90array3dcreatescalar_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array3d *PETSC_F90_2PTR_PROTO_NOVAR);
358 PETSC_EXTERN void f90array3daccessscalar_(F90Array3d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
359 PETSC_EXTERN void f90array3ddestroyscalar_(F90Array3d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
360 PETSC_EXTERN void f90array3dcreatereal_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array3d *PETSC_F90_2PTR_PROTO_NOVAR);
361 PETSC_EXTERN void f90array3daccessreal_(F90Array3d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
362 PETSC_EXTERN void f90array3ddestroyreal_(F90Array3d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
363 PETSC_EXTERN void f90array3dcreateint_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array3d *PETSC_F90_2PTR_PROTO_NOVAR);
364 PETSC_EXTERN void f90array3daccessint_(F90Array3d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
365 PETSC_EXTERN void f90array3ddestroyint_(F90Array3d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
366 PETSC_EXTERN void f90array3dcreatefortranaddr_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array3d *PETSC_F90_2PTR_PROTO_NOVAR);
367 PETSC_EXTERN void f90array3daccessfortranaddr_(F90Array3d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
368 PETSC_EXTERN void f90array3ddestroyfortranaddr_(F90Array3d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
369
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))370 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))
371 {
372 PetscFunctionBegin;
373 if (type == MPIU_SCALAR) {
374 f90array3dcreatescalar_(array, &start1, &len1, &start2, &len2, &start3, &len3, ptr PETSC_F90_2PTR_PARAM(ptrd));
375 } else if (type == MPIU_REAL) {
376 f90array3dcreatereal_(array, &start1, &len1, &start2, &len2, &start3, &len3, ptr PETSC_F90_2PTR_PARAM(ptrd));
377 } else if (type == MPIU_INT) {
378 f90array3dcreateint_(array, &start1, &len1, &start2, &len2, &start3, &len3, ptr PETSC_F90_2PTR_PARAM(ptrd));
379 } else if (type == MPIU_FORTRANADDR) {
380 f90array3dcreatefortranaddr_(array, &start1, &len1, &start2, &len2, &start3, &len3, ptr PETSC_F90_2PTR_PARAM(ptrd));
381 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
382 PetscFunctionReturn(PETSC_SUCCESS);
383 }
384
F90Array3dAccess(F90Array3d * ptr,MPI_Datatype type,void ** array PETSC_F90_2PTR_PROTO (ptrd))385 PetscErrorCode F90Array3dAccess(F90Array3d *ptr, MPI_Datatype type, void **array PETSC_F90_2PTR_PROTO(ptrd))
386 {
387 PetscFunctionBegin;
388 if (type == MPIU_SCALAR) {
389 f90array3daccessscalar_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
390 } else if (type == MPIU_REAL) {
391 f90array3daccessreal_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
392 } else if (type == MPIU_INT) {
393 f90array3daccessint_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
394 } else if (type == MPIU_FORTRANADDR) {
395 f90array3daccessfortranaddr_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
396 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
397 PetscFunctionReturn(PETSC_SUCCESS);
398 }
399
F90Array3dDestroy(F90Array3d * ptr,MPI_Datatype type PETSC_F90_2PTR_PROTO (ptrd))400 PetscErrorCode F90Array3dDestroy(F90Array3d *ptr, MPI_Datatype type PETSC_F90_2PTR_PROTO(ptrd))
401 {
402 PetscFunctionBegin;
403 if (type == MPIU_SCALAR) {
404 f90array3ddestroyscalar_(ptr PETSC_F90_2PTR_PARAM(ptrd));
405 } else if (type == MPIU_REAL) {
406 f90array3ddestroyreal_(ptr PETSC_F90_2PTR_PARAM(ptrd));
407 } else if (type == MPIU_INT) {
408 f90array3ddestroyint_(ptr PETSC_F90_2PTR_PARAM(ptrd));
409 } else if (type == MPIU_FORTRANADDR) {
410 f90array3ddestroyfortranaddr_(ptr PETSC_F90_2PTR_PARAM(ptrd));
411 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
412 PetscFunctionReturn(PETSC_SUCCESS);
413 }
414
415 #if defined(PETSC_HAVE_FORTRAN_CAPS)
416 #define f90array4dcreatescalar_ F90ARRAY4DCREATESCALAR
417 #define f90array4daccessscalar_ F90ARRAY4DACCESSSCALAR
418 #define f90array4ddestroyscalar_ F90ARRAY4DDESTROYSCALAR
419 #define f90array4dcreatereal_ F90ARRAY4DCREATEREAL
420 #define f90array4daccessreal_ F90ARRAY4DACCESSREAL
421 #define f90array4ddestroyreal_ F90ARRAY4DDESTROYREAL
422 #define f90array4dcreateint_ F90ARRAY4DCREATEINT
423 #define f90array4daccessint_ F90ARRAY4DACCESSINT
424 #define f90array4ddestroyint_ F90ARRAY4DDESTROYINT
425 #define f90array4dcreatefortranaddr_ F90ARRAY4DCREATEFORTRANADDR
426 #define f90array4daccessfortranaddr_ F90ARRAY4DACCESSFORTRANADDR
427 #define f90array4ddestroyfortranaddr_ F90ARRAY4DDESTROYFORTRANADDR
428 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
429 #define f90array4dcreatescalar_ f90array4dcreatescalar
430 #define f90array4daccessscalar_ f90array4daccessscalar
431 #define f90array4ddestroyscalar_ f90array4ddestroyscalar
432 #define f90array4dcreatereal_ f90array4dcreatereal
433 #define f90array4daccessreal_ f90array4daccessreal
434 #define f90array4ddestroyreal_ f90array4ddestroyreal
435 #define f90array4dcreateint_ f90array4dcreateint
436 #define f90array4daccessint_ f90array4daccessint
437 #define f90array4ddestroyint_ f90array4ddestroyint
438 #define f90array4dcreatefortranaddr_ f90array4dcreatefortranaddr
439 #define f90array4daccessfortranaddr_ f90array4daccessfortranaddr
440 #define f90array4ddestroyfortranaddr_ f90array4ddestroyfortranaddr
441 #endif
442
443 PETSC_EXTERN void f90array4dcreatescalar_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array4d *PETSC_F90_2PTR_PROTO_NOVAR);
444 PETSC_EXTERN void f90array4daccessscalar_(F90Array4d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
445 PETSC_EXTERN void f90array4ddestroyscalar_(F90Array4d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
446 PETSC_EXTERN void f90array4dcreatereal_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array4d *PETSC_F90_2PTR_PROTO_NOVAR);
447 PETSC_EXTERN void f90array4daccessreal_(F90Array4d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
448 PETSC_EXTERN void f90array4ddestroyreal_(F90Array4d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
449 PETSC_EXTERN void f90array4dcreateint_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array4d *PETSC_F90_2PTR_PROTO_NOVAR);
450 PETSC_EXTERN void f90array4daccessint_(F90Array4d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
451 PETSC_EXTERN void f90array4ddestroyint_(F90Array4d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
452 PETSC_EXTERN void f90array4dcreatefortranaddr_(void *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, F90Array4d *PETSC_F90_2PTR_PROTO_NOVAR);
453 PETSC_EXTERN void f90array4daccessfortranaddr_(F90Array4d *, void **PETSC_F90_2PTR_PROTO_NOVAR);
454 PETSC_EXTERN void f90array4ddestroyfortranaddr_(F90Array4d *ptr PETSC_F90_2PTR_PROTO_NOVAR);
455
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))456 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))
457 {
458 PetscFunctionBegin;
459 PetscCheck(type == MPIU_SCALAR, PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
460 f90array4dcreatescalar_(array, &start1, &len1, &start2, &len2, &start3, &len3, &start4, &len4, ptr PETSC_F90_2PTR_PARAM(ptrd));
461 PetscFunctionReturn(PETSC_SUCCESS);
462 }
463
F90Array4dAccess(F90Array4d * ptr,MPI_Datatype type,void ** array PETSC_F90_2PTR_PROTO (ptrd))464 PetscErrorCode F90Array4dAccess(F90Array4d *ptr, MPI_Datatype type, void **array PETSC_F90_2PTR_PROTO(ptrd))
465 {
466 PetscFunctionBegin;
467 if (type == MPIU_SCALAR) {
468 f90array4daccessscalar_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
469 } else if (type == MPIU_REAL) {
470 f90array4daccessreal_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
471 } else if (type == MPIU_INT) {
472 f90array4daccessint_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
473 } else if (type == MPIU_FORTRANADDR) {
474 f90array4daccessfortranaddr_(ptr, array PETSC_F90_2PTR_PARAM(ptrd));
475 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
476 PetscFunctionReturn(PETSC_SUCCESS);
477 }
478
F90Array4dDestroy(F90Array4d * ptr,MPI_Datatype type PETSC_F90_2PTR_PROTO (ptrd))479 PetscErrorCode F90Array4dDestroy(F90Array4d *ptr, MPI_Datatype type PETSC_F90_2PTR_PROTO(ptrd))
480 {
481 PetscFunctionBegin;
482 PetscCheck(type == MPIU_SCALAR, PETSC_COMM_SELF, PETSC_ERR_SUP, "Unsupported MPI_Datatype");
483 f90array4ddestroyscalar_(ptr PETSC_F90_2PTR_PARAM(ptrd));
484 PetscFunctionReturn(PETSC_SUCCESS);
485 }
486
487 #if defined(PETSC_HAVE_FORTRAN_CAPS)
488 #define f90array1dgetaddrscalar_ F90ARRAY1DGETADDRSCALAR
489 #define f90array1dgetaddrreal_ F90ARRAY1DGETADDRREAL
490 #define f90array1dgetaddrint_ F90ARRAY1DGETADDRINT
491 #define f90array1dgetaddrmpiint_ F90ARRAY1DGETADDRMPIINT
492 #define f90array1dgetaddrfortranaddr_ F90ARRAY1DGETADDRFORTRANADDR
493 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
494 #define f90array1dgetaddrscalar_ f90array1dgetaddrscalar
495 #define f90array1dgetaddrreal_ f90array1dgetaddrreal
496 #define f90array1dgetaddrint_ f90array1dgetaddrint
497 #define f90array1dgetaddrmpiint_ f90array1dgetaddrmpiint
498 #define f90array1dgetaddrfortranaddr_ f90array1dgetaddrfortranaddr
499 #endif
500
f90array1dgetaddrscalar_(void * array,PetscFortranAddr * address)501 PETSC_EXTERN void f90array1dgetaddrscalar_(void *array, PetscFortranAddr *address)
502 {
503 *address = (PetscFortranAddr)array;
504 }
f90array1dgetaddrreal_(void * array,PetscFortranAddr * address)505 PETSC_EXTERN void f90array1dgetaddrreal_(void *array, PetscFortranAddr *address)
506 {
507 *address = (PetscFortranAddr)array;
508 }
f90array1dgetaddrint_(void * array,PetscFortranAddr * address)509 PETSC_EXTERN void f90array1dgetaddrint_(void *array, PetscFortranAddr *address)
510 {
511 *address = (PetscFortranAddr)array;
512 }
f90array1dgetaddrmpiint_(void * array,PetscFortranAddr * address)513 PETSC_EXTERN void f90array1dgetaddrmpiint_(void *array, PetscFortranAddr *address)
514 {
515 *address = (PetscFortranAddr)array;
516 }
f90array1dgetaddrfortranaddr_(void * array,PetscFortranAddr * address)517 PETSC_EXTERN void f90array1dgetaddrfortranaddr_(void *array, PetscFortranAddr *address)
518 {
519 *address = (PetscFortranAddr)array;
520 }
521
522 #if defined(PETSC_HAVE_FORTRAN_CAPS)
523 #define f90array2dgetaddrscalar_ F90ARRAY2DGETADDRSCALAR
524 #define f90array2dgetaddrreal_ F90ARRAY2DGETADDRREAL
525 #define f90array2dgetaddrint_ F90ARRAY2DGETADDRINT
526 #define f90array2dgetaddrfortranaddr_ F90ARRAY2DGETADDRFORTRANADDR
527 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
528 #define f90array2dgetaddrscalar_ f90array2dgetaddrscalar
529 #define f90array2dgetaddrreal_ f90array2dgetaddrreal
530 #define f90array2dgetaddrint_ f90array2dgetaddrint
531 #define f90array2dgetaddrfortranaddr_ f90array2dgetaddrfortranaddr
532 #endif
533
f90array2dgetaddrscalar_(void * array,PetscFortranAddr * address)534 PETSC_EXTERN void f90array2dgetaddrscalar_(void *array, PetscFortranAddr *address)
535 {
536 *address = (PetscFortranAddr)array;
537 }
f90array2dgetaddrreal_(void * array,PetscFortranAddr * address)538 PETSC_EXTERN void f90array2dgetaddrreal_(void *array, PetscFortranAddr *address)
539 {
540 *address = (PetscFortranAddr)array;
541 }
f90array2dgetaddrint_(void * array,PetscFortranAddr * address)542 PETSC_EXTERN void f90array2dgetaddrint_(void *array, PetscFortranAddr *address)
543 {
544 *address = (PetscFortranAddr)array;
545 }
f90array2dgetaddrfortranaddr_(void * array,PetscFortranAddr * address)546 PETSC_EXTERN void f90array2dgetaddrfortranaddr_(void *array, PetscFortranAddr *address)
547 {
548 *address = (PetscFortranAddr)array;
549 }
550
551 #if defined(PETSC_HAVE_FORTRAN_CAPS)
552 #define f90array3dgetaddrscalar_ F90ARRAY3DGETADDRSCALAR
553 #define f90array3dgetaddrreal_ F90ARRAY3DGETADDRREAL
554 #define f90array3dgetaddrint_ F90ARRAY3DGETADDRINT
555 #define f90array3dgetaddrfortranaddr_ F90ARRAY3DGETADDRFORTRANADDR
556 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
557 #define f90array3dgetaddrscalar_ f90array3dgetaddrscalar
558 #define f90array3dgetaddrreal_ f90array3dgetaddrreal
559 #define f90array3dgetaddrint_ f90array3dgetaddrint
560 #define f90array3dgetaddrfortranaddr_ f90array3dgetaddrfortranaddr
561 #endif
562
f90array3dgetaddrscalar_(void * array,PetscFortranAddr * address)563 PETSC_EXTERN void f90array3dgetaddrscalar_(void *array, PetscFortranAddr *address)
564 {
565 *address = (PetscFortranAddr)array;
566 }
f90array3dgetaddrreal_(void * array,PetscFortranAddr * address)567 PETSC_EXTERN void f90array3dgetaddrreal_(void *array, PetscFortranAddr *address)
568 {
569 *address = (PetscFortranAddr)array;
570 }
f90array3dgetaddrint_(void * array,PetscFortranAddr * address)571 PETSC_EXTERN void f90array3dgetaddrint_(void *array, PetscFortranAddr *address)
572 {
573 *address = (PetscFortranAddr)array;
574 }
f90array3dgetaddrfortranaddr_(void * array,PetscFortranAddr * address)575 PETSC_EXTERN void f90array3dgetaddrfortranaddr_(void *array, PetscFortranAddr *address)
576 {
577 *address = (PetscFortranAddr)array;
578 }
579
580 #if defined(PETSC_HAVE_FORTRAN_CAPS)
581 #define f90array4dgetaddrscalar_ F90ARRAY4DGETADDRSCALAR
582 #define f90array4dgetaddrreal_ F90ARRAY4DGETADDRREAL
583 #define f90array4dgetaddrint_ F90ARRAY4DGETADDRINT
584 #define f90array4dgetaddrfortranaddr_ F90ARRAY4DGETADDRFORTRANADDR
585 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
586 #define f90array4dgetaddrscalar_ f90array4dgetaddrscalar
587 #define f90array4dgetaddrreal_ f90array4dgetaddrreal
588 #define f90array4dgetaddrint_ f90array4dgetaddrint
589 #define f90array4dgetaddrfortranaddr_ f90array4dgetaddrfortranaddr
590 #endif
591
f90array4dgetaddrscalar_(void * array,PetscFortranAddr * address)592 PETSC_EXTERN void f90array4dgetaddrscalar_(void *array, PetscFortranAddr *address)
593 {
594 *address = (PetscFortranAddr)array;
595 }
f90array4dgetaddrreal_(void * array,PetscFortranAddr * address)596 PETSC_EXTERN void f90array4dgetaddrreal_(void *array, PetscFortranAddr *address)
597 {
598 *address = (PetscFortranAddr)array;
599 }
f90array4dgetaddrint_(void * array,PetscFortranAddr * address)600 PETSC_EXTERN void f90array4dgetaddrint_(void *array, PetscFortranAddr *address)
601 {
602 *address = (PetscFortranAddr)array;
603 }
f90array4dgetaddrfortranaddr_(void * array,PetscFortranAddr * address)604 PETSC_EXTERN void f90array4dgetaddrfortranaddr_(void *array, PetscFortranAddr *address)
605 {
606 *address = (PetscFortranAddr)array;
607 }
608