xref: /petsc/src/dm/impls/da/dacreate.c (revision 2c9ec6dfe874b911fa49ef7e759d29a8430d6aff)
1 
2 #include <petsc-private/dmdaimpl.h>    /*I   "petscdmda.h"   I*/
3 
4 #undef __FUNCT__
5 #define __FUNCT__ "DMSetFromOptions_DA"
6 PetscErrorCode  DMSetFromOptions_DA(DM da)
7 {
8   PetscErrorCode ierr;
9   DM_DA          *dd         = (DM_DA*)da->data;
10   PetscInt       refine      = 0,maxnlevels = 100,refx[100],refy[100],refz[100],n,i;
11   PetscBool      negativeMNP = PETSC_FALSE,bM = PETSC_FALSE,bN = PETSC_FALSE, bP = PETSC_FALSE,flg;
12 
13   PetscFunctionBegin;
14   PetscValidHeaderSpecific(da,DM_CLASSID,1);
15 
16   if (dd->M < 0) {
17     dd->M       = -dd->M;
18     bM          = PETSC_TRUE;
19     negativeMNP = PETSC_TRUE;
20   }
21   if (dd->dim > 1 && dd->N < 0) {
22     dd->N       = -dd->N;
23     bN          = PETSC_TRUE;
24     negativeMNP = PETSC_TRUE;
25   }
26   if (dd->dim > 2 && dd->P < 0) {
27     dd->P       = -dd->P;
28     bP          = PETSC_TRUE;
29     negativeMNP = PETSC_TRUE;
30   }
31 
32   ierr = PetscOptionsHead("DMDA Options");CHKERRQ(ierr);
33   if (bM) {ierr = PetscOptionsInt("-da_grid_x","Number of grid points in x direction","DMDASetSizes",dd->M,&dd->M,NULL);CHKERRQ(ierr);}
34   if (bN) {ierr = PetscOptionsInt("-da_grid_y","Number of grid points in y direction","DMDASetSizes",dd->N,&dd->N,NULL);CHKERRQ(ierr);}
35   if (bP) {ierr = PetscOptionsInt("-da_grid_z","Number of grid points in z direction","DMDASetSizes",dd->P,&dd->P,NULL);CHKERRQ(ierr);}
36 
37   ierr = PetscOptionsInt("-da_overlap","Decomposition overlap in all directions","DMDASetOverlap",dd->xol,&dd->xol,&flg);CHKERRQ(ierr);
38   if (flg) {ierr = DMDASetOverlap(da,dd->xol,dd->xol,dd->xol);CHKERRQ(ierr);}
39   ierr = PetscOptionsInt("-da_overlap_x","Decomposition overlap in x direction","DMDASetOverlap",dd->xol,&dd->xol,NULL);CHKERRQ(ierr);
40   if (dd->dim > 1) {ierr = PetscOptionsInt("-da_overlap_y","Decomposition overlap in y direction","DMDASetOverlap",dd->yol,&dd->yol,NULL);CHKERRQ(ierr);}
41   if (dd->dim > 2) {ierr = PetscOptionsInt("-da_overlap_z","Decomposition overlap in z direction","DMDASetOverlap",dd->zol,&dd->zol,NULL);CHKERRQ(ierr);}
42 
43   ierr = PetscOptionsInt("-da_local_subdomains","","DMDASetNumLocalSubdomains",dd->Nsub,&dd->Nsub,&flg);CHKERRQ(ierr);
44   if (flg) {ierr = DMDASetNumLocalSubDomains(da,dd->Nsub);CHKERRQ(ierr);}
45 
46   /* Handle DMDA parallel distibution */
47   ierr = PetscOptionsInt("-da_processors_x","Number of processors in x direction","DMDASetNumProcs",dd->m,&dd->m,NULL);CHKERRQ(ierr);
48   if (dd->dim > 1) {ierr = PetscOptionsInt("-da_processors_y","Number of processors in y direction","DMDASetNumProcs",dd->n,&dd->n,NULL);CHKERRQ(ierr);}
49   if (dd->dim > 2) {ierr = PetscOptionsInt("-da_processors_z","Number of processors in z direction","DMDASetNumProcs",dd->p,&dd->p,NULL);CHKERRQ(ierr);}
50   /* Handle DMDA refinement */
51   ierr = PetscOptionsInt("-da_refine_x","Refinement ratio in x direction","DMDASetRefinementFactor",dd->refine_x,&dd->refine_x,NULL);CHKERRQ(ierr);
52   if (dd->dim > 1) {ierr = PetscOptionsInt("-da_refine_y","Refinement ratio in y direction","DMDASetRefinementFactor",dd->refine_y,&dd->refine_y,NULL);CHKERRQ(ierr);}
53   if (dd->dim > 2) {ierr = PetscOptionsInt("-da_refine_z","Refinement ratio in z direction","DMDASetRefinementFactor",dd->refine_z,&dd->refine_z,NULL);CHKERRQ(ierr);}
54   dd->coarsen_x = dd->refine_x; dd->coarsen_y = dd->refine_y; dd->coarsen_z = dd->refine_z;
55 
56   /* Get refinement factors, defaults taken from the coarse DMDA */
57   ierr = DMDAGetRefinementFactor(da,&refx[0],&refy[0],&refz[0]);CHKERRQ(ierr);
58   for (i=1; i<maxnlevels; i++) {
59     refx[i] = refx[0];
60     refy[i] = refy[0];
61     refz[i] = refz[0];
62   }
63   n    = maxnlevels;
64   ierr = PetscOptionsGetIntArray(((PetscObject)da)->prefix,"-da_refine_hierarchy_x",refx,&n,NULL);CHKERRQ(ierr);
65   if (da->levelup - da->leveldown >= 0) dd->refine_x = refx[da->levelup - da->leveldown];
66   if (da->levelup - da->leveldown >= 1) dd->coarsen_x = refx[da->levelup - da->leveldown - 1];
67   if (dd->dim > 1) {
68     n    = maxnlevels;
69     ierr = PetscOptionsGetIntArray(((PetscObject)da)->prefix,"-da_refine_hierarchy_y",refy,&n,NULL);CHKERRQ(ierr);
70     if (da->levelup - da->leveldown >= 0) dd->refine_y = refy[da->levelup - da->leveldown];
71     if (da->levelup - da->leveldown >= 1) dd->coarsen_y = refy[da->levelup - da->leveldown - 1];
72   }
73   if (dd->dim > 2) {
74     n    = maxnlevels;
75     ierr = PetscOptionsGetIntArray(((PetscObject)da)->prefix,"-da_refine_hierarchy_z",refz,&n,NULL);CHKERRQ(ierr);
76     if (da->levelup - da->leveldown >= 0) dd->refine_z = refz[da->levelup - da->leveldown];
77     if (da->levelup - da->leveldown >= 1) dd->coarsen_z = refz[da->levelup - da->leveldown - 1];
78   }
79 
80   if (negativeMNP) {ierr = PetscOptionsInt("-da_refine","Uniformly refine DA one or more times","None",refine,&refine,NULL);CHKERRQ(ierr);}
81   ierr = PetscOptionsTail();CHKERRQ(ierr);
82 
83   while (refine--) {
84     if (dd->bx == DM_BOUNDARY_PERIODIC || dd->interptype == DMDA_Q0) {
85       dd->M = dd->refine_x*dd->M;
86     } else {
87       dd->M = 1 + dd->refine_x*(dd->M - 1);
88     }
89     if (dd->by == DM_BOUNDARY_PERIODIC || dd->interptype == DMDA_Q0) {
90       dd->N = dd->refine_y*dd->N;
91     } else {
92       dd->N = 1 + dd->refine_y*(dd->N - 1);
93     }
94     if (dd->bz == DM_BOUNDARY_PERIODIC || dd->interptype == DMDA_Q0) {
95       dd->P = dd->refine_z*dd->P;
96     } else {
97       dd->P = 1 + dd->refine_z*(dd->P - 1);
98     }
99     da->levelup++;
100     if (da->levelup - da->leveldown >= 0) {
101       dd->refine_x = refx[da->levelup - da->leveldown];
102       dd->refine_y = refy[da->levelup - da->leveldown];
103       dd->refine_z = refz[da->levelup - da->leveldown];
104     }
105     if (da->levelup - da->leveldown >= 1) {
106       dd->coarsen_x = refx[da->levelup - da->leveldown - 1];
107       dd->coarsen_y = refy[da->levelup - da->leveldown - 1];
108       dd->coarsen_z = refz[da->levelup - da->leveldown - 1];
109     }
110   }
111   PetscFunctionReturn(0);
112 }
113 
114 extern PetscErrorCode  DMCreateGlobalVector_DA(DM,Vec*);
115 extern PetscErrorCode  DMCreateLocalVector_DA(DM,Vec*);
116 extern PetscErrorCode  DMGlobalToLocalBegin_DA(DM,Vec,InsertMode,Vec);
117 extern PetscErrorCode  DMGlobalToLocalEnd_DA(DM,Vec,InsertMode,Vec);
118 extern PetscErrorCode  DMLocalToGlobalBegin_DA(DM,Vec,InsertMode,Vec);
119 extern PetscErrorCode  DMLocalToGlobalEnd_DA(DM,Vec,InsertMode,Vec);
120 extern PetscErrorCode  DMLocalToLocalBegin_DA(DM,Vec,InsertMode,Vec);
121 extern PetscErrorCode  DMLocalToLocalEnd_DA(DM,Vec,InsertMode,Vec);
122 extern PetscErrorCode  DMCreateInterpolation_DA(DM,DM,Mat*,Vec*);
123 extern PetscErrorCode  DMCreateColoring_DA(DM,ISColoringType,ISColoring*);
124 extern PetscErrorCode  DMCreateMatrix_DA(DM,Mat*);
125 extern PetscErrorCode  DMCreateCoordinateDM_DA(DM,DM*);
126 extern PetscErrorCode  DMRefine_DA(DM,MPI_Comm,DM*);
127 extern PetscErrorCode  DMCoarsen_DA(DM,MPI_Comm,DM*);
128 extern PetscErrorCode  DMRefineHierarchy_DA(DM,PetscInt,DM[]);
129 extern PetscErrorCode  DMCoarsenHierarchy_DA(DM,PetscInt,DM[]);
130 extern PetscErrorCode  DMCreateInjection_DA(DM,DM,VecScatter*);
131 extern PetscErrorCode  DMCreateAggregates_DA(DM,DM,Mat*);
132 extern PetscErrorCode  DMView_DA(DM,PetscViewer);
133 extern PetscErrorCode  DMSetUp_DA(DM);
134 extern PetscErrorCode  DMDestroy_DA(DM);
135 extern PetscErrorCode  DMCreateDomainDecomposition_DA(DM,PetscInt*,char***,IS**,IS**,DM**);
136 extern PetscErrorCode  DMCreateDomainDecompositionScatters_DA(DM,PetscInt,DM*,VecScatter**,VecScatter**,VecScatter**);
137 
138 #undef __FUNCT__
139 #define __FUNCT__ "DMLoad_DA"
140 PetscErrorCode DMLoad_DA(DM da,PetscViewer viewer)
141 {
142   PetscErrorCode   ierr;
143   PetscInt         dim,m,n,p,dof,swidth;
144   DMDAStencilType  stencil;
145   DMBoundaryType   bx,by,bz;
146   PetscBool        coors;
147   DM               dac;
148   Vec              c;
149 
150   PetscFunctionBegin;
151   ierr = PetscViewerBinaryRead(viewer,&dim,1,PETSC_INT);CHKERRQ(ierr);
152   ierr = PetscViewerBinaryRead(viewer,&m,1,PETSC_INT);CHKERRQ(ierr);
153   ierr = PetscViewerBinaryRead(viewer,&n,1,PETSC_INT);CHKERRQ(ierr);
154   ierr = PetscViewerBinaryRead(viewer,&p,1,PETSC_INT);CHKERRQ(ierr);
155   ierr = PetscViewerBinaryRead(viewer,&dof,1,PETSC_INT);CHKERRQ(ierr);
156   ierr = PetscViewerBinaryRead(viewer,&swidth,1,PETSC_INT);CHKERRQ(ierr);
157   ierr = PetscViewerBinaryRead(viewer,&bx,1,PETSC_ENUM);CHKERRQ(ierr);
158   ierr = PetscViewerBinaryRead(viewer,&by,1,PETSC_ENUM);CHKERRQ(ierr);
159   ierr = PetscViewerBinaryRead(viewer,&bz,1,PETSC_ENUM);CHKERRQ(ierr);
160   ierr = PetscViewerBinaryRead(viewer,&stencil,1,PETSC_ENUM);CHKERRQ(ierr);
161 
162   ierr = DMDASetDim(da, dim);CHKERRQ(ierr);
163   ierr = DMDASetSizes(da, m,n,p);CHKERRQ(ierr);
164   ierr = DMDASetBoundaryType(da, bx, by, bz);CHKERRQ(ierr);
165   ierr = DMDASetDof(da, dof);CHKERRQ(ierr);
166   ierr = DMDASetStencilType(da, stencil);CHKERRQ(ierr);
167   ierr = DMDASetStencilWidth(da, swidth);CHKERRQ(ierr);
168   ierr = DMSetUp(da);CHKERRQ(ierr);
169   ierr = PetscViewerBinaryRead(viewer,&coors,1,PETSC_ENUM);CHKERRQ(ierr);
170   if (coors) {
171     ierr = DMGetCoordinateDM(da,&dac);CHKERRQ(ierr);
172     ierr = DMCreateGlobalVector(dac,&c);CHKERRQ(ierr);
173     ierr = VecLoad(c,viewer);CHKERRQ(ierr);
174     ierr = DMSetCoordinates(da,c);CHKERRQ(ierr);
175     ierr = VecDestroy(&c);CHKERRQ(ierr);
176   }
177   PetscFunctionReturn(0);
178 }
179 
180 #undef __FUNCT__
181 #define __FUNCT__ "DMCreateSubDM_DA"
182 PetscErrorCode DMCreateSubDM_DA(DM dm, PetscInt numFields, PetscInt fields[], IS *is, DM *subdm)
183 {
184   DM_DA         *da = (DM_DA*) dm->data;
185   PetscSection   section;
186   PetscErrorCode ierr;
187 
188   PetscFunctionBegin;
189   if (subdm) {
190     PetscSF sf;
191     Vec     coords;
192     void   *ctx;
193     /* Cannot use DMClone since the dof stuff is mixed in. Ugh
194     ierr = DMClone(dm, subdm);CHKERRQ(ierr); */
195     ierr = DMCreate(PetscObjectComm((PetscObject)dm), subdm);CHKERRQ(ierr);
196     ierr = DMGetPointSF(dm, &sf);CHKERRQ(ierr);
197     ierr = DMSetPointSF(*subdm, sf);CHKERRQ(ierr);
198     ierr = DMGetApplicationContext(dm, &ctx);CHKERRQ(ierr);
199     ierr = DMSetApplicationContext(*subdm, ctx);CHKERRQ(ierr);
200     ierr = DMGetCoordinatesLocal(dm, &coords);CHKERRQ(ierr);
201     if (coords) {
202       ierr = DMSetCoordinatesLocal(*subdm, coords);CHKERRQ(ierr);
203     } else {
204       ierr = DMGetCoordinates(dm, &coords);CHKERRQ(ierr);
205       if (coords) {ierr = DMSetCoordinates(*subdm, coords);CHKERRQ(ierr);}
206     }
207 
208     ierr = DMSetType(*subdm, DMDA);CHKERRQ(ierr);
209     ierr = DMDASetDim(*subdm, da->dim);CHKERRQ(ierr);
210     ierr = DMDASetSizes(*subdm, da->M, da->N, da->P);CHKERRQ(ierr);
211     ierr = DMDASetNumProcs(*subdm, da->m, da->n, da->p);CHKERRQ(ierr);
212     ierr = DMDASetBoundaryType(*subdm, da->bx, da->by, da->bz);CHKERRQ(ierr);
213     ierr = DMDASetDof(*subdm, numFields);CHKERRQ(ierr);
214     ierr = DMDASetStencilType(*subdm, da->stencil_type);CHKERRQ(ierr);
215     ierr = DMDASetStencilWidth(*subdm, da->s);CHKERRQ(ierr);
216     ierr = DMDASetOwnershipRanges(*subdm, da->lx, da->ly, da->lz);CHKERRQ(ierr);
217   }
218   ierr = DMGetDefaultSection(dm, &section);CHKERRQ(ierr);
219   if (section) {
220     ierr = DMCreateSubDM_Section_Private(dm, numFields, fields, is, subdm);CHKERRQ(ierr);
221   } else {
222     if (is) {
223       PetscInt *indices, cnt = 0, dof = da->w, i, j;
224 
225       ierr = PetscMalloc1(da->Nlocal*numFields/dof, &indices);CHKERRQ(ierr);
226       for (i = da->base/dof; i < (da->base+da->Nlocal)/dof; ++i) {
227         for (j = 0; j < numFields; ++j) {
228           indices[cnt++] = dof*i + fields[j];
229         }
230       }
231       if (cnt != da->Nlocal*numFields/dof) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Count %d does not equal expected value %d", cnt, da->Nlocal*numFields/dof);
232       ierr = ISCreateGeneral(PetscObjectComm((PetscObject) dm), cnt, indices, PETSC_OWN_POINTER, is);CHKERRQ(ierr);
233     }
234   }
235   PetscFunctionReturn(0);
236 }
237 
238 #undef __FUNCT__
239 #define __FUNCT__ "DMCreateFieldDecomposition_DA"
240 PetscErrorCode DMCreateFieldDecomposition_DA(DM dm, PetscInt *len,char ***namelist, IS **islist, DM **dmlist)
241 {
242   PetscInt       i;
243   PetscErrorCode ierr;
244   DM_DA          *dd = (DM_DA*)dm->data;
245   PetscInt       dof = dd->w;
246 
247   PetscFunctionBegin;
248   if (len) *len = dof;
249   if (islist) {
250     Vec      v;
251     PetscInt rstart,n;
252 
253     ierr = DMGetGlobalVector(dm,&v);CHKERRQ(ierr);
254     ierr = VecGetOwnershipRange(v,&rstart,NULL);CHKERRQ(ierr);
255     ierr = VecGetLocalSize(v,&n);CHKERRQ(ierr);
256     ierr = DMRestoreGlobalVector(dm,&v);CHKERRQ(ierr);
257     ierr = PetscMalloc1(dof,islist);CHKERRQ(ierr);
258     for (i=0; i<dof; i++) {
259       ierr = ISCreateStride(PetscObjectComm((PetscObject)dm),n/dof,rstart+i,dof,&(*islist)[i]);CHKERRQ(ierr);
260     }
261   }
262   if (namelist) {
263     ierr = PetscMalloc1(dof, namelist);CHKERRQ(ierr);
264     if (dd->fieldname) {
265       for (i=0; i<dof; i++) {
266         ierr = PetscStrallocpy(dd->fieldname[i],&(*namelist)[i]);CHKERRQ(ierr);
267       }
268     } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Currently DMDA must have fieldnames");
269   }
270   if (dmlist) {
271     DM da;
272 
273     ierr = DMDACreate(PetscObjectComm((PetscObject)dm), &da);CHKERRQ(ierr);
274     ierr = DMDASetDim(da, dd->dim);CHKERRQ(ierr);
275     ierr = DMDASetSizes(da, dd->M, dd->N, dd->P);CHKERRQ(ierr);
276     ierr = DMDASetNumProcs(da, dd->m, dd->n, dd->p);CHKERRQ(ierr);
277     ierr = DMDASetBoundaryType(da, dd->bx, dd->by, dd->bz);CHKERRQ(ierr);
278     ierr = DMDASetDof(da, 1);CHKERRQ(ierr);
279     ierr = DMDASetStencilType(da, dd->stencil_type);CHKERRQ(ierr);
280     ierr = DMDASetStencilWidth(da, dd->s);CHKERRQ(ierr);
281     ierr = DMSetUp(da);CHKERRQ(ierr);
282     ierr = PetscMalloc1(dof,dmlist);CHKERRQ(ierr);
283     for (i=0; i<dof-1; i++) {ierr = PetscObjectReference((PetscObject)da);CHKERRQ(ierr);}
284     for (i=0; i<dof; i++) (*dmlist)[i] = da;
285   }
286   PetscFunctionReturn(0);
287 }
288 
289 #undef __FUNCT__
290 #define __FUNCT__ "DMClone_DA"
291 PetscErrorCode DMClone_DA(DM dm, DM *newdm)
292 {
293   DM_DA         *da = (DM_DA *) dm->data;
294   PetscErrorCode ierr;
295 
296   PetscFunctionBegin;
297   ierr = DMSetType(*newdm, DMDA);CHKERRQ(ierr);
298   ierr = DMDASetDim(*newdm, da->dim);CHKERRQ(ierr);
299   ierr = DMDASetSizes(*newdm, da->M, da->N, da->P);CHKERRQ(ierr);
300   ierr = DMDASetNumProcs(*newdm, da->m, da->n, da->p);CHKERRQ(ierr);
301   ierr = DMDASetBoundaryType(*newdm, da->bx, da->by, da->bz);CHKERRQ(ierr);
302   ierr = DMDASetDof(*newdm, da->w);CHKERRQ(ierr);
303   ierr = DMDASetStencilType(*newdm, da->stencil_type);CHKERRQ(ierr);
304   ierr = DMDASetStencilWidth(*newdm, da->s);CHKERRQ(ierr);
305   ierr = DMDASetOwnershipRanges(*newdm, da->lx, da->ly, da->lz);CHKERRQ(ierr);
306   ierr = DMSetUp(*newdm);CHKERRQ(ierr);
307   PetscFunctionReturn(0);
308 }
309 
310 /*MC
311    DMDA = "da" - A DM object that is used to manage data for a structured grid in 1, 2, or 3 dimensions.
312          In the global representation of the vector each process stores a non-overlapping rectangular (or slab in 3d) portion of the grid points.
313          In the local representation these rectangular regions (slabs) are extended in all directions by a stencil width.
314 
315          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
316          vertex centered.
317 
318   Level: intermediate
319 
320 .seealso: DMType, DMCOMPOSITE, DMDACreate(), DMCreate(), DMSetType()
321 M*/
322 
323 
324 #undef __FUNCT__
325 #define __FUNCT__ "DMCreate_DA"
326 PETSC_EXTERN PetscErrorCode DMCreate_DA(DM da)
327 {
328   PetscErrorCode ierr;
329   DM_DA          *dd;
330 
331   PetscFunctionBegin;
332   PetscValidPointer(da,1);
333   ierr     = PetscNewLog(da,&dd);CHKERRQ(ierr);
334   da->data = dd;
335 
336   dd->dim        = -1;
337   dd->interptype = DMDA_Q1;
338   dd->refine_x   = 2;
339   dd->refine_y   = 2;
340   dd->refine_z   = 2;
341   dd->coarsen_x  = 2;
342   dd->coarsen_y  = 2;
343   dd->coarsen_z  = 2;
344   dd->fieldname  = NULL;
345   dd->nlocal     = -1;
346   dd->Nlocal     = -1;
347   dd->M          = -1;
348   dd->N          = -1;
349   dd->P          = -1;
350   dd->m          = -1;
351   dd->n          = -1;
352   dd->p          = -1;
353   dd->w          = -1;
354   dd->s          = -1;
355 
356   dd->xs = -1; dd->xe = -1; dd->ys = -1; dd->ye = -1; dd->zs = -1; dd->ze = -1;
357   dd->Xs = -1; dd->Xe = -1; dd->Ys = -1; dd->Ye = -1; dd->Zs = -1; dd->Ze = -1;
358 
359   dd->Nsub            = 1;
360   dd->xol             = 0;
361   dd->yol             = 0;
362   dd->zol             = 0;
363   dd->xo              = 0;
364   dd->yo              = 0;
365   dd->zo              = 0;
366   dd->Mo              = -1;
367   dd->No              = -1;
368   dd->Po              = -1;
369 
370   dd->gtol         = NULL;
371   dd->ltol         = NULL;
372   dd->ao           = NULL;
373   dd->base         = -1;
374   dd->bx           = DM_BOUNDARY_NONE;
375   dd->by           = DM_BOUNDARY_NONE;
376   dd->bz           = DM_BOUNDARY_NONE;
377   dd->stencil_type = DMDA_STENCIL_BOX;
378   dd->interptype   = DMDA_Q1;
379   dd->lx           = NULL;
380   dd->ly           = NULL;
381   dd->lz           = NULL;
382 
383   dd->elementtype = DMDA_ELEMENT_Q1;
384 
385   da->ops->globaltolocalbegin          = DMGlobalToLocalBegin_DA;
386   da->ops->globaltolocalend            = DMGlobalToLocalEnd_DA;
387   da->ops->localtoglobalbegin          = DMLocalToGlobalBegin_DA;
388   da->ops->localtoglobalend            = DMLocalToGlobalEnd_DA;
389   da->ops->localtolocalbegin           = DMLocalToLocalBegin_DA;
390   da->ops->localtolocalend             = DMLocalToLocalEnd_DA;
391   da->ops->createglobalvector          = DMCreateGlobalVector_DA;
392   da->ops->createlocalvector           = DMCreateLocalVector_DA;
393   da->ops->createinterpolation         = DMCreateInterpolation_DA;
394   da->ops->getcoloring                 = DMCreateColoring_DA;
395   da->ops->creatematrix                = DMCreateMatrix_DA;
396   da->ops->refine                      = DMRefine_DA;
397   da->ops->coarsen                     = DMCoarsen_DA;
398   da->ops->refinehierarchy             = DMRefineHierarchy_DA;
399   da->ops->coarsenhierarchy            = DMCoarsenHierarchy_DA;
400   da->ops->getinjection                = DMCreateInjection_DA;
401   da->ops->getaggregates               = DMCreateAggregates_DA;
402   da->ops->destroy                     = DMDestroy_DA;
403   da->ops->view                        = 0;
404   da->ops->setfromoptions              = DMSetFromOptions_DA;
405   da->ops->setup                       = DMSetUp_DA;
406   da->ops->clone                       = DMClone_DA;
407   da->ops->load                        = DMLoad_DA;
408   da->ops->createcoordinatedm          = DMCreateCoordinateDM_DA;
409   da->ops->createsubdm                 = DMCreateSubDM_DA;
410   da->ops->createfielddecomposition    = DMCreateFieldDecomposition_DA;
411   da->ops->createdomaindecomposition   = DMCreateDomainDecomposition_DA;
412   da->ops->createddscatters            = DMCreateDomainDecompositionScatters_DA;
413   PetscFunctionReturn(0);
414 }
415 
416 #undef __FUNCT__
417 #define __FUNCT__ "DMDACreate"
418 /*@
419   DMDACreate - Creates a DMDA object.
420 
421   Collective on MPI_Comm
422 
423   Input Parameter:
424 . comm - The communicator for the DMDA object
425 
426   Output Parameter:
427 . da  - The DMDA object
428 
429   Level: advanced
430 
431   Developers Note: Since there exists DMDACreate1/2/3d() should this routine even exist?
432 
433 .keywords: DMDA, create
434 .seealso:  DMDASetSizes(), DMDADuplicate(),  DMDACreate1d(), DMDACreate2d(), DMDACreate3d()
435 @*/
436 PetscErrorCode  DMDACreate(MPI_Comm comm, DM *da)
437 {
438   PetscErrorCode ierr;
439 
440   PetscFunctionBegin;
441   PetscValidPointer(da,2);
442   ierr = DMCreate(comm,da);CHKERRQ(ierr);
443   ierr = DMSetType(*da,DMDA);CHKERRQ(ierr);
444   PetscFunctionReturn(0);
445 }
446 
447 
448