Actual source code: dmdaimpl.h
1: /*
2: Distributed arrays - communication tools for parallel, rectangular grids.
3: */
5: #pragma once
7: #include <petscdmda.h>
8: #include <petsc/private/dmimpl.h>
10: typedef struct {
11: PetscInt M, N, P; /* array dimensions */
12: PetscInt m, n, p; /* processor layout */
13: PetscInt w; /* degrees of freedom per node */
14: PetscInt s; /* stencil width */
15: PetscInt xs, xe, ys, ye, zs, ze; /* range of local values */
16: PetscInt Xs, Xe, Ys, Ye, Zs, Ze; /* range including ghost values
17: values above already scaled by w */
18: PetscInt base; /* global number of 1st local node, includes the * w term */
19: DMBoundaryType bx, by, bz; /* indicates type of ghost nodes at boundary */
20: VecScatter gtol, ltol; /* scatters, see below for details */
21: DMDAStencilType stencil_type; /* stencil, either box or star */
22: DMDAInterpolationType interptype;
24: PetscInt nlocal, Nlocal; /* local size of local vector and global vector, includes the * w term */
26: PetscInt xol, yol, zol; /* overlap of local subdomains */
27: PetscInt xo, yo, zo; /* offsets for the indices in x y and z */
28: PetscInt Mo, No, Po; /* the size of the problem the offset is in to */
29: PetscInt Nsub; /* number of local subdomains to decompose into */
30: PetscInt nonxs, nonys, nonzs; /* the nonoverlapping starts in the case of a subdomain da */
31: PetscInt nonxm, nonym, nonzm; /* the nonoverlapping sizes in the case of a subdomain da */
33: AO ao; /* application ordering context */
34: AOType aotype; /* type of application ordering */
36: char **fieldname; /* names of individual components in vectors */
37: char **coordinatename; /* names of coordinate directions, for example, x, y, z */
39: PetscInt *lx, *ly, *lz; /* number of nodes in each partition block along 3 axis */
40: Vec natural; /* global vector for storing items in natural order */
41: VecScatter gton; /* vector scatter from global to natural */
42: PetscMPIInt *neighbors; /* ranks of all neighbors and self */
44: ISColoring localcoloring; /* set by DMCreateColoring() */
45: ISColoring ghostedcoloring;
47: DMDAElementType elementtype;
48: PetscInt ne; /* number of elements */
49: PetscInt nen; /* number of nodes per element */
50: PetscInt *e; /* the elements */
51: IS ecorners; /* corners of the subdomain */
53: PetscInt refine_x, refine_y, refine_z; /* ratio used in refining */
54: PetscInt coarsen_x, coarsen_y, coarsen_z; /* ratio used for coarsening */
55: /* if the refinement is done differently on different levels */
56: PetscInt refine_x_hier_n, *refine_x_hier, refine_y_hier_n, *refine_y_hier, refine_z_hier_n, *refine_z_hier;
58: #define DMDA_MAX_WORK_ARRAYS 2 /* work arrays for holding work via DMDAGetArray() */
59: void *arrayin[DMDA_MAX_WORK_ARRAYS], *arrayout[DMDA_MAX_WORK_ARRAYS];
60: void *arrayghostedin[DMDA_MAX_WORK_ARRAYS], *arrayghostedout[DMDA_MAX_WORK_ARRAYS];
61: void *startin[DMDA_MAX_WORK_ARRAYS], *startout[DMDA_MAX_WORK_ARRAYS];
62: void *startghostedin[DMDA_MAX_WORK_ARRAYS], *startghostedout[DMDA_MAX_WORK_ARRAYS];
64: PetscErrorCode (*lf)(DM, Vec, Vec, void *);
65: PetscErrorCode (*lj)(DM, Vec, Vec, void *);
67: /* used by DMDASetBlockFills() */
68: PetscInt *ofill, *dfill;
69: PetscInt *ofillcols;
71: /* used by DMDASetMatPreallocateOnly() */
72: PetscBool prealloc_only;
73: PetscInt preallocCenterDim; /* Dimension of the points which connect adjacent points for preallocation */
74: } DM_DA;
76: /*
77: Vectors:
78: Global has on each processor the interior degrees of freedom and
79: no ghost points. This vector is what the solvers usually see.
80: Local has on each processor the ghost points as well. This is
81: what code to calculate Jacobians, etc. usually sees.
82: Vector scatters:
83: gtol - Global representation to local
84: ltol - Local representation to local representation, updates the
85: ghostpoint values in the second vector from (correct) interior
86: values in the first vector. This is good for explicit
87: nearest neighbor timestepping.
88: */
90: PETSC_INTERN PetscErrorCode VecView_MPI_DA(Vec, PetscViewer);
91: PETSC_INTERN PetscErrorCode VecLoad_Default_DA(Vec, PetscViewer);
92: PETSC_INTERN PetscErrorCode DMView_DA_Matlab(DM, PetscViewer);
93: PETSC_INTERN PetscErrorCode DMView_DA_Binary(DM, PetscViewer);
94: PETSC_INTERN PetscErrorCode DMView_DA_VTK(DM, PetscViewer);
95: PETSC_INTERN PetscErrorCode DMView_DA_GLVis(DM, PetscViewer);
96: PETSC_INTERN PetscErrorCode DMDASelectFields(DM, PetscInt *, PetscInt **);
98: PETSC_INTERN PetscErrorCode DMCreateGlobalVector_DA(DM, Vec *);
99: PETSC_INTERN PetscErrorCode DMCreateLocalVector_DA(DM, Vec *);
100: PETSC_INTERN PetscErrorCode DMGlobalToLocalBegin_DA(DM, Vec, InsertMode, Vec);
101: PETSC_INTERN PetscErrorCode DMGlobalToLocalEnd_DA(DM, Vec, InsertMode, Vec);
102: PETSC_INTERN PetscErrorCode DMLocalToGlobalBegin_DA(DM, Vec, InsertMode, Vec);
103: PETSC_INTERN PetscErrorCode DMLocalToGlobalEnd_DA(DM, Vec, InsertMode, Vec);
104: PETSC_INTERN PetscErrorCode DMLocalToLocalBegin_DA(DM, Vec, InsertMode, Vec);
105: PETSC_INTERN PetscErrorCode DMLocalToLocalEnd_DA(DM, Vec, InsertMode, Vec);
106: PETSC_INTERN PetscErrorCode DMCreateInterpolation_DA(DM, DM, Mat *, Vec *);
107: PETSC_INTERN PetscErrorCode DMCreateColoring_DA(DM, ISColoringType, ISColoring *);
108: PETSC_INTERN PetscErrorCode DMCreateMatrix_DA(DM, Mat *);
109: PETSC_INTERN PetscErrorCode DMCreateCoordinateDM_DA(DM, DM *);
110: PETSC_INTERN PetscErrorCode DMRefine_DA(DM, MPI_Comm, DM *);
111: PETSC_INTERN PetscErrorCode DMCoarsen_DA(DM, MPI_Comm, DM *);
112: PETSC_INTERN PetscErrorCode DMRefineHierarchy_DA(DM, PetscInt, DM[]);
113: PETSC_INTERN PetscErrorCode DMCoarsenHierarchy_DA(DM, PetscInt, DM[]);
114: PETSC_INTERN PetscErrorCode DMCreateInjection_DA(DM, DM, Mat *);
115: PETSC_INTERN PetscErrorCode DMView_DA(DM, PetscViewer);
116: PETSC_INTERN PetscErrorCode DMSetUp_DA(DM);
117: PETSC_INTERN PetscErrorCode DMDestroy_DA(DM);
118: PETSC_INTERN PetscErrorCode DMCreateDomainDecomposition_DA(DM, PetscInt *, char ***, IS **, IS **, DM **);
119: PETSC_INTERN PetscErrorCode DMCreateDomainDecompositionScatters_DA(DM, PetscInt, DM *, VecScatter **, VecScatter **, VecScatter **);
120: PETSC_INTERN PetscErrorCode DMGetCompatibility_DA(DM, DM, PetscBool *, PetscBool *);
121: PETSC_INTERN PetscErrorCode DMLocatePoints_DA_Regular(DM, Vec, DMPointLocationType, PetscSF);
122: PETSC_INTERN PetscErrorCode DMGetLocalBoundingBox_DA(DM, PetscReal[], PetscReal[], PetscInt[], PetscInt[]);
123: PETSC_INTERN PetscErrorCode DMSetUpGLVisViewer_DMDA(PetscObject, PetscViewer);
124: PETSC_INTERN PetscErrorCode DMLocalToLocalCreate_DA(DM);
125: PETSC_INTERN PetscErrorCode DMDAGetNatural_Private(DM, PetscInt *, IS *);
126: PETSC_INTERN PetscErrorCode DMSetUp_DA_1D(DM);
127: PETSC_INTERN PetscErrorCode DMSetUp_DA_2D(DM);
128: PETSC_INTERN PetscErrorCode DMSetUp_DA_3D(DM);