Actual source code: petscdmda.h
petsc-3.3-p7 2013-05-11
4: #include petscdm.h
5: #include petscpf.h
6: #include petscao.h
8: /*E
9: DMDAStencilType - Determines if the stencil extends only along the coordinate directions, or also
10: to the northeast, northwest etc
12: Level: beginner
14: .seealso: DMDACreate1d(), DMDACreate2d(), DMDACreate3d(), DMDACreate()
15: E*/
16: typedef enum { DMDA_STENCIL_STAR,DMDA_STENCIL_BOX } DMDAStencilType;
18: /*MC
19: DMDA_STENCIL_STAR - "Star"-type stencil. In logical grid coordinates, only (i,j,k), (i+s,j,k), (i,j+s,k),
20: (i,j,k+s) are in the stencil NOT, for example, (i+s,j+s,k)
22: Level: beginner
24: .seealso: DMDA_STENCIL_BOX, DMDAStencilType
25: M*/
27: /*MC
28: DMDA_STENCIL_BOX - "Box"-type stencil. In logical grid coordinates, any of (i,j,k), (i+s,j+r,k+t) may
29: be in the stencil.
31: Level: beginner
33: .seealso: DMDA_STENCIL_STAR, DMDAStencilType
34: M*/
36: /*E
37: DMDABoundaryType - Describes the choice for fill of ghost cells on physical domain boundaries.
39: Level: beginner
41: A boundary may be of type DMDA_BOUNDARY_NONE (no ghost nodes), DMDA_BOUNDARY_GHOST (ghost nodes
42: exist but aren't filled, you can put values into them and then apply a stencil that uses those ghost locations),
43: DMDA_BOUNDARY_MIRROR (not yet implemented), or DMDA_BOUNDARY_PERIODIC
44: (ghost nodes filled by the opposite edge of the domain).
46: Note: This is information for the boundary of the __PHYSICAL__ domain. It has nothing to do with boundaries between
47: processes, that width is always determined by the stencil width, see DMDASetStencilWidth().
49: .seealso: DMDASetBoundaryType(), DMDACreate1d(), DMDACreate2d(), DMDACreate3d(), DMDACreate()
50: E*/
51: typedef enum { DMDA_BOUNDARY_NONE, DMDA_BOUNDARY_GHOSTED, DMDA_BOUNDARY_MIRROR, DMDA_BOUNDARY_PERIODIC } DMDABoundaryType;
53: PETSC_EXTERN const char *DMDABoundaryTypes[];
55: /*E
56: DMDAInterpolationType - Defines the type of interpolation that will be returned by
57: DMCreateInterpolation.
59: Level: beginner
61: .seealso: DMDACreate1d(), DMDACreate2d(), DMDACreate3d(), DMCreateInterpolation(), DMDASetInterpolationType(), DMDACreate()
62: E*/
63: typedef enum { DMDA_Q0, DMDA_Q1 } DMDAInterpolationType;
65: PETSC_EXTERN PetscErrorCode DMDASetInterpolationType(DM,DMDAInterpolationType);
66: PETSC_EXTERN PetscErrorCode DMDAGetInterpolationType(DM,DMDAInterpolationType*);
68: /*E
69: DMDAElementType - Defines the type of elements that will be returned by
70: DMDAGetElements()
72: Level: beginner
74: .seealso: DMDACreate1d(), DMDACreate2d(), DMDACreate3d(), DMCreateInterpolation(), DMDASetInterpolationType(),
75: DMDASetElementType(), DMDAGetElements(), DMDARestoreElements(), DMDACreate()
76: E*/
77: typedef enum { DMDA_ELEMENT_P1, DMDA_ELEMENT_Q1 } DMDAElementType;
79: PETSC_EXTERN PetscErrorCode DMDASetElementType(DM,DMDAElementType);
80: PETSC_EXTERN PetscErrorCode DMDAGetElementType(DM,DMDAElementType*);
81: PETSC_EXTERN PetscErrorCode DMDAGetElements(DM,PetscInt *,PetscInt *,const PetscInt*[]);
82: PETSC_EXTERN PetscErrorCode DMDARestoreElements(DM,PetscInt *,PetscInt *,const PetscInt*[]);
84: typedef enum { DMDA_X,DMDA_Y,DMDA_Z } DMDADirection;
86: #define MATSEQUSFFT "sequsfft"
88: PETSC_EXTERN PetscErrorCode DMDACreate(MPI_Comm,DM*);
89: PETSC_EXTERN PetscErrorCode DMDASetDim(DM,PetscInt);
90: PETSC_EXTERN PetscErrorCode DMDASetSizes(DM,PetscInt,PetscInt,PetscInt);
91: PETSC_EXTERN PetscErrorCode DMDACreate1d(MPI_Comm,DMDABoundaryType,PetscInt,PetscInt,PetscInt,const PetscInt[],DM *);
92: PETSC_EXTERN PetscErrorCode DMDACreate2d(MPI_Comm,DMDABoundaryType,DMDABoundaryType,DMDAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,const PetscInt[],const PetscInt[],DM*);
93: PETSC_EXTERN PetscErrorCode DMDACreate3d(MPI_Comm,DMDABoundaryType,DMDABoundaryType,DMDABoundaryType,DMDAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,const PetscInt[],const PetscInt[],const PetscInt[],DM*);
95: PETSC_EXTERN PetscErrorCode DMDAGlobalToNaturalBegin(DM,Vec,InsertMode,Vec);
96: PETSC_EXTERN PetscErrorCode DMDAGlobalToNaturalEnd(DM,Vec,InsertMode,Vec);
97: PETSC_EXTERN PetscErrorCode DMDANaturalToGlobalBegin(DM,Vec,InsertMode,Vec);
98: PETSC_EXTERN PetscErrorCode DMDANaturalToGlobalEnd(DM,Vec,InsertMode,Vec);
99: PETSC_EXTERN PetscErrorCode DMDALocalToLocalBegin(DM,Vec,InsertMode,Vec);
100: PETSC_EXTERN PetscErrorCode DMDALocalToLocalEnd(DM,Vec,InsertMode,Vec);
101: PETSC_EXTERN PetscErrorCode DMDACreateNaturalVector(DM,Vec *);
103: PETSC_EXTERN PetscErrorCode DMDAGetCorners(DM,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
104: PETSC_EXTERN PetscErrorCode DMDAGetGhostCorners(DM,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
105: PETSC_EXTERN PetscErrorCode DMDAGetInfo(DM,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,DMDABoundaryType*,DMDABoundaryType*,DMDABoundaryType*,DMDAStencilType*);
106: PETSC_EXTERN PetscErrorCode DMDAGetProcessorSubset(DM,DMDADirection,PetscInt,MPI_Comm*);
107: PETSC_EXTERN PetscErrorCode DMDAGetProcessorSubsets(DM,DMDADirection,MPI_Comm*);
109: PETSC_EXTERN PetscErrorCode DMDAGlobalToNaturalAllCreate(DM,VecScatter*);
110: PETSC_EXTERN PetscErrorCode DMDANaturalAllToGlobalCreate(DM,VecScatter*);
112: PETSC_EXTERN PetscErrorCode DMDAGetGlobalIndices(DM,PetscInt*,PetscInt**);
114: PETSC_EXTERN PetscErrorCode DMDAGetScatter(DM,VecScatter*,VecScatter*,VecScatter*);
115: PETSC_EXTERN PetscErrorCode DMDAGetNeighbors(DM,const PetscMPIInt**);
117: PETSC_EXTERN PetscErrorCode DMDAGetAO(DM,AO*);
118: PETSC_EXTERN PetscErrorCode DMDASetCoordinates(DM,Vec);
119: PETSC_EXTERN PetscErrorCode DMDASetGhostedCoordinates(DM,Vec);
120: PETSC_EXTERN PetscErrorCode DMDAGetCoordinates(DM,Vec *);
121: PETSC_EXTERN PetscErrorCode DMDAGetGhostedCoordinates(DM,Vec *);
122: PETSC_EXTERN PetscErrorCode DMDAGetCoordinateDA(DM,DM *);
123: PETSC_EXTERN PetscErrorCode DMDASetUniformCoordinates(DM,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal);
124: PETSC_EXTERN PetscErrorCode DMDAGetBoundingBox(DM,PetscReal[],PetscReal[]);
125: PETSC_EXTERN PetscErrorCode DMDAGetLocalBoundingBox(DM,PetscReal[],PetscReal[]);
126: /* function to wrap coordinates around boundary */
127: PETSC_EXTERN PetscErrorCode DMDAMapCoordsToPeriodicDomain(DM,PetscScalar*,PetscScalar*);
129: PETSC_EXTERN PetscErrorCode DMDAGetReducedDA(DM,PetscInt,DM*);
131: PETSC_EXTERN PetscErrorCode DMDASetFieldName(DM,PetscInt,const char[]);
132: PETSC_EXTERN PetscErrorCode DMDAGetFieldName(DM,PetscInt,const char**);
134: PETSC_EXTERN PetscErrorCode DMDASetBoundaryType(DM,DMDABoundaryType,DMDABoundaryType,DMDABoundaryType);
135: PETSC_EXTERN PetscErrorCode DMDASetDof(DM, PetscInt);
136: PETSC_EXTERN PetscErrorCode DMDASetStencilWidth(DM, PetscInt);
137: PETSC_EXTERN PetscErrorCode DMDASetOwnershipRanges(DM,const PetscInt[],const PetscInt[],const PetscInt[]);
138: PETSC_EXTERN PetscErrorCode DMDAGetOwnershipRanges(DM,const PetscInt**,const PetscInt**,const PetscInt**);
139: PETSC_EXTERN PetscErrorCode DMDASetNumProcs(DM, PetscInt, PetscInt, PetscInt);
140: PETSC_EXTERN PetscErrorCode DMDASetStencilType(DM, DMDAStencilType);
142: PETSC_EXTERN PetscErrorCode DMDAVecGetArray(DM,Vec,void *);
143: PETSC_EXTERN PetscErrorCode DMDAVecRestoreArray(DM,Vec,void *);
145: PETSC_EXTERN PetscErrorCode DMDAVecGetArrayDOF(DM,Vec,void *);
146: PETSC_EXTERN PetscErrorCode DMDAVecRestoreArrayDOF(DM,Vec,void *);
148: PETSC_EXTERN PetscErrorCode DMDASplitComm2d(MPI_Comm,PetscInt,PetscInt,PetscInt,MPI_Comm*);
150: /*S
151: DMDALocalInfo - C struct that contains information about a structured grid and a processors logical
152: location in it.
154: Level: beginner
156: Concepts: distributed array
158: Developer note: Then entries in this struct are int instead of PetscInt so that the elements may
159: be extracted in Fortran as if from an integer array
161: .seealso: DMDACreate1d(), DMDACreate2d(), DMDACreate3d(), DMDestroy(), DM, DMDAGetLocalInfo(), DMDAGetInfo()
162: S*/
163: typedef struct {
164: PetscInt dim,dof,sw;
165: PetscInt mx,my,mz; /* global number of grid points in each direction */
166: PetscInt xs,ys,zs; /* starting point of this processor, excluding ghosts */
167: PetscInt xm,ym,zm; /* number of grid points on this processor, excluding ghosts */
168: PetscInt gxs,gys,gzs; /* starting point of this processor including ghosts */
169: PetscInt gxm,gym,gzm; /* number of grid points on this processor including ghosts */
170: DMDABoundaryType bx,by,bz; /* type of ghost nodes at boundary */
171: DMDAStencilType st;
172: DM da;
173: } DMDALocalInfo;
175: /*MC
176: DMDAForEachPointBegin2d - Starts a loop over the local part of a two dimensional DMDA
178: Synopsis:
179: void DMDAForEachPointBegin2d(DALocalInfo *info,PetscInt i,PetscInt j);
180:
181: Not Collective
183: Level: intermediate
185: .seealso: DMDAForEachPointEnd2d(), DMDAVecGetArray()
186: M*/
187: #define DMDAForEachPointBegin2d(info,i,j) {\
188: PetscInt _xints = info->xs,_xinte = info->xs+info->xm,_yints = info->ys,_yinte = info->ys+info->ym;\
189: for (j=_yints; j<_yinte; j++) {\
190: for (i=_xints; i<_xinte; i++) {\
192: /*MC
193: DMDAForEachPointEnd2d - Ends a loop over the local part of a two dimensional DMDA
195: Synopsis:
196: void DMDAForEachPointEnd2d;
197:
198: Not Collective
200: Level: intermediate
202: .seealso: DMDAForEachPointBegin2d(), DMDAVecGetArray()
203: M*/
204: #define DMDAForEachPointEnd2d }}}
206: /*MC
207: DMDACoor2d - Structure for holding 2d (x and y) coordinates.
209: Level: intermediate
211: Sample Usage:
212: DMDACoor2d **coors;
213: Vec vcoors;
214: DM cda;
216: DMDAGetCoordinates(da,&vcoors);
217: DMDAGetCoordinateDA(da,&cda);
218: DMDAVecGetArray(cda,vcoors,&coors);
219: DMDAGetCorners(cda,&mstart,&nstart,0,&m,&n,0)
220: for (i=mstart; i<mstart+m; i++) {
221: for (j=nstart; j<nstart+n; j++) {
222: x = coors[j][i].x;
223: y = coors[j][i].y;
224: ......
225: }
226: }
227: DMDAVecRestoreArray(dac,vcoors,&coors);
229: .seealso: DMDACoor3d, DMDAForEachPointBegin(), DMDAGetCoordinateDA(), DMDAGetCoordinates(), DMDAGetGhostCoordinates()
230: M*/
231: typedef struct {PetscScalar x,y;} DMDACoor2d;
233: /*MC
234: DMDACoor3d - Structure for holding 3d (x, y and z) coordinates.
236: Level: intermediate
238: Sample Usage:
239: DMDACoor3d ***coors;
240: Vec vcoors;
241: DM cda;
243: DMDAGetCoordinates(da,&vcoors);
244: DMDAGetCoordinateDA(da,&cda);
245: DMDAVecGetArray(cda,vcoors,&coors);
246: DMDAGetCorners(cda,&mstart,&nstart,&pstart,&m,&n,&p)
247: for (i=mstart; i<mstart+m; i++) {
248: for (j=nstart; j<nstart+n; j++) {
249: for (k=pstart; k<pstart+p; k++) {
250: x = coors[k][j][i].x;
251: y = coors[k][j][i].y;
252: z = coors[k][j][i].z;
253: ......
254: }
255: }
256: DMDAVecRestoreArray(dac,vcoors,&coors);
258: .seealso: DMDACoor2d, DMDAForEachPointBegin(), DMDAGetCoordinateDA(), DMDAGetCoordinates(), DMDAGetGhostCoordinates()
259: M*/
260: typedef struct {PetscScalar x,y,z;} DMDACoor3d;
261:
262: PETSC_EXTERN PetscErrorCode DMDAGetLocalInfo(DM,DMDALocalInfo*);
263: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*DMDALocalFunction1)(DMDALocalInfo*,void*,void*,void*);
264: PETSC_EXTERN PetscErrorCode DMDAComputeFunctionLocal(DM, DMDALocalFunction1, Vec, Vec, void *);
265: PETSC_EXTERN PetscErrorCode DMDAComputeFunctionLocalGhost(DM, DMDALocalFunction1, Vec, Vec, void *);
266: PETSC_EXTERN PetscErrorCode DMDAFormJacobianLocal(DM, DMDALocalFunction1, Vec, Mat, void *);
267: PETSC_EXTERN PetscErrorCode DMDAComputeFunction1(DM,Vec,Vec,void*);
268: PETSC_EXTERN PetscErrorCode DMDAComputeFunction(DM,PetscErrorCode (*)(void),Vec,Vec,void*);
269: PETSC_EXTERN PetscErrorCode DMDAComputeFunctioni1(DM,PetscInt,Vec,PetscScalar*,void*);
270: PETSC_EXTERN PetscErrorCode DMDAComputeFunctionib1(DM,PetscInt,Vec,PetscScalar*,void*);
271: PETSC_EXTERN PetscErrorCode DMDAComputeJacobian1WithAdic(DM,Vec,Mat,void*);
272: PETSC_EXTERN PetscErrorCode DMDAComputeJacobian1WithAdifor(DM,Vec,Mat,void*);
273: PETSC_EXTERN PetscErrorCode DMDAMultiplyByJacobian1WithAdic(DM,Vec,Vec,Vec,void*);
274: PETSC_EXTERN PetscErrorCode DMDAMultiplyByJacobian1WithAdifor(DM,Vec,Vec,Vec,void*);
275: PETSC_EXTERN PetscErrorCode DMDAMultiplyByJacobian1WithAD(DM,Vec,Vec,Vec,void*);
276: PETSC_EXTERN PetscErrorCode DMDAComputeJacobian1(DM,Vec,Mat,void*);
277: PETSC_EXTERN PetscErrorCode DMDAGetLocalFunction(DM,DMDALocalFunction1*);
278: PETSC_EXTERN PetscErrorCode DMDASetLocalFunction(DM,DMDALocalFunction1);
279: PETSC_EXTERN PetscErrorCode DMDASetLocalFunctioni(DM,PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,PetscScalar*,void*));
280: PETSC_EXTERN PetscErrorCode DMDASetLocalFunctionib(DM,PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,PetscScalar*,void*));
281: PETSC_EXTERN PetscErrorCode DMDAGetLocalJacobian(DM,DMDALocalFunction1*);
282: PETSC_EXTERN PetscErrorCode DMDASetLocalJacobian(DM,DMDALocalFunction1);
283: PETSC_EXTERN PetscErrorCode DMDASetLocalAdicFunction_Private(DM,DMDALocalFunction1);
285: PETSC_EXTERN PetscErrorCode MatSetDM(Mat,DM);
286: PETSC_EXTERN PetscErrorCode MatRegisterDAAD(void);
287: PETSC_EXTERN PetscErrorCode MatCreateDAAD(DM,Mat*);
288: PETSC_EXTERN PetscErrorCode MatCreateSeqUSFFT(Vec,DM,Mat*);
290: /*MC
291: DMDASetLocalAdicFunction - Caches in a DM a local function computed by ADIC/ADIFOR
293: Synopsis:
294: PetscErrorCode DMDASetLocalAdicFunction(DM da,DMDALocalFunction1 ad_lf)
295:
296: Logically Collective on DM
298: Input Parameter:
299: + da - initial distributed array
300: - ad_lf - the local function as computed by ADIC/ADIFOR
302: Level: intermediate
304: .keywords: distributed array, refine
306: .seealso: DMDACreate1d(), DMDACreate2d(), DMDACreate3d(), DMDestroy(), DMDAGetLocalFunction(), DMDASetLocalFunction(),
307: DMDASetLocalJacobian()
308: M*/
309: #if defined(PETSC_HAVE_ADIC)
310: # define DMDASetLocalAdicFunction(a,d) DMDASetLocalAdicFunction_Private(a,(DMDALocalFunction1)d)
311: #else
312: # define DMDASetLocalAdicFunction(a,d) DMDASetLocalAdicFunction_Private(a,0)
313: #endif
315: PETSC_EXTERN PetscErrorCode DMDASetLocalAdicMFFunction_Private(DM,DMDALocalFunction1);
316: #if defined(PETSC_HAVE_ADIC)
317: # define DMDASetLocalAdicMFFunction(a,d) DMDASetLocalAdicMFFunction_Private(a,(DMDALocalFunction1)d)
318: #else
319: # define DMDASetLocalAdicMFFunction(a,d) DMDASetLocalAdicMFFunction_Private(a,0)
320: #endif
321: PETSC_EXTERN PetscErrorCode DMDASetLocalAdicFunctioni_Private(DM,PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*));
322: #if defined(PETSC_HAVE_ADIC)
323: # define DMDASetLocalAdicFunctioni(a,d) DMDASetLocalAdicFunctioni_Private(a,(PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*))d)
324: #else
325: # define DMDASetLocalAdicFunctioni(a,d) DMDASetLocalAdicFunctioni_Private(a,0)
326: #endif
327: PETSC_EXTERN PetscErrorCode DMDASetLocalAdicMFFunctioni_Private(DM,PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*));
328: #if defined(PETSC_HAVE_ADIC)
329: # define DMDASetLocalAdicMFFunctioni(a,d) DMDASetLocalAdicMFFunctioni_Private(a,(PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*))d)
330: #else
331: # define DMDASetLocalAdicMFFunctioni(a,d) DMDASetLocalAdicMFFunctioni_Private(a,0)
332: #endif
334: PETSC_EXTERN PetscErrorCode DMDASetLocalAdicFunctionib_Private(DM,PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*));
335: #if defined(PETSC_HAVE_ADIC)
336: # define DMDASetLocalAdicFunctionib(a,d) DMDASetLocalAdicFunctionib_Private(a,(PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*))d)
337: #else
338: # define DMDASetLocalAdicFunctionib(a,d) DMDASetLocalAdicFunctionib_Private(a,0)
339: #endif
340: PETSC_EXTERN PetscErrorCode DMDASetLocalAdicMFFunctionib_Private(DM,PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*));
341: #if defined(PETSC_HAVE_ADIC)
342: # define DMDASetLocalAdicMFFunctionib(a,d) DMDASetLocalAdicMFFunctionib_Private(a,(PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*))d)
343: #else
344: # define DMDASetLocalAdicMFFunctionib(a,d) DMDASetLocalAdicMFFunctionib_Private(a,0)
345: #endif
347: PETSC_EXTERN PetscErrorCode DMDAComputeFunctioniTest1(DM,void*);
348: PETSC_EXTERN PetscErrorCode DMDASetGetMatrix(DM,PetscErrorCode (*)(DM, const MatType,Mat *));
349: PETSC_EXTERN PetscErrorCode DMDASetBlockFills(DM,PetscInt*,PetscInt*);
350: PETSC_EXTERN PetscErrorCode DMDASetRefinementFactor(DM,PetscInt,PetscInt,PetscInt);
351: PETSC_EXTERN PetscErrorCode DMDAGetRefinementFactor(DM,PetscInt*,PetscInt*,PetscInt*);
353: PETSC_EXTERN PetscErrorCode DMDAGetAdicArray(DM,PetscBool ,void*,void*,PetscInt*);
354: PETSC_EXTERN PetscErrorCode DMDARestoreAdicArray(DM,PetscBool ,void*,void*,PetscInt*);
355: PETSC_EXTERN PetscErrorCode DMDAGetAdicMFArray(DM,PetscBool ,void*,void*,PetscInt*);
356: PETSC_EXTERN PetscErrorCode DMDAGetAdicMFArray4(DM,PetscBool ,void*,void*,PetscInt*);
357: PETSC_EXTERN PetscErrorCode DMDAGetAdicMFArray9(DM,PetscBool ,void*,void*,PetscInt*);
358: PETSC_EXTERN PetscErrorCode DMDAGetAdicMFArrayb(DM,PetscBool ,void*,void*,PetscInt*);
359: PETSC_EXTERN PetscErrorCode DMDARestoreAdicMFArray(DM,PetscBool ,void*,void*,PetscInt*);
360: PETSC_EXTERN PetscErrorCode DMDAGetArray(DM,PetscBool ,void*);
361: PETSC_EXTERN PetscErrorCode DMDARestoreArray(DM,PetscBool ,void*);
362: PETSC_EXTERN PetscErrorCode ad_DAGetArray(DM,PetscBool ,void*);
363: PETSC_EXTERN PetscErrorCode ad_DARestoreArray(DM,PetscBool ,void*);
364: PETSC_EXTERN PetscErrorCode admf_DAGetArray(DM,PetscBool ,void*);
365: PETSC_EXTERN PetscErrorCode admf_DARestoreArray(DM,PetscBool ,void*);
367: PETSC_EXTERN PetscErrorCode DMDACreatePF(DM,PF*);
369: PETSC_EXTERN PetscErrorCode DMDACreateSection(DM, PetscInt[], PetscInt[], PetscInt[], PetscInt[]);
371: #define DMDA_FILE_CLASSID 1211220
372: #endif