Actual source code: daimpl.h
1: /*
2: Distributed arrays - communication tools for parallel, rectangular grids.
3: */
5: #if !defined(_DAIMPL_H)
6: #define _DAIMPL_H
8: #include petscda.h
10: /*
11: The DM interface is shared by DA, VecPack, and any other object that may
12: be used with the DMMG class. If you change this MAKE SURE you change
13: struct _DAOps and struct _VecPackOps!
14: */
15: typedef struct _DMOps *DMOps;
16: struct _DMOps {
17: PetscErrorCode (*view)(DM,PetscViewer);
18: PetscErrorCode (*createglobalvector)(DM,Vec*);
19: PetscErrorCode (*getcoloring)(DM,ISColoringType,ISColoring*);
20: PetscErrorCode (*getmatrix)(DM, MatType,Mat*);
21: PetscErrorCode (*getinterpolation)(DM,DM,Mat*,Vec*);
22: PetscErrorCode (*getinjection)(DM,DM,VecScatter*);
23: PetscErrorCode (*refine)(DM,MPI_Comm,DM*);
24: PetscErrorCode (*coarsen)(DM,MPI_Comm,DM*);
25: PetscErrorCode (*refinehierarchy)(DM,PetscInt,DM**);
26: PetscErrorCode (*coarsenhierarchy)(DM,PetscInt,DM**);
27: PetscErrorCode (*forminitialguess)(DM,PetscErrorCode (*)(void),Vec,void*);
28: PetscErrorCode (*formfunction)(DM,Vec*);
29: };
31: struct _p_DM {
32: PETSCHEADER(struct _DMOps);
33: };
35: typedef struct _DAOps *DAOps;
36: struct _DAOps {
37: PetscErrorCode (*view)(DA,PetscViewer);
38: PetscErrorCode (*createglobalvector)(DA,Vec*);
39: PetscErrorCode (*getcoloring)(DA,ISColoringType,ISColoring*);
40: PetscErrorCode (*getmatrix)(DA, MatType,Mat*);
41: PetscErrorCode (*getinterpolation)(DA,DA,Mat*,Vec*);
42: PetscErrorCode (*getinjection)(DA,DA,VecScatter*);
43: PetscErrorCode (*refine)(DA,MPI_Comm,DA*);
44: PetscErrorCode (*getelements)(DA,PetscInt*,const PetscInt*[]);
45: PetscErrorCode (*restoreelements)(DA,PetscInt*,const PetscInt*[]);
46: };
48: struct _p_DA {
49: PETSCHEADER(struct _DAOps);
50: PetscInt M,N,P; /* array dimensions */
51: PetscInt m,n,p; /* processor layout */
52: PetscInt w; /* degrees of freedom per node */
53: PetscInt s; /* stencil width */
54: PetscInt xs,xe,ys,ye,zs,ze; /* range of local values */
55: PetscInt Xs,Xe,Ys,Ye,Zs,Ze; /* range including ghost values
56: values above already scaled by w */
57: PetscInt *idx,Nl; /* local to global map */
58: PetscInt base; /* global number of 1st local node */
59: DAPeriodicType wrap; /* indicates type of periodic boundaries */
60: VecScatter gtol,ltog,ltol; /* scatters, see below for details */
61: DAStencilType stencil_type; /* stencil, either box or star */
62: PetscInt dim; /* DA dimension (1,2, or 3) */
63: DAInterpolationType interptype;
65: PetscInt nlocal,Nlocal; /* local size of local vector and global vector */
67: AO ao; /* application ordering context */
69: ISLocalToGlobalMapping ltogmap,ltogmapb; /* local to global mapping for associated vectors */
70: Vec coordinates; /* coordinates (x,y,z) of local nodes, not including ghosts*/
71: DA da_coordinates; /* da for getting ghost values of coordinates */
72: Vec ghosted_coordinates;/* coordinates with ghost nodes */
73: char **fieldname; /* names of individual components in vectors */
75: PetscInt *lx,*ly,*lz; /* number of nodes in each partition block along 3 axis */
76: Vec natural; /* global vector for storing items in natural order */
77: VecScatter gton; /* vector scatter from global to natural */
79: ISColoring localcoloring; /* set by DAGetColoring() */
80: ISColoring ghostedcoloring;
82: DAElementType elementtype;
83: PetscInt ne; /* number of elements */
84: PetscInt *e; /* the elements */
86: #define DA_MAX_WORK_VECTORS 10 /* work vectors available to users via DAVecGetArray() */
87: Vec localin[DA_MAX_WORK_VECTORS],localout[DA_MAX_WORK_VECTORS];
88: Vec globalin[DA_MAX_WORK_VECTORS],globalout[DA_MAX_WORK_VECTORS];
90: PetscInt refine_x,refine_y,refine_z; /* ratio used in refining */
92: #define DA_MAX_AD_ARRAYS 2 /* work arrays for holding derivative type data, via DAGetAdicArray() */
93: void *adarrayin[DA_MAX_AD_ARRAYS],*adarrayout[DA_MAX_AD_ARRAYS];
94: void *adarrayghostedin[DA_MAX_AD_ARRAYS],*adarrayghostedout[DA_MAX_AD_ARRAYS];
95: void *adstartin[DA_MAX_AD_ARRAYS],*adstartout[DA_MAX_AD_ARRAYS];
96: void *adstartghostedin[DA_MAX_AD_ARRAYS],*adstartghostedout[DA_MAX_AD_ARRAYS];
97: PetscInt tdof,ghostedtdof;
99: /* work arrays for holding derivative type data, via DAGetAdicMFArray() */
100: void *admfarrayin[DA_MAX_AD_ARRAYS],*admfarrayout[DA_MAX_AD_ARRAYS];
101: void *admfarrayghostedin[DA_MAX_AD_ARRAYS],*admfarrayghostedout[DA_MAX_AD_ARRAYS];
102: void *admfstartin[DA_MAX_AD_ARRAYS],*admfstartout[DA_MAX_AD_ARRAYS];
103: void *admfstartghostedin[DA_MAX_AD_ARRAYS],*admfstartghostedout[DA_MAX_AD_ARRAYS];
105: #define DA_MAX_WORK_ARRAYS 2 /* work arrays for holding work via DAGetArray() */
106: void *arrayin[DA_MAX_WORK_ARRAYS],*arrayout[DA_MAX_WORK_ARRAYS];
107: void *arrayghostedin[DA_MAX_WORK_ARRAYS],*arrayghostedout[DA_MAX_WORK_ARRAYS];
108: void *startin[DA_MAX_WORK_ARRAYS],*startout[DA_MAX_WORK_ARRAYS];
109: void *startghostedin[DA_MAX_WORK_ARRAYS],*startghostedout[DA_MAX_WORK_ARRAYS];
111: DALocalFunction1 lf;
112: DALocalFunction1 lj;
113: DALocalFunction1 adic_lf;
114: DALocalFunction1 adicmf_lf;
115: DALocalFunction1 adifor_lf;
116: DALocalFunction1 adiformf_lf;
118: PetscErrorCode (*lfi)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*);
119: PetscErrorCode (*adic_lfi)(DALocalInfo*,MatStencil*,void*,void*,void*);
120: PetscErrorCode (*adicmf_lfi)(DALocalInfo*,MatStencil*,void*,void*,void*);
121: PetscErrorCode (*lfib)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*);
122: PetscErrorCode (*adic_lfib)(DALocalInfo*,MatStencil*,void*,void*,void*);
123: PetscErrorCode (*adicmf_lfib)(DALocalInfo*,MatStencil*,void*,void*,void*);
125: /* used by DASetBlockFills() */
126: PetscInt *ofill,*dfill;
128: /* used by DASetMatPreallocateOnly() */
129: PetscTruth prealloc_only;
130: };
132: /*
133: Vectors:
134: Global has on each processor the interior degrees of freedom and
135: no ghost points. This vector is what the solvers usually see.
136: Local has on each processor the ghost points as well. This is
137: what code to calculate Jacobians, etc. usually sees.
138: Vector scatters:
139: gtol - Global representation to local
140: ltog - Local representation to global (involves no communication)
141: ltol - Local representation to local representation, updates the
142: ghostpoint values in the second vector from (correct) interior
143: values in the first vector. This is good for explicit
144: nearest neighbor timestepping.
145: */
148: EXTERN PetscErrorCode VecView_MPI_DA(Vec,PetscViewer);
149: EXTERN PetscErrorCode VecLoadIntoVector_Binary_DA(PetscViewer,Vec);
151: EXTERN PetscErrorCode DAView_Private(DA);
155: #endif