Actual source code: dacorn.c
1: /*
2: Code for manipulating distributed regular arrays in parallel.
3: */
5: #include <petsc/private/dmdaimpl.h>
6: #include <petscdmfield.h>
8: PetscErrorCode DMCreateCoordinateDM_DA(DM dm, DM *cdm)
9: {
10: const char *prefix;
12: PetscFunctionBegin;
13: PetscCall(DMDACreateCompatibleDMDA(dm, dm->dim, cdm));
14: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)dm, &prefix));
15: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)*cdm, prefix));
16: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)*cdm, "cdm_"));
17: PetscFunctionReturn(PETSC_SUCCESS);
18: }
20: /*@
21: DMDASetFieldName - Sets the names of individual field components in multicomponent
22: vectors associated with a `DMDA`.
24: Logically Collective; name must contain a common value
26: Input Parameters:
27: + da - the `DMDA`
28: . nf - field number for the `DMDA` (0, 1, ... dof-1), where dof indicates the
29: number of degrees of freedom per node within the `DMDA`
30: - name - the name of the field (component)
32: Level: intermediate
34: Note:
35: It must be called after having called `DMSetUp()`.
37: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAGetFieldName()`, `DMDASetCoordinateName()`, `DMDAGetCoordinateName()`, `DMDASetFieldNames()`, `DMSetUp()`
38: @*/
39: PetscErrorCode DMDASetFieldName(DM da, PetscInt nf, const char name[])
40: {
41: DM_DA *dd = (DM_DA *)da->data;
43: PetscFunctionBegin;
45: PetscCheck(nf >= 0 && nf < dd->w, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Invalid field number: %" PetscInt_FMT, nf);
46: PetscCheck(dd->fieldname, PetscObjectComm((PetscObject)da), PETSC_ERR_ORDER, "You should call DMSetUp() first");
47: PetscCall(PetscFree(dd->fieldname[nf]));
48: PetscCall(PetscStrallocpy(name, &dd->fieldname[nf]));
49: PetscFunctionReturn(PETSC_SUCCESS);
50: }
52: /*@C
53: DMDAGetFieldNames - Gets the name of all the components in the vector associated with the `DMDA`
55: Not Collective; names will contain a common value; No Fortran Support
57: Input Parameter:
58: . da - the `DMDA` object
60: Output Parameter:
61: . names - the names of the components, final string is `NULL`, will have the same number of entries as the dof used in creating the `DMDA`
63: Level: intermediate
65: Fortran Note:
66: Use `DMDAGetFieldName()`
68: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAGetFieldName()`, `DMDASetCoordinateName()`, `DMDAGetCoordinateName()`, `DMDASetFieldName()`, `DMDASetFieldNames()`
69: @*/
70: PetscErrorCode DMDAGetFieldNames(DM da, const char *const **names)
71: {
72: DM_DA *dd = (DM_DA *)da->data;
74: PetscFunctionBegin;
75: *names = (const char *const *)dd->fieldname;
76: PetscFunctionReturn(PETSC_SUCCESS);
77: }
79: /*@C
80: DMDASetFieldNames - Sets the name of each component in the vector associated with the `DMDA`
82: Logically Collective; names must contain a common value; No Fortran Support
84: Input Parameters:
85: + da - the `DMDA` object
86: - names - the names of the components, final string must be `NULL`, must have the same number of entries as the dof used in creating the `DMDA`
88: Level: intermediate
90: Note:
91: It must be called after having called `DMSetUp()`.
93: Fortran Note:
94: Use `DMDASetFieldName()`
96: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAGetFieldName()`, `DMDASetCoordinateName()`, `DMDAGetCoordinateName()`, `DMDASetFieldName()`, `DMSetUp()`
97: @*/
98: PetscErrorCode DMDASetFieldNames(DM da, const char *const *names)
99: {
100: DM_DA *dd = (DM_DA *)da->data;
101: char **fieldname;
102: PetscInt nf = 0;
104: PetscFunctionBegin;
105: PetscCheck(dd->fieldname, PetscObjectComm((PetscObject)da), PETSC_ERR_ORDER, "You should call DMSetUp() first");
106: while (names[nf++]) { }
107: PetscCheck(nf == dd->w + 1, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Invalid number of fields %" PetscInt_FMT, nf - 1);
108: PetscCall(PetscStrArrayallocpy(names, &fieldname));
109: PetscCall(PetscStrArrayDestroy(&dd->fieldname));
110: dd->fieldname = fieldname;
111: PetscFunctionReturn(PETSC_SUCCESS);
112: }
114: /*@
115: DMDAGetFieldName - Gets the names of individual field components in multicomponent
116: vectors associated with a `DMDA`.
118: Not Collective; name will contain a common value
120: Input Parameters:
121: + da - the `DMDA`
122: - nf - field number for the `DMDA` (0, 1, ... dof-1), where dof indicates the
123: number of degrees of freedom per node within the `DMDA`
125: Output Parameter:
126: . name - the name of the field (component)
128: Level: intermediate
130: Note:
131: It must be called after having called `DMSetUp()`.
133: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDASetFieldName()`, `DMDASetCoordinateName()`, `DMDAGetCoordinateName()`, `DMSetUp()`
134: @*/
135: PetscErrorCode DMDAGetFieldName(DM da, PetscInt nf, const char *name[])
136: {
137: DM_DA *dd = (DM_DA *)da->data;
139: PetscFunctionBegin;
141: PetscAssertPointer(name, 3);
142: PetscCheck(nf >= 0 && nf < dd->w, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Invalid field number: %" PetscInt_FMT, nf);
143: PetscCheck(dd->fieldname, PetscObjectComm((PetscObject)da), PETSC_ERR_ORDER, "You should call DMSetUp() first");
144: *name = dd->fieldname[nf];
145: PetscFunctionReturn(PETSC_SUCCESS);
146: }
148: /*@
149: DMDASetCoordinateName - Sets the name of the coordinate directions associated with a `DMDA`, for example "x" or "y"
151: Logically Collective; name must contain a common value; No Fortran Support
153: Input Parameters:
154: + dm - the `DMDA`
155: . nf - coordinate number for the `DMDA` (0, 1, ... dim-1),
156: - name - the name of the coordinate
158: Level: intermediate
160: Note:
161: Must be called after having called `DMSetUp()`.
163: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAGetCoordinateName()`, `DMDASetFieldName()`, `DMDAGetFieldName()`, `DMSetUp()`
164: @*/
165: PetscErrorCode DMDASetCoordinateName(DM dm, PetscInt nf, const char name[])
166: {
167: DM_DA *dd = (DM_DA *)dm->data;
169: PetscFunctionBegin;
171: PetscCheck(nf >= 0 && nf < dm->dim, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Invalid coordinate number: %" PetscInt_FMT, nf);
172: PetscCheck(dd->coordinatename, PetscObjectComm((PetscObject)dm), PETSC_ERR_ORDER, "You should call DMSetUp() first");
173: PetscCall(PetscFree(dd->coordinatename[nf]));
174: PetscCall(PetscStrallocpy(name, &dd->coordinatename[nf]));
175: PetscFunctionReturn(PETSC_SUCCESS);
176: }
178: /*@
179: DMDAGetCoordinateName - Gets the name of a coordinate direction associated with a `DMDA`.
181: Not Collective; name will contain a common value; No Fortran Support
183: Input Parameters:
184: + dm - the `DMDA`
185: - nf - number for the `DMDA` (0, 1, ... dim-1)
187: Output Parameter:
188: . name - the name of the coordinate direction
190: Level: intermediate
192: Note:
193: It must be called after having called `DMSetUp()`.
195: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDASetCoordinateName()`, `DMDASetFieldName()`, `DMDAGetFieldName()`, `DMSetUp()`
196: @*/
197: PetscErrorCode DMDAGetCoordinateName(DM dm, PetscInt nf, const char *name[])
198: {
199: DM_DA *dd = (DM_DA *)dm->data;
201: PetscFunctionBegin;
203: PetscAssertPointer(name, 3);
204: PetscCheck(nf >= 0 && nf < dm->dim, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Invalid coordinate number: %" PetscInt_FMT, nf);
205: PetscCheck(dd->coordinatename, PetscObjectComm((PetscObject)dm), PETSC_ERR_ORDER, "You should call DMSetUp() first");
206: *name = dd->coordinatename[nf];
207: PetscFunctionReturn(PETSC_SUCCESS);
208: }
210: /*@
211: DMDAGetCorners - Returns the global (`x`,`y`,`z`) indices of the lower left
212: corner and size of the local region, excluding ghost points.
214: Not Collective
216: Input Parameter:
217: . da - the `DMDA`
219: Output Parameters:
220: + x - the corner index for the first dimension
221: . y - the corner index for the second dimension (only used in 2D and 3D problems)
222: . z - the corner index for the third dimension (only used in 3D problems)
223: . m - the width in the first dimension
224: . n - the width in the second dimension (only used in 2D and 3D problems)
225: - p - the width in the third dimension (only used in 3D problems)
227: Level: beginner
229: Notes:
230: Any of `y`, `z`, `n`, and `p` can be passed in as `NULL` if not needed.
232: The corner information is independent of the number of degrees of
233: freedom per node set with the `DMDACreateXX()` routine. Thus the `x`, `y`, and `z`
234: can be thought of as the lower left coordinates of the patch of values on process on a logical grid and `m`, `n`, and `p` as the
235: extent of the patch, where each grid point has (potentially) several degrees of freedom.
237: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetOwnershipRanges()`, `DMStagGetCorners()`, `DMSTAG`
238: @*/
239: PetscErrorCode DMDAGetCorners(DM da, PetscInt *x, PetscInt *y, PetscInt *z, PetscInt *m, PetscInt *n, PetscInt *p)
240: {
241: PetscInt w;
242: DM_DA *dd = (DM_DA *)da->data;
244: PetscFunctionBegin;
246: /* since the xs, xe ... have all been multiplied by the number of degrees
247: of freedom per cell, w = dd->w, we divide that out before returning.*/
248: w = dd->w;
249: if (x) *x = dd->xs / w + dd->xo;
250: /* the y and z have NOT been multiplied by w */
251: if (y) *y = dd->ys + dd->yo;
252: if (z) *z = dd->zs + dd->zo;
253: if (m) *m = (dd->xe - dd->xs) / w;
254: if (n) *n = (dd->ye - dd->ys);
255: if (p) *p = (dd->ze - dd->zs);
256: PetscFunctionReturn(PETSC_SUCCESS);
257: }
259: PetscErrorCode DMGetLocalBoundingIndices_DMDA(DM dm, PetscReal lmin[], PetscReal lmax[])
260: {
261: DMDALocalInfo info;
263: PetscFunctionBegin;
264: PetscCall(DMDAGetLocalInfo(dm, &info));
265: lmin[0] = info.xs;
266: lmin[1] = info.ys;
267: lmin[2] = info.zs;
268: lmax[0] = info.xs + info.xm - 1;
269: lmax[1] = info.ys + info.ym - 1;
270: lmax[2] = info.zs + info.zm - 1;
271: PetscFunctionReturn(PETSC_SUCCESS);
272: }
274: // PetscClangLinter pragma ignore: -fdoc-*
275: /*@
276: DMDAGetReducedDMDA - Deprecated; use DMDACreateCompatibleDMDA()
278: Level: deprecated
279: @*/
280: PetscErrorCode DMDAGetReducedDMDA(DM da, PetscInt nfields, DM *nda)
281: {
282: PetscFunctionBegin;
283: PetscCall(DMDACreateCompatibleDMDA(da, nfields, nda));
284: PetscFunctionReturn(PETSC_SUCCESS);
285: }
287: /*@
288: DMDACreateCompatibleDMDA - Creates a `DMDA` with the same layout as given `DMDA` but with fewer or more fields
290: Collective
292: Input Parameters:
293: + da - the `DMDA`
294: - nfields - number of fields in new `DMDA`
296: Output Parameter:
297: . nda - the new `DMDA`
299: Level: intermediate
301: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMSetCoordinates()`, `DMDASetUniformCoordinates()`, `DMGetCoordinates()`, `DMDAGetGhostedCoordinates()`,
302: `DMStagCreateCompatibleDMStag()`
303: @*/
304: PetscErrorCode DMDACreateCompatibleDMDA(DM da, PetscInt nfields, DM *nda)
305: {
306: DM_DA *dd = (DM_DA *)da->data;
307: PetscInt s, m, n, p, M, N, P, dim, Mo, No, Po;
308: const PetscInt *lx, *ly, *lz;
309: DMBoundaryType bx, by, bz;
310: DMDAStencilType stencil_type;
311: Vec coords;
312: PetscInt ox, oy, oz;
313: PetscInt cl, rl;
315: PetscFunctionBegin;
316: dim = da->dim;
317: M = dd->M;
318: N = dd->N;
319: P = dd->P;
320: m = dd->m;
321: n = dd->n;
322: p = dd->p;
323: s = dd->s;
324: bx = dd->bx;
325: by = dd->by;
326: bz = dd->bz;
328: stencil_type = dd->stencil_type;
330: PetscCall(DMDAGetOwnershipRanges(da, &lx, &ly, &lz));
331: if (dim == 1) {
332: PetscCall(DMDACreate1d(PetscObjectComm((PetscObject)da), bx, M, nfields, s, dd->lx, nda));
333: } else if (dim == 2) {
334: PetscCall(DMDACreate2d(PetscObjectComm((PetscObject)da), bx, by, stencil_type, M, N, m, n, nfields, s, lx, ly, nda));
335: } else if (dim == 3) {
336: PetscCall(DMDACreate3d(PetscObjectComm((PetscObject)da), bx, by, bz, stencil_type, M, N, P, m, n, p, nfields, s, lx, ly, lz, nda));
337: }
338: PetscCall(DMSetUp(*nda));
339: PetscCall(DMGetCoordinates(da, &coords));
340: PetscCall(DMSetCoordinates(*nda, coords));
342: /* allow for getting a reduced DA corresponding to a domain decomposition */
343: PetscCall(DMDAGetOffset(da, &ox, &oy, &oz, &Mo, &No, &Po));
344: PetscCall(DMDASetOffset(*nda, ox, oy, oz, Mo, No, Po));
346: /* allow for getting a reduced DA corresponding to a coarsened DA */
347: PetscCall(DMGetCoarsenLevel(da, &cl));
348: PetscCall(DMGetRefineLevel(da, &rl));
350: (*nda)->levelup = rl;
351: (*nda)->leveldown = cl;
352: PetscFunctionReturn(PETSC_SUCCESS);
353: }
355: /*@C
356: DMDAGetCoordinateArray - Gets an array containing the coordinates of the `DMDA`
358: Not Collective; No Fortran Support
360: Input Parameter:
361: . dm - the `DMDA`
363: Output Parameter:
364: . xc - the coordinates
366: Level: intermediate
368: Note:
369: Use `DMDARestoreCoordinateArray()` to return the array
371: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDASetCoordinateName()`, `DMDASetFieldName()`, `DMDAGetFieldName()`, `DMDARestoreCoordinateArray()`
372: @*/
373: PetscErrorCode DMDAGetCoordinateArray(DM dm, void *xc)
374: {
375: DM cdm;
376: Vec x;
378: PetscFunctionBegin;
380: PetscCall(DMGetCoordinates(dm, &x));
381: PetscCall(DMGetCoordinateDM(dm, &cdm));
382: PetscCall(DMDAVecGetArray(cdm, x, xc));
383: PetscFunctionReturn(PETSC_SUCCESS);
384: }
386: /*@C
387: DMDARestoreCoordinateArray - Returns an array containing the coordinates of the `DMDA` obtained with `DMDAGetCoordinateArray()`
389: Not Collective; No Fortran Support
391: Input Parameters:
392: + dm - the `DMDA`
393: - xc - the coordinates
395: Level: intermediate
397: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDASetCoordinateName()`, `DMDASetFieldName()`, `DMDAGetFieldName()`, `DMDAGetCoordinateArray()`
398: @*/
399: PetscErrorCode DMDARestoreCoordinateArray(DM dm, void *xc)
400: {
401: DM cdm;
402: Vec x;
404: PetscFunctionBegin;
406: PetscCall(DMGetCoordinates(dm, &x));
407: PetscCall(DMGetCoordinateDM(dm, &cdm));
408: PetscCall(DMDAVecRestoreArray(cdm, x, xc));
409: PetscFunctionReturn(PETSC_SUCCESS);
410: }