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 *idx,Nl; /* local to global map */
20: PetscInt base; /* global number of 1st local node */
21: DMDABoundaryType bx,by,bz; /* indicates type of ghost nodes at boundary */
22: VecScatter gtol,ltog,ltol; /* scatters, see below for details */
23: DMDAStencilType stencil_type; /* stencil, either box or star */
24: PetscInt dim; /* DMDA dimension (1,2, or 3) */
25: DMDAInterpolationType interptype;
27: PetscInt nlocal,Nlocal; /* local size of local vector and global vector */
29: AO ao; /* application ordering context */
31: Vec coordinates; /* coordinates (x,y,z) of local nodes, not including ghosts*/
32: DM da_coordinates; /* da for getting ghost values of coordinates */
33: Vec ghosted_coordinates;/* coordinates with ghost nodes */
34: char **fieldname; /* names of individual components in vectors */
36: PetscInt *lx,*ly,*lz; /* number of nodes in each partition block along 3 axis */
37: Vec natural; /* global vector for storing items in natural order */
38: VecScatter gton; /* vector scatter from global to natural */
39: PetscMPIInt *neighbors; /* ranks of all neighbors and self */
41: ISColoring localcoloring; /* set by DMCreateColoring() */
42: ISColoring ghostedcoloring;
44: DMDAElementType elementtype;
45: PetscInt ne; /* number of elements */
46: PetscInt *e; /* the elements */
48: PetscInt refine_x,refine_y,refine_z; /* ratio used in refining */
49: PetscInt coarsen_x,coarsen_y,coarsen_z; /* ratio used for coarsening */
51: #define DMDA_MAX_AD_ARRAYS 2 /* work arrays for holding derivative type data, via DMDAGetAdicArray() */ 52: void *adarrayin[DMDA_MAX_AD_ARRAYS],*adarrayout[DMDA_MAX_AD_ARRAYS];
53: void *adarrayghostedin[DMDA_MAX_AD_ARRAYS],*adarrayghostedout[DMDA_MAX_AD_ARRAYS];
54: void *adstartin[DMDA_MAX_AD_ARRAYS],*adstartout[DMDA_MAX_AD_ARRAYS];
55: void *adstartghostedin[DMDA_MAX_AD_ARRAYS],*adstartghostedout[DMDA_MAX_AD_ARRAYS];
56: PetscInt tdof,ghostedtdof;
58: /* work arrays for holding derivative type data, via DMDAGetAdicMFArray() */
59: void *admfarrayin[DMDA_MAX_AD_ARRAYS],*admfarrayout[DMDA_MAX_AD_ARRAYS];
60: void *admfarrayghostedin[DMDA_MAX_AD_ARRAYS],*admfarrayghostedout[DMDA_MAX_AD_ARRAYS];
61: void *admfstartin[DMDA_MAX_AD_ARRAYS],*admfstartout[DMDA_MAX_AD_ARRAYS];
62: void *admfstartghostedin[DMDA_MAX_AD_ARRAYS],*admfstartghostedout[DMDA_MAX_AD_ARRAYS];
64: #define DMDA_MAX_WORK_ARRAYS 2 /* work arrays for holding work via DMDAGetArray() */ 65: void *arrayin[DMDA_MAX_WORK_ARRAYS],*arrayout[DMDA_MAX_WORK_ARRAYS];
66: void *arrayghostedin[DMDA_MAX_WORK_ARRAYS],*arrayghostedout[DMDA_MAX_WORK_ARRAYS];
67: void *startin[DMDA_MAX_WORK_ARRAYS],*startout[DMDA_MAX_WORK_ARRAYS];
68: void *startghostedin[DMDA_MAX_WORK_ARRAYS],*startghostedout[DMDA_MAX_WORK_ARRAYS];
70: DMDALocalFunction1 lf;
71: DMDALocalFunction1 lj;
72: DMDALocalFunction1 adic_lf;
73: DMDALocalFunction1 adicmf_lf;
74: DMDALocalFunction1 adifor_lf;
75: DMDALocalFunction1 adiformf_lf;
77: PetscErrorCode (*lfi)(DMDALocalInfo*,MatStencil*,void*,PetscScalar*,void*);
78: PetscErrorCode (*adic_lfi)(DMDALocalInfo*,MatStencil*,void*,void*,void*);
79: PetscErrorCode (*adicmf_lfi)(DMDALocalInfo*,MatStencil*,void*,void*,void*);
80: PetscErrorCode (*lfib)(DMDALocalInfo*,MatStencil*,void*,PetscScalar*,void*);
81: PetscErrorCode (*adic_lfib)(DMDALocalInfo*,MatStencil*,void*,void*,void*);
82: PetscErrorCode (*adicmf_lfib)(DMDALocalInfo*,MatStencil*,void*,void*,void*);
84: /* used by DMDASetBlockFills() */
85: PetscInt *ofill,*dfill;
87: /* used by DMDASetMatPreallocateOnly() */
88: PetscBool prealloc_only;
89: } DM_DA;
91: /*
92: Vectors:
93: Global has on each processor the interior degrees of freedom and
94: no ghost points. This vector is what the solvers usually see.
95: Local has on each processor the ghost points as well. This is
96: what code to calculate Jacobians, etc. usually sees.
97: Vector scatters:
98: gtol - Global representation to local
99: ltog - Local representation to global (involves no communication)
100: ltol - Local representation to local representation, updates the
101: ghostpoint values in the second vector from (correct) interior
102: values in the first vector. This is good for explicit
103: nearest neighbor timestepping.
104: */
106: EXTERN_C_BEGIN107: PETSC_EXTERN PetscErrorCode VecView_MPI_DA(Vec,PetscViewer);
108: PETSC_EXTERN PetscErrorCode VecLoad_Default_DA(Vec, PetscViewer);
109: EXTERN_C_END110: PETSC_EXTERN PetscErrorCode DMView_DA_Private(DM);
111: PETSC_EXTERN PetscErrorCode DMView_DA_Matlab(DM,PetscViewer);
112: PETSC_EXTERN PetscErrorCode DMView_DA_Binary(DM,PetscViewer);
113: PETSC_EXTERN PetscErrorCode DMView_DA_VTK(DM,PetscViewer);
114: PETSC_EXTERN PetscErrorCodeDMDAVTKWriteAll(PetscObject,PetscViewer);
116: PETSC_EXTERN PetscLogEvent DMDA_LocalADFunction;
118: #endif