Actual source code: petscdmda.h
4: #include petscdm.h
5: #include petscpf.h
6: #include petscao.h
7: PETSC_EXTERN_CXX_BEGIN
9: /*E
10: DMDAStencilType - Determines if the stencil extends only along the coordinate directions, or also
11: to the northeast, northwest etc
13: Level: beginner
15: .seealso: DMDACreate1d(), DMDACreate2d(), DMDACreate3d(), DMDACreate()
16: E*/
17: typedef enum { DMDA_STENCIL_STAR,DMDA_STENCIL_BOX } DMDAStencilType;
19: /*MC
20: DMDA_STENCIL_STAR - "Star"-type stencil. In logical grid coordinates, only (i,j,k), (i+s,j,k), (i,j+s,k),
21: (i,j,k+s) are in the stencil NOT, for example, (i+s,j+s,k)
23: Level: beginner
25: .seealso: DMDA_STENCIL_BOX, DMDAStencilType
26: M*/
28: /*MC
29: DMDA_STENCIL_BOX - "Box"-type stencil. In logical grid coordinates, any of (i,j,k), (i+s,j+r,k+t) may
30: be in the stencil.
32: Level: beginner
34: .seealso: DMDA_STENCIL_STAR, DMDAStencilType
35: M*/
37: /*E
38: DMDABoundaryType - Describes the choice for fill of ghost cells on domain boundaries.
40: Level: beginner
42: A boundary may be of type DMDA_BOUNDARY_NONE (no ghost nodes), DMDA_BOUNDARY_GHOST (ghost nodes
43: exist but aren't filled), DMDA_BOUNDARY_MIRROR (not yet implemented), or DMDA_BOUNDARY_PERIODIC
44: (ghost nodes filled by the opposite edge of the domain).
46: .seealso: DMDASetBoundaryType(), DMDACreate1d(), DMDACreate2d(), DMDACreate3d(), DMDACreate()
47: E*/
48: typedef enum { DMDA_BOUNDARY_NONE, DMDA_BOUNDARY_GHOSTED, DMDA_BOUNDARY_MIRROR, DMDA_BOUNDARY_PERIODIC } DMDABoundaryType;
50: extern const char *DMDABoundaryTypes[];
52: /*E
53: DMDAInterpolationType - Defines the type of interpolation that will be returned by
54: DMGetInterpolation.
56: Level: beginner
58: .seealso: DMDACreate1d(), DMDACreate2d(), DMDACreate3d(), DMGetInterpolation(), DMDASetInterpolationType(), DMDACreate()
59: E*/
60: typedef enum { DMDA_Q0, DMDA_Q1 } DMDAInterpolationType;
62: extern PetscErrorCode DMDASetInterpolationType(DM,DMDAInterpolationType);
63: extern PetscErrorCode DMDAGetInterpolationType(DM,DMDAInterpolationType*);
65: /*E
66: DMDAElementType - Defines the type of elements that will be returned by
67: DMDAGetElements()
69: Level: beginner
71: .seealso: DMDACreate1d(), DMDACreate2d(), DMDACreate3d(), DMGetInterpolation(), DMDASetInterpolationType(),
72: DMDASetElementType(), DMDAGetElements(), DMDARestoreElements(), DMDACreate()
73: E*/
74: typedef enum { DMDA_ELEMENT_P1, DMDA_ELEMENT_Q1 } DMDAElementType;
76: extern PetscErrorCode DMDASetElementType(DM,DMDAElementType);
77: extern PetscErrorCode DMDAGetElementType(DM,DMDAElementType*);
78: extern PetscErrorCode DMDAGetElements(DM,PetscInt *,PetscInt *,const PetscInt*[]);
79: extern PetscErrorCode DMDARestoreElements(DM,PetscInt *,PetscInt *,const PetscInt*[]);
81: typedef enum { DMDA_X,DMDA_Y,DMDA_Z } DMDADirection;
83: #define MATSEQUSFFT "sequsfft"
85: extern PetscErrorCode DMDACreate(MPI_Comm,DM*);
86: extern PetscErrorCode DMDASetDim(DM,PetscInt);
87: extern PetscErrorCode DMDASetSizes(DM,PetscInt,PetscInt,PetscInt);
88: extern PetscErrorCode DMDACreate1d(MPI_Comm,DMDABoundaryType,PetscInt,PetscInt,PetscInt,const PetscInt[],DM *);
89: extern PetscErrorCode DMDACreate2d(MPI_Comm,DMDABoundaryType,DMDABoundaryType,DMDAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,const PetscInt[],const PetscInt[],DM*);
90: extern PetscErrorCode DMDACreate3d(MPI_Comm,DMDABoundaryType,DMDABoundaryType,DMDABoundaryType,DMDAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,const PetscInt[],const PetscInt[],const PetscInt[],DM*);
92: extern PetscErrorCode DMDAGlobalToNaturalBegin(DM,Vec,InsertMode,Vec);
93: extern PetscErrorCode DMDAGlobalToNaturalEnd(DM,Vec,InsertMode,Vec);
94: extern PetscErrorCode DMDANaturalToGlobalBegin(DM,Vec,InsertMode,Vec);
95: extern PetscErrorCode DMDANaturalToGlobalEnd(DM,Vec,InsertMode,Vec);
96: extern PetscErrorCode DMDALocalToLocalBegin(DM,Vec,InsertMode,Vec);
97: extern PetscErrorCode DMDALocalToLocalEnd(DM,Vec,InsertMode,Vec);
98: extern PetscErrorCode DMDACreateNaturalVector(DM,Vec *);
100: extern PetscErrorCode DMDAGetCorners(DM,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
101: extern PetscErrorCode DMDAGetGhostCorners(DM,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
102: extern PetscErrorCode DMDAGetInfo(DM,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,DMDABoundaryType*,DMDABoundaryType*,DMDABoundaryType*,DMDAStencilType*);
103: extern PetscErrorCode DMDAGetProcessorSubset(DM,DMDADirection,PetscInt,MPI_Comm*);
104: extern PetscErrorCode DMDAGetProcessorSubsets(DM,DMDADirection,MPI_Comm*);
106: extern PetscErrorCode DMDAGlobalToNaturalAllCreate(DM,VecScatter*);
107: extern PetscErrorCode DMDANaturalAllToGlobalCreate(DM,VecScatter*);
109: extern PetscErrorCode DMDAGetGlobalIndices(DM,PetscInt*,PetscInt**);
111: extern PetscErrorCode DMDAGetScatter(DM,VecScatter*,VecScatter*,VecScatter*);
112: extern PetscErrorCode DMDAGetNeighbors(DM,const PetscMPIInt**);
114: extern PetscErrorCode DMDAGetAO(DM,AO*);
115: extern PetscErrorCode DMDASetCoordinates(DM,Vec);
116: extern PetscErrorCode DMDASetGhostedCoordinates(DM,Vec);
117: extern PetscErrorCode DMDAGetCoordinates(DM,Vec *);
118: extern PetscErrorCode DMDAGetGhostedCoordinates(DM,Vec *);
119: extern PetscErrorCode DMDAGetCoordinateDA(DM,DM *);
120: extern PetscErrorCode DMDASetUniformCoordinates(DM,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal);
121: extern PetscErrorCode DMDAGetBoundingBox(DM,PetscReal[],PetscReal[]);
122: extern PetscErrorCode DMDAGetLocalBoundingBox(DM,PetscReal[],PetscReal[]);
123: /* function to wrap coordinates around boundary */
124: extern PetscErrorCode DMDAMapCoordsToPeriodicDomain(DM,PetscScalar*,PetscScalar*);
126: extern PetscErrorCode DMDAGetReducedDA(DM,PetscInt,DM*);
128: extern PetscErrorCode DMDASetFieldName(DM,PetscInt,const char[]);
129: extern PetscErrorCode DMDAGetFieldName(DM,PetscInt,const char**);
131: extern PetscErrorCode DMDASetBoundaryType(DM,DMDABoundaryType,DMDABoundaryType,DMDABoundaryType);
132: extern PetscErrorCode DMDASetDof(DM, PetscInt);
133: extern PetscErrorCode DMDASetStencilWidth(DM, PetscInt);
134: extern PetscErrorCode DMDASetOwnershipRanges(DM,const PetscInt[],const PetscInt[],const PetscInt[]);
135: extern PetscErrorCode DMDAGetOwnershipRanges(DM,const PetscInt**,const PetscInt**,const PetscInt**);
136: extern PetscErrorCode DMDASetNumProcs(DM, PetscInt, PetscInt, PetscInt);
137: extern PetscErrorCode DMDASetStencilType(DM, DMDAStencilType);
139: extern PetscErrorCode DMDAVecGetArray(DM,Vec,void *);
140: extern PetscErrorCode DMDAVecRestoreArray(DM,Vec,void *);
142: extern PetscErrorCode DMDAVecGetArrayDOF(DM,Vec,void *);
143: extern PetscErrorCode DMDAVecRestoreArrayDOF(DM,Vec,void *);
145: extern PetscErrorCode DMDASplitComm2d(MPI_Comm,PetscInt,PetscInt,PetscInt,MPI_Comm*);
147: /*S
148: DMDALocalInfo - C struct that contains information about a structured grid and a processors logical
149: location in it.
151: Level: beginner
153: Concepts: distributed array
155: Developer note: Then entries in this struct are int instead of PetscInt so that the elements may
156: be extracted in Fortran as if from an integer array
158: .seealso: DMDACreate1d(), DMDACreate2d(), DMDACreate3d(), DMDestroy(), DM, DMDAGetLocalInfo(), DMDAGetInfo()
159: S*/
160: typedef struct {
161: PetscInt dim,dof,sw;
162: PetscInt mx,my,mz; /* global number of grid points in each direction */
163: PetscInt xs,ys,zs; /* starting point of this processor, excluding ghosts */
164: PetscInt xm,ym,zm; /* number of grid points on this processor, excluding ghosts */
165: PetscInt gxs,gys,gzs; /* starting point of this processor including ghosts */
166: PetscInt gxm,gym,gzm; /* number of grid points on this processor including ghosts */
167: DMDABoundaryType bx,by,bz; /* type of ghost nodes at boundary */
168: DMDAStencilType st;
169: DM da;
170: } DMDALocalInfo;
172: /*MC
173: DMDAForEachPointBegin2d - Starts a loop over the local part of a two dimensional DMDA
175: Synopsis:
176: void DMDAForEachPointBegin2d(DALocalInfo *info,PetscInt i,PetscInt j);
177:
178: Not Collective
180: Level: intermediate
182: .seealso: DMDAForEachPointEnd2d(), DMDAVecGetArray()
183: M*/
184: #define DMDAForEachPointBegin2d(info,i,j) {\
185: PetscInt _xints = info->xs,_xinte = info->xs+info->xm,_yints = info->ys,_yinte = info->ys+info->ym;\
186: for (j=_yints; j<_yinte; j++) {\
187: for (i=_xints; i<_xinte; i++) {\
189: /*MC
190: DMDAForEachPointEnd2d - Ends a loop over the local part of a two dimensional DMDA
192: Synopsis:
193: void DMDAForEachPointEnd2d;
194:
195: Not Collective
197: Level: intermediate
199: .seealso: DMDAForEachPointBegin2d(), DMDAVecGetArray()
200: M*/
201: #define DMDAForEachPointEnd2d }}}
203: /*MC
204: DMDACoor2d - Structure for holding 2d (x and y) coordinates.
206: Level: intermediate
208: Sample Usage:
209: DMDACoor2d **coors;
210: Vec vcoors;
211: DM cda;
213: DMDAGetCoordinates(da,&vcoors);
214: DMDAGetCoordinateDA(da,&cda);
215: DMDAVecGetArray(cda,vcoors,&coors);
216: DMDAGetCorners(cda,&mstart,&nstart,0,&m,&n,0)
217: for (i=mstart; i<mstart+m; i++) {
218: for (j=nstart; j<nstart+n; j++) {
219: x = coors[j][i].x;
220: y = coors[j][i].y;
221: ......
222: }
223: }
224: DMDAVecRestoreArray(dac,vcoors,&coors);
226: .seealso: DMDACoor3d, DMDAForEachPointBegin(), DMDAGetCoordinateDA(), DMDAGetCoordinates(), DMDAGetGhostCoordinates()
227: M*/
228: typedef struct {PetscScalar x,y;} DMDACoor2d;
230: /*MC
231: DMDACoor3d - Structure for holding 3d (x, y and z) coordinates.
233: Level: intermediate
235: Sample Usage:
236: DMDACoor3d ***coors;
237: Vec vcoors;
238: DM cda;
240: DMDAGetCoordinates(da,&vcoors);
241: DMDAGetCoordinateDA(da,&cda);
242: DMDAVecGetArray(cda,vcoors,&coors);
243: DMDAGetCorners(cda,&mstart,&nstart,&pstart,&m,&n,&p)
244: for (i=mstart; i<mstart+m; i++) {
245: for (j=nstart; j<nstart+n; j++) {
246: for (k=pstart; k<pstart+p; k++) {
247: x = coors[k][j][i].x;
248: y = coors[k][j][i].y;
249: z = coors[k][j][i].z;
250: ......
251: }
252: }
253: DMDAVecRestoreArray(dac,vcoors,&coors);
255: .seealso: DMDACoor2d, DMDAForEachPointBegin(), DMDAGetCoordinateDA(), DMDAGetCoordinates(), DMDAGetGhostCoordinates()
256: M*/
257: typedef struct {PetscScalar x,y,z;} DMDACoor3d;
258:
259: extern PetscErrorCode DMDAGetLocalInfo(DM,DMDALocalInfo*);
260: typedef PetscErrorCode (*DMDALocalFunction1)(DMDALocalInfo*,void*,void*,void*);
261: extern PetscErrorCode DMDAFormFunctionLocal(DM, DMDALocalFunction1, Vec, Vec, void *);
262: extern PetscErrorCode DMDAFormFunctionLocalGhost(DM, DMDALocalFunction1, Vec, Vec, void *);
263: extern PetscErrorCode DMDAFormJacobianLocal(DM, DMDALocalFunction1, Vec, Mat, void *);
264: extern PetscErrorCode DMDAFormFunction1(DM,Vec,Vec,void*);
265: extern PetscErrorCode DMDAFormFunction(DM,PetscErrorCode (*)(void),Vec,Vec,void*);
266: extern PetscErrorCode DMDAFormFunctioni1(DM,PetscInt,Vec,PetscScalar*,void*);
267: extern PetscErrorCode DMDAFormFunctionib1(DM,PetscInt,Vec,PetscScalar*,void*);
268: extern PetscErrorCode DMDAComputeJacobian1WithAdic(DM,Vec,Mat,void*);
269: extern PetscErrorCode DMDAComputeJacobian1WithAdifor(DM,Vec,Mat,void*);
270: extern PetscErrorCode DMDAMultiplyByJacobian1WithAdic(DM,Vec,Vec,Vec,void*);
271: extern PetscErrorCode DMDAMultiplyByJacobian1WithAdifor(DM,Vec,Vec,Vec,void*);
272: extern PetscErrorCode DMDAMultiplyByJacobian1WithAD(DM,Vec,Vec,Vec,void*);
273: extern PetscErrorCode DMDAComputeJacobian1(DM,Vec,Mat,void*);
274: extern PetscErrorCode DMDAGetLocalFunction(DM,DMDALocalFunction1*);
275: extern PetscErrorCode DMDASetLocalFunction(DM,DMDALocalFunction1);
276: extern PetscErrorCode DMDASetLocalFunctioni(DM,PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,PetscScalar*,void*));
277: extern PetscErrorCode DMDASetLocalFunctionib(DM,PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,PetscScalar*,void*));
278: extern PetscErrorCode DMDAGetLocalJacobian(DM,DMDALocalFunction1*);
279: extern PetscErrorCode DMDASetLocalJacobian(DM,DMDALocalFunction1);
280: extern PetscErrorCode DMDASetLocalAdicFunction_Private(DM,DMDALocalFunction1);
282: extern PetscErrorCode MatSetDM(Mat,DM);
283: extern PetscErrorCode MatRegisterDAAD(void);
284: extern PetscErrorCode MatCreateDAAD(DM,Mat*);
285: extern PetscErrorCode MatCreateSeqUSFFT(Vec,DM,Mat*);
287: /*MC
288: DMDASetLocalAdicFunction - Caches in a DM a local function computed by ADIC/ADIFOR
290: Synopsis:
291: PetscErrorCode DMDASetLocalAdicFunction(DM da,DMDALocalFunction1 ad_lf)
292:
293: Logically Collective on DM
295: Input Parameter:
296: + da - initial distributed array
297: - ad_lf - the local function as computed by ADIC/ADIFOR
299: Level: intermediate
301: .keywords: distributed array, refine
303: .seealso: DMDACreate1d(), DMDACreate2d(), DMDACreate3d(), DMDestroy(), DMDAGetLocalFunction(), DMDASetLocalFunction(),
304: DMDASetLocalJacobian()
305: M*/
306: #if defined(PETSC_HAVE_ADIC)
307: # define DMDASetLocalAdicFunction(a,d) DMDASetLocalAdicFunction_Private(a,(DMDALocalFunction1)d)
308: #else
309: # define DMDASetLocalAdicFunction(a,d) DMDASetLocalAdicFunction_Private(a,0)
310: #endif
312: extern PetscErrorCode DMDASetLocalAdicMFFunction_Private(DM,DMDALocalFunction1);
313: #if defined(PETSC_HAVE_ADIC)
314: # define DMDASetLocalAdicMFFunction(a,d) DMDASetLocalAdicMFFunction_Private(a,(DMDALocalFunction1)d)
315: #else
316: # define DMDASetLocalAdicMFFunction(a,d) DMDASetLocalAdicMFFunction_Private(a,0)
317: #endif
318: extern PetscErrorCode DMDASetLocalAdicFunctioni_Private(DM,PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*));
319: #if defined(PETSC_HAVE_ADIC)
320: # define DMDASetLocalAdicFunctioni(a,d) DMDASetLocalAdicFunctioni_Private(a,(PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*))d)
321: #else
322: # define DMDASetLocalAdicFunctioni(a,d) DMDASetLocalAdicFunctioni_Private(a,0)
323: #endif
324: extern PetscErrorCode DMDASetLocalAdicMFFunctioni_Private(DM,PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*));
325: #if defined(PETSC_HAVE_ADIC)
326: # define DMDASetLocalAdicMFFunctioni(a,d) DMDASetLocalAdicMFFunctioni_Private(a,(PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*))d)
327: #else
328: # define DMDASetLocalAdicMFFunctioni(a,d) DMDASetLocalAdicMFFunctioni_Private(a,0)
329: #endif
331: extern PetscErrorCode DMDASetLocalAdicFunctionib_Private(DM,PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*));
332: #if defined(PETSC_HAVE_ADIC)
333: # define DMDASetLocalAdicFunctionib(a,d) DMDASetLocalAdicFunctionib_Private(a,(PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*))d)
334: #else
335: # define DMDASetLocalAdicFunctionib(a,d) DMDASetLocalAdicFunctionib_Private(a,0)
336: #endif
337: extern PetscErrorCode DMDASetLocalAdicMFFunctionib_Private(DM,PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*));
338: #if defined(PETSC_HAVE_ADIC)
339: # define DMDASetLocalAdicMFFunctionib(a,d) DMDASetLocalAdicMFFunctionib_Private(a,(PetscErrorCode (*)(DMDALocalInfo*,MatStencil*,void*,void*,void*))d)
340: #else
341: # define DMDASetLocalAdicMFFunctionib(a,d) DMDASetLocalAdicMFFunctionib_Private(a,0)
342: #endif
344: extern PetscErrorCode DMDAFormFunctioniTest1(DM,void*);
345: extern PetscErrorCode DMDASetGetMatrix(DM,PetscErrorCode (*)(DM, const MatType,Mat *));
346: extern PetscErrorCode DMDASetBlockFills(DM,PetscInt*,PetscInt*);
347: extern PetscErrorCode DMDASetRefinementFactor(DM,PetscInt,PetscInt,PetscInt);
348: extern PetscErrorCode DMDAGetRefinementFactor(DM,PetscInt*,PetscInt*,PetscInt*);
350: extern PetscErrorCode DMDAGetAdicArray(DM,PetscBool ,void*,void*,PetscInt*);
351: extern PetscErrorCode DMDARestoreAdicArray(DM,PetscBool ,void*,void*,PetscInt*);
352: extern PetscErrorCode DMDAGetAdicMFArray(DM,PetscBool ,void*,void*,PetscInt*);
353: extern PetscErrorCode DMDAGetAdicMFArray4(DM,PetscBool ,void*,void*,PetscInt*);
354: extern PetscErrorCode DMDAGetAdicMFArray9(DM,PetscBool ,void*,void*,PetscInt*);
355: extern PetscErrorCode DMDAGetAdicMFArrayb(DM,PetscBool ,void*,void*,PetscInt*);
356: extern PetscErrorCode DMDARestoreAdicMFArray(DM,PetscBool ,void*,void*,PetscInt*);
357: extern PetscErrorCode DMDAGetArray(DM,PetscBool ,void*);
358: extern PetscErrorCode DMDARestoreArray(DM,PetscBool ,void*);
359: extern PetscErrorCode ad_DAGetArray(DM,PetscBool ,void*);
360: extern PetscErrorCode ad_DARestoreArray(DM,PetscBool ,void*);
361: extern PetscErrorCode admf_DAGetArray(DM,PetscBool ,void*);
362: extern PetscErrorCode admf_DARestoreArray(DM,PetscBool ,void*);
364: extern PetscErrorCode DMDACreatePF(DM,PF*);
366: #define DMDA_FILE_CLASSID 1211220
367: PETSC_EXTERN_CXX_END
368: #endif