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