1*6dd63270SBarry Smith #include <petsc/private/ftnimpl.h> 2*6dd63270SBarry Smith #include <petscdmda.h> 3*6dd63270SBarry Smith 4*6dd63270SBarry Smith #if defined(PETSC_HAVE_FORTRAN_CAPS) 5*6dd63270SBarry Smith #define dmdavecgetarrayf901_ DMDAVECGETARRAYF901 6*6dd63270SBarry Smith #define dmdavecrestorearrayf901_ DMDAVECRESTOREARRAYF901 7*6dd63270SBarry Smith #define dmdavecgetarrayf902_ DMDAVECGETARRAYF902 8*6dd63270SBarry Smith #define dmdavecrestorearrayf902_ DMDAVECRESTOREARRAYF902 9*6dd63270SBarry Smith #define dmdavecgetarrayf903_ DMDAVECGETARRAYF903 10*6dd63270SBarry Smith #define dmdavecrestorearrayf903_ DMDAVECRESTOREARRAYF903 11*6dd63270SBarry Smith #define dmdavecgetarrayf904_ DMDAVECGETARRAYF904 12*6dd63270SBarry Smith #define dmdavecrestorearrayf904_ DMDAVECRESTOREARRAYF904 13*6dd63270SBarry Smith #define dmdavecgetarrayreadf901_ DMDAVECGETARRAYREADF901 14*6dd63270SBarry Smith #define dmdavecrestorearrayreadf901_ DMDAVECRESTOREARRAYREADF901 15*6dd63270SBarry Smith #define dmdavecgetarrayreadf902_ DMDAVECGETARRAYREADF902 16*6dd63270SBarry Smith #define dmdavecrestorearrayreadf902_ DMDAVECRESTOREARRAYREADF902 17*6dd63270SBarry Smith #define dmdavecgetarrayreadf903_ DMDAVECGETARRAYREADF903 18*6dd63270SBarry Smith #define dmdavecrestorearrayreadf903_ DMDAVECRESTOREARRAYREADF903 19*6dd63270SBarry Smith #define dmdavecgetarrayreadf904_ DMDAVECGETARRAYREADF904 20*6dd63270SBarry Smith #define dmdavecrestorearrayreadf904_ DMDAVECRESTOREARRAYREADF904 21*6dd63270SBarry Smith #define dmdagetelements_ DMDAGETELEMENTS 22*6dd63270SBarry Smith #define dmdarestoreelements_ DMDARESTOREELEMENTS 23*6dd63270SBarry Smith #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 24*6dd63270SBarry Smith #define dmdavecgetarrayf901_ dmdavecgetarrayf901 25*6dd63270SBarry Smith #define dmdavecrestorearrayf901_ dmdavecrestorearrayf901 26*6dd63270SBarry Smith #define dmdavecgetarrayf902_ dmdavecgetarrayf902 27*6dd63270SBarry Smith #define dmdavecrestorearrayf902_ dmdavecrestorearrayf902 28*6dd63270SBarry Smith #define dmdavecgetarrayf903_ dmdavecgetarrayf903 29*6dd63270SBarry Smith #define dmdavecrestorearrayf903_ dmdavecrestorearrayf903 30*6dd63270SBarry Smith #define dmdavecgetarrayf904_ dmdavecgetarrayf904 31*6dd63270SBarry Smith #define dmdavecrestorearrayf904_ dmdavecrestorearrayf904 32*6dd63270SBarry Smith #define dmdavecgetarrayreadf901_ dmdavecgetarrayreadf901 33*6dd63270SBarry Smith #define dmdavecrestorearrayreadf901_ dmdavecrestorearrayreadf901 34*6dd63270SBarry Smith #define dmdavecgetarrayreadf902_ dmdavecgetarrayreadf902 35*6dd63270SBarry Smith #define dmdavecrestorearrayreadf902_ dmdavecrestorearrayreadf902 36*6dd63270SBarry Smith #define dmdavecgetarrayreadf903_ dmdavecgetarrayreadf903 37*6dd63270SBarry Smith #define dmdavecrestorearrayreadf903_ dmdavecrestorearrayreadf903 38*6dd63270SBarry Smith #define dmdavecgetarrayreadf904_ dmdavecgetarrayreadf904 39*6dd63270SBarry Smith #define dmdavecrestorearrayreadf904_ dmdavecrestorearrayreadf904 40*6dd63270SBarry Smith #define dmdagetelements_ dmdagetelements 41*6dd63270SBarry Smith #define dmdarestoreelements_ dmdarestoreelements 42*6dd63270SBarry Smith #endif 43*6dd63270SBarry Smith 44*6dd63270SBarry Smith PETSC_EXTERN void dmdagetelements_(DM *dm, PetscInt *nel, PetscInt *nen, F90Array1d *e, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 45*6dd63270SBarry Smith { 46*6dd63270SBarry Smith const PetscInt *fa; 47*6dd63270SBarry Smith 48*6dd63270SBarry Smith if (!e) { 49*6dd63270SBarry Smith *ierr = PetscError(((PetscObject)e)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "e==NULL, maybe #include <petsc/finclude/petscvec.h> is missing?"); 50*6dd63270SBarry Smith return; 51*6dd63270SBarry Smith } 52*6dd63270SBarry Smith *ierr = DMDAGetElements(*dm, nel, nen, &fa); 53*6dd63270SBarry Smith if (*ierr) return; 54*6dd63270SBarry Smith *ierr = F90Array1dCreate((PetscInt *)fa, MPIU_INT, 1, (*nel) * (*nen), e PETSC_F90_2PTR_PARAM(ptrd)); 55*6dd63270SBarry Smith } 56*6dd63270SBarry Smith 57*6dd63270SBarry Smith PETSC_EXTERN void dmdarestoreelements_(DM *dm, PetscInt *nel, PetscInt *nen, F90Array1d *e, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 58*6dd63270SBarry Smith { 59*6dd63270SBarry Smith if (!e) { 60*6dd63270SBarry Smith *ierr = PetscError(((PetscObject)e)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "e==NULL, maybe #include <petsc/finclude/petscvec.h> is missing?"); 61*6dd63270SBarry Smith return; 62*6dd63270SBarry Smith } 63*6dd63270SBarry Smith *ierr = F90Array1dDestroy(e, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 64*6dd63270SBarry Smith } 65*6dd63270SBarry Smith 66*6dd63270SBarry Smith PETSC_EXTERN void dmdavecgetarrayf901_(DM *da, Vec *v, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 67*6dd63270SBarry Smith { 68*6dd63270SBarry Smith PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 69*6dd63270SBarry Smith PetscScalar *aa; 70*6dd63270SBarry Smith 71*6dd63270SBarry Smith *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 72*6dd63270SBarry Smith if (*ierr) return; 73*6dd63270SBarry Smith *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 74*6dd63270SBarry Smith if (*ierr) return; 75*6dd63270SBarry Smith *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 76*6dd63270SBarry Smith if (*ierr) return; 77*6dd63270SBarry Smith 78*6dd63270SBarry Smith /* Handle case where user passes in global vector as opposed to local */ 79*6dd63270SBarry Smith *ierr = VecGetLocalSize(*v, &N); 80*6dd63270SBarry Smith if (*ierr) return; 81*6dd63270SBarry Smith if (N == xm * ym * zm * dof) { 82*6dd63270SBarry Smith gxm = xm; 83*6dd63270SBarry Smith gym = ym; 84*6dd63270SBarry Smith gzm = zm; 85*6dd63270SBarry Smith gxs = xs; 86*6dd63270SBarry Smith gys = ys; 87*6dd63270SBarry Smith gzs = zs; 88*6dd63270SBarry Smith } else if (N != gxm * gym * gzm * dof) { 89*6dd63270SBarry Smith *ierr = PETSC_ERR_ARG_INCOMP; 90*6dd63270SBarry Smith return; 91*6dd63270SBarry Smith } 92*6dd63270SBarry Smith *ierr = VecGetArray(*v, &aa); 93*6dd63270SBarry Smith if (*ierr) return; 94*6dd63270SBarry Smith *ierr = F90Array1dCreate(aa, MPIU_SCALAR, gxs, gxm, a PETSC_F90_2PTR_PARAM(ptrd)); 95*6dd63270SBarry Smith if (*ierr) return; 96*6dd63270SBarry Smith } 97*6dd63270SBarry Smith 98*6dd63270SBarry Smith PETSC_EXTERN void dmdavecrestorearrayf901_(DM *da, Vec *v, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 99*6dd63270SBarry Smith { 100*6dd63270SBarry Smith PetscScalar *fa; 101*6dd63270SBarry Smith *ierr = F90Array1dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 102*6dd63270SBarry Smith *ierr = VecRestoreArray(*v, &fa); 103*6dd63270SBarry Smith if (*ierr) return; 104*6dd63270SBarry Smith *ierr = F90Array1dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 105*6dd63270SBarry Smith } 106*6dd63270SBarry Smith 107*6dd63270SBarry Smith PETSC_EXTERN void dmdavecgetarrayf902_(DM *da, Vec *v, F90Array2d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 108*6dd63270SBarry Smith { 109*6dd63270SBarry Smith PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 110*6dd63270SBarry Smith PetscScalar *aa; 111*6dd63270SBarry Smith 112*6dd63270SBarry Smith *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 113*6dd63270SBarry Smith if (*ierr) return; 114*6dd63270SBarry Smith *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 115*6dd63270SBarry Smith if (*ierr) return; 116*6dd63270SBarry Smith *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 117*6dd63270SBarry Smith if (*ierr) return; 118*6dd63270SBarry Smith 119*6dd63270SBarry Smith /* Handle case where user passes in global vector as opposed to local */ 120*6dd63270SBarry Smith *ierr = VecGetLocalSize(*v, &N); 121*6dd63270SBarry Smith if (*ierr) return; 122*6dd63270SBarry Smith if (N == xm * ym * zm * dof) { 123*6dd63270SBarry Smith gxm = xm; 124*6dd63270SBarry Smith gym = ym; 125*6dd63270SBarry Smith gzm = zm; 126*6dd63270SBarry Smith gxs = xs; 127*6dd63270SBarry Smith gys = ys; 128*6dd63270SBarry Smith gzs = zs; 129*6dd63270SBarry Smith } else if (N != gxm * gym * gzm * dof) { 130*6dd63270SBarry Smith *ierr = PETSC_ERR_ARG_INCOMP; 131*6dd63270SBarry Smith return; 132*6dd63270SBarry Smith } 133*6dd63270SBarry Smith if (dim == 1) { 134*6dd63270SBarry Smith gys = gxs; 135*6dd63270SBarry Smith gym = gxm; 136*6dd63270SBarry Smith gxs = 0; 137*6dd63270SBarry Smith gxm = dof; 138*6dd63270SBarry Smith } 139*6dd63270SBarry Smith *ierr = VecGetArray(*v, &aa); 140*6dd63270SBarry Smith if (*ierr) return; 141*6dd63270SBarry Smith *ierr = F90Array2dCreate(aa, MPIU_SCALAR, gxs, gxm, gys, gym, a PETSC_F90_2PTR_PARAM(ptrd)); 142*6dd63270SBarry Smith if (*ierr) return; 143*6dd63270SBarry Smith } 144*6dd63270SBarry Smith 145*6dd63270SBarry Smith PETSC_EXTERN void dmdavecrestorearrayf902_(DM *da, Vec *v, F90Array2d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 146*6dd63270SBarry Smith { 147*6dd63270SBarry Smith PetscScalar *fa; 148*6dd63270SBarry Smith *ierr = F90Array2dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 149*6dd63270SBarry Smith *ierr = VecRestoreArray(*v, &fa); 150*6dd63270SBarry Smith if (*ierr) return; 151*6dd63270SBarry Smith *ierr = F90Array2dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 152*6dd63270SBarry Smith } 153*6dd63270SBarry Smith 154*6dd63270SBarry Smith PETSC_EXTERN void dmdavecgetarrayf903_(DM *da, Vec *v, F90Array3d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 155*6dd63270SBarry Smith { 156*6dd63270SBarry Smith PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 157*6dd63270SBarry Smith PetscScalar *aa; 158*6dd63270SBarry Smith 159*6dd63270SBarry Smith *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 160*6dd63270SBarry Smith if (*ierr) return; 161*6dd63270SBarry Smith *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 162*6dd63270SBarry Smith if (*ierr) return; 163*6dd63270SBarry Smith *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 164*6dd63270SBarry Smith if (*ierr) return; 165*6dd63270SBarry Smith 166*6dd63270SBarry Smith /* Handle case where user passes in global vector as opposed to local */ 167*6dd63270SBarry Smith *ierr = VecGetLocalSize(*v, &N); 168*6dd63270SBarry Smith if (*ierr) return; 169*6dd63270SBarry Smith if (N == xm * ym * zm * dof) { 170*6dd63270SBarry Smith gxm = xm; 171*6dd63270SBarry Smith gym = ym; 172*6dd63270SBarry Smith gzm = zm; 173*6dd63270SBarry Smith gxs = xs; 174*6dd63270SBarry Smith gys = ys; 175*6dd63270SBarry Smith gzs = zs; 176*6dd63270SBarry Smith } else if (N != gxm * gym * gzm * dof) { 177*6dd63270SBarry Smith *ierr = PETSC_ERR_ARG_INCOMP; 178*6dd63270SBarry Smith return; 179*6dd63270SBarry Smith } 180*6dd63270SBarry Smith if (dim == 2) { 181*6dd63270SBarry Smith gzs = gys; 182*6dd63270SBarry Smith gzm = gym; 183*6dd63270SBarry Smith gys = gxs; 184*6dd63270SBarry Smith gym = gxm; 185*6dd63270SBarry Smith gxs = 0; 186*6dd63270SBarry Smith gxm = dof; 187*6dd63270SBarry Smith } 188*6dd63270SBarry Smith *ierr = VecGetArray(*v, &aa); 189*6dd63270SBarry Smith if (*ierr) return; 190*6dd63270SBarry Smith *ierr = F90Array3dCreate(aa, MPIU_SCALAR, gxs, gxm, gys, gym, gzs, gzm, a PETSC_F90_2PTR_PARAM(ptrd)); 191*6dd63270SBarry Smith if (*ierr) return; 192*6dd63270SBarry Smith } 193*6dd63270SBarry Smith 194*6dd63270SBarry Smith PETSC_EXTERN void dmdavecrestorearrayf903_(DM *da, Vec *v, F90Array3d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 195*6dd63270SBarry Smith { 196*6dd63270SBarry Smith PetscScalar *fa; 197*6dd63270SBarry Smith *ierr = F90Array3dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 198*6dd63270SBarry Smith *ierr = VecRestoreArray(*v, &fa); 199*6dd63270SBarry Smith if (*ierr) return; 200*6dd63270SBarry Smith *ierr = F90Array3dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 201*6dd63270SBarry Smith } 202*6dd63270SBarry Smith 203*6dd63270SBarry Smith PETSC_EXTERN void dmdavecgetarrayf904_(DM *da, Vec *v, F90Array4d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 204*6dd63270SBarry Smith { 205*6dd63270SBarry Smith PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof, zero = 0; 206*6dd63270SBarry Smith PetscScalar *aa; 207*6dd63270SBarry Smith 208*6dd63270SBarry Smith *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 209*6dd63270SBarry Smith if (*ierr) return; 210*6dd63270SBarry Smith *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 211*6dd63270SBarry Smith if (*ierr) return; 212*6dd63270SBarry Smith *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 213*6dd63270SBarry Smith if (*ierr) return; 214*6dd63270SBarry Smith 215*6dd63270SBarry Smith /* Handle case where user passes in global vector as opposed to local */ 216*6dd63270SBarry Smith *ierr = VecGetLocalSize(*v, &N); 217*6dd63270SBarry Smith if (*ierr) return; 218*6dd63270SBarry Smith if (N == xm * ym * zm * dof) { 219*6dd63270SBarry Smith gxm = xm; 220*6dd63270SBarry Smith gym = ym; 221*6dd63270SBarry Smith gzm = zm; 222*6dd63270SBarry Smith gxs = xs; 223*6dd63270SBarry Smith gys = ys; 224*6dd63270SBarry Smith gzs = zs; 225*6dd63270SBarry Smith } else if (N != gxm * gym * gzm * dof) { 226*6dd63270SBarry Smith *ierr = PETSC_ERR_ARG_INCOMP; 227*6dd63270SBarry Smith return; 228*6dd63270SBarry Smith } 229*6dd63270SBarry Smith *ierr = VecGetArray(*v, &aa); 230*6dd63270SBarry Smith if (*ierr) return; 231*6dd63270SBarry Smith *ierr = F90Array4dCreate(aa, MPIU_SCALAR, zero, dof, gxs, gxm, gys, gym, gzs, gzm, a PETSC_F90_2PTR_PARAM(ptrd)); 232*6dd63270SBarry Smith if (*ierr) return; 233*6dd63270SBarry Smith } 234*6dd63270SBarry Smith 235*6dd63270SBarry Smith PETSC_EXTERN void dmdavecrestorearrayf904_(DM *da, Vec *v, F90Array4d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 236*6dd63270SBarry Smith { 237*6dd63270SBarry Smith PetscScalar *fa; 238*6dd63270SBarry Smith /* 239*6dd63270SBarry Smith F90Array4dAccess is not implemented, so the following call would fail 240*6dd63270SBarry Smith */ 241*6dd63270SBarry Smith *ierr = F90Array4dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 242*6dd63270SBarry Smith *ierr = VecRestoreArray(*v, &fa); 243*6dd63270SBarry Smith if (*ierr) return; 244*6dd63270SBarry Smith *ierr = F90Array4dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 245*6dd63270SBarry Smith } 246*6dd63270SBarry Smith 247*6dd63270SBarry Smith PETSC_EXTERN void dmdavecgetarrayreadf901_(DM *da, Vec *v, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 248*6dd63270SBarry Smith { 249*6dd63270SBarry Smith PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 250*6dd63270SBarry Smith const PetscScalar *aa; 251*6dd63270SBarry Smith 252*6dd63270SBarry Smith *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 253*6dd63270SBarry Smith if (*ierr) return; 254*6dd63270SBarry Smith *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 255*6dd63270SBarry Smith if (*ierr) return; 256*6dd63270SBarry Smith *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 257*6dd63270SBarry Smith if (*ierr) return; 258*6dd63270SBarry Smith 259*6dd63270SBarry Smith /* Handle case where user passes in global vector as opposed to local */ 260*6dd63270SBarry Smith *ierr = VecGetLocalSize(*v, &N); 261*6dd63270SBarry Smith if (*ierr) return; 262*6dd63270SBarry Smith if (N == xm * ym * zm * dof) { 263*6dd63270SBarry Smith gxm = xm; 264*6dd63270SBarry Smith gym = ym; 265*6dd63270SBarry Smith gzm = zm; 266*6dd63270SBarry Smith gxs = xs; 267*6dd63270SBarry Smith gys = ys; 268*6dd63270SBarry Smith gzs = zs; 269*6dd63270SBarry Smith } else if (N != gxm * gym * gzm * dof) { 270*6dd63270SBarry Smith *ierr = PETSC_ERR_ARG_INCOMP; 271*6dd63270SBarry Smith return; 272*6dd63270SBarry Smith } 273*6dd63270SBarry Smith *ierr = VecGetArrayRead(*v, &aa); 274*6dd63270SBarry Smith if (*ierr) return; 275*6dd63270SBarry Smith *ierr = F90Array1dCreate((void *)aa, MPIU_SCALAR, gxs, gxm, a PETSC_F90_2PTR_PARAM(ptrd)); 276*6dd63270SBarry Smith if (*ierr) return; 277*6dd63270SBarry Smith } 278*6dd63270SBarry Smith 279*6dd63270SBarry Smith PETSC_EXTERN void dmdavecrestorearrayreadf901_(DM *da, Vec *v, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 280*6dd63270SBarry Smith { 281*6dd63270SBarry Smith const PetscScalar *fa; 282*6dd63270SBarry Smith *ierr = F90Array1dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 283*6dd63270SBarry Smith *ierr = VecRestoreArrayRead(*v, &fa); 284*6dd63270SBarry Smith if (*ierr) return; 285*6dd63270SBarry Smith *ierr = F90Array1dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 286*6dd63270SBarry Smith } 287*6dd63270SBarry Smith 288*6dd63270SBarry Smith PETSC_EXTERN void dmdavecgetarrayreadf902_(DM *da, Vec *v, F90Array2d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 289*6dd63270SBarry Smith { 290*6dd63270SBarry Smith PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 291*6dd63270SBarry Smith const PetscScalar *aa; 292*6dd63270SBarry Smith 293*6dd63270SBarry Smith *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 294*6dd63270SBarry Smith if (*ierr) return; 295*6dd63270SBarry Smith *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 296*6dd63270SBarry Smith if (*ierr) return; 297*6dd63270SBarry Smith *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 298*6dd63270SBarry Smith if (*ierr) return; 299*6dd63270SBarry Smith 300*6dd63270SBarry Smith /* Handle case where user passes in global vector as opposed to local */ 301*6dd63270SBarry Smith *ierr = VecGetLocalSize(*v, &N); 302*6dd63270SBarry Smith if (*ierr) return; 303*6dd63270SBarry Smith if (N == xm * ym * zm * dof) { 304*6dd63270SBarry Smith gxm = xm; 305*6dd63270SBarry Smith gym = ym; 306*6dd63270SBarry Smith gzm = zm; 307*6dd63270SBarry Smith gxs = xs; 308*6dd63270SBarry Smith gys = ys; 309*6dd63270SBarry Smith gzs = zs; 310*6dd63270SBarry Smith } else if (N != gxm * gym * gzm * dof) { 311*6dd63270SBarry Smith *ierr = PETSC_ERR_ARG_INCOMP; 312*6dd63270SBarry Smith return; 313*6dd63270SBarry Smith } 314*6dd63270SBarry Smith if (dim == 1) { 315*6dd63270SBarry Smith gys = gxs; 316*6dd63270SBarry Smith gym = gxm; 317*6dd63270SBarry Smith gxs = 0; 318*6dd63270SBarry Smith gxm = dof; 319*6dd63270SBarry Smith } 320*6dd63270SBarry Smith *ierr = VecGetArrayRead(*v, &aa); 321*6dd63270SBarry Smith if (*ierr) return; 322*6dd63270SBarry Smith *ierr = F90Array2dCreate((void *)aa, MPIU_SCALAR, gxs, gxm, gys, gym, a PETSC_F90_2PTR_PARAM(ptrd)); 323*6dd63270SBarry Smith if (*ierr) return; 324*6dd63270SBarry Smith } 325*6dd63270SBarry Smith 326*6dd63270SBarry Smith PETSC_EXTERN void dmdavecrestorearrayreadf902_(DM *da, Vec *v, F90Array2d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 327*6dd63270SBarry Smith { 328*6dd63270SBarry Smith const PetscScalar *fa; 329*6dd63270SBarry Smith *ierr = F90Array2dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 330*6dd63270SBarry Smith *ierr = VecRestoreArrayRead(*v, &fa); 331*6dd63270SBarry Smith if (*ierr) return; 332*6dd63270SBarry Smith *ierr = F90Array2dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 333*6dd63270SBarry Smith } 334*6dd63270SBarry Smith 335*6dd63270SBarry Smith PETSC_EXTERN void dmdavecgetarrayreadf903_(DM *da, Vec *v, F90Array3d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 336*6dd63270SBarry Smith { 337*6dd63270SBarry Smith PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 338*6dd63270SBarry Smith const PetscScalar *aa; 339*6dd63270SBarry Smith 340*6dd63270SBarry Smith *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 341*6dd63270SBarry Smith if (*ierr) return; 342*6dd63270SBarry Smith *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 343*6dd63270SBarry Smith if (*ierr) return; 344*6dd63270SBarry Smith *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 345*6dd63270SBarry Smith if (*ierr) return; 346*6dd63270SBarry Smith 347*6dd63270SBarry Smith /* Handle case where user passes in global vector as opposed to local */ 348*6dd63270SBarry Smith *ierr = VecGetLocalSize(*v, &N); 349*6dd63270SBarry Smith if (*ierr) return; 350*6dd63270SBarry Smith if (N == xm * ym * zm * dof) { 351*6dd63270SBarry Smith gxm = xm; 352*6dd63270SBarry Smith gym = ym; 353*6dd63270SBarry Smith gzm = zm; 354*6dd63270SBarry Smith gxs = xs; 355*6dd63270SBarry Smith gys = ys; 356*6dd63270SBarry Smith gzs = zs; 357*6dd63270SBarry Smith } else if (N != gxm * gym * gzm * dof) { 358*6dd63270SBarry Smith *ierr = PETSC_ERR_ARG_INCOMP; 359*6dd63270SBarry Smith return; 360*6dd63270SBarry Smith } 361*6dd63270SBarry Smith if (dim == 2) { 362*6dd63270SBarry Smith gzs = gys; 363*6dd63270SBarry Smith gzm = gym; 364*6dd63270SBarry Smith gys = gxs; 365*6dd63270SBarry Smith gym = gxm; 366*6dd63270SBarry Smith gxs = 0; 367*6dd63270SBarry Smith gxm = dof; 368*6dd63270SBarry Smith } 369*6dd63270SBarry Smith *ierr = VecGetArrayRead(*v, &aa); 370*6dd63270SBarry Smith if (*ierr) return; 371*6dd63270SBarry Smith *ierr = F90Array3dCreate((void *)aa, MPIU_SCALAR, gxs, gxm, gys, gym, gzs, gzm, a PETSC_F90_2PTR_PARAM(ptrd)); 372*6dd63270SBarry Smith if (*ierr) return; 373*6dd63270SBarry Smith } 374*6dd63270SBarry Smith 375*6dd63270SBarry Smith PETSC_EXTERN void dmdavecrestorearrayreadf903_(DM *da, Vec *v, F90Array3d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 376*6dd63270SBarry Smith { 377*6dd63270SBarry Smith const PetscScalar *fa; 378*6dd63270SBarry Smith *ierr = F90Array3dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 379*6dd63270SBarry Smith *ierr = VecRestoreArrayRead(*v, &fa); 380*6dd63270SBarry Smith if (*ierr) return; 381*6dd63270SBarry Smith *ierr = F90Array3dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 382*6dd63270SBarry Smith } 383*6dd63270SBarry Smith 384*6dd63270SBarry Smith PETSC_EXTERN void dmdavecgetarrayreadf904_(DM *da, Vec *v, F90Array4d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 385*6dd63270SBarry Smith { 386*6dd63270SBarry Smith PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof, zero = 0; 387*6dd63270SBarry Smith const PetscScalar *aa; 388*6dd63270SBarry Smith 389*6dd63270SBarry Smith *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm); 390*6dd63270SBarry Smith if (*ierr) return; 391*6dd63270SBarry Smith *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm); 392*6dd63270SBarry Smith if (*ierr) return; 393*6dd63270SBarry Smith *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL); 394*6dd63270SBarry Smith if (*ierr) return; 395*6dd63270SBarry Smith 396*6dd63270SBarry Smith /* Handle case where user passes in global vector as opposed to local */ 397*6dd63270SBarry Smith *ierr = VecGetLocalSize(*v, &N); 398*6dd63270SBarry Smith if (*ierr) return; 399*6dd63270SBarry Smith if (N == xm * ym * zm * dof) { 400*6dd63270SBarry Smith gxm = xm; 401*6dd63270SBarry Smith gym = ym; 402*6dd63270SBarry Smith gzm = zm; 403*6dd63270SBarry Smith gxs = xs; 404*6dd63270SBarry Smith gys = ys; 405*6dd63270SBarry Smith gzs = zs; 406*6dd63270SBarry Smith } else if (N != gxm * gym * gzm * dof) { 407*6dd63270SBarry Smith *ierr = PETSC_ERR_ARG_INCOMP; 408*6dd63270SBarry Smith return; 409*6dd63270SBarry Smith } 410*6dd63270SBarry Smith *ierr = VecGetArrayRead(*v, &aa); 411*6dd63270SBarry Smith if (*ierr) return; 412*6dd63270SBarry Smith *ierr = F90Array4dCreate((void *)aa, MPIU_SCALAR, zero, dof, gxs, gxm, gys, gym, gzs, gzm, a PETSC_F90_2PTR_PARAM(ptrd)); 413*6dd63270SBarry Smith if (*ierr) return; 414*6dd63270SBarry Smith } 415*6dd63270SBarry Smith 416*6dd63270SBarry Smith PETSC_EXTERN void dmdavecrestorearrayreadf904_(DM *da, Vec *v, F90Array4d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd)) 417*6dd63270SBarry Smith { 418*6dd63270SBarry Smith const PetscScalar *fa; 419*6dd63270SBarry Smith /* 420*6dd63270SBarry Smith F90Array4dAccess is not implemented, so the following call would fail 421*6dd63270SBarry Smith */ 422*6dd63270SBarry Smith *ierr = F90Array4dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 423*6dd63270SBarry Smith *ierr = VecRestoreArrayRead(*v, &fa); 424*6dd63270SBarry Smith if (*ierr) return; 425*6dd63270SBarry Smith *ierr = F90Array4dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 426*6dd63270SBarry Smith } 427