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 src/dm/dmimpl.h

 10: typedef struct _DAOps *DAOps;
 11: struct _DAOps {
 12:   DMOPS(DA)
 13:   PetscErrorCode (*getelements)(DA,PetscInt*,const PetscInt*[]);
 14:   PetscErrorCode (*restoreelements)(DA,PetscInt*,const PetscInt*[]);
 15: };

 17: struct _p_DA {
 18:   PETSCHEADER(struct _DAOps);
 19:   PetscInt            M,N,P;                 /* array dimensions */
 20:   PetscInt            m,n,p;                 /* processor layout */
 21:   PetscInt            w;                     /* degrees of freedom per node */
 22:   PetscInt            s;                     /* stencil width */
 23:   PetscInt            xs,xe,ys,ye,zs,ze;     /* range of local values */
 24:   PetscInt            Xs,Xe,Ys,Ye,Zs,Ze;     /* range including ghost values
 25:                                                    values above already scaled by w */
 26:   PetscInt            *idx,Nl;               /* local to global map */
 27:   PetscInt            base;                  /* global number of 1st local node */
 28:   DAPeriodicType      wrap;                  /* indicates type of periodic boundaries */
 29:   VecScatter          gtol,ltog,ltol;        /* scatters, see below for details */
 30:   DAStencilType       stencil_type;          /* stencil, either box or star */
 31:   PetscInt            dim;                   /* DA dimension (1,2, or 3) */
 32:   DAInterpolationType interptype;

 34:   PetscInt            nlocal,Nlocal;         /* local size of local vector and global vector */

 36:   AO                  ao;                    /* application ordering context */

 38:   ISLocalToGlobalMapping ltogmap,ltogmapb;   /* local to global mapping for associated vectors */
 39:   Vec                    coordinates;        /* coordinates (x,y,z) of local nodes, not including ghosts*/
 40:   DA                     da_coordinates;     /* da for getting ghost values of coordinates */
 41:   Vec                    ghosted_coordinates;/* coordinates with ghost nodes */
 42:   char                   **fieldname;        /* names of individual components in vectors */

 44:   PetscInt               *lx,*ly,*lz;        /* number of nodes in each partition block along 3 axis */
 45:   Vec                    natural;            /* global vector for storing items in natural order */
 46:   VecScatter             gton;               /* vector scatter from global to natural */

 48:   ISColoring             localcoloring;       /* set by DAGetColoring() */
 49:   ISColoring             ghostedcoloring;

 51:   DAElementType          elementtype;
 52:   PetscInt               ne;                  /* number of elements */
 53:   PetscInt               *e;                  /* the elements */

 55: #define DA_MAX_WORK_VECTORS 10 /* work vectors available to users  via DAVecGetArray() */
 56:   Vec                    localin[DA_MAX_WORK_VECTORS],localout[DA_MAX_WORK_VECTORS];
 57:   Vec                    globalin[DA_MAX_WORK_VECTORS],globalout[DA_MAX_WORK_VECTORS];

 59:   PetscInt               refine_x,refine_y,refine_z; /* ratio used in refining */

 61: #define DA_MAX_AD_ARRAYS 2 /* work arrays for holding derivative type data, via DAGetAdicArray() */
 62:   void                   *adarrayin[DA_MAX_AD_ARRAYS],*adarrayout[DA_MAX_AD_ARRAYS];
 63:   void                   *adarrayghostedin[DA_MAX_AD_ARRAYS],*adarrayghostedout[DA_MAX_AD_ARRAYS];
 64:   void                   *adstartin[DA_MAX_AD_ARRAYS],*adstartout[DA_MAX_AD_ARRAYS];
 65:   void                   *adstartghostedin[DA_MAX_AD_ARRAYS],*adstartghostedout[DA_MAX_AD_ARRAYS];
 66:   PetscInt                    tdof,ghostedtdof;

 68:                             /* work arrays for holding derivative type data, via DAGetAdicMFArray() */
 69:   void                   *admfarrayin[DA_MAX_AD_ARRAYS],*admfarrayout[DA_MAX_AD_ARRAYS];
 70:   void                   *admfarrayghostedin[DA_MAX_AD_ARRAYS],*admfarrayghostedout[DA_MAX_AD_ARRAYS];
 71:   void                   *admfstartin[DA_MAX_AD_ARRAYS],*admfstartout[DA_MAX_AD_ARRAYS];
 72:   void                   *admfstartghostedin[DA_MAX_AD_ARRAYS],*admfstartghostedout[DA_MAX_AD_ARRAYS];

 74: #define DA_MAX_WORK_ARRAYS 2 /* work arrays for holding work via DAGetArray() */
 75:   void                   *arrayin[DA_MAX_WORK_ARRAYS],*arrayout[DA_MAX_WORK_ARRAYS];
 76:   void                   *arrayghostedin[DA_MAX_WORK_ARRAYS],*arrayghostedout[DA_MAX_WORK_ARRAYS];
 77:   void                   *startin[DA_MAX_WORK_ARRAYS],*startout[DA_MAX_WORK_ARRAYS];
 78:   void                   *startghostedin[DA_MAX_WORK_ARRAYS],*startghostedout[DA_MAX_WORK_ARRAYS];

 80:   DALocalFunction1       lf;
 81:   DALocalFunction1       lj;
 82:   DALocalFunction1       adic_lf;
 83:   DALocalFunction1       adicmf_lf;
 84:   DALocalFunction1       adifor_lf;
 85:   DALocalFunction1       adiformf_lf;

 87:   PetscErrorCode (*lfi)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*);
 88:   PetscErrorCode (*adic_lfi)(DALocalInfo*,MatStencil*,void*,void*,void*);
 89:   PetscErrorCode (*adicmf_lfi)(DALocalInfo*,MatStencil*,void*,void*,void*);
 90:   PetscErrorCode (*lfib)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*);
 91:   PetscErrorCode (*adic_lfib)(DALocalInfo*,MatStencil*,void*,void*,void*);
 92:   PetscErrorCode (*adicmf_lfib)(DALocalInfo*,MatStencil*,void*,void*,void*);

 94:   /* used by DASetBlockFills() */
 95:   PetscInt               *ofill,*dfill;

 97:   /* used by DASetMatPreallocateOnly() */
 98:   PetscTruth             prealloc_only;
 99: };

101: /*
102:   Vectors:
103:      Global has on each processor the interior degrees of freedom and
104:          no ghost points. This vector is what the solvers usually see.
105:      Local has on each processor the ghost points as well. This is 
106:           what code to calculate Jacobians, etc. usually sees.
107:   Vector scatters:
108:      gtol - Global representation to local
109:      ltog - Local representation to global (involves no communication)
110:      ltol - Local representation to local representation, updates the
111:             ghostpoint values in the second vector from (correct) interior
112:             values in the first vector.  This is good for explicit
113:             nearest neighbor timestepping.
114: */

117: EXTERN PetscErrorCode  VecView_MPI_DA(Vec,PetscViewer);
118: EXTERN PetscErrorCode  VecLoadIntoVector_Binary_DA(PetscViewer,Vec);
120: EXTERN PetscErrorCode DAView_Private(DA);


124: #endif