Actual source code: petscdmda.h

4: #include petscdm.h 5: #include petscpf.h 6: #include petscao.h 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; 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; 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; 81: typedef enum { DMDA_X,DMDA_Y,DMDA_Z } DMDADirection; 83: #define MATSEQUSFFT "sequsfft" 123: /* function to wrap coordinates around boundary */ 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: 260: typedef PetscErrorCode (*DMDALocalFunction1)(DMDALocalInfo*,void*,void*,void*); 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 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 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 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 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 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 366: #define DMDA_FILE_CLASSID 1211220 368: #endif