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: PetscCheck(dd->fieldname, PETSC_COMM_SELF, PETSC_ERR_SUP, "Currently DMDA must have fieldnames");
229: for (i = 0; i < dof; i++) PetscCall(PetscStrallocpy(dd->fieldname[i], &(*namelist)[i]));
230: }
231: if (dmlist) {
232: DM da;
234: PetscCall(DMDACreate(PetscObjectComm((PetscObject)dm), &da));
235: PetscCall(DMSetDimension(da, dm->dim));
236: PetscCall(DMDASetSizes(da, dd->M, dd->N, dd->P));
237: PetscCall(DMDASetNumProcs(da, dd->m, dd->n, dd->p));
238: PetscCall(DMDASetBoundaryType(da, dd->bx, dd->by, dd->bz));
239: PetscCall(DMDASetDof(da, 1));
240: PetscCall(DMDASetStencilType(da, dd->stencil_type));
241: PetscCall(DMDASetStencilWidth(da, dd->s));
242: PetscCall(DMSetUp(da));
243: PetscCall(PetscMalloc1(dof, dmlist));
244: for (i = 0; i < dof - 1; i++) PetscCall(PetscObjectReference((PetscObject)da));
245: for (i = 0; i < dof; i++) (*dmlist)[i] = da;
246: }
247: PetscFunctionReturn(PETSC_SUCCESS);
248: }
250: static PetscErrorCode DMClone_DA(DM dm, DM *newdm)
251: {
252: DM_DA *da = (DM_DA *)dm->data;
254: PetscFunctionBegin;
255: PetscCall(DMSetType(*newdm, DMDA));
256: PetscCall(DMSetDimension(*newdm, dm->dim));
257: PetscCall(DMDASetSizes(*newdm, da->M, da->N, da->P));
258: PetscCall(DMDASetNumProcs(*newdm, da->m, da->n, da->p));
259: PetscCall(DMDASetBoundaryType(*newdm, da->bx, da->by, da->bz));
260: PetscCall(DMDASetDof(*newdm, da->w));
261: PetscCall(DMDASetStencilType(*newdm, da->stencil_type));
262: PetscCall(DMDASetStencilWidth(*newdm, da->s));
263: PetscCall(DMDASetOwnershipRanges(*newdm, da->lx, da->ly, da->lz));
264: PetscCall(DMSetUp(*newdm));
265: PetscFunctionReturn(PETSC_SUCCESS);
266: }
268: static PetscErrorCode DMHasCreateInjection_DA(DM dm, PetscBool *flg)
269: {
270: DM_DA *da = (DM_DA *)dm->data;
272: PetscFunctionBegin;
274: PetscAssertPointer(flg, 2);
275: *flg = da->interptype == DMDA_Q1 ? PETSC_TRUE : PETSC_FALSE;
276: PetscFunctionReturn(PETSC_SUCCESS);
277: }
279: static PetscErrorCode DMGetDimPoints_DA(DM dm, PetscInt dim, PetscInt *pStart, PetscInt *pEnd)
280: {
281: PetscFunctionBegin;
282: PetscCall(DMDAGetDepthStratum(dm, dim, pStart, pEnd));
283: PetscFunctionReturn(PETSC_SUCCESS);
284: }
286: static PetscErrorCode DMGetNeighbors_DA(DM dm, PetscInt *nranks, const PetscMPIInt *ranks[])
287: {
288: PetscInt dim;
289: DMDAStencilType st;
291: PetscFunctionBegin;
292: PetscCall(DMDAGetNeighbors(dm, ranks));
293: PetscCall(DMDAGetInfo(dm, &dim, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, &st));
295: switch (dim) {
296: case 1:
297: *nranks = 3;
298: /* if (st == DMDA_STENCIL_STAR) *nranks = 3; */
299: break;
300: case 2:
301: *nranks = 9;
302: /* if (st == DMDA_STENCIL_STAR) *nranks = 5; */
303: break;
304: case 3:
305: *nranks = 27;
306: /* if (st == DMDA_STENCIL_STAR) *nranks = 7; */
307: break;
308: default:
309: break;
310: }
311: PetscFunctionReturn(PETSC_SUCCESS);
312: }
314: /*MC
315: DMDA = "da" - A `DM` object that is used to help solve PDEs on a structured grid (or mesh) in 1, 2, or 3 dimensions.
317: Level: intermediate
319: Notes:
320: In the global representation of the vectors each process stores a non-overlapping rectangular (or slab in 3d) portion of the grid points.
321: In the local representation these rectangular regions (slabs) are extended in all directions by a stencil width set with `DMDASetStencilWidth()`.
323: The vectors can be thought of as either cell centered or vertex centered on the grid (or mesh). But some variables cannot be cell centered and others
324: vertex centered; see the documentation for `DMSTAG`, a similar `DM` implementation which supports more general staggered grids.
326: Periodic boundary conditions can be handled by using a `DMBoundaryType` of `DM_BOUNDARY_PERIODIC` provided with `DMDASetBoundaryType()`.
327: Other `DMBoundaryType`values allow for different handling of terms along the boundary of the grid (or mesh).
329: .seealso: [](sec_struct), `DMType`, `DMCOMPOSITE`, `DMSTAG`, `DMDACreate()`, `DMCreate()`, `DMSetType()`, `DMDASetStencilWidth()`, `DMDASetStencilType()`,
330: `DMDAStencilType`
331: M*/
333: PETSC_EXTERN PetscErrorCode DMCreate_DA(DM da)
334: {
335: DM_DA *dd;
337: PetscFunctionBegin;
338: PetscAssertPointer(da, 1);
339: PetscCall(PetscNew(&dd));
340: da->data = dd;
342: da->dim = -1;
343: dd->interptype = DMDA_Q1;
344: dd->refine_x = 2;
345: dd->refine_y = 2;
346: dd->refine_z = 2;
347: dd->coarsen_x = 2;
348: dd->coarsen_y = 2;
349: dd->coarsen_z = 2;
350: dd->fieldname = NULL;
351: dd->nlocal = -1;
352: dd->Nlocal = -1;
353: dd->M = -1;
354: dd->N = -1;
355: dd->P = -1;
356: dd->m = -1;
357: dd->n = -1;
358: dd->p = -1;
359: dd->w = -1;
360: dd->s = -1;
362: dd->xs = -1;
363: dd->xe = -1;
364: dd->ys = -1;
365: dd->ye = -1;
366: dd->zs = -1;
367: dd->ze = -1;
368: dd->Xs = -1;
369: dd->Xe = -1;
370: dd->Ys = -1;
371: dd->Ye = -1;
372: dd->Zs = -1;
373: dd->Ze = -1;
375: dd->Nsub = 1;
376: dd->xol = 0;
377: dd->yol = 0;
378: dd->zol = 0;
379: dd->xo = 0;
380: dd->yo = 0;
381: dd->zo = 0;
382: dd->Mo = -1;
383: dd->No = -1;
384: dd->Po = -1;
386: dd->gtol = NULL;
387: dd->ltol = NULL;
388: dd->ao = NULL;
389: PetscCall(PetscStrallocpy(AOBASIC, (char **)&dd->aotype));
390: dd->base = -1;
391: dd->bx = DM_BOUNDARY_NONE;
392: dd->by = DM_BOUNDARY_NONE;
393: dd->bz = DM_BOUNDARY_NONE;
394: dd->stencil_type = DMDA_STENCIL_BOX;
395: dd->interptype = DMDA_Q1;
396: dd->lx = NULL;
397: dd->ly = NULL;
398: dd->lz = NULL;
400: dd->elementtype = DMDA_ELEMENT_Q1;
402: da->ops->globaltolocalbegin = DMGlobalToLocalBegin_DA;
403: da->ops->globaltolocalend = DMGlobalToLocalEnd_DA;
404: da->ops->localtoglobalbegin = DMLocalToGlobalBegin_DA;
405: da->ops->localtoglobalend = DMLocalToGlobalEnd_DA;
406: da->ops->localtolocalbegin = DMLocalToLocalBegin_DA;
407: da->ops->localtolocalend = DMLocalToLocalEnd_DA;
408: da->ops->createglobalvector = DMCreateGlobalVector_DA;
409: da->ops->createlocalvector = DMCreateLocalVector_DA;
410: da->ops->createinterpolation = DMCreateInterpolation_DA;
411: da->ops->getcoloring = DMCreateColoring_DA;
412: da->ops->creatematrix = DMCreateMatrix_DA;
413: da->ops->refine = DMRefine_DA;
414: da->ops->coarsen = DMCoarsen_DA;
415: da->ops->refinehierarchy = DMRefineHierarchy_DA;
416: da->ops->coarsenhierarchy = DMCoarsenHierarchy_DA;
417: da->ops->createinjection = DMCreateInjection_DA;
418: da->ops->hascreateinjection = DMHasCreateInjection_DA;
419: da->ops->destroy = DMDestroy_DA;
420: da->ops->view = NULL;
421: da->ops->setfromoptions = DMSetFromOptions_DA;
422: da->ops->setup = DMSetUp_DA;
423: da->ops->clone = DMClone_DA;
424: da->ops->load = DMLoad_DA;
425: da->ops->createcoordinatedm = DMCreateCoordinateDM_DA;
426: da->ops->createcellcoordinatedm = NULL;
427: da->ops->createsubdm = DMCreateSubDM_DA;
428: da->ops->createfielddecomposition = DMCreateFieldDecomposition_DA;
429: da->ops->createdomaindecomposition = DMCreateDomainDecomposition_DA;
430: da->ops->createddscatters = DMCreateDomainDecompositionScatters_DA;
431: da->ops->getdimpoints = DMGetDimPoints_DA;
432: da->ops->getneighbors = DMGetNeighbors_DA;
433: da->ops->getlocalboundingbox = DMGetLocalBoundingBox_DA;
434: da->ops->locatepoints = DMLocatePoints_DA_Regular;
435: da->ops->getcompatibility = DMGetCompatibility_DA;
436: PetscCall(PetscObjectComposeFunction((PetscObject)da, "DMSetUpGLVisViewer_C", DMSetUpGLVisViewer_DMDA));
437: PetscFunctionReturn(PETSC_SUCCESS);
438: }
440: /*@
441: DMDACreate - Creates a `DMDA` object for managing structured grids.
443: Collective
445: Input Parameter:
446: . comm - The communicator for the `DMDA` object
448: Output Parameter:
449: . da - the `DMDA` object
451: Level: advanced
453: Notes:
454: See [](sec_struct) for details on the construction of a `DMDA`
456: `DMDACreate1d()`, `DMDACreate2d()`, and `DMDACreate3d()` are convenience routines to quickly completely create a `DMDA`
458: .seealso: [](sec_struct), `DM`, `DMDA`, `DMSetUp()`, `DMDASetSizes()`, `DMClone()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
459: @*/
460: PetscErrorCode DMDACreate(MPI_Comm comm, DM *da)
461: {
462: PetscFunctionBegin;
463: PetscAssertPointer(da, 2);
464: PetscCall(DMCreate(comm, da));
465: PetscCall(DMSetType(*da, DMDA));
466: PetscFunctionReturn(PETSC_SUCCESS);
467: }