Actual source code: petscda.h

  1: /* $Id: petscda.h,v 1.66 2001/04/10 19:37:22 bsmith Exp $ */

  3: /*
  4:       Regular array object, for easy parallelism of simple grid 
  5:    problems on regular distributed arrays.
  6: */
  9: #include "petscvec.h"
 10: #include "petscao.h"

 12: #define DA_COOKIE PETSC_COOKIE+14

 14: /*S
 15:      DA - Abstract PETSc object that manages distributed field data for a single structured grid

 17:    Level: beginner

 19:   Concepts: distributed array

 21: .seealso:  DACreate1d(), DACreate2d(), DACreate3d(), DADestroy(), VecScatter
 22: S*/
 23: typedef struct _p_DA* DA;

 25: /*E
 26:     DAStencilType - Determines if the stencil extends only along the coordinate directions, or also
 27:       to the northest, northwest etc

 29:    Level: beginner

 31: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DA
 32: E*/
 33: typedef enum { DA_STENCIL_STAR,DA_STENCIL_BOX } DAStencilType;

 35: /*E
 36:     DAPeriodicType - Is the domain periodic in one or more directions

 38:    Level: beginner

 40: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DA
 41: E*/
 42: typedef enum { DA_NONPERIODIC,DA_XPERIODIC,DA_YPERIODIC,DA_XYPERIODIC,
 43:                DA_XYZPERIODIC,DA_XZPERIODIC,DA_YZPERIODIC,DA_ZPERIODIC}
 44:                DAPeriodicType;

 46: #define DAXPeriodic(pt) ((pt)==DA_XPERIODIC||(pt)==DA_XYPERIODIC||(pt)==DA_XZPERIODIC||(pt)==DA_XYZPERIODIC)
 47: #define DAYPeriodic(pt) ((pt)==DA_YPERIODIC||(pt)==DA_XYPERIODIC||(pt)==DA_YZPERIODIC||(pt)==DA_XYZPERIODIC)
 48: #define DAZPeriodic(pt) ((pt)==DA_ZPERIODIC||(pt)==DA_XZPERIODIC||(pt)==DA_YZPERIODIC||(pt)==DA_XYZPERIODIC)

 50: typedef enum { DA_X,DA_Y,DA_Z } DADirection;

 52: EXTERN int   DACreate1d(MPI_Comm,DAPeriodicType,int,int,int,int*,DA *);
 53: EXTERN int   DACreate2d(MPI_Comm,DAPeriodicType,DAStencilType,int,int,int,int,int,int,int*,int*,DA *);
 54: EXTERN int   DACreate3d(MPI_Comm,DAPeriodicType,DAStencilType,
 55:                         int,int,int,int,int,int,int,int,int *,int *,int *,DA *);
 56: EXTERN int   DADestroy(DA);
 57: EXTERN int   DAView(DA,PetscViewer);

 59: EXTERN int   DAPrintHelp(DA);

 61: EXTERN int   DAGlobalToLocalBegin(DA,Vec,InsertMode,Vec);
 62: EXTERN int   DAGlobalToLocalEnd(DA,Vec,InsertMode,Vec);
 63: EXTERN int   DAGlobalToNaturalBegin(DA,Vec,InsertMode,Vec);
 64: EXTERN int   DAGlobalToNaturalEnd(DA,Vec,InsertMode,Vec);
 65: EXTERN int   DANaturalToGlobalBegin(DA,Vec,InsertMode,Vec);
 66: EXTERN int   DANaturalToGlobalEnd(DA,Vec,InsertMode,Vec);
 67: EXTERN int   DALocalToLocalBegin(DA,Vec,InsertMode,Vec);
 68: EXTERN int   DALocalToLocalEnd(DA,Vec,InsertMode,Vec);
 69: EXTERN int   DALocalToGlobal(DA,Vec,InsertMode,Vec);
 70: EXTERN int   DAGetOwnershipRange(DA,int **,int **,int **);
 71: EXTERN int   DACreateGlobalVector(DA,Vec *);
 72: EXTERN int   DACreateNaturalVector(DA,Vec *);
 73: EXTERN int   DACreateLocalVector(DA,Vec *);
 74: EXTERN int   DAGetLocalVector(DA,Vec *);
 75: EXTERN int   DARestoreLocalVector(DA,Vec *);
 76: EXTERN int   DAGetGlobalVector(DA,Vec *);
 77: EXTERN int   DARestoreGlobalVector(DA,Vec *);
 78: EXTERN int   DALoad(PetscViewer,int,int,int,DA *);
 79: EXTERN int   DAGetCorners(DA,int*,int*,int*,int*,int*,int*);
 80: EXTERN int   DAGetGhostCorners(DA,int*,int*,int*,int*,int*,int*);
 81: EXTERN int   DAGetInfo(DA,int*,int*,int*,int*,int*,int*,int*,int*,int*,DAPeriodicType*,DAStencilType*);
 82: EXTERN int   DAGetProcessorSubset(DA,DADirection,int,MPI_Comm*);
 83: EXTERN int   DARefine(DA,MPI_Comm,DA*);

 85: EXTERN int   DAGlobalToNaturalAllCreate(DA,VecScatter*);
 86: EXTERN int   DANaturalAllToGlobalCreate(DA,VecScatter*);

 88: EXTERN int   DAGetGlobalIndices(DA,int*,int**);
 89: EXTERN int   DAGetISLocalToGlobalMapping(DA,ISLocalToGlobalMapping*);
 90: EXTERN int   DAGetISLocalToGlobalMappingBlck(DA,ISLocalToGlobalMapping*);

 92: EXTERN int   DAGetScatter(DA,VecScatter*,VecScatter*,VecScatter*);

 94: EXTERN int   DAGetAO(DA,AO*);
 95: EXTERN int   DASetCoordinates(DA,Vec);
 96: EXTERN int   DAGetCoordinates(DA,Vec *);
 97: EXTERN int   DASetUniformCoordinates(DA,double,double,double,double,double,double);
 98: EXTERN int   DASetFieldName(DA,int,const char[]);
 99: EXTERN int   DAGetFieldName(DA,int,char **);

101: EXTERN int   DAVecGetArray(DA,Vec,void **);
102: EXTERN int   DAVecRestoreArray(DA,Vec,void **);

104: EXTERN int   DASplitComm2d(MPI_Comm,int,int,int,MPI_Comm*);

106: /*S
107:      DALocalInfo - C struct that contains information about a structured grid and a processors logical
108:               location in it.

110:    Level: beginner

112:   Concepts: distributed array

114: .seealso:  DACreate1d(), DACreate2d(), DACreate3d(), DADestroy(), DA, DAGetLocalInfo(), DAGetInfo()
115: S*/
116: typedef struct {
117:   int            dim,dof,sw;
118:   DAPeriodicType pt;
119:   DAStencilType  st;
120:   int            mx,my,mz;    /* global number of grid points in each direction */
121:   int            xs,ys,zs;    /* starting point of this processor */
122:   int            xm,ym,zm;    /* number of grid points on this processor */
123: } DALocalInfo;

125: EXTERN int DAGetLocalInfo(DA,DALocalInfo*);

127: #include "petscmat.h"
128: EXTERN int   DAGetColoring(DA,ISColoringType,MatType,ISColoring *,Mat *);
129: EXTERN int   DAGetInterpolation(DA,DA,Mat*,Vec*);

131: #include "petscpf.h"
132: EXTERN int DACreatePF(DA,PF*);

134: /*S
135:      VecPack - Abstract PETSc object that manages treating several distinct vectors as if they
136:         were one.   The VecPack routines allow one to manage a nonlinear solver that works on a
137:         vector that consists of several distinct parts. This is mostly used for LNKS solvers, 
138:         that is design optimization problems that are written as a nonlinear system

140:    Level: beginner

142:   Concepts: multi-component, LNKS solvers

144: .seealso:  VecPackCreate(), VecPackDestroy()
145: S*/
146: typedef struct _p_VecPack *VecPack;

148: EXTERN int VecPackCreate(MPI_Comm,VecPack*);
149: EXTERN int VecPackDestroy(VecPack);
150: EXTERN int VecPackAddArray(VecPack,int);
151: EXTERN int VecPackAddDA(VecPack,DA);
152: EXTERN int VecPackAddVecScatter(VecPack,VecScatter);
153: EXTERN int VecPackScatter(VecPack,Vec,...);
154: EXTERN int VecPackGather(VecPack,Vec,...);
155: EXTERN int VecPackGetAccess(VecPack,Vec,...);
156: EXTERN int VecPackRestoreAccess(VecPack,Vec,...);
157: EXTERN int VecPackGetLocalVectors(VecPack,...);
158: EXTERN int VecPackGetEntries(VecPack,...);
159: EXTERN int VecPackRestoreLocalVectors(VecPack,...);
160: EXTERN int VecPackCreateGlobalVector(VecPack,Vec*);
161: EXTERN int VecPackGetGlobalIndices(VecPack,...);
162: EXTERN int VecPackRefine(VecPack,MPI_Comm,VecPack*);
163: EXTERN int VecPackGetInterpolation(VecPack,VecPack,Mat*,Vec*);

165: #include "petscsnes.h"

167: /*S
168:      DM - Abstract PETSc object that manages an abstract grid object
169:           
170:    Level: intermediate

172:   Concepts: grids, grid refinement

174:    Notes: The DA object and the VecPack object are examples of DMs

176: .seealso:  VecPackCreate(), DA, VecPack
177: S*/
178: typedef struct _p_DM* DM;

180: EXTERN int DMView(DM,PetscViewer);
181: EXTERN int DMDestroy(DM);
182: EXTERN int DMCreateGlobalVector(DM,Vec*);
183: EXTERN int DMGetColoring(DM,ISColoringType,MatType,ISColoring*,Mat*);
184: EXTERN int DMGetInterpolation(DM,DM,Mat*,Vec*);
185: EXTERN int DMRefine(DM,MPI_Comm,DM*);
186: EXTERN int DMGetInterpolationScale(DM,DM,Mat,Vec*);

188: /*S
189:      DMMG -  Data structure to easily manage multi-level non-linear solvers on grids managed by DM
190:           
191:    Level: intermediate

193:   Concepts: multigrid, Newton-multigrid

195: .seealso:  VecPackCreate(), DA, VecPack, DM, DMMGCreate()
196: S*/
197: typedef struct _p_DMMG *DMMG;
198: struct _p_DMMG {
199:   DM         dm;                   /* grid information for this level */
200:   Vec        x,b,r;                /* global vectors used in multigrid preconditioner for this level*/
201:   Mat        J;                    /* matrix on this level */
202:   Mat        R;                    /* restriction to next coarser level (not defined on level 0) */
203:   int        nlevels;              /* number of levels above this one (total number of levels on level 0)*/
204:   MPI_Comm   comm;
205:   int        (*solve)(DMMG*,int);
206:   void       *user;

208:   /* SLES only */
209:   SLES       sles;
210:   int        (*rhs)(DMMG,Vec);

212:   /* SNES only */
213:   PetscTruth    matrixfree;
214:   Mat           B;
215:   Vec           Rscale;                /* scaling to restriction before computing Jacobian */
216:   int           (*computejacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
217:   int           (*computefunction)(SNES,Vec,Vec,void*);
218:   int           (*computefunctionlocal)(Scalar**,Scalar**,DALocalInfo*,void*);
219:   int           (*ad_computefunctionlocal)(Scalar**,Scalar**,DALocalInfo*,void*);
220:   ISColoring    iscoloring;            /* used with AD or FD coloring for Jacobian */
221:   MatFDColoring fdcoloring;            /* only used with FD coloring for Jacobian */
222:   SNES          snes;
223:   int           (*initialguess)(SNES,Vec,void*);
224:   Vec           work1,work2;
225: };

227: EXTERN int DMMGCreate(MPI_Comm,int,void*,DMMG**);
228: EXTERN int DMMGDestroy(DMMG*);
229: EXTERN int DMMGSetUp(DMMG*);
230: EXTERN int DMMGSetSLES(DMMG*,int (*)(DMMG,Vec),int (*)(DMMG,Mat));
231: EXTERN int DMMGSetSNES(DMMG*,int (*)(SNES,Vec,Vec,void*),int (*)(SNES,Vec,Mat*,Mat*,MatStructure*,void*));
232: EXTERN int DMMGSetSNESLocal(DMMG*,int(*)(Scalar**,Scalar**,DALocalInfo*,void*),int(*jacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*));
233: EXTERN int DMMGSetInitialGuess(DMMG*,int (*)(SNES,Vec,void*));
234: EXTERN int DMMGView(DMMG*,PetscViewer);
235: EXTERN int DMMGSolve(DMMG*);
236: EXTERN int DMMGSetUseMatrixFree(DMMG*);
237: EXTERN int DMMGSetDM(DMMG*,DM);

239: #define DMMGGetb(ctx)    (ctx)[(ctx)[0]->nlevels-1]->b
240: #define DMMGGetr(ctx)    (ctx)[(ctx)[0]->nlevels-1]->r
241: #define DMMGGetx(ctx)    (ctx)[(ctx)[0]->nlevels-1]->x
242: #define DMMGGetJ(ctx)    (ctx)[(ctx)[0]->nlevels-1]->J
243: #define DMMGGetB(ctx)    (ctx)[(ctx)[0]->nlevels-1]->B
244: #define DMMGGetFine(ctx) (ctx)[(ctx)[0]->nlevels-1]
245: #define DMMGGetSLES(ctx) (ctx)[(ctx)[0]->nlevels-1]->sles
246: #define DMMGGetSNES(ctx) (ctx)[(ctx)[0]->nlevels-1]->snes
247: #define DMMGGetDA(ctx)   (DA)((ctx)[(ctx)[0]->nlevels-1]->dm)
248: #define DMMGGetVecPack(ctx)   (VecPack)((ctx)[(ctx)[0]->nlevels-1]->dm)
249: #define DMMGGetUser(ctx)      (VecPack)((ctx)[(ctx)[0]->nlevels-1]->user)
250: #define DMMGGetLevels(ctx)  (ctx)[0]->nlevels

252: #endif