1 #include <petsc/private/dmdaimpl.h> /*I "petscdmda.h" I*/
2
DMDAGetElements_1D(DM dm,PetscInt * nel,PetscInt * nen,const PetscInt * e[])3 static PetscErrorCode DMDAGetElements_1D(DM dm, PetscInt *nel, PetscInt *nen, const PetscInt *e[])
4 {
5 DM_DA *da = (DM_DA *)dm->data;
6 PetscInt i, xs, xe, Xs, Xe;
7 PetscInt cnt = 0;
8
9 PetscFunctionBegin;
10 if (!da->e) {
11 PetscInt corners[2];
12
13 PetscCheck(da->s, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Cannot get elements for DMDA with zero stencil width");
14 PetscCall(DMDAGetCorners(dm, &xs, NULL, NULL, &xe, NULL, NULL));
15 PetscCall(DMDAGetGhostCorners(dm, &Xs, NULL, NULL, &Xe, NULL, NULL));
16 xe += xs;
17 Xe += Xs;
18 if (xs != Xs) xs -= 1;
19 da->ne = 1 * (xe - xs - 1);
20 PetscCall(PetscMalloc1(1 + 2 * da->ne, &da->e));
21 for (i = xs; i < xe - 1; i++) {
22 da->e[cnt++] = (i - Xs);
23 da->e[cnt++] = (i - Xs + 1);
24 }
25 da->nen = 2;
26
27 corners[0] = (xs - Xs);
28 corners[1] = (xe - 1 - Xs);
29 PetscCall(ISCreateGeneral(PETSC_COMM_SELF, 2, corners, PETSC_COPY_VALUES, &da->ecorners));
30 }
31 *nel = da->ne;
32 *nen = da->nen;
33 *e = da->e;
34 PetscFunctionReturn(PETSC_SUCCESS);
35 }
36
DMDAGetElements_2D(DM dm,PetscInt * nel,PetscInt * nen,const PetscInt * e[])37 static PetscErrorCode DMDAGetElements_2D(DM dm, PetscInt *nel, PetscInt *nen, const PetscInt *e[])
38 {
39 DM_DA *da = (DM_DA *)dm->data;
40 PetscInt i, xs, xe, Xs, Xe;
41 PetscInt j, ys, ye, Ys, Ye;
42 PetscInt cnt = 0, cell[4], ns = 2;
43 PetscInt c, split[] = {0, 1, 3, 2, 3, 1};
44
45 PetscFunctionBegin;
46 if (!da->e) {
47 PetscInt corners[4], nn = 0;
48
49 PetscCheck(da->s, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Cannot get elements for DMDA with zero stencil width");
50
51 switch (da->elementtype) {
52 case DMDA_ELEMENT_Q1:
53 da->nen = 4;
54 break;
55 case DMDA_ELEMENT_P1:
56 da->nen = 3;
57 break;
58 default:
59 SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Unknown element type %d", da->elementtype);
60 }
61 nn = da->nen;
62
63 if (da->elementtype == DMDA_ELEMENT_P1) ns = 2;
64 if (da->elementtype == DMDA_ELEMENT_Q1) ns = 1;
65 PetscCall(DMDAGetCorners(dm, &xs, &ys, NULL, &xe, &ye, NULL));
66 PetscCall(DMDAGetGhostCorners(dm, &Xs, &Ys, NULL, &Xe, &Ye, NULL));
67 xe += xs;
68 Xe += Xs;
69 if (xs != Xs) xs -= 1;
70 ye += ys;
71 Ye += Ys;
72 if (ys != Ys) ys -= 1;
73 da->ne = ns * (xe - xs - 1) * (ye - ys - 1);
74 PetscCall(PetscMalloc1(1 + nn * da->ne, &da->e));
75 for (j = ys; j < ye - 1; j++) {
76 for (i = xs; i < xe - 1; i++) {
77 cell[0] = (i - Xs) + (j - Ys) * (Xe - Xs);
78 cell[1] = (i - Xs + 1) + (j - Ys) * (Xe - Xs);
79 cell[2] = (i - Xs + 1) + (j - Ys + 1) * (Xe - Xs);
80 cell[3] = (i - Xs) + (j - Ys + 1) * (Xe - Xs);
81 if (da->elementtype == DMDA_ELEMENT_P1) {
82 for (c = 0; c < ns * nn; c++) da->e[cnt++] = cell[split[c]];
83 }
84 if (da->elementtype == DMDA_ELEMENT_Q1) {
85 for (c = 0; c < ns * nn; c++) da->e[cnt++] = cell[c];
86 }
87 }
88 }
89
90 corners[0] = (xs - Xs) + (ys - Ys) * (Xe - Xs);
91 corners[1] = (xe - 1 - Xs) + (ys - Ys) * (Xe - Xs);
92 corners[2] = (xs - Xs) + (ye - 1 - Ys) * (Xe - Xs);
93 corners[3] = (xe - 1 - Xs) + (ye - 1 - Ys) * (Xe - Xs);
94 PetscCall(ISCreateGeneral(PETSC_COMM_SELF, 4, corners, PETSC_COPY_VALUES, &da->ecorners));
95 }
96 *nel = da->ne;
97 *nen = da->nen;
98 *e = da->e;
99 PetscFunctionReturn(PETSC_SUCCESS);
100 }
101
DMDAGetElements_3D(DM dm,PetscInt * nel,PetscInt * nen,const PetscInt * e[])102 static PetscErrorCode DMDAGetElements_3D(DM dm, PetscInt *nel, PetscInt *nen, const PetscInt *e[])
103 {
104 DM_DA *da = (DM_DA *)dm->data;
105 PetscInt i, xs, xe, Xs, Xe;
106 PetscInt j, ys, ye, Ys, Ye;
107 PetscInt k, zs, ze, Zs, Ze;
108 PetscInt cnt = 0, cell[8], ns = 6;
109 PetscInt c, split[] = {0, 1, 3, 7, 0, 1, 7, 4, 1, 2, 3, 7, 1, 2, 7, 6, 1, 4, 5, 7, 1, 5, 6, 7};
110
111 PetscFunctionBegin;
112 if (!da->e) {
113 PetscInt corners[8], nn = 0;
114
115 PetscCheck(da->s, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Cannot get elements for DMDA with zero stencil width");
116
117 switch (da->elementtype) {
118 case DMDA_ELEMENT_Q1:
119 da->nen = 8;
120 break;
121 case DMDA_ELEMENT_P1:
122 da->nen = 4;
123 break;
124 default:
125 SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Unknown element type %d", da->elementtype);
126 }
127 nn = da->nen;
128
129 if (da->elementtype == DMDA_ELEMENT_P1) ns = 6;
130 if (da->elementtype == DMDA_ELEMENT_Q1) ns = 1;
131 PetscCall(DMDAGetCorners(dm, &xs, &ys, &zs, &xe, &ye, &ze));
132 PetscCall(DMDAGetGhostCorners(dm, &Xs, &Ys, &Zs, &Xe, &Ye, &Ze));
133 xe += xs;
134 Xe += Xs;
135 if (xs != Xs) xs -= 1;
136 ye += ys;
137 Ye += Ys;
138 if (ys != Ys) ys -= 1;
139 ze += zs;
140 Ze += Zs;
141 if (zs != Zs) zs -= 1;
142 da->ne = ns * (xe - xs - 1) * (ye - ys - 1) * (ze - zs - 1);
143 PetscCall(PetscMalloc1(1 + nn * da->ne, &da->e));
144 for (k = zs; k < ze - 1; k++) {
145 for (j = ys; j < ye - 1; j++) {
146 for (i = xs; i < xe - 1; i++) {
147 cell[0] = (i - Xs) + (j - Ys) * (Xe - Xs) + (k - Zs) * (Xe - Xs) * (Ye - Ys);
148 cell[1] = (i + 1 - Xs) + (j - Ys) * (Xe - Xs) + (k - Zs) * (Xe - Xs) * (Ye - Ys);
149 cell[2] = (i + 1 - Xs) + (j + 1 - Ys) * (Xe - Xs) + (k - Zs) * (Xe - Xs) * (Ye - Ys);
150 cell[3] = (i - Xs) + (j + 1 - Ys) * (Xe - Xs) + (k - Zs) * (Xe - Xs) * (Ye - Ys);
151 cell[4] = (i - Xs) + (j - Ys) * (Xe - Xs) + (k + 1 - Zs) * (Xe - Xs) * (Ye - Ys);
152 cell[5] = (i + 1 - Xs) + (j - Ys) * (Xe - Xs) + (k + 1 - Zs) * (Xe - Xs) * (Ye - Ys);
153 cell[6] = (i + 1 - Xs) + (j + 1 - Ys) * (Xe - Xs) + (k + 1 - Zs) * (Xe - Xs) * (Ye - Ys);
154 cell[7] = (i - Xs) + (j + 1 - Ys) * (Xe - Xs) + (k + 1 - Zs) * (Xe - Xs) * (Ye - Ys);
155 if (da->elementtype == DMDA_ELEMENT_P1) {
156 for (c = 0; c < ns * nn; c++) da->e[cnt++] = cell[split[c]];
157 }
158 if (da->elementtype == DMDA_ELEMENT_Q1) {
159 for (c = 0; c < ns * nn; c++) da->e[cnt++] = cell[c];
160 }
161 }
162 }
163 }
164
165 corners[0] = (xs - Xs) + (ys - Ys) * (Xe - Xs) + (zs - Zs) * (Xe - Xs) * (Ye - Ys);
166 corners[1] = (xe - 1 - Xs) + (ys - Ys) * (Xe - Xs) + (zs - Zs) * (Xe - Xs) * (Ye - Ys);
167 corners[2] = (xs - Xs) + (ye - 1 - Ys) * (Xe - Xs) + (zs - Zs) * (Xe - Xs) * (Ye - Ys);
168 corners[3] = (xe - 1 - Xs) + (ye - 1 - Ys) * (Xe - Xs) + (zs - Zs) * (Xe - Xs) * (Ye - Ys);
169 corners[4] = (xs - Xs) + (ys - Ys) * (Xe - Xs) + (ze - 1 - Zs) * (Xe - Xs) * (Ye - Ys);
170 corners[5] = (xe - 1 - Xs) + (ys - Ys) * (Xe - Xs) + (ze - 1 - Zs) * (Xe - Xs) * (Ye - Ys);
171 corners[6] = (xs - Xs) + (ye - 1 - Ys) * (Xe - Xs) + (ze - 1 - Zs) * (Xe - Xs) * (Ye - Ys);
172 corners[7] = (xe - 1 - Xs) + (ye - 1 - Ys) * (Xe - Xs) + (ze - 1 - Zs) * (Xe - Xs) * (Ye - Ys);
173 PetscCall(ISCreateGeneral(PETSC_COMM_SELF, 8, corners, PETSC_COPY_VALUES, &da->ecorners));
174 }
175 *nel = da->ne;
176 *nen = da->nen;
177 *e = da->e;
178 PetscFunctionReturn(PETSC_SUCCESS);
179 }
180
181 /*@
182 DMDAGetElementsCorners - Returns the global (i,j,k) indices of the lower left
183 corner of the non-overlapping decomposition of elements identified by `DMDAGetElements()`
184
185 Not Collective
186
187 Input Parameter:
188 . da - the `DMDA` object
189
190 Output Parameters:
191 + gx - the i index
192 . gy - the j index
193 - gz - the k index
194
195 Level: intermediate
196
197 .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAElementType`, `DMDASetElementType()`, `DMDAGetElements()`, `DMDAGetCorners()`, `DMDAGetGhostCorners()`, `DMDAGetElementsSizes()`,
198 `DMDAGetElementsCornersIS()`, `DMDARestoreElementsCornersIS()`
199 @*/
DMDAGetElementsCorners(DM da,PeOp PetscInt * gx,PeOp PetscInt * gy,PeOp PetscInt * gz)200 PetscErrorCode DMDAGetElementsCorners(DM da, PeOp PetscInt *gx, PeOp PetscInt *gy, PeOp PetscInt *gz)
201 {
202 PetscInt xs, Xs;
203 PetscInt ys, Ys;
204 PetscInt zs, Zs;
205 PetscBool isda;
206
207 PetscFunctionBegin;
208 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA);
209 if (gx) PetscAssertPointer(gx, 2);
210 if (gy) PetscAssertPointer(gy, 3);
211 if (gz) PetscAssertPointer(gz, 4);
212 PetscCall(PetscObjectTypeCompare((PetscObject)da, DMDA, &isda));
213 PetscCheck(isda, PetscObjectComm((PetscObject)da), PETSC_ERR_USER, "Not for DM type %s", ((PetscObject)da)->type_name);
214 PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, NULL, NULL, NULL));
215 PetscCall(DMDAGetGhostCorners(da, &Xs, &Ys, &Zs, NULL, NULL, NULL));
216 if (xs != Xs) xs -= 1;
217 if (ys != Ys) ys -= 1;
218 if (zs != Zs) zs -= 1;
219 if (gx) *gx = xs;
220 if (gy) *gy = ys;
221 if (gz) *gz = zs;
222 PetscFunctionReturn(PETSC_SUCCESS);
223 }
224
225 /*@
226 DMDAGetElementsSizes - Gets the local number of elements per coordinate direction for the non-overlapping decomposition identified by `DMDAGetElements()`
227
228 Not Collective
229
230 Input Parameter:
231 . da - the `DMDA` object
232
233 Output Parameters:
234 + mx - number of local elements in x-direction
235 . my - number of local elements in y-direction
236 - mz - number of local elements in z-direction
237
238 Level: intermediate
239
240 Note:
241 Returns the same number of elements, irrespective of the `DMDAElementType`
242
243 .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAElementType`, `DMDASetElementType()`, `DMDAGetElements()`, `DMDAGetElementsCorners()`
244 @*/
DMDAGetElementsSizes(DM da,PeOp PetscInt * mx,PeOp PetscInt * my,PeOp PetscInt * mz)245 PetscErrorCode DMDAGetElementsSizes(DM da, PeOp PetscInt *mx, PeOp PetscInt *my, PeOp PetscInt *mz)
246 {
247 PetscInt xs, xe, Xs;
248 PetscInt ys, ye, Ys;
249 PetscInt zs, ze, Zs;
250 PetscInt dim;
251 PetscBool isda;
252
253 PetscFunctionBegin;
254 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA);
255 if (mx) PetscAssertPointer(mx, 2);
256 if (my) PetscAssertPointer(my, 3);
257 if (mz) PetscAssertPointer(mz, 4);
258 PetscCall(PetscObjectTypeCompare((PetscObject)da, DMDA, &isda));
259 PetscCheck(isda, PetscObjectComm((PetscObject)da), PETSC_ERR_USER, "Not for DM type %s", ((PetscObject)da)->type_name);
260 PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xe, &ye, &ze));
261 PetscCall(DMDAGetGhostCorners(da, &Xs, &Ys, &Zs, NULL, NULL, NULL));
262 xe += xs;
263 if (xs != Xs) xs -= 1;
264 ye += ys;
265 if (ys != Ys) ys -= 1;
266 ze += zs;
267 if (zs != Zs) zs -= 1;
268 if (mx) *mx = 0;
269 if (my) *my = 0;
270 if (mz) *mz = 0;
271 PetscCall(DMGetDimension(da, &dim));
272 switch (dim) {
273 case 3:
274 if (mz) *mz = ze - zs - 1; /* fall through */
275 case 2:
276 if (my) *my = ye - ys - 1; /* fall through */
277 case 1:
278 if (mx) *mx = xe - xs - 1;
279 break;
280 }
281 PetscFunctionReturn(PETSC_SUCCESS);
282 }
283
284 /*@
285 DMDASetElementType - Sets the element type to be returned by `DMDAGetElements()`
286
287 Not Collective
288
289 Input Parameter:
290 . da - the `DMDA` object
291
292 Output Parameter:
293 . etype - the element type, currently either `DMDA_ELEMENT_P1` or `DMDA_ELEMENT_Q1`
294
295 Level: intermediate
296
297 .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAElementType`, `DMDAGetElementType()`, `DMDAGetElements()`, `DMDARestoreElements()`,
298 `DMDA_ELEMENT_P1`, `DMDA_ELEMENT_Q1`
299 @*/
DMDASetElementType(DM da,DMDAElementType etype)300 PetscErrorCode DMDASetElementType(DM da, DMDAElementType etype)
301 {
302 DM_DA *dd = (DM_DA *)da->data;
303 PetscBool isda;
304
305 PetscFunctionBegin;
306 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA);
307 PetscValidLogicalCollectiveEnum(da, etype, 2);
308 PetscCall(PetscObjectTypeCompare((PetscObject)da, DMDA, &isda));
309 if (!isda) PetscFunctionReturn(PETSC_SUCCESS);
310 if (dd->elementtype != etype) {
311 PetscCall(PetscFree(dd->e));
312 PetscCall(ISDestroy(&dd->ecorners));
313
314 dd->elementtype = etype;
315 dd->ne = 0;
316 dd->nen = 0;
317 dd->e = NULL;
318 }
319 PetscFunctionReturn(PETSC_SUCCESS);
320 }
321
322 /*@
323 DMDAGetElementType - Gets the element type to be returned by `DMDAGetElements()`
324
325 Not Collective
326
327 Input Parameter:
328 . da - the `DMDA` object
329
330 Output Parameter:
331 . etype - the element type, currently either `DMDA_ELEMENT_P1` or `DMDA_ELEMENT_Q1`
332
333 Level: intermediate
334
335 .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAElementType`, `DMDASetElementType()`, `DMDAGetElements()`, `DMDARestoreElements()`,
336 `DMDA_ELEMENT_P1`, `DMDA_ELEMENT_Q1`
337 @*/
DMDAGetElementType(DM da,DMDAElementType * etype)338 PetscErrorCode DMDAGetElementType(DM da, DMDAElementType *etype)
339 {
340 DM_DA *dd = (DM_DA *)da->data;
341 PetscBool isda;
342
343 PetscFunctionBegin;
344 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA);
345 PetscAssertPointer(etype, 2);
346 PetscCall(PetscObjectTypeCompare((PetscObject)da, DMDA, &isda));
347 PetscCheck(isda, PetscObjectComm((PetscObject)da), PETSC_ERR_USER, "Not for DM type %s", ((PetscObject)da)->type_name);
348 *etype = dd->elementtype;
349 PetscFunctionReturn(PETSC_SUCCESS);
350 }
351
352 /*@C
353 DMDAGetElements - Gets an array containing the indices (in local indexing)
354 of all the local elements
355
356 Not Collective
357
358 Input Parameter:
359 . dm - the `DMDA` object
360
361 Output Parameters:
362 + nel - number of local elements
363 . nen - number of nodes in each element (for example in one dimension it is 2, in two dimensions it is 3 (for `DMDA_ELEMENT_P1`) and 4
364 (for `DMDA_ELEMENT_Q1`)
365 - e - the local indices of the elements' vertices, of length `nel` * `nen`
366
367 Level: intermediate
368
369 Notes:
370 Call `DMDARestoreElements()` once you have finished accessing the elements.
371
372 Each process uniquely owns a subset of the elements. That is no element is owned by two or more processes.
373
374 If on each process you integrate over its owned elements and use `ADD_VALUES` in `Vec`/`MatSetValuesLocal()` then you'll obtain the correct result.
375
376 Fortran Note:
377 Use
378 .vb
379 PetscScalar, pointer :: e(:)
380 .ve
381 to declare the element array
382
383 .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAElementType`, `DMDASetElementType()`, `VecSetValuesLocal()`, `MatSetValuesLocal()`,
384 `DMGlobalToLocalBegin()`, `DMLocalToGlobalBegin()`, `DMDARestoreElements()`, `DMDA_ELEMENT_P1`, `DMDA_ELEMENT_Q1`, `DMDAGetElementsSizes()`,
385 `DMDAGetElementsCorners()`
386 @*/
DMDAGetElements(DM dm,PetscInt * nel,PetscInt * nen,const PetscInt * e[])387 PetscErrorCode DMDAGetElements(DM dm, PetscInt *nel, PetscInt *nen, const PetscInt *e[])
388 {
389 PetscInt dim;
390 DM_DA *dd = (DM_DA *)dm->data;
391 PetscBool isda;
392
393 PetscFunctionBegin;
394 PetscValidHeaderSpecificType(dm, DM_CLASSID, 1, DMDA);
395 PetscAssertPointer(nel, 2);
396 PetscAssertPointer(nen, 3);
397 PetscAssertPointer(e, 4);
398 PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMDA, &isda));
399 PetscCheck(isda, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for DM type %s", ((PetscObject)dm)->type_name);
400 PetscCheck(dd->stencil_type != DMDA_STENCIL_STAR, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "DMDAGetElements() requires you use a stencil type of DMDA_STENCIL_BOX");
401 PetscCall(DMGetDimension(dm, &dim));
402 if (dd->e) {
403 *nel = dd->ne;
404 *nen = dd->nen;
405 *e = dd->e;
406 PetscFunctionReturn(PETSC_SUCCESS);
407 }
408 if (dim == -1) {
409 *nel = 0;
410 *nen = 0;
411 *e = NULL;
412 } else if (dim == 1) {
413 PetscCall(DMDAGetElements_1D(dm, nel, nen, e));
414 } else if (dim == 2) {
415 PetscCall(DMDAGetElements_2D(dm, nel, nen, e));
416 } else if (dim == 3) {
417 PetscCall(DMDAGetElements_3D(dm, nel, nen, e));
418 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
419 PetscFunctionReturn(PETSC_SUCCESS);
420 }
421
422 /*@
423 DMDAGetSubdomainCornersIS - Gets an index set containing the corner indices (in local indexing)
424 of the non-overlapping decomposition identified by `DMDAGetElements()`
425
426 Not Collective
427
428 Input Parameter:
429 . dm - the `DMDA` object
430
431 Output Parameter:
432 . is - the index set
433
434 Level: intermediate
435
436 Note:
437 Call `DMDARestoreSubdomainCornersIS()` once you have finished accessing the index set.
438
439 .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAElementType`, `DMDASetElementType()`, `DMDAGetElements()`, `DMDARestoreElementsCornersIS()`,
440 `DMDAGetElementsSizes()`, `DMDAGetElementsCorners()`
441 @*/
DMDAGetSubdomainCornersIS(DM dm,IS * is)442 PetscErrorCode DMDAGetSubdomainCornersIS(DM dm, IS *is)
443 {
444 DM_DA *dd = (DM_DA *)dm->data;
445 PetscBool isda;
446
447 PetscFunctionBegin;
448 PetscValidHeaderSpecificType(dm, DM_CLASSID, 1, DMDA);
449 PetscAssertPointer(is, 2);
450 PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMDA, &isda));
451 PetscCheck(isda, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for DM type %s", ((PetscObject)dm)->type_name);
452 PetscCheck(dd->stencil_type != DMDA_STENCIL_STAR, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "DMDAGetElement() requires you use a stencil type of DMDA_STENCIL_BOX");
453 if (!dd->ecorners) { /* compute elements if not yet done */
454 const PetscInt *e;
455 PetscInt nel, nen;
456
457 PetscCall(DMDAGetElements(dm, &nel, &nen, &e));
458 PetscCall(DMDARestoreElements(dm, &nel, &nen, &e));
459 }
460 *is = dd->ecorners;
461 PetscFunctionReturn(PETSC_SUCCESS);
462 }
463
464 /*@C
465 DMDARestoreElements - Restores the array obtained with `DMDAGetElements()`
466
467 Not Collective
468
469 Input Parameters:
470 + dm - the `DM` object
471 . nel - number of local elements
472 . nen - number of nodes in each element
473 - e - the local indices of the elements' vertices
474
475 Level: intermediate
476
477 Note:
478 This restore signals the `DMDA` object that you no longer need access to the array information.
479
480 .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAElementType`, `DMDASetElementType()`, `DMDAGetElements()`
481 @*/
DMDARestoreElements(DM dm,PetscInt * nel,PetscInt * nen,const PetscInt * e[])482 PetscErrorCode DMDARestoreElements(DM dm, PetscInt *nel, PetscInt *nen, const PetscInt *e[])
483 {
484 PetscFunctionBegin;
485 PetscValidHeaderSpecificType(dm, DM_CLASSID, 1, DMDA);
486 PetscAssertPointer(nel, 2);
487 PetscAssertPointer(nen, 3);
488 PetscAssertPointer(e, 4);
489 if (nel) *nel = 0;
490 if (nen) *nen = -1;
491 if (e) *e = NULL;
492 PetscFunctionReturn(PETSC_SUCCESS);
493 }
494
495 /*@
496 DMDARestoreSubdomainCornersIS - Restores the `IS` obtained with `DMDAGetSubdomainCornersIS()`
497
498 Not Collective
499
500 Input Parameters:
501 + dm - the `DM` object
502 - is - the index set
503
504 Level: intermediate
505
506 .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAElementType`, `DMDASetElementType()`, `DMDAGetSubdomainCornersIS()`
507 @*/
DMDARestoreSubdomainCornersIS(DM dm,IS * is)508 PetscErrorCode DMDARestoreSubdomainCornersIS(DM dm, IS *is)
509 {
510 PetscFunctionBegin;
511 PetscValidHeaderSpecificType(dm, DM_CLASSID, 1, DMDA);
512 PetscValidHeaderSpecific(*is, IS_CLASSID, 2);
513 *is = NULL;
514 PetscFunctionReturn(PETSC_SUCCESS);
515 }
516