xref: /petsc/src/dm/impls/da/dacreate.c (revision c688d0420b4e513ff34944d1e1ad7d4e50aafa8d)
1 
2 #include <petsc/private/dmdaimpl.h>    /*I   "petscdmda.h"   I*/
3 
4 #undef __FUNCT__
5 #define __FUNCT__ "DMSetFromOptions_DA"
6 PetscErrorCode  DMSetFromOptions_DA(PetscOptionItems *PetscOptionsObject,DM da)
7 {
8   PetscErrorCode ierr;
9   DM_DA          *dd    = (DM_DA*)da->data;
10   PetscInt       refine = 0,dim = da->dim,maxnlevels = 100,refx[100],refy[100],refz[100],n,i;
11   PetscBool      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     dd->negativeMNP = PETSC_TRUE;
20   }
21   if (da->dim > 1 && dd->N < 0) {
22     dd->N           = -dd->N;
23     bN              = PETSC_TRUE;
24     dd->negativeMNP = PETSC_TRUE;
25   }
26   if (da->dim > 2 && dd->P < 0) {
27     dd->P           = -dd->P;
28     bP              = PETSC_TRUE;
29     dd->negativeMNP = PETSC_TRUE;
30   }
31 
32   ierr = PetscOptionsHead(PetscOptionsObject,"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 (dim > 1) {ierr = PetscOptionsInt("-da_overlap_y","Decomposition overlap in y direction","DMDASetOverlap",dd->yol,&dd->yol,NULL);CHKERRQ(ierr);}
41   if (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 (dim > 1) {ierr = PetscOptionsInt("-da_processors_y","Number of processors in y direction","DMDASetNumProcs",dd->n,&dd->n,NULL);CHKERRQ(ierr);}
49   if (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 (dim > 1) {ierr = PetscOptionsInt("-da_refine_y","Refinement ratio in y direction","DMDASetRefinementFactor",dd->refine_y,&dd->refine_y,NULL);CHKERRQ(ierr);}
53   if (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 = PetscOptionsIntArray("-da_refine_hierarchy_x","Refinement factor for each level","None",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 (dim > 1) {
68     n    = maxnlevels;
69     ierr = PetscOptionsIntArray("-da_refine_hierarchy_y","Refinement factor for each level","None",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 (dim > 2) {
74     n    = maxnlevels;
75     ierr = PetscOptionsIntArray("-da_refine_hierarchy_z","Refinement factor for each level","None",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 (dd->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,Mat*);
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,NULL,PETSC_INT);CHKERRQ(ierr);
152   ierr = PetscViewerBinaryRead(viewer,&m,1,NULL,PETSC_INT);CHKERRQ(ierr);
153   ierr = PetscViewerBinaryRead(viewer,&n,1,NULL,PETSC_INT);CHKERRQ(ierr);
154   ierr = PetscViewerBinaryRead(viewer,&p,1,NULL,PETSC_INT);CHKERRQ(ierr);
155   ierr = PetscViewerBinaryRead(viewer,&dof,1,NULL,PETSC_INT);CHKERRQ(ierr);
156   ierr = PetscViewerBinaryRead(viewer,&swidth,1,NULL,PETSC_INT);CHKERRQ(ierr);
157   ierr = PetscViewerBinaryRead(viewer,&bx,1,NULL,PETSC_ENUM);CHKERRQ(ierr);
158   ierr = PetscViewerBinaryRead(viewer,&by,1,NULL,PETSC_ENUM);CHKERRQ(ierr);
159   ierr = PetscViewerBinaryRead(viewer,&bz,1,NULL,PETSC_ENUM);CHKERRQ(ierr);
160   ierr = PetscViewerBinaryRead(viewer,&stencil,1,NULL,PETSC_ENUM);CHKERRQ(ierr);
161 
162   ierr = DMSetDimension(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,NULL,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 = DMSetDimension(*subdm, dm->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 = DMSetDimension(da, dm->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 = DMSetDimension(*newdm, dm->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 #undef __FUNCT__
311 #define __FUNCT__ "DMGetDimPoints_DA"
312 static PetscErrorCode DMGetDimPoints_DA(DM dm, PetscInt dim, PetscInt *pStart, PetscInt *pEnd)
313 {
314   PetscErrorCode ierr;
315 
316   PetscFunctionBegin;
317   ierr = DMDAGetDepthStratum(dm, dim, pStart, pEnd);CHKERRQ(ierr);
318   PetscFunctionReturn(0);
319 }
320 
321 #undef __FUNCT__
322 #define __FUNCT__ "DMGetNeighbors_DA"
323 static PetscErrorCode DMGetNeighbors_DA(DM dm, PetscInt *nranks, const PetscMPIInt *ranks[])
324 {
325   PetscErrorCode ierr;
326   PetscInt dim;
327   DMDAStencilType st;
328 
329   PetscFunctionBegin;
330   ierr = DMDAGetNeighbors(dm,ranks);CHKERRQ(ierr);
331   ierr = DMDAGetInfo(dm,&dim,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,&st);CHKERRQ(ierr);
332 
333   switch (dim) {
334     case 1:
335       *nranks = 3;
336       /* if (st == DMDA_STENCIL_STAR) { *nranks = 3; } */
337       break;
338     case 2:
339       *nranks = 9;
340       /* if (st == DMDA_STENCIL_STAR) { *nranks = 5; } */
341       break;
342     case 3:
343       *nranks = 27;
344       /* if (st == DMDA_STENCIL_STAR) { *nranks = 7; } */
345       break;
346     default:
347       break;
348   }
349   PetscFunctionReturn(0);
350 }
351 
352 /*MC
353    DMDA = "da" - A DM object that is used to manage data for a structured grid in 1, 2, or 3 dimensions.
354          In the global representation of the vector each process stores a non-overlapping rectangular (or slab in 3d) portion of the grid points.
355          In the local representation these rectangular regions (slabs) are extended in all directions by a stencil width.
356 
357          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
358          vertex centered.
359 
360   Level: intermediate
361 
362 .seealso: DMType, DMCOMPOSITE, DMDACreate(), DMCreate(), DMSetType()
363 M*/
364 
365 extern PetscErrorCode DMProjectFunctionLocal_DA(DM, PetscReal, PetscErrorCode (**)(PetscInt, PetscReal, const PetscReal [], PetscInt, PetscScalar *, void *), void **, InsertMode, Vec);
366 extern PetscErrorCode DMComputeL2Diff_DA(DM, PetscReal, PetscErrorCode (**)(PetscInt, PetscReal, const PetscReal [], PetscInt, PetscScalar *, void *), void **, Vec, PetscReal *);
367 extern PetscErrorCode DMComputeL2GradientDiff_DA(DM, PetscReal, PetscErrorCode (**)(PetscInt, PetscReal, const PetscReal [], const PetscReal [],PetscInt, PetscScalar *, void *), void **, Vec,const PetscReal [], PetscReal *);
368 
369 
370 #undef __FUNCT__
371 #define __FUNCT__ "DMCreate_DA"
372 PETSC_EXTERN PetscErrorCode DMCreate_DA(DM da)
373 {
374   PetscErrorCode ierr;
375   DM_DA          *dd;
376 
377   PetscFunctionBegin;
378   PetscValidPointer(da,1);
379   ierr     = PetscNewLog(da,&dd);CHKERRQ(ierr);
380   da->data = dd;
381 
382   da->dim        = -1;
383   dd->interptype = DMDA_Q1;
384   dd->refine_x   = 2;
385   dd->refine_y   = 2;
386   dd->refine_z   = 2;
387   dd->coarsen_x  = 2;
388   dd->coarsen_y  = 2;
389   dd->coarsen_z  = 2;
390   dd->fieldname  = NULL;
391   dd->nlocal     = -1;
392   dd->Nlocal     = -1;
393   dd->M          = -1;
394   dd->N          = -1;
395   dd->P          = -1;
396   dd->m          = -1;
397   dd->n          = -1;
398   dd->p          = -1;
399   dd->w          = -1;
400   dd->s          = -1;
401 
402   dd->xs = -1; dd->xe = -1; dd->ys = -1; dd->ye = -1; dd->zs = -1; dd->ze = -1;
403   dd->Xs = -1; dd->Xe = -1; dd->Ys = -1; dd->Ye = -1; dd->Zs = -1; dd->Ze = -1;
404 
405   dd->Nsub            = 1;
406   dd->xol             = 0;
407   dd->yol             = 0;
408   dd->zol             = 0;
409   dd->xo              = 0;
410   dd->yo              = 0;
411   dd->zo              = 0;
412   dd->Mo              = -1;
413   dd->No              = -1;
414   dd->Po              = -1;
415 
416   dd->gtol         = NULL;
417   dd->ltol         = NULL;
418   dd->ao           = NULL;
419   PetscStrallocpy(AOBASIC,(char**)&dd->aotype);
420   dd->base         = -1;
421   dd->bx           = DM_BOUNDARY_NONE;
422   dd->by           = DM_BOUNDARY_NONE;
423   dd->bz           = DM_BOUNDARY_NONE;
424   dd->stencil_type = DMDA_STENCIL_BOX;
425   dd->interptype   = DMDA_Q1;
426   dd->lx           = NULL;
427   dd->ly           = NULL;
428   dd->lz           = NULL;
429 
430   dd->elementtype = DMDA_ELEMENT_Q1;
431 
432   da->ops->globaltolocalbegin          = DMGlobalToLocalBegin_DA;
433   da->ops->globaltolocalend            = DMGlobalToLocalEnd_DA;
434   da->ops->localtoglobalbegin          = DMLocalToGlobalBegin_DA;
435   da->ops->localtoglobalend            = DMLocalToGlobalEnd_DA;
436   da->ops->localtolocalbegin           = DMLocalToLocalBegin_DA;
437   da->ops->localtolocalend             = DMLocalToLocalEnd_DA;
438   da->ops->createglobalvector          = DMCreateGlobalVector_DA;
439   da->ops->createlocalvector           = DMCreateLocalVector_DA;
440   da->ops->createinterpolation         = DMCreateInterpolation_DA;
441   da->ops->getcoloring                 = DMCreateColoring_DA;
442   da->ops->creatematrix                = DMCreateMatrix_DA;
443   da->ops->refine                      = DMRefine_DA;
444   da->ops->coarsen                     = DMCoarsen_DA;
445   da->ops->refinehierarchy             = DMRefineHierarchy_DA;
446   da->ops->coarsenhierarchy            = DMCoarsenHierarchy_DA;
447   da->ops->getinjection                = DMCreateInjection_DA;
448   da->ops->getaggregates               = DMCreateAggregates_DA;
449   da->ops->destroy                     = DMDestroy_DA;
450   da->ops->view                        = 0;
451   da->ops->setfromoptions              = DMSetFromOptions_DA;
452   da->ops->setup                       = DMSetUp_DA;
453   da->ops->clone                       = DMClone_DA;
454   da->ops->load                        = DMLoad_DA;
455   da->ops->createcoordinatedm          = DMCreateCoordinateDM_DA;
456   da->ops->createsubdm                 = DMCreateSubDM_DA;
457   da->ops->createfielddecomposition    = DMCreateFieldDecomposition_DA;
458   da->ops->createdomaindecomposition   = DMCreateDomainDecomposition_DA;
459   da->ops->createddscatters            = DMCreateDomainDecompositionScatters_DA;
460   da->ops->getdimpoints                = DMGetDimPoints_DA;
461   da->ops->projectfunctionlocal        = DMProjectFunctionLocal_DA;
462   da->ops->computel2diff               = DMComputeL2Diff_DA;
463   da->ops->computel2gradientdiff       = DMComputeL2GradientDiff_DA;
464   da->ops->getneighbors                = DMGetNeighbors_DA;
465   PetscFunctionReturn(0);
466 }
467 
468 #undef __FUNCT__
469 #define __FUNCT__ "DMDACreate"
470 /*@
471   DMDACreate - Creates a DMDA object.
472 
473   Collective on MPI_Comm
474 
475   Input Parameter:
476 . comm - The communicator for the DMDA object
477 
478   Output Parameter:
479 . da  - The DMDA object
480 
481   Level: advanced
482 
483   Developers Note: Since there exists DMDACreate1/2/3d() should this routine even exist?
484 
485 .keywords: DMDA, create
486 .seealso:  DMDASetSizes(), DMDADuplicate(),  DMDACreate1d(), DMDACreate2d(), DMDACreate3d()
487 @*/
488 PetscErrorCode  DMDACreate(MPI_Comm comm, DM *da)
489 {
490   PetscErrorCode ierr;
491 
492   PetscFunctionBegin;
493   PetscValidPointer(da,2);
494   ierr = DMCreate(comm,da);CHKERRQ(ierr);
495   ierr = DMSetType(*da,DMDA);CHKERRQ(ierr);
496   PetscFunctionReturn(0);
497 }
498 
499 
500