1: /*
2: Distributed arrays - communication tools for parallel, rectangular grids.
3: */
5: #if !defined(_DAIMPL_H)
6: #define _DAIMPL_H 8: #include <petscdmda.h>
9: #include <petsc-private/dmimpl.h>
11: typedef struct {
12: PetscInt M,N,P; /* array dimensions */
13: PetscInt m,n,p; /* processor layout */
14: PetscInt w; /* degrees of freedom per node */
15: PetscInt s; /* stencil width */
16: PetscInt xs,xe,ys,ye,zs,ze; /* range of local values */
17: PetscInt Xs,Xe,Ys,Ye,Zs,Ze; /* range including ghost values
18: values above already scaled by w */
19: PetscInt base; /* global number of 1st local node, includes the * w term */
20: DMBoundaryType bx,by,bz; /* indicates type of ghost nodes at boundary */
21: VecScatter gtol,ltol; /* scatters, see below for details */
22: DMDAStencilType stencil_type; /* stencil, either box or star */
23: PetscInt dim; /* DMDA dimension (1,2, or 3) */
24: DMDAInterpolationType interptype;
26: PetscInt nlocal,Nlocal; /* local size of local vector and global vector, includes the * w term */
28: PetscInt xol,yol,zol; /* overlap of local subdomains */
29: PetscInt xo,yo,zo; /* offsets for the indices in x y and z */
30: PetscInt Mo,No,Po; /* the size of the problem the offset is in to */
31: PetscInt Nsub; /* number of local subdomains to decompose into */
32: PetscInt nonxs,nonys,nonzs; /* the nonoverlapping starts in the case of a subdomain da */
33: PetscInt nonxm,nonym,nonzm; /* the nonoverlapping sizes in the case of a subdomain da */
35: AO ao; /* application ordering context */
37: char **fieldname; /* names of individual components in vectors */
38: char **coordinatename; /* names of coordinate directions, for example, x, y, z */
40: PetscInt *lx,*ly,*lz; /* number of nodes in each partition block along 3 axis */
41: Vec natural; /* global vector for storing items in natural order */
42: VecScatter gton; /* vector scatter from global to natural */
43: PetscMPIInt *neighbors; /* ranks of all neighbors and self */
45: ISColoring localcoloring; /* set by DMCreateColoring() */
46: ISColoring ghostedcoloring;
48: DMDAElementType elementtype;
49: PetscInt ne; /* number of elements */
50: PetscInt *e; /* the elements */
52: PetscInt refine_x,refine_y,refine_z; /* ratio used in refining */
53: PetscInt coarsen_x,coarsen_y,coarsen_z; /* ratio used for coarsening */
55: #define DMDA_MAX_WORK_ARRAYS 2 /* work arrays for holding work via DMDAGetArray() */ 56: void *arrayin[DMDA_MAX_WORK_ARRAYS],*arrayout[DMDA_MAX_WORK_ARRAYS];
57: void *arrayghostedin[DMDA_MAX_WORK_ARRAYS],*arrayghostedout[DMDA_MAX_WORK_ARRAYS];
58: void *startin[DMDA_MAX_WORK_ARRAYS],*startout[DMDA_MAX_WORK_ARRAYS];
59: void *startghostedin[DMDA_MAX_WORK_ARRAYS],*startghostedout[DMDA_MAX_WORK_ARRAYS];
61: PetscErrorCode (*lf)(DM, Vec, Vec, void *);
62: PetscErrorCode (*lj)(DM, Vec, Vec, void *);
64: /* used by DMDASetBlockFills() */
65: PetscInt *ofill,*dfill;
66: PetscInt *ofillcols;
68: /* used by DMDASetMatPreallocateOnly() */
69: PetscBool prealloc_only;
70: PetscInt preallocCenterDim; /* Dimension of the points which connect adjacent points for preallocation */
71: } DM_DA;
73: /*
74: Vectors:
75: Global has on each processor the interior degrees of freedom and
76: no ghost points. This vector is what the solvers usually see.
77: Local has on each processor the ghost points as well. This is
78: what code to calculate Jacobians, etc. usually sees.
79: Vector scatters:
80: gtol - Global representation to local
81: ltog - Local representation to global (involves no communication)
82: ltol - Local representation to local representation, updates the
83: ghostpoint values in the second vector from (correct) interior
84: values in the first vector. This is good for explicit
85: nearest neighbor timestepping.
86: */
88: PETSC_EXTERN PetscErrorCode VecView_MPI_DA(Vec,PetscViewer);
89: PETSC_EXTERN PetscErrorCode VecLoad_Default_DA(Vec, PetscViewer);
90: PETSC_EXTERN PetscErrorCode DMView_DA_Matlab(DM,PetscViewer);
91: PETSC_EXTERN PetscErrorCode DMView_DA_Binary(DM,PetscViewer);
92: PETSC_EXTERN PetscErrorCode DMView_DA_VTK(DM,PetscViewer);
93: PETSC_EXTERN PetscErrorCodeDMDAVTKWriteAll(PetscObject,PetscViewer);
94: PETSC_EXTERN PetscErrorCode DMDASelectFields(DM,PetscInt*,PetscInt**);
96: PETSC_EXTERN PetscLogEvent DMDA_LocalADFunction;
98: #endif