Actual source code: dmdaimpl.h
1: /*
2: Distributed arrays - communication tools for parallel, rectangular grids.
3: */
5: #ifndef _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: DMDAInterpolationType interptype;
25: PetscInt nlocal, Nlocal; /* local size of local vector and global vector, includes the * w term */
27: PetscInt xol, yol, zol; /* overlap of local subdomains */
28: PetscInt xo, yo, zo; /* offsets for the indices in x y and z */
29: PetscInt Mo, No, Po; /* the size of the problem the offset is in to */
30: PetscInt Nsub; /* number of local subdomains to decompose into */
31: PetscInt nonxs, nonys, nonzs; /* the nonoverlapping starts in the case of a subdomain da */
32: PetscInt nonxm, nonym, nonzm; /* the nonoverlapping sizes in the case of a subdomain da */
34: AO ao; /* application ordering context */
35: AOType aotype; /* type of application ordering */
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 nen; /* number of nodes per element */
51: PetscInt *e; /* the elements */
52: IS ecorners; /* corners of the subdomain */
54: PetscInt refine_x, refine_y, refine_z; /* ratio used in refining */
55: PetscInt coarsen_x, coarsen_y, coarsen_z; /* ratio used for coarsening */
56: /* if the refinement is done differently on different levels */
57: PetscInt refine_x_hier_n, *refine_x_hier, refine_y_hier_n, *refine_y_hier, refine_z_hier_n, *refine_z_hier;
59: #define DMDA_MAX_WORK_ARRAYS 2 /* work arrays for holding work via DMDAGetArray() */
60: void *arrayin[DMDA_MAX_WORK_ARRAYS], *arrayout[DMDA_MAX_WORK_ARRAYS];
61: void *arrayghostedin[DMDA_MAX_WORK_ARRAYS], *arrayghostedout[DMDA_MAX_WORK_ARRAYS];
62: void *startin[DMDA_MAX_WORK_ARRAYS], *startout[DMDA_MAX_WORK_ARRAYS];
63: void *startghostedin[DMDA_MAX_WORK_ARRAYS], *startghostedout[DMDA_MAX_WORK_ARRAYS];
65: PetscErrorCode (*lf)(DM, Vec, Vec, void *);
66: PetscErrorCode (*lj)(DM, Vec, Vec, void *);
68: /* used by DMDASetBlockFills() */
69: PetscInt *ofill, *dfill;
70: PetscInt *ofillcols;
72: /* used by DMDASetMatPreallocateOnly() */
73: PetscBool prealloc_only;
74: PetscInt preallocCenterDim; /* Dimension of the points which connect adjacent points for preallocation */
75: } DM_DA;
77: /*
78: Vectors:
79: Global has on each processor the interior degrees of freedom and
80: no ghost points. This vector is what the solvers usually see.
81: Local has on each processor the ghost points as well. This is
82: what code to calculate Jacobians, etc. usually sees.
83: Vector scatters:
84: gtol - Global representation to local
85: ltol - Local representation to local representation, updates the
86: ghostpoint values in the second vector from (correct) interior
87: values in the first vector. This is good for explicit
88: nearest neighbor timestepping.
89: */
91: PETSC_INTERN PetscErrorCode VecView_MPI_DA(Vec, PetscViewer);
92: PETSC_INTERN PetscErrorCode VecLoad_Default_DA(Vec, PetscViewer);
93: PETSC_INTERN PetscErrorCode DMView_DA_Matlab(DM, PetscViewer);
94: PETSC_INTERN PetscErrorCode DMView_DA_Binary(DM, PetscViewer);
95: PETSC_INTERN PetscErrorCode DMView_DA_VTK(DM, PetscViewer);
96: PETSC_INTERN PetscErrorCode DMView_DA_GLVis(DM, PetscViewer);
97: PETSC_EXTERN PetscErrorCode DMDAVTKWriteAll(PetscObject, PetscViewer);
98: PETSC_EXTERN PetscErrorCode DMDASelectFields(DM, PetscInt *, PetscInt **);
100: #endif