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