Actual source code: zda1f90.c

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

  4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
  5:   #define dmdagetlocalinfof90_         DMDAGETLOCALINFOF90
  6:   #define dmdavecgetarrayf901_         DMDAVECGETARRAYF901
  7:   #define dmdavecrestorearrayf901_     DMDAVECRESTOREARRAYF901
  8:   #define dmdavecgetarrayf902_         DMDAVECGETARRAYF902
  9:   #define dmdavecrestorearrayf902_     DMDAVECRESTOREARRAYF902
 10:   #define dmdavecgetarrayf903_         DMDAVECGETARRAYF903
 11:   #define dmdavecrestorearrayf903_     DMDAVECRESTOREARRAYF903
 12:   #define dmdavecgetarrayf904_         DMDAVECGETARRAYF904
 13:   #define dmdavecrestorearrayf904_     DMDAVECRESTOREARRAYF904
 14:   #define dmdavecgetarrayreadf901_     DMDAVECGETARRAYREADF901
 15:   #define dmdavecrestorearrayreadf901_ DMDAVECRESTOREARRAYREADF901
 16:   #define dmdavecgetarrayreadf902_     DMDAVECGETARRAYREADF902
 17:   #define dmdavecrestorearrayreadf902_ DMDAVECRESTOREARRAYREADF902
 18:   #define dmdavecgetarrayreadf903_     DMDAVECGETARRAYREADF903
 19:   #define dmdavecrestorearrayreadf903_ DMDAVECRESTOREARRAYREADF903
 20:   #define dmdavecgetarrayreadf904_     DMDAVECGETARRAYREADF904
 21:   #define dmdavecrestorearrayreadf904_ DMDAVECRESTOREARRAYREADF904
 22:   #define dmdagetelements_             DMDAGETELEMENTS
 23:   #define dmdarestoreelements_         DMDARESTOREELEMENTS
 24: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 25:   #define dmdagetlocalinfof90_         dmdagetlocalinfof90
 26:   #define dmdavecgetarrayf901_         dmdavecgetarrayf901
 27:   #define dmdavecrestorearrayf901_     dmdavecrestorearrayf901
 28:   #define dmdavecgetarrayf902_         dmdavecgetarrayf902
 29:   #define dmdavecrestorearrayf902_     dmdavecrestorearrayf902
 30:   #define dmdavecgetarrayf903_         dmdavecgetarrayf903
 31:   #define dmdavecrestorearrayf903_     dmdavecrestorearrayf903
 32:   #define dmdavecgetarrayf904_         dmdavecgetarrayf904
 33:   #define dmdavecrestorearrayf904_     dmdavecrestorearrayf904
 34:   #define dmdavecgetarrayreadf901_     dmdavecgetarrayreadf901
 35:   #define dmdavecrestorearrayreadf901_ dmdavecrestorearrayreadf901
 36:   #define dmdavecgetarrayreadf902_     dmdavecgetarrayreadf902
 37:   #define dmdavecrestorearrayreadf902_ dmdavecrestorearrayreadf902
 38:   #define dmdavecgetarrayreadf903_     dmdavecgetarrayreadf903
 39:   #define dmdavecrestorearrayreadf903_ dmdavecrestorearrayreadf903
 40:   #define dmdavecgetarrayreadf904_     dmdavecgetarrayreadf904
 41:   #define dmdavecrestorearrayreadf904_ dmdavecrestorearrayreadf904
 42:   #define dmdagetelements_             dmdagetelements
 43:   #define dmdarestoreelements_         dmdarestoreelements
 44: #endif

 46: PETSC_EXTERN void dmdagetelements_(DM *dm, PetscInt *nel, PetscInt *nen, F90Array1d *e, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 47: {
 48:   const PetscInt *fa;

 50:   if (!e) {
 51:     *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?");
 52:     return;
 53:   }
 54:   *ierr = DMDAGetElements(*dm, nel, nen, &fa);
 55:   if (*ierr) return;
 56:   *ierr = F90Array1dCreate((PetscInt *)fa, MPIU_INT, 1, (*nel) * (*nen), e PETSC_F90_2PTR_PARAM(ptrd));
 57: }

 59: PETSC_EXTERN void dmdarestoreelements_(DM *dm, PetscInt *nel, PetscInt *nen, F90Array1d *e, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 60: {
 61:   if (!e) {
 62:     *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?");
 63:     return;
 64:   }
 65:   *ierr = F90Array1dDestroy(e, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
 66: }

 68: PETSC_EXTERN void dmdagetlocalinfof90_(DM *da, DMDALocalInfo *info, PetscErrorCode *ierr)
 69: {
 70:   *ierr = DMDAGetLocalInfo(*da, info);
 71: }

 73: PETSC_EXTERN void dmdavecgetarrayf901_(DM *da, Vec *v, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
 74: {
 75:   PetscInt     xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
 76:   PetscScalar *aa;

 78:   *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm);
 79:   if (*ierr) return;
 80:   *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
 81:   if (*ierr) return;
 82:   *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);
 83:   if (*ierr) return;

 85:   /* Handle case where user passes in global vector as opposed to local */
 86:   *ierr = VecGetLocalSize(*v, &N);
 87:   if (*ierr) return;
 88:   if (N == xm * ym * zm * dof) {
 89:     gxm = xm;
 90:     gym = ym;
 91:     gzm = zm;
 92:     gxs = xs;
 93:     gys = ys;
 94:     gzs = zs;
 95:   } else if (N != gxm * gym * gzm * dof) {
 96:     *ierr = PETSC_ERR_ARG_INCOMP;
 97:     return;
 98:   }
 99:   *ierr = VecGetArray(*v, &aa);
100:   if (*ierr) return;
101:   *ierr = F90Array1dCreate(aa, MPIU_SCALAR, gxs, gxm, a PETSC_F90_2PTR_PARAM(ptrd));
102:   if (*ierr) return;
103: }

105: PETSC_EXTERN void dmdavecrestorearrayf901_(DM *da, Vec *v, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
106: {
107:   PetscScalar *fa;
108:   *ierr = F90Array1dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
109:   *ierr = VecRestoreArray(*v, &fa);
110:   if (*ierr) return;
111:   *ierr = F90Array1dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
112: }

114: PETSC_EXTERN void dmdavecgetarrayf902_(DM *da, Vec *v, F90Array2d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
115: {
116:   PetscInt     xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
117:   PetscScalar *aa;

119:   *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm);
120:   if (*ierr) return;
121:   *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
122:   if (*ierr) return;
123:   *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);
124:   if (*ierr) return;

126:   /* Handle case where user passes in global vector as opposed to local */
127:   *ierr = VecGetLocalSize(*v, &N);
128:   if (*ierr) return;
129:   if (N == xm * ym * zm * dof) {
130:     gxm = xm;
131:     gym = ym;
132:     gzm = zm;
133:     gxs = xs;
134:     gys = ys;
135:     gzs = zs;
136:   } else if (N != gxm * gym * gzm * dof) {
137:     *ierr = PETSC_ERR_ARG_INCOMP;
138:     return;
139:   }
140:   if (dim == 1) {
141:     gys = gxs;
142:     gym = gxm;
143:     gxs = 0;
144:     gxm = dof;
145:   }
146:   *ierr = VecGetArray(*v, &aa);
147:   if (*ierr) return;
148:   *ierr = F90Array2dCreate(aa, MPIU_SCALAR, gxs, gxm, gys, gym, a PETSC_F90_2PTR_PARAM(ptrd));
149:   if (*ierr) return;
150: }

152: PETSC_EXTERN void dmdavecrestorearrayf902_(DM *da, Vec *v, F90Array2d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
153: {
154:   PetscScalar *fa;
155:   *ierr = F90Array2dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
156:   *ierr = VecRestoreArray(*v, &fa);
157:   if (*ierr) return;
158:   *ierr = F90Array2dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
159: }

161: PETSC_EXTERN void dmdavecgetarrayf903_(DM *da, Vec *v, F90Array3d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
162: {
163:   PetscInt     xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
164:   PetscScalar *aa;

166:   *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm);
167:   if (*ierr) return;
168:   *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
169:   if (*ierr) return;
170:   *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);
171:   if (*ierr) return;

173:   /* Handle case where user passes in global vector as opposed to local */
174:   *ierr = VecGetLocalSize(*v, &N);
175:   if (*ierr) return;
176:   if (N == xm * ym * zm * dof) {
177:     gxm = xm;
178:     gym = ym;
179:     gzm = zm;
180:     gxs = xs;
181:     gys = ys;
182:     gzs = zs;
183:   } else if (N != gxm * gym * gzm * dof) {
184:     *ierr = PETSC_ERR_ARG_INCOMP;
185:     return;
186:   }
187:   if (dim == 2) {
188:     gzs = gys;
189:     gzm = gym;
190:     gys = gxs;
191:     gym = gxm;
192:     gxs = 0;
193:     gxm = dof;
194:   }
195:   *ierr = VecGetArray(*v, &aa);
196:   if (*ierr) return;
197:   *ierr = F90Array3dCreate(aa, MPIU_SCALAR, gxs, gxm, gys, gym, gzs, gzm, a PETSC_F90_2PTR_PARAM(ptrd));
198:   if (*ierr) return;
199: }

201: PETSC_EXTERN void dmdavecrestorearrayf903_(DM *da, Vec *v, F90Array3d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
202: {
203:   PetscScalar *fa;
204:   *ierr = F90Array3dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
205:   *ierr = VecRestoreArray(*v, &fa);
206:   if (*ierr) return;
207:   *ierr = F90Array3dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
208: }

210: PETSC_EXTERN void dmdavecgetarrayf904_(DM *da, Vec *v, F90Array4d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
211: {
212:   PetscInt     xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof, zero = 0;
213:   PetscScalar *aa;

215:   *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm);
216:   if (*ierr) return;
217:   *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
218:   if (*ierr) return;
219:   *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);
220:   if (*ierr) return;

222:   /* Handle case where user passes in global vector as opposed to local */
223:   *ierr = VecGetLocalSize(*v, &N);
224:   if (*ierr) return;
225:   if (N == xm * ym * zm * dof) {
226:     gxm = xm;
227:     gym = ym;
228:     gzm = zm;
229:     gxs = xs;
230:     gys = ys;
231:     gzs = zs;
232:   } else if (N != gxm * gym * gzm * dof) {
233:     *ierr = PETSC_ERR_ARG_INCOMP;
234:     return;
235:   }
236:   *ierr = VecGetArray(*v, &aa);
237:   if (*ierr) return;
238:   *ierr = F90Array4dCreate(aa, MPIU_SCALAR, zero, dof, gxs, gxm, gys, gym, gzs, gzm, a PETSC_F90_2PTR_PARAM(ptrd));
239:   if (*ierr) return;
240: }

242: PETSC_EXTERN void dmdavecrestorearrayf904_(DM *da, Vec *v, F90Array4d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
243: {
244:   PetscScalar *fa;
245:   /*
246:     F90Array4dAccess is not implemented, so the following call would fail
247:   */
248:   *ierr = F90Array4dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
249:   *ierr = VecRestoreArray(*v, &fa);
250:   if (*ierr) return;
251:   *ierr = F90Array4dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
252: }

254: PETSC_EXTERN void dmdavecgetarrayreadf901_(DM *da, Vec *v, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
255: {
256:   PetscInt           xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
257:   const PetscScalar *aa;

259:   *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm);
260:   if (*ierr) return;
261:   *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
262:   if (*ierr) return;
263:   *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);
264:   if (*ierr) return;

266:   /* Handle case where user passes in global vector as opposed to local */
267:   *ierr = VecGetLocalSize(*v, &N);
268:   if (*ierr) return;
269:   if (N == xm * ym * zm * dof) {
270:     gxm = xm;
271:     gym = ym;
272:     gzm = zm;
273:     gxs = xs;
274:     gys = ys;
275:     gzs = zs;
276:   } else if (N != gxm * gym * gzm * dof) {
277:     *ierr = PETSC_ERR_ARG_INCOMP;
278:     return;
279:   }
280:   *ierr = VecGetArrayRead(*v, &aa);
281:   if (*ierr) return;
282:   *ierr = F90Array1dCreate((void *)aa, MPIU_SCALAR, gxs, gxm, a PETSC_F90_2PTR_PARAM(ptrd));
283:   if (*ierr) return;
284: }

286: PETSC_EXTERN void dmdavecrestorearrayreadf901_(DM *da, Vec *v, F90Array1d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
287: {
288:   const PetscScalar *fa;
289:   *ierr = F90Array1dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
290:   *ierr = VecRestoreArrayRead(*v, &fa);
291:   if (*ierr) return;
292:   *ierr = F90Array1dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
293: }

295: PETSC_EXTERN void dmdavecgetarrayreadf902_(DM *da, Vec *v, F90Array2d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
296: {
297:   PetscInt           xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
298:   const PetscScalar *aa;

300:   *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm);
301:   if (*ierr) return;
302:   *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
303:   if (*ierr) return;
304:   *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);
305:   if (*ierr) return;

307:   /* Handle case where user passes in global vector as opposed to local */
308:   *ierr = VecGetLocalSize(*v, &N);
309:   if (*ierr) return;
310:   if (N == xm * ym * zm * dof) {
311:     gxm = xm;
312:     gym = ym;
313:     gzm = zm;
314:     gxs = xs;
315:     gys = ys;
316:     gzs = zs;
317:   } else if (N != gxm * gym * gzm * dof) {
318:     *ierr = PETSC_ERR_ARG_INCOMP;
319:     return;
320:   }
321:   if (dim == 1) {
322:     gys = gxs;
323:     gym = gxm;
324:     gxs = 0;
325:     gxm = dof;
326:   }
327:   *ierr = VecGetArrayRead(*v, &aa);
328:   if (*ierr) return;
329:   *ierr = F90Array2dCreate((void *)aa, MPIU_SCALAR, gxs, gxm, gys, gym, a PETSC_F90_2PTR_PARAM(ptrd));
330:   if (*ierr) return;
331: }

333: PETSC_EXTERN void dmdavecrestorearrayreadf902_(DM *da, Vec *v, F90Array2d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
334: {
335:   const PetscScalar *fa;
336:   *ierr = F90Array2dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
337:   *ierr = VecRestoreArrayRead(*v, &fa);
338:   if (*ierr) return;
339:   *ierr = F90Array2dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
340: }

342: PETSC_EXTERN void dmdavecgetarrayreadf903_(DM *da, Vec *v, F90Array3d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
343: {
344:   PetscInt           xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
345:   const PetscScalar *aa;

347:   *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm);
348:   if (*ierr) return;
349:   *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
350:   if (*ierr) return;
351:   *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);
352:   if (*ierr) return;

354:   /* Handle case where user passes in global vector as opposed to local */
355:   *ierr = VecGetLocalSize(*v, &N);
356:   if (*ierr) return;
357:   if (N == xm * ym * zm * dof) {
358:     gxm = xm;
359:     gym = ym;
360:     gzm = zm;
361:     gxs = xs;
362:     gys = ys;
363:     gzs = zs;
364:   } else if (N != gxm * gym * gzm * dof) {
365:     *ierr = PETSC_ERR_ARG_INCOMP;
366:     return;
367:   }
368:   if (dim == 2) {
369:     gzs = gys;
370:     gzm = gym;
371:     gys = gxs;
372:     gym = gxm;
373:     gxs = 0;
374:     gxm = dof;
375:   }
376:   *ierr = VecGetArrayRead(*v, &aa);
377:   if (*ierr) return;
378:   *ierr = F90Array3dCreate((void *)aa, MPIU_SCALAR, gxs, gxm, gys, gym, gzs, gzm, a PETSC_F90_2PTR_PARAM(ptrd));
379:   if (*ierr) return;
380: }

382: PETSC_EXTERN void dmdavecrestorearrayreadf903_(DM *da, Vec *v, F90Array3d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
383: {
384:   const PetscScalar *fa;
385:   *ierr = F90Array3dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
386:   *ierr = VecRestoreArrayRead(*v, &fa);
387:   if (*ierr) return;
388:   *ierr = F90Array3dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
389: }

391: PETSC_EXTERN void dmdavecgetarrayreadf904_(DM *da, Vec *v, F90Array4d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
392: {
393:   PetscInt           xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof, zero = 0;
394:   const PetscScalar *aa;

396:   *ierr = DMDAGetCorners(*da, &xs, &ys, &zs, &xm, &ym, &zm);
397:   if (*ierr) return;
398:   *ierr = DMDAGetGhostCorners(*da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
399:   if (*ierr) return;
400:   *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);
401:   if (*ierr) return;

403:   /* Handle case where user passes in global vector as opposed to local */
404:   *ierr = VecGetLocalSize(*v, &N);
405:   if (*ierr) return;
406:   if (N == xm * ym * zm * dof) {
407:     gxm = xm;
408:     gym = ym;
409:     gzm = zm;
410:     gxs = xs;
411:     gys = ys;
412:     gzs = zs;
413:   } else if (N != gxm * gym * gzm * dof) {
414:     *ierr = PETSC_ERR_ARG_INCOMP;
415:     return;
416:   }
417:   *ierr = VecGetArrayRead(*v, &aa);
418:   if (*ierr) return;
419:   *ierr = F90Array4dCreate((void *)aa, MPIU_SCALAR, zero, dof, gxs, gxm, gys, gym, gzs, gzm, a PETSC_F90_2PTR_PARAM(ptrd));
420:   if (*ierr) return;
421: }

423: PETSC_EXTERN void dmdavecrestorearrayreadf904_(DM *da, Vec *v, F90Array4d *a, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
424: {
425:   const PetscScalar *fa;
426:   /*
427:     F90Array4dAccess is not implemented, so the following call would fail
428:   */
429:   *ierr = F90Array4dAccess(a, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd));
430:   *ierr = VecRestoreArrayRead(*v, &fa);
431:   if (*ierr) return;
432:   *ierr = F90Array4dDestroy(a, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd));
433: }