Actual source code: dacreate.c
1: #include <petsc/private/dmdaimpl.h>
3: static PetscErrorCode DMSetFromOptions_DA(DM da, PetscOptionItems *PetscOptionsObject)
4: {
5: DM_DA *dd = (DM_DA *)da->data;
6: PetscInt refine = 0, dim = da->dim, maxnlevels = 100, refx[100], refy[100], refz[100], n, i;
7: DMBoundaryType bt = DM_BOUNDARY_NONE;
8: PetscBool flg;
10: PetscFunctionBegin;
11: PetscCheck(!da->setupcalled, PetscObjectComm((PetscObject)da), PETSC_ERR_ARG_WRONGSTATE, "Cannot call DMSetFromOptions() after DMSetUp() for DMDA");
12: PetscCheck(dd->M >= 0, PetscObjectComm((PetscObject)da), PETSC_ERR_ARG_OUTOFRANGE, "Dimension must be non-negative, call DMSetFromOptions() if you want to change the value at runtime");
13: PetscCheck(dd->N >= 0, PetscObjectComm((PetscObject)da), PETSC_ERR_ARG_OUTOFRANGE, "Dimension must be non-negative, call DMSetFromOptions() if you want to change the value at runtime");
14: PetscCheck(dd->P >= 0, PetscObjectComm((PetscObject)da), PETSC_ERR_ARG_OUTOFRANGE, "Dimension must be non-negative, call DMSetFromOptions() if you want to change the value at runtime");
16: PetscOptionsHeadBegin(PetscOptionsObject, "DMDA Options");
17: PetscCall(PetscOptionsBoundedInt("-da_grid_x", "Number of grid points in x direction", "DMDASetSizes", dd->M, &dd->M, NULL, 1));
18: if (dim > 1) PetscCall(PetscOptionsBoundedInt("-da_grid_y", "Number of grid points in y direction", "DMDASetSizes", dd->N, &dd->N, NULL, 1));
19: if (dim > 2) PetscCall(PetscOptionsBoundedInt("-da_grid_z", "Number of grid points in z direction", "DMDASetSizes", dd->P, &dd->P, NULL, 1));
21: PetscCall(PetscOptionsBoundedInt("-da_overlap", "Decomposition overlap in all directions", "DMDASetOverlap", dd->xol, &dd->xol, &flg, 0));
22: if (flg) PetscCall(DMDASetOverlap(da, dd->xol, dd->xol, dd->xol));
23: PetscCall(PetscOptionsBoundedInt("-da_overlap_x", "Decomposition overlap in x direction", "DMDASetOverlap", dd->xol, &dd->xol, NULL, 0));
24: if (dim > 1) PetscCall(PetscOptionsBoundedInt("-da_overlap_y", "Decomposition overlap in y direction", "DMDASetOverlap", dd->yol, &dd->yol, NULL, 0));
25: if (dim > 2) PetscCall(PetscOptionsBoundedInt("-da_overlap_z", "Decomposition overlap in z direction", "DMDASetOverlap", dd->zol, &dd->zol, NULL, 0));
27: PetscCall(PetscOptionsBoundedInt("-da_local_subdomains", "", "DMDASetNumLocalSubdomains", dd->Nsub, &dd->Nsub, &flg, PETSC_DECIDE));
28: if (flg) PetscCall(DMDASetNumLocalSubDomains(da, dd->Nsub));
30: /* Handle DMDA parallel distribution */
31: PetscCall(PetscOptionsBoundedInt("-da_processors_x", "Number of processors in x direction", "DMDASetNumProcs", dd->m, &dd->m, NULL, PETSC_DECIDE));
32: if (dim > 1) PetscCall(PetscOptionsBoundedInt("-da_processors_y", "Number of processors in y direction", "DMDASetNumProcs", dd->n, &dd->n, NULL, PETSC_DECIDE));
33: if (dim > 2) PetscCall(PetscOptionsBoundedInt("-da_processors_z", "Number of processors in z direction", "DMDASetNumProcs", dd->p, &dd->p, NULL, PETSC_DECIDE));
34: // Handle boundaries
35: PetscCall(PetscOptionsEnum("-da_bd_x", "Boundary type for x direction", "DMDASetBoundaryType", DMBoundaryTypes, (PetscEnum)dd->bx, (PetscEnum *)&dd->bx, NULL));
36: if (dim > 1) PetscCall(PetscOptionsEnum("-da_bd_y", "Boundary type for y direction", "DMDASetBoundaryType", DMBoundaryTypes, (PetscEnum)dd->by, (PetscEnum *)&dd->by, NULL));
37: if (dim > 2) PetscCall(PetscOptionsEnum("-da_bd_z", "Boundary type for z direction", "DMDASetBoundaryType", DMBoundaryTypes, (PetscEnum)dd->bz, (PetscEnum *)&dd->bz, NULL));
38: PetscCall(PetscOptionsEnum("-da_bd_all", "Boundary type for every direction", "DMDASetBoundaryType", DMBoundaryTypes, (PetscEnum)bt, (PetscEnum *)&bt, &flg));
39: if (flg) PetscCall(DMDASetBoundaryType(da, bt, bt, bt));
40: /* Handle DMDA refinement */
41: PetscCall(PetscOptionsBoundedInt("-da_refine_x", "Refinement ratio in x direction", "DMDASetRefinementFactor", dd->refine_x, &dd->refine_x, NULL, 1));
42: if (dim > 1) PetscCall(PetscOptionsBoundedInt("-da_refine_y", "Refinement ratio in y direction", "DMDASetRefinementFactor", dd->refine_y, &dd->refine_y, NULL, 1));
43: if (dim > 2) PetscCall(PetscOptionsBoundedInt("-da_refine_z", "Refinement ratio in z direction", "DMDASetRefinementFactor", dd->refine_z, &dd->refine_z, NULL, 1));
44: dd->coarsen_x = dd->refine_x;
45: dd->coarsen_y = dd->refine_y;
46: dd->coarsen_z = dd->refine_z;
48: /* Get refinement factors, defaults taken from the coarse DMDA */
49: PetscCall(DMDAGetRefinementFactor(da, &refx[0], &refy[0], &refz[0]));
50: for (i = 1; i < maxnlevels; i++) {
51: refx[i] = refx[0];
52: refy[i] = refy[0];
53: refz[i] = refz[0];
54: }
55: n = maxnlevels;
56: PetscCall(PetscOptionsIntArray("-da_refine_hierarchy_x", "Refinement factor for each level", "None", refx, &n, &flg));
57: if (flg) {
58: dd->refine_x = refx[0];
59: dd->refine_x_hier_n = n;
60: PetscCall(PetscMalloc1(n, &dd->refine_x_hier));
61: PetscCall(PetscArraycpy(dd->refine_x_hier, refx, n));
62: }
63: if (dim > 1) {
64: n = maxnlevels;
65: PetscCall(PetscOptionsIntArray("-da_refine_hierarchy_y", "Refinement factor for each level", "None", refy, &n, &flg));
66: if (flg) {
67: dd->refine_y = refy[0];
68: dd->refine_y_hier_n = n;
69: PetscCall(PetscMalloc1(n, &dd->refine_y_hier));
70: PetscCall(PetscArraycpy(dd->refine_y_hier, refy, n));
71: }
72: }
73: if (dim > 2) {
74: n = maxnlevels;
75: PetscCall(PetscOptionsIntArray("-da_refine_hierarchy_z", "Refinement factor for each level", "None", refz, &n, &flg));
76: if (flg) {
77: dd->refine_z = refz[0];
78: dd->refine_z_hier_n = n;
79: PetscCall(PetscMalloc1(n, &dd->refine_z_hier));
80: PetscCall(PetscArraycpy(dd->refine_z_hier, refz, n));
81: }
82: }
84: PetscCall(PetscOptionsBoundedInt("-da_refine", "Uniformly refine DA one or more times", "None", refine, &refine, NULL, 0));
85: PetscOptionsHeadEnd();
87: while (refine--) {
88: if (dd->bx == DM_BOUNDARY_PERIODIC || dd->interptype == DMDA_Q0) {
89: PetscCall(PetscIntMultError(dd->refine_x, dd->M, &dd->M));
90: } else {
91: PetscCall(PetscIntMultError(dd->refine_x, dd->M - 1, &dd->M));
92: dd->M += 1;
93: }
94: if (dim > 1 && (dd->by == DM_BOUNDARY_PERIODIC || dd->interptype == DMDA_Q0)) {
95: PetscCall(PetscIntMultError(dd->refine_y, dd->N, &dd->N));
96: } else if (dim > 1) {
97: PetscCall(PetscIntMultError(dd->refine_y, dd->N - 1, &dd->N));
98: dd->N += 1;
99: }
100: if (dim > 2 && (dd->bz == DM_BOUNDARY_PERIODIC || dd->interptype == DMDA_Q0)) {
101: PetscCall(PetscIntMultError(dd->refine_z, dd->P, &dd->P));
102: } else if (dim > 2) {
103: PetscCall(PetscIntMultError(dd->refine_z, dd->P - 1, &dd->P));
104: dd->P += 1;
105: }
106: da->levelup++;
107: if (da->levelup - da->leveldown >= 0) {
108: dd->refine_x = refx[da->levelup - da->leveldown];
109: dd->refine_y = refy[da->levelup - da->leveldown];
110: dd->refine_z = refz[da->levelup - da->leveldown];
111: }
112: if (da->levelup - da->leveldown >= 1) {
113: dd->coarsen_x = refx[da->levelup - da->leveldown - 1];
114: dd->coarsen_y = refy[da->levelup - da->leveldown - 1];
115: dd->coarsen_z = refz[da->levelup - da->leveldown - 1];
116: }
117: }
118: PetscFunctionReturn(PETSC_SUCCESS);
119: }
121: static PetscErrorCode DMLoad_DA(DM da, PetscViewer viewer)
122: {
123: PetscInt dim, m, n, p, dof, swidth;
124: DMDAStencilType stencil;
125: DMBoundaryType bx, by, bz;
126: PetscBool coors;
127: DM dac;
128: Vec c;
130: PetscFunctionBegin;
131: PetscCall(PetscViewerBinaryRead(viewer, &dim, 1, NULL, PETSC_INT));
132: PetscCall(PetscViewerBinaryRead(viewer, &m, 1, NULL, PETSC_INT));
133: PetscCall(PetscViewerBinaryRead(viewer, &n, 1, NULL, PETSC_INT));
134: PetscCall(PetscViewerBinaryRead(viewer, &p, 1, NULL, PETSC_INT));
135: PetscCall(PetscViewerBinaryRead(viewer, &dof, 1, NULL, PETSC_INT));
136: PetscCall(PetscViewerBinaryRead(viewer, &swidth, 1, NULL, PETSC_INT));
137: PetscCall(PetscViewerBinaryRead(viewer, &bx, 1, NULL, PETSC_ENUM));
138: PetscCall(PetscViewerBinaryRead(viewer, &by, 1, NULL, PETSC_ENUM));
139: PetscCall(PetscViewerBinaryRead(viewer, &bz, 1, NULL, PETSC_ENUM));
140: PetscCall(PetscViewerBinaryRead(viewer, &stencil, 1, NULL, PETSC_ENUM));
142: PetscCall(DMSetDimension(da, dim));
143: PetscCall(DMDASetSizes(da, m, n, p));
144: PetscCall(DMDASetBoundaryType(da, bx, by, bz));
145: PetscCall(DMDASetDof(da, dof));
146: PetscCall(DMDASetStencilType(da, stencil));
147: PetscCall(DMDASetStencilWidth(da, swidth));
148: PetscCall(DMSetUp(da));
149: PetscCall(PetscViewerBinaryRead(viewer, &coors, 1, NULL, PETSC_ENUM));
150: if (coors) {
151: PetscCall(DMGetCoordinateDM(da, &dac));
152: PetscCall(DMCreateGlobalVector(dac, &c));
153: PetscCall(VecLoad(c, viewer));
154: PetscCall(DMSetCoordinates(da, c));
155: PetscCall(VecDestroy(&c));
156: }
157: PetscFunctionReturn(PETSC_SUCCESS);
158: }
160: static PetscErrorCode DMCreateSubDM_DA(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm)
161: {
162: DM_DA *da = (DM_DA *)dm->data;
164: PetscFunctionBegin;
165: if (subdm) {
166: PetscSF sf;
167: Vec coords;
168: void *ctx;
169: /* Cannot use DMClone since the dof stuff is mixed in. Ugh
170: PetscCall(DMClone(dm, subdm)); */
171: PetscCall(DMCreate(PetscObjectComm((PetscObject)dm), subdm));
172: PetscCall(DMGetPointSF(dm, &sf));
173: PetscCall(DMSetPointSF(*subdm, sf));
174: PetscCall(DMGetApplicationContext(dm, &ctx));
175: PetscCall(DMSetApplicationContext(*subdm, ctx));
176: PetscCall(DMGetCoordinatesLocal(dm, &coords));
177: if (coords) {
178: PetscCall(DMSetCoordinatesLocal(*subdm, coords));
179: } else {
180: PetscCall(DMGetCoordinates(dm, &coords));
181: if (coords) PetscCall(DMSetCoordinates(*subdm, coords));
182: }
184: PetscCall(DMSetType(*subdm, DMDA));
185: PetscCall(DMSetDimension(*subdm, dm->dim));
186: PetscCall(DMDASetSizes(*subdm, da->M, da->N, da->P));
187: PetscCall(DMDASetNumProcs(*subdm, da->m, da->n, da->p));
188: PetscCall(DMDASetBoundaryType(*subdm, da->bx, da->by, da->bz));
189: PetscCall(DMDASetDof(*subdm, numFields));
190: PetscCall(DMDASetStencilType(*subdm, da->stencil_type));
191: PetscCall(DMDASetStencilWidth(*subdm, da->s));
192: PetscCall(DMDASetOwnershipRanges(*subdm, da->lx, da->ly, da->lz));
193: }
194: if (is) {
195: PetscInt *indices, cnt = 0, dof = da->w, i, j;
197: PetscCall(PetscMalloc1(da->Nlocal * numFields / dof, &indices));
198: for (i = da->base / dof; i < (da->base + da->Nlocal) / dof; ++i) {
199: for (j = 0; j < numFields; ++j) indices[cnt++] = dof * i + fields[j];
200: }
201: PetscCheck(cnt == da->Nlocal * numFields / dof, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Count %" PetscInt_FMT " does not equal expected value %" PetscInt_FMT, cnt, da->Nlocal * numFields / dof);
202: PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)dm), cnt, indices, PETSC_OWN_POINTER, is));
203: }
204: PetscFunctionReturn(PETSC_SUCCESS);
205: }
207: static PetscErrorCode DMCreateFieldDecomposition_DA(DM dm, PetscInt *len, char ***namelist, IS **islist, DM **dmlist)
208: {
209: PetscInt i;
210: DM_DA *dd = (DM_DA *)dm->data;
211: PetscInt dof = dd->w;
213: PetscFunctionBegin;
214: if (len) *len = dof;
215: if (islist) {
216: Vec v;
217: PetscInt rstart, n;
219: PetscCall(DMGetGlobalVector(dm, &v));
220: PetscCall(VecGetOwnershipRange(v, &rstart, NULL));
221: PetscCall(VecGetLocalSize(v, &n));
222: PetscCall(DMRestoreGlobalVector(dm, &v));
223: PetscCall(PetscMalloc1(dof, islist));
224: for (i = 0; i < dof; i++) PetscCall(ISCreateStride(PetscObjectComm((PetscObject)dm), n / dof, rstart + i, dof, &(*islist)[i]));
225: }
226: if (namelist) {
227: PetscCall(PetscMalloc1(dof, namelist));
228: if (dd->fieldname) {
229: for (i = 0; i < dof; i++) PetscCall(PetscStrallocpy(dd->fieldname[i], &(*namelist)[i]));
230: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Currently DMDA must have fieldnames");
231: }
232: if (dmlist) {
233: DM da;
235: PetscCall(DMDACreate(PetscObjectComm((PetscObject)dm), &da));
236: PetscCall(DMSetDimension(da, dm->dim));
237: PetscCall(DMDASetSizes(da, dd->M, dd->N, dd->P));
238: PetscCall(DMDASetNumProcs(da, dd->m, dd->n, dd->p));
239: PetscCall(DMDASetBoundaryType(da, dd->bx, dd->by, dd->bz));
240: PetscCall(DMDASetDof(da, 1));
241: PetscCall(DMDASetStencilType(da, dd->stencil_type));
242: PetscCall(DMDASetStencilWidth(da, dd->s));
243: PetscCall(DMSetUp(da));
244: PetscCall(PetscMalloc1(dof, dmlist));
245: for (i = 0; i < dof - 1; i++) PetscCall(PetscObjectReference((PetscObject)da));
246: for (i = 0; i < dof; i++) (*dmlist)[i] = da;
247: }
248: PetscFunctionReturn(PETSC_SUCCESS);
249: }
251: static PetscErrorCode DMClone_DA(DM dm, DM *newdm)
252: {
253: DM_DA *da = (DM_DA *)dm->data;
255: PetscFunctionBegin;
256: PetscCall(DMSetType(*newdm, DMDA));
257: PetscCall(DMSetDimension(*newdm, dm->dim));
258: PetscCall(DMDASetSizes(*newdm, da->M, da->N, da->P));
259: PetscCall(DMDASetNumProcs(*newdm, da->m, da->n, da->p));
260: PetscCall(DMDASetBoundaryType(*newdm, da->bx, da->by, da->bz));
261: PetscCall(DMDASetDof(*newdm, da->w));
262: PetscCall(DMDASetStencilType(*newdm, da->stencil_type));
263: PetscCall(DMDASetStencilWidth(*newdm, da->s));
264: PetscCall(DMDASetOwnershipRanges(*newdm, da->lx, da->ly, da->lz));
265: PetscCall(DMSetUp(*newdm));
266: PetscFunctionReturn(PETSC_SUCCESS);
267: }
269: static PetscErrorCode DMHasCreateInjection_DA(DM dm, PetscBool *flg)
270: {
271: DM_DA *da = (DM_DA *)dm->data;
273: PetscFunctionBegin;
275: PetscAssertPointer(flg, 2);
276: *flg = da->interptype == DMDA_Q1 ? PETSC_TRUE : PETSC_FALSE;
277: PetscFunctionReturn(PETSC_SUCCESS);
278: }
280: static PetscErrorCode DMGetDimPoints_DA(DM dm, PetscInt dim, PetscInt *pStart, PetscInt *pEnd)
281: {
282: PetscFunctionBegin;
283: PetscCall(DMDAGetDepthStratum(dm, dim, pStart, pEnd));
284: PetscFunctionReturn(PETSC_SUCCESS);
285: }
287: static PetscErrorCode DMGetNeighbors_DA(DM dm, PetscInt *nranks, const PetscMPIInt *ranks[])
288: {
289: PetscInt dim;
290: DMDAStencilType st;
292: PetscFunctionBegin;
293: PetscCall(DMDAGetNeighbors(dm, ranks));
294: PetscCall(DMDAGetInfo(dm, &dim, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, &st));
296: switch (dim) {
297: case 1:
298: *nranks = 3;
299: /* if (st == DMDA_STENCIL_STAR) *nranks = 3; */
300: break;
301: case 2:
302: *nranks = 9;
303: /* if (st == DMDA_STENCIL_STAR) *nranks = 5; */
304: break;
305: case 3:
306: *nranks = 27;
307: /* if (st == DMDA_STENCIL_STAR) *nranks = 7; */
308: break;
309: default:
310: break;
311: }
312: PetscFunctionReturn(PETSC_SUCCESS);
313: }
315: /*MC
316: DMDA = "da" - A `DM` object that is used to manage data for a structured grid in 1, 2, or 3 dimensions.
317: In the global representation of the vector each process stores a non-overlapping rectangular (or slab in 3d) portion of the grid points.
318: In the local representation these rectangular regions (slabs) are extended in all directions by a stencil width set with `DMDASetStencilWidth()`.
320: The vectors can be thought of as either cell centered or vertex centered on the mesh. But some variables cannot be cell centered and others
321: vertex centered; see the documentation for `DMSTAG`, a similar `DM` implementation which supports more general staggered grids.
323: Level: intermediate
325: .seealso: [](sec_struct), `DMType`, `DMCOMPOSITE`, `DMSTAG`, `DMDACreate()`, `DMCreate()`, `DMSetType()`, `DMDASetStencilWidth()`, `DMDASetStencilType()`,
326: `DMDAStencilType`
327: M*/
329: PETSC_EXTERN PetscErrorCode DMCreate_DA(DM da)
330: {
331: DM_DA *dd;
333: PetscFunctionBegin;
334: PetscAssertPointer(da, 1);
335: PetscCall(PetscNew(&dd));
336: da->data = dd;
338: da->dim = -1;
339: dd->interptype = DMDA_Q1;
340: dd->refine_x = 2;
341: dd->refine_y = 2;
342: dd->refine_z = 2;
343: dd->coarsen_x = 2;
344: dd->coarsen_y = 2;
345: dd->coarsen_z = 2;
346: dd->fieldname = NULL;
347: dd->nlocal = -1;
348: dd->Nlocal = -1;
349: dd->M = -1;
350: dd->N = -1;
351: dd->P = -1;
352: dd->m = -1;
353: dd->n = -1;
354: dd->p = -1;
355: dd->w = -1;
356: dd->s = -1;
358: dd->xs = -1;
359: dd->xe = -1;
360: dd->ys = -1;
361: dd->ye = -1;
362: dd->zs = -1;
363: dd->ze = -1;
364: dd->Xs = -1;
365: dd->Xe = -1;
366: dd->Ys = -1;
367: dd->Ye = -1;
368: dd->Zs = -1;
369: dd->Ze = -1;
371: dd->Nsub = 1;
372: dd->xol = 0;
373: dd->yol = 0;
374: dd->zol = 0;
375: dd->xo = 0;
376: dd->yo = 0;
377: dd->zo = 0;
378: dd->Mo = -1;
379: dd->No = -1;
380: dd->Po = -1;
382: dd->gtol = NULL;
383: dd->ltol = NULL;
384: dd->ao = NULL;
385: PetscCall(PetscStrallocpy(AOBASIC, (char **)&dd->aotype));
386: dd->base = -1;
387: dd->bx = DM_BOUNDARY_NONE;
388: dd->by = DM_BOUNDARY_NONE;
389: dd->bz = DM_BOUNDARY_NONE;
390: dd->stencil_type = DMDA_STENCIL_BOX;
391: dd->interptype = DMDA_Q1;
392: dd->lx = NULL;
393: dd->ly = NULL;
394: dd->lz = NULL;
396: dd->elementtype = DMDA_ELEMENT_Q1;
398: da->ops->globaltolocalbegin = DMGlobalToLocalBegin_DA;
399: da->ops->globaltolocalend = DMGlobalToLocalEnd_DA;
400: da->ops->localtoglobalbegin = DMLocalToGlobalBegin_DA;
401: da->ops->localtoglobalend = DMLocalToGlobalEnd_DA;
402: da->ops->localtolocalbegin = DMLocalToLocalBegin_DA;
403: da->ops->localtolocalend = DMLocalToLocalEnd_DA;
404: da->ops->createglobalvector = DMCreateGlobalVector_DA;
405: da->ops->createlocalvector = DMCreateLocalVector_DA;
406: da->ops->createinterpolation = DMCreateInterpolation_DA;
407: da->ops->getcoloring = DMCreateColoring_DA;
408: da->ops->creatematrix = DMCreateMatrix_DA;
409: da->ops->refine = DMRefine_DA;
410: da->ops->coarsen = DMCoarsen_DA;
411: da->ops->refinehierarchy = DMRefineHierarchy_DA;
412: da->ops->coarsenhierarchy = DMCoarsenHierarchy_DA;
413: da->ops->createinjection = DMCreateInjection_DA;
414: da->ops->hascreateinjection = DMHasCreateInjection_DA;
415: da->ops->destroy = DMDestroy_DA;
416: da->ops->view = NULL;
417: da->ops->setfromoptions = DMSetFromOptions_DA;
418: da->ops->setup = DMSetUp_DA;
419: da->ops->clone = DMClone_DA;
420: da->ops->load = DMLoad_DA;
421: da->ops->createcoordinatedm = DMCreateCoordinateDM_DA;
422: da->ops->createsubdm = DMCreateSubDM_DA;
423: da->ops->createfielddecomposition = DMCreateFieldDecomposition_DA;
424: da->ops->createdomaindecomposition = DMCreateDomainDecomposition_DA;
425: da->ops->createddscatters = DMCreateDomainDecompositionScatters_DA;
426: da->ops->getdimpoints = DMGetDimPoints_DA;
427: da->ops->getneighbors = DMGetNeighbors_DA;
428: da->ops->getlocalboundingbox = DMGetLocalBoundingBox_DA;
429: da->ops->locatepoints = DMLocatePoints_DA_Regular;
430: da->ops->getcompatibility = DMGetCompatibility_DA;
431: PetscCall(PetscObjectComposeFunction((PetscObject)da, "DMSetUpGLVisViewer_C", DMSetUpGLVisViewer_DMDA));
432: PetscFunctionReturn(PETSC_SUCCESS);
433: }
435: /*@
436: DMDACreate - Creates a `DMDA` object for managing structured grids.
438: Collective
440: Input Parameter:
441: . comm - The communicator for the `DMDA` object
443: Output Parameter:
444: . da - the `DMDA` object
446: Level: advanced
448: Notes:
449: See [](sec_struct) for details on the construction of a `DMDA`
451: `DMDACreate1d()`, `DMDACreate2d()`, and `DMDACreate3d()` are convenience routines to quickly completely create a `DMDA`
453: .seealso: [](sec_struct), `DM`, `DMDA`, `DMSetUp()`, `DMDASetSizes()`, `DMClone()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
454: @*/
455: PetscErrorCode DMDACreate(MPI_Comm comm, DM *da)
456: {
457: PetscFunctionBegin;
458: PetscAssertPointer(da, 2);
459: PetscCall(DMCreate(comm, da));
460: PetscCall(DMSetType(*da, DMDA));
461: PetscFunctionReturn(PETSC_SUCCESS);
462: }