Actual source code: dalocal.c
1: /*
2: Code for manipulating distributed regular arrays in parallel.
3: */
5: #include <petsc/private/dmdaimpl.h>
6: #include <petscbt.h>
7: #include <petscsf.h>
8: #include <petscds.h>
9: #include <petscfe.h>
11: /*
12: This allows the DMDA vectors to properly tell MATLAB their dimensions
13: */
14: #if PetscDefined(HAVE_MATLAB)
15: #include <engine.h> /* MATLAB include file */
16: #include <mex.h> /* MATLAB include file */
17: static PetscErrorCode VecMatlabEnginePut_DA2d(PetscObject obj, void *mengine)
18: {
19: PetscInt n, m;
20: Vec vec = (Vec)obj;
21: PetscScalar *array;
22: mxArray *mat;
23: DM da;
25: PetscFunctionBegin;
26: PetscCall(VecGetDM(vec, &da));
27: PetscCheck(da, PetscObjectComm((PetscObject)vec), PETSC_ERR_ARG_WRONGSTATE, "Vector not associated with a DMDA");
28: PetscCall(DMDAGetGhostCorners(da, 0, 0, 0, &m, &n, 0));
30: PetscCall(VecGetArray(vec, &array));
31: #if !PetscDefined(USE_COMPLEX)
32: mat = mxCreateDoubleMatrix(m, n, mxREAL);
33: #else
34: mat = mxCreateDoubleMatrix(m, n, mxCOMPLEX);
35: #endif
36: PetscCall(PetscArraycpy(mxGetPr(mat), array, n * m));
37: PetscCall(PetscObjectName(obj));
38: engPutVariable((Engine *)mengine, obj->name, mat);
40: PetscCall(VecRestoreArray(vec, &array));
41: PetscFunctionReturn(PETSC_SUCCESS);
42: }
43: #endif
45: PetscErrorCode DMCreateLocalVector_DA(DM da, Vec *g)
46: {
47: DM_DA *dd = (DM_DA *)da->data;
49: PetscFunctionBegin;
51: PetscAssertPointer(g, 2);
52: PetscCall(VecCreate(PETSC_COMM_SELF, g));
53: PetscCall(VecSetSizes(*g, dd->nlocal, PETSC_DETERMINE));
54: PetscCall(VecSetBlockSize(*g, dd->w));
55: PetscCall(VecSetType(*g, da->vectype));
56: if (dd->nlocal < da->bind_below) {
57: PetscCall(VecSetBindingPropagates(*g, PETSC_TRUE));
58: PetscCall(VecBindToCPU(*g, PETSC_TRUE));
59: }
60: PetscCall(VecSetDM(*g, da));
61: #if PetscDefined(HAVE_MATLAB)
62: if (dd->w == 1 && da->dim == 2) PetscCall(PetscObjectComposeFunction((PetscObject)*g, "PetscMatlabEnginePut_C", VecMatlabEnginePut_DA2d));
63: #endif
64: PetscFunctionReturn(PETSC_SUCCESS);
65: }
67: /*@
68: DMDAGetNumCells - Get the number of cells (or vertices) in the local piece of the `DMDA`. This includes ghost cells.
70: Input Parameter:
71: . dm - The `DMDA` object
73: Output Parameters:
74: + numCellsX - The number of local cells in the x-direction
75: . numCellsY - The number of local cells in the y-direction
76: . numCellsZ - The number of local cells in the z-direction
77: - numCells - The number of local cells
79: Level: developer
81: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAGetCellPoint()`
82: @*/
83: PetscErrorCode DMDAGetNumCells(DM dm, PeOp PetscInt *numCellsX, PeOp PetscInt *numCellsY, PeOp PetscInt *numCellsZ, PeOp PetscInt *numCells)
84: {
85: DM_DA *da = (DM_DA *)dm->data;
86: const PetscInt dim = dm->dim;
87: const PetscInt mx = (da->Xe - da->Xs) / da->w, my = da->Ye - da->Ys, mz = da->Ze - da->Zs;
88: const PetscInt nC = (mx) * (dim > 1 ? (my) * (dim > 2 ? (mz) : 1) : 1);
90: PetscFunctionBegin;
92: if (numCellsX) {
93: PetscAssertPointer(numCellsX, 2);
94: *numCellsX = mx;
95: }
96: if (numCellsY) {
97: PetscAssertPointer(numCellsY, 3);
98: *numCellsY = my;
99: }
100: if (numCellsZ) {
101: PetscAssertPointer(numCellsZ, 4);
102: *numCellsZ = mz;
103: }
104: if (numCells) {
105: PetscAssertPointer(numCells, 5);
106: *numCells = nC;
107: }
108: PetscFunctionReturn(PETSC_SUCCESS);
109: }
111: /*@
112: DMDAGetCellPoint - Get the `DM` point corresponding to the tuple (i, j, k) in the `DMDA`
114: Input Parameters:
115: + dm - The `DMDA` object
116: . i - The global x index for the cell
117: . j - The global y index for the cell
118: - k - The global z index for the cell
120: Output Parameter:
121: . point - The local `DM` point
123: Level: developer
125: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAGetNumCells()`
126: @*/
127: PetscErrorCode DMDAGetCellPoint(DM dm, PetscInt i, PetscInt j, PetscInt k, PetscInt *point)
128: {
129: const PetscInt dim = dm->dim;
130: DMDALocalInfo info;
132: PetscFunctionBegin;
134: PetscAssertPointer(point, 5);
135: PetscCall(DMDAGetLocalInfo(dm, &info));
136: if (dim > 0) PetscCheck(!(i < info.gxs) && !(i >= info.gxs + info.gxm), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "X index %" PetscInt_FMT " not in [%" PetscInt_FMT ", %" PetscInt_FMT ")", i, info.gxs, info.gxs + info.gxm);
137: if (dim > 1) PetscCheck(!(j < info.gys) && !(j >= info.gys + info.gym), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Y index %" PetscInt_FMT " not in [%" PetscInt_FMT ", %" PetscInt_FMT ")", j, info.gys, info.gys + info.gym);
138: if (dim > 2) PetscCheck(!(k < info.gzs) && !(k >= info.gzs + info.gzm), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Z index %" PetscInt_FMT " not in [%" PetscInt_FMT ", %" PetscInt_FMT ")", k, info.gzs, info.gzs + info.gzm);
139: *point = i + (dim > 1 ? (j + (dim > 2 ? k * info.gym : 0)) * info.gxm : 0);
140: PetscFunctionReturn(PETSC_SUCCESS);
141: }
143: /*@
144: DMDAGetNumVertices - Return the number of local vertices (including ghost vertices) of a `DMDA` in each dimension and in total.
146: Not Collective
148: Input Parameter:
149: . dm - the `DMDA`
151: Output Parameters:
152: + numVerticesX - number of vertices in the x direction, or `NULL` if not needed
153: . numVerticesY - number of vertices in the y direction (`1` if `dim < 2`), or `NULL` if not needed
154: . numVerticesZ - number of vertices in the z direction (`1` if `dim < 3`), or `NULL` if not needed
155: - numVertices - total number of vertices, or `NULL` if not needed
157: Level: developer
159: .seealso: `DM`, `DMDA`, `DMDAGetNumCells()`, `DMDAGetNumFaces()`
160: @*/
161: PetscErrorCode DMDAGetNumVertices(DM dm, PetscInt *numVerticesX, PetscInt *numVerticesY, PetscInt *numVerticesZ, PetscInt *numVertices)
162: {
163: DM_DA *da = (DM_DA *)dm->data;
164: const PetscInt dim = dm->dim;
165: const PetscInt mx = (da->Xe - da->Xs) / da->w, my = da->Ye - da->Ys, mz = da->Ze - da->Zs;
166: const PetscInt nVx = mx + 1;
167: const PetscInt nVy = dim > 1 ? (my + 1) : 1;
168: const PetscInt nVz = dim > 2 ? (mz + 1) : 1;
169: const PetscInt nV = nVx * nVy * nVz;
171: PetscFunctionBegin;
172: if (numVerticesX) {
173: PetscAssertPointer(numVerticesX, 2);
174: *numVerticesX = nVx;
175: }
176: if (numVerticesY) {
177: PetscAssertPointer(numVerticesY, 3);
178: *numVerticesY = nVy;
179: }
180: if (numVerticesZ) {
181: PetscAssertPointer(numVerticesZ, 4);
182: *numVerticesZ = nVz;
183: }
184: if (numVertices) {
185: PetscAssertPointer(numVertices, 5);
186: *numVertices = nV;
187: }
188: PetscFunctionReturn(PETSC_SUCCESS);
189: }
191: /*@
192: DMDAGetNumFaces - Return the number of local mesh faces of each orientation (including ghost faces) for a `DMDA`.
194: Not Collective
196: Input Parameter:
197: . dm - the `DMDA`
199: Output Parameters:
200: + numXFacesX - number of X-normal faces along the x direction, or `NULL` if not needed
201: . numXFaces - total number of X-normal faces, or `NULL` if not needed
202: . numYFacesY - number of Y-normal faces along the y direction, or `NULL` if not needed
203: . numYFaces - total number of Y-normal faces (`0` for 1D), or `NULL` if not needed
204: . numZFacesZ - number of Z-normal faces along the z direction, or `NULL` if not needed
205: - numZFaces - total number of Z-normal faces (`0` for 1D/2D), or `NULL` if not needed
207: Level: developer
209: .seealso: `DM`, `DMDA`, `DMDAGetNumVertices()`, `DMDAGetNumCells()`
210: @*/
211: PetscErrorCode DMDAGetNumFaces(DM dm, PetscInt *numXFacesX, PetscInt *numXFaces, PetscInt *numYFacesY, PetscInt *numYFaces, PetscInt *numZFacesZ, PetscInt *numZFaces)
212: {
213: DM_DA *da = (DM_DA *)dm->data;
214: const PetscInt dim = dm->dim;
215: const PetscInt mx = (da->Xe - da->Xs) / da->w, my = da->Ye - da->Ys, mz = da->Ze - da->Zs;
216: const PetscInt nxF = (dim > 1 ? (my) * (dim > 2 ? (mz) : 1) : 1);
217: const PetscInt nXF = (mx + 1) * nxF;
218: const PetscInt nyF = mx * (dim > 2 ? mz : 1);
219: const PetscInt nYF = dim > 1 ? (my + 1) * nyF : 0;
220: const PetscInt nzF = mx * (dim > 1 ? my : 0);
221: const PetscInt nZF = dim > 2 ? (mz + 1) * nzF : 0;
223: PetscFunctionBegin;
224: if (numXFacesX) {
225: PetscAssertPointer(numXFacesX, 2);
226: *numXFacesX = nxF;
227: }
228: if (numXFaces) {
229: PetscAssertPointer(numXFaces, 3);
230: *numXFaces = nXF;
231: }
232: if (numYFacesY) {
233: PetscAssertPointer(numYFacesY, 4);
234: *numYFacesY = nyF;
235: }
236: if (numYFaces) {
237: PetscAssertPointer(numYFaces, 5);
238: *numYFaces = nYF;
239: }
240: if (numZFacesZ) {
241: PetscAssertPointer(numZFacesZ, 6);
242: *numZFacesZ = nzF;
243: }
244: if (numZFaces) {
245: PetscAssertPointer(numZFaces, 7);
246: *numZFaces = nZF;
247: }
248: PetscFunctionReturn(PETSC_SUCCESS);
249: }
251: /*@
252: DMDAGetHeightStratum - Get the bounds [`start`, `end`) for all points at a certain height.
254: Not Collective
256: Input Parameters:
257: + dm - The `DMDA` object
258: - height - The requested height
260: Output Parameters:
261: + pStart - The first point at this `height`
262: - pEnd - One beyond the last point at this `height`
264: Level: developer
266: Note:
267: See `DMPlexGetHeightStratum()` for the meaning of these values
269: .seealso: [](ch_unstructured), `DM`, `DMDA`, `DMPlexGetDepthStratum()`, `DMPlexGetHeightStratum()`, `DMPlexGetCellTypeStratum()`, `DMPlexGetDepth()`,
270: `DMPlexGetDepthLabel()`, `DMPlexGetPointDepth()`, `DMPlexSymmetrize()`, `DMPlexInterpolate()`, `DMDAGetDepthStratum()`
271: @*/
272: PetscErrorCode DMDAGetHeightStratum(DM dm, PetscInt height, PeOp PetscInt *pStart, PeOp PetscInt *pEnd)
273: {
274: const PetscInt dim = dm->dim;
275: PetscInt nC, nV, nXF, nYF, nZF;
277: PetscFunctionBegin;
278: if (pStart) PetscAssertPointer(pStart, 3);
279: if (pEnd) PetscAssertPointer(pEnd, 4);
280: PetscCall(DMDAGetNumCells(dm, NULL, NULL, NULL, &nC));
281: PetscCall(DMDAGetNumVertices(dm, NULL, NULL, NULL, &nV));
282: PetscCall(DMDAGetNumFaces(dm, NULL, &nXF, NULL, &nYF, NULL, &nZF));
283: if (height == 0) {
284: /* Cells */
285: if (pStart) *pStart = 0;
286: if (pEnd) *pEnd = nC;
287: } else if (height == 1) {
288: /* Faces */
289: if (pStart) *pStart = nC + nV;
290: if (pEnd) *pEnd = nC + nV + nXF + nYF + nZF;
291: } else if (height == dim) {
292: /* Vertices */
293: if (pStart) *pStart = nC;
294: if (pEnd) *pEnd = nC + nV;
295: } else if (height < 0) {
296: /* All points */
297: if (pStart) *pStart = 0;
298: if (pEnd) *pEnd = nC + nV + nXF + nYF + nZF;
299: } else SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "No points of height %" PetscInt_FMT " in the DA", height);
300: PetscFunctionReturn(PETSC_SUCCESS);
301: }
303: /*@
304: DMDAGetDepthStratum - Get the bounds [`start`, `end`) for all points at a certain depth.
306: Not Collective
308: Input Parameters:
309: + dm - The `DMDA` object
310: - depth - The requested depth
312: Output Parameters:
313: + pStart - The first point at this `depth`
314: - pEnd - One beyond the last point at this `depth`
316: Level: developer
318: Note:
319: See `DMPlexGetDepthStratum()` for the meaning of these values
321: .seealso: [](ch_unstructured), `DM`, `DMDA`, `DMPlexGetDepthStratum()`, `DMPlexGetHeightStratum()`, `DMPlexGetCellTypeStratum()`, `DMPlexGetDepth()`,
322: `DMPlexGetDepthLabel()`, `DMPlexGetPointDepth()`, `DMPlexSymmetrize()`, `DMPlexInterpolate()`, `DMDAGetHeightStratum()`
323: @*/
324: PetscErrorCode DMDAGetDepthStratum(DM dm, PetscInt depth, PeOp PetscInt *pStart, PeOp PetscInt *pEnd)
325: {
326: const PetscInt dim = dm->dim;
327: PetscInt nC, nV, nXF, nYF, nZF;
329: PetscFunctionBegin;
330: if (pStart) PetscAssertPointer(pStart, 3);
331: if (pEnd) PetscAssertPointer(pEnd, 4);
332: PetscCall(DMDAGetNumCells(dm, NULL, NULL, NULL, &nC));
333: PetscCall(DMDAGetNumVertices(dm, NULL, NULL, NULL, &nV));
334: PetscCall(DMDAGetNumFaces(dm, NULL, &nXF, NULL, &nYF, NULL, &nZF));
335: if (depth == dim) {
336: /* Cells */
337: if (pStart) *pStart = 0;
338: if (pEnd) *pEnd = nC;
339: } else if (depth == dim - 1) {
340: /* Faces */
341: if (pStart) *pStart = nC + nV;
342: if (pEnd) *pEnd = nC + nV + nXF + nYF + nZF;
343: } else if (depth == 0) {
344: /* Vertices */
345: if (pStart) *pStart = nC;
346: if (pEnd) *pEnd = nC + nV;
347: } else if (depth < 0) {
348: /* All points */
349: if (pStart) *pStart = 0;
350: if (pEnd) *pEnd = nC + nV + nXF + nYF + nZF;
351: } else SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "No points of depth %" PetscInt_FMT " in the DA", depth);
352: PetscFunctionReturn(PETSC_SUCCESS);
353: }
355: /*@
356: DMDASetVertexCoordinates - Sets the lower and upper coordinates for a `DMDA`
358: Logically Collective
360: Input Parameters:
361: + dm - The `DMDA` object
362: . xl - the lower x coordinate
363: . xu - the upper x coordinate
364: . yl - the lower y coordinate
365: . yu - the upper y coordinate
366: . zl - the lower z coordinate
367: - zu - the upper z coordinate
369: Level: intermediate
371: .seealso: [](ch_unstructured), `DM`, `DMDA`
372: @*/
373: PetscErrorCode DMDASetVertexCoordinates(DM dm, PetscReal xl, PetscReal xu, PetscReal yl, PetscReal yu, PetscReal zl, PetscReal zu)
374: {
375: DM_DA *da = (DM_DA *)dm->data;
376: Vec coordinates;
377: PetscSection section;
378: PetscScalar *coords;
379: PetscReal h[3];
380: PetscInt dim, size, M, N, P, nVx, nVy, nVz, nV, vStart, vEnd, v, i, j, k;
382: PetscFunctionBegin;
384: PetscCall(DMDAGetInfo(dm, &dim, &M, &N, &P, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL));
385: PetscCheck(dim <= 3, PetscObjectComm((PetscObject)dm), PETSC_ERR_PLIB, "The following code only works for dim <= 3");
386: h[0] = (xu - xl) / M;
387: h[1] = (yu - yl) / N;
388: h[2] = (zu - zl) / P;
389: PetscCall(DMDAGetDepthStratum(dm, 0, &vStart, &vEnd));
390: PetscCall(DMDAGetNumVertices(dm, &nVx, &nVy, &nVz, &nV));
391: PetscCall(PetscSectionCreate(PetscObjectComm((PetscObject)dm), §ion));
392: PetscCall(PetscSectionSetNumFields(section, 1));
393: PetscCall(PetscSectionSetFieldComponents(section, 0, dim));
394: PetscCall(PetscSectionSetChart(section, vStart, vEnd));
395: for (v = vStart; v < vEnd; ++v) PetscCall(PetscSectionSetDof(section, v, dim));
396: PetscCall(PetscSectionSetUp(section));
397: PetscCall(PetscSectionGetStorageSize(section, &size));
398: PetscCall(VecCreateSeq(PETSC_COMM_SELF, size, &coordinates));
399: PetscCall(PetscObjectSetName((PetscObject)coordinates, "coordinates"));
400: PetscCall(VecGetArray(coordinates, &coords));
401: for (k = 0; k < nVz; ++k) {
402: PetscInt ind[3], d, off;
404: ind[0] = 0;
405: ind[1] = 0;
406: ind[2] = k + da->zs;
407: for (j = 0; j < nVy; ++j) {
408: ind[1] = j + da->ys;
409: for (i = 0; i < nVx; ++i) {
410: const PetscInt vertex = (k * nVy + j) * nVx + i + vStart;
412: PetscCall(PetscSectionGetOffset(section, vertex, &off));
413: ind[0] = i + da->xs;
414: for (d = 0; d < dim; ++d) coords[off + d] = h[d] * ind[d];
415: }
416: }
417: }
418: PetscCall(VecRestoreArray(coordinates, &coords));
419: PetscCall(DMSetCoordinateSection(dm, PETSC_DETERMINE, section));
420: PetscCall(DMSetCoordinatesLocal(dm, coordinates));
421: PetscCall(PetscSectionDestroy(§ion));
422: PetscCall(VecDestroy(&coordinates));
423: PetscFunctionReturn(PETSC_SUCCESS);
424: }
426: /*@C
427: DMDAGetArray - Gets a work array for a `DMDA`
429: Input Parameters:
430: + da - a `DMDA`
431: - ghosted - do you want arrays for the ghosted or nonghosted patch
433: Output Parameter:
434: . vptr - array data structured
436: Level: advanced
438: Notes:
439: The vector values are NOT initialized and may have garbage in them, so you may need
440: to zero them.
442: Use `DMDARestoreArray()` to return the array
444: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDARestoreArray()`
445: @*/
446: PetscErrorCode DMDAGetArray(DM da, PetscBool ghosted, void *vptr)
447: {
448: PetscInt j, i, xs, ys, xm, ym, zs, zm;
449: char *iarray_start;
450: void **iptr = (void **)vptr;
451: DM_DA *dd = (DM_DA *)da->data;
453: PetscFunctionBegin;
455: if (ghosted) {
456: for (i = 0; i < DMDA_MAX_WORK_ARRAYS; i++) {
457: if (dd->arrayghostedin[i]) {
458: *iptr = dd->arrayghostedin[i];
459: iarray_start = (char *)dd->startghostedin[i];
460: dd->arrayghostedin[i] = NULL;
461: dd->startghostedin[i] = NULL;
463: goto done;
464: }
465: }
466: xs = dd->Xs;
467: ys = dd->Ys;
468: zs = dd->Zs;
469: xm = dd->Xe - dd->Xs;
470: ym = dd->Ye - dd->Ys;
471: zm = dd->Ze - dd->Zs;
472: } else {
473: for (i = 0; i < DMDA_MAX_WORK_ARRAYS; i++) {
474: if (dd->arrayin[i]) {
475: *iptr = dd->arrayin[i];
476: iarray_start = (char *)dd->startin[i];
477: dd->arrayin[i] = NULL;
478: dd->startin[i] = NULL;
480: goto done;
481: }
482: }
483: xs = dd->xs;
484: ys = dd->ys;
485: zs = dd->zs;
486: xm = dd->xe - dd->xs;
487: ym = dd->ye - dd->ys;
488: zm = dd->ze - dd->zs;
489: }
491: switch (da->dim) {
492: case 1: {
493: void *ptr;
495: PetscCall(PetscMalloc(xm * sizeof(PetscScalar), &iarray_start));
497: ptr = (void *)((PetscScalar *)iarray_start - xs);
498: *iptr = ptr;
499: break;
500: }
501: case 2: {
502: void **ptr;
504: PetscCall(PetscMalloc((ym + 1) * sizeof(void *) + xm * ym * sizeof(PetscScalar), &iarray_start));
506: ptr = (void **)(iarray_start + xm * ym * sizeof(PetscScalar) - ys * sizeof(void *));
507: for (j = ys; j < ys + ym; j++) ptr[j] = iarray_start + sizeof(PetscScalar) * (xm * (j - ys) - xs);
508: *iptr = (void *)ptr;
509: break;
510: }
511: case 3: {
512: void ***ptr, **bptr;
514: PetscCall(PetscMalloc((zm + 1) * sizeof(void **) + (ym * zm + 1) * sizeof(void *) + xm * ym * zm * sizeof(PetscScalar), &iarray_start));
516: ptr = (void ***)(iarray_start + xm * ym * zm * sizeof(PetscScalar) - zs * sizeof(void *));
517: bptr = (void **)(iarray_start + xm * ym * zm * sizeof(PetscScalar) + zm * sizeof(void **));
518: for (i = zs; i < zs + zm; i++) ptr[i] = bptr + ((i - zs) * ym - ys);
519: for (i = zs; i < zs + zm; i++) {
520: for (j = ys; j < ys + ym; j++) ptr[i][j] = iarray_start + sizeof(PetscScalar) * (xm * ym * (i - zs) + xm * (j - ys) - xs);
521: }
522: *iptr = (void *)ptr;
523: break;
524: }
525: default:
526: SETERRQ(PetscObjectComm((PetscObject)da), PETSC_ERR_SUP, "Dimension %" PetscInt_FMT " not supported", da->dim);
527: }
529: done:
530: /* add arrays to the checked out list */
531: if (ghosted) {
532: for (i = 0; i < DMDA_MAX_WORK_ARRAYS; i++) {
533: if (!dd->arrayghostedout[i]) {
534: dd->arrayghostedout[i] = *iptr;
535: dd->startghostedout[i] = iarray_start;
536: break;
537: }
538: }
539: } else {
540: for (i = 0; i < DMDA_MAX_WORK_ARRAYS; i++) {
541: if (!dd->arrayout[i]) {
542: dd->arrayout[i] = *iptr;
543: dd->startout[i] = iarray_start;
544: break;
545: }
546: }
547: }
548: PetscFunctionReturn(PETSC_SUCCESS);
549: }
551: /*@C
552: DMDARestoreArray - Restores an array for a `DMDA` obtained with `DMDAGetArray()`
554: Input Parameters:
555: + da - information about my local patch
556: . ghosted - do you want arrays for the ghosted or nonghosted patch
557: - vptr - array data structured
559: Level: advanced
561: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAGetArray()`
562: @*/
563: PetscErrorCode DMDARestoreArray(DM da, PetscBool ghosted, void *vptr)
564: {
565: void **iptr = (void **)vptr, *iarray_start = NULL;
566: DM_DA *dd = (DM_DA *)da->data;
568: PetscFunctionBegin;
570: if (ghosted) {
571: for (PetscInt i = 0; i < DMDA_MAX_WORK_ARRAYS; i++) {
572: if (dd->arrayghostedout[i] == *iptr) {
573: iarray_start = dd->startghostedout[i];
574: dd->arrayghostedout[i] = NULL;
575: dd->startghostedout[i] = NULL;
576: break;
577: }
578: }
579: for (PetscInt i = 0; i < DMDA_MAX_WORK_ARRAYS; i++) {
580: if (!dd->arrayghostedin[i]) {
581: dd->arrayghostedin[i] = *iptr;
582: dd->startghostedin[i] = iarray_start;
583: break;
584: }
585: }
586: } else {
587: for (PetscInt i = 0; i < DMDA_MAX_WORK_ARRAYS; i++) {
588: if (dd->arrayout[i] == *iptr) {
589: iarray_start = dd->startout[i];
590: dd->arrayout[i] = NULL;
591: dd->startout[i] = NULL;
592: break;
593: }
594: }
595: for (PetscInt i = 0; i < DMDA_MAX_WORK_ARRAYS; i++) {
596: if (!dd->arrayin[i]) {
597: dd->arrayin[i] = *iptr;
598: dd->startin[i] = iarray_start;
599: break;
600: }
601: }
602: }
603: PetscFunctionReturn(PETSC_SUCCESS);
604: }