Actual source code: petscda.h

  1: /* $Id: petscda.h,v 1.77 2001/09/11 16:34:35 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: /*S
 13:      DA - Abstract PETSc object that manages distributed field data for a single structured grid

 15:    Level: beginner

 17:   Concepts: distributed array

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

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

 27:    Level: beginner

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

 33: /*E
 34:     DAPeriodicType - Is the domain periodic in one or more directions

 36:    Level: beginner

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

 44: /*E
 45:     DAInterpolationType - Defines the type of interpolation that will be returned by 
 46:        DAGetInterpolation.

 48:    Level: beginner

 50: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DA, DAGetInterpolation(), DASetInterpolationType()
 51: E*/
 52: typedef enum { DA_Q0, DA_Q1 } DAInterpolationType;

 54: EXTERN int DASetInterpolationType(DA,DAInterpolationType);

 56: #define DAXPeriodic(pt) ((pt)==DA_XPERIODIC||(pt)==DA_XYPERIODIC||(pt)==DA_XZPERIODIC||(pt)==DA_XYZPERIODIC)
 57: #define DAYPeriodic(pt) ((pt)==DA_YPERIODIC||(pt)==DA_XYPERIODIC||(pt)==DA_YZPERIODIC||(pt)==DA_XYZPERIODIC)
 58: #define DAZPeriodic(pt) ((pt)==DA_ZPERIODIC||(pt)==DA_XZPERIODIC||(pt)==DA_YZPERIODIC||(pt)==DA_XYZPERIODIC)

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

 62: /* Logging support */
 63: extern int DA_COOKIE;
 64: enum {DA_GlobalToLocal, DA_LocalToGlobal, DA_MAX_EVENTS};
 65: extern int DAEvents[DA_MAX_EVENTS];
 66: #define DALogEventBegin(e,o1,o2,o3,o4) PetscLogEventBegin(DAEvents[e],o1,o2,o3,o4)
 67: #define DALogEventEnd(e,o1,o2,o3,o4)   PetscLogEventEnd(DAEvents[e],o1,o2,o3,o4)

 69: EXTERN int   DACreate1d(MPI_Comm,DAPeriodicType,int,int,int,int*,DA *);
 70: EXTERN int   DACreate2d(MPI_Comm,DAPeriodicType,DAStencilType,int,int,int,int,int,int,int*,int*,DA *);
 71: EXTERN int   DACreate3d(MPI_Comm,DAPeriodicType,DAStencilType,
 72:                         int,int,int,int,int,int,int,int,int *,int *,int *,DA *);
 73: EXTERN int   DADestroy(DA);
 74: EXTERN int   DAView(DA,PetscViewer);

 76: EXTERN int   DAPrintHelp(DA);

 78: EXTERN int   DAGlobalToLocalBegin(DA,Vec,InsertMode,Vec);
 79: EXTERN int   DAGlobalToLocalEnd(DA,Vec,InsertMode,Vec);
 80: EXTERN int   DAGlobalToNaturalBegin(DA,Vec,InsertMode,Vec);
 81: EXTERN int   DAGlobalToNaturalEnd(DA,Vec,InsertMode,Vec);
 82: EXTERN int   DANaturalToGlobalBegin(DA,Vec,InsertMode,Vec);
 83: EXTERN int   DANaturalToGlobalEnd(DA,Vec,InsertMode,Vec);
 84: EXTERN int   DALocalToLocalBegin(DA,Vec,InsertMode,Vec);
 85: EXTERN int   DALocalToLocalEnd(DA,Vec,InsertMode,Vec);
 86: EXTERN int   DALocalToGlobal(DA,Vec,InsertMode,Vec);
 87: EXTERN int   DALocalToGlobalBegin(DA,Vec,Vec);
 88: EXTERN int   DALocalToGlobalEnd(DA,Vec,Vec);
 89: EXTERN int   DAGetOwnershipRange(DA,int **,int **,int **);
 90: EXTERN int   DACreateGlobalVector(DA,Vec *);
 91: EXTERN int   DACreateNaturalVector(DA,Vec *);
 92: EXTERN int   DACreateLocalVector(DA,Vec *);
 93: EXTERN int   DAGetLocalVector(DA,Vec *);
 94: EXTERN int   DARestoreLocalVector(DA,Vec *);
 95: EXTERN int   DAGetGlobalVector(DA,Vec *);
 96: EXTERN int   DARestoreGlobalVector(DA,Vec *);
 97: EXTERN int   DALoad(PetscViewer,int,int,int,DA *);
 98: EXTERN int   DAGetCorners(DA,int*,int*,int*,int*,int*,int*);
 99: EXTERN int   DAGetGhostCorners(DA,int*,int*,int*,int*,int*,int*);
100: EXTERN int   DAGetInfo(DA,int*,int*,int*,int*,int*,int*,int*,int*,int*,DAPeriodicType*,DAStencilType*);
101: EXTERN int   DAGetProcessorSubset(DA,DADirection,int,MPI_Comm*);
102: EXTERN int   DARefine(DA,MPI_Comm,DA*);

104: EXTERN int   DAGlobalToNaturalAllCreate(DA,VecScatter*);
105: EXTERN int   DANaturalAllToGlobalCreate(DA,VecScatter*);

107: EXTERN int   DAGetGlobalIndices(DA,int*,int**);
108: EXTERN int   DAGetISLocalToGlobalMapping(DA,ISLocalToGlobalMapping*);
109: EXTERN int   DAGetISLocalToGlobalMappingBlck(DA,ISLocalToGlobalMapping*);

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

113: EXTERN int   DAGetAO(DA,AO*);
114: EXTERN int   DASetCoordinates(DA,Vec);
115: EXTERN int   DAGetCoordinates(DA,Vec *);
116: EXTERN int   DASetUniformCoordinates(DA,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal);
117: EXTERN int   DASetFieldName(DA,int,const char[]);
118: EXTERN int   DAGetFieldName(DA,int,char **);

120: EXTERN int   DAVecGetArray(DA,Vec,void **);
121: EXTERN int   DAVecRestoreArray(DA,Vec,void **);

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

125: EXTERN int   MatRegisterDAAD(void);
126: EXTERN int   MatCreateDAAD(DA,Mat*);

128: /*S
129:      DALocalInfo - C struct that contains information about a structured grid and a processors logical
130:               location in it.

132:    Level: beginner

134:   Concepts: distributed array

136: .seealso:  DACreate1d(), DACreate2d(), DACreate3d(), DADestroy(), DA, DAGetLocalInfo(), DAGetInfo()
137: S*/
138: typedef struct {
139:   int            dim,dof,sw;
140:   DAPeriodicType pt;
141:   DAStencilType  st;
142:   int            mx,my,mz;    /* global number of grid points in each direction */
143:   int            xs,ys,zs;    /* starting point of this processor, excluding ghosts */
144:   int            xm,ym,zm;    /* number of grid points on this processor, excluding ghosts */
145:   int            gxs,gys,gzs;    /* starting point of this processor including ghosts */
146:   int            gxm,gym,gzm;    /* number of grid points on this processor including ghosts */
147:   DA             da;
148: } DALocalInfo;

150: EXTERN int DAGetLocalInfo(DA,DALocalInfo*);
151: typedef int (*DALocalFunction1)(DALocalInfo*,void*,void*,void*);
152: EXTERN int DAFormFunction1(DA,Vec,Vec,void*);
153: EXTERN int DAFormFunctioni1(DA,int,Vec,PetscScalar*,void*);
154: EXTERN int DAComputeJacobian1WithAdic(DA,Vec,Mat,void*);
155: EXTERN int DAComputeJacobian1WithAdifor(DA,Vec,Mat,void*);
156: EXTERN int DAMultiplyByJacobian1WithAdic(DA,Vec,Vec,Vec,void*);
157: EXTERN int DAMultiplyByJacobian1WithAdifor(DA,Vec,Vec,Vec,void*);
158: EXTERN int DAMultiplyByJacobian1WithAD(DA,Vec,Vec,Vec,void*);
159: EXTERN int DAComputeJacobian1(DA,Vec,Mat,void*);
160: EXTERN int DAGetLocalFunction(DA,DALocalFunction1*);
161: EXTERN int DASetLocalFunction(DA,DALocalFunction1);
162: EXTERN int DASetLocalFunctioni(DA,int (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*));
163: EXTERN int DASetLocalJacobian(DA,DALocalFunction1);
164: EXTERN int DASetLocalAdicFunction_Private(DA,DALocalFunction1);
165: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
166: #  define DASetLocalAdicFunction(a,d) DASetLocalAdicFunction_Private(a,(DALocalFunction1)d)
167: #else
168: #  define DASetLocalAdicFunction(a,d) DASetLocalAdicFunction_Private(a,0)
169: #endif
170: EXTERN int DASetLocalAdicMFFunction_Private(DA,DALocalFunction1);
171: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
172: #  define DASetLocalAdicMFFunction(a,d) DASetLocalAdicMFFunction_Private(a,(DALocalFunction1)d)
173: #else
174: #  define DASetLocalAdicMFFunction(a,d) DASetLocalAdicMFFunction_Private(a,0)
175: #endif
176: EXTERN int DASetLocalAdicFunctioni_Private(DA,int (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
177: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
178: #  define DASetLocalAdicFunctioni(a,d) DASetLocalAdicFunctioni_Private(a,(int (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
179: #else
180: #  define DASetLocalAdicFunctioni(a,d) DASetLocalAdicFunctioni_Private(a,0)
181: #endif
182: EXTERN int DASetLocalAdicMFFunctioni_Private(DA,int (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
183: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
184: #  define DASetLocalAdicMFFunctioni(a,d) DASetLocalAdicMFFunctioni_Private(a,(int (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
185: #else
186: #  define DASetLocalAdicMFFunctioni(a,d) DASetLocalAdicMFFunctioni_Private(a,0)
187: #endif
188: EXTERN int DAFormFunctioniTest1(DA,void*);


191:  #include petscmat.h
192: EXTERN int DAGetColoring(DA,ISColoringType,ISColoring *);
193: EXTERN int DAGetMatrix(DA,MatType,Mat *);
194: EXTERN int DAGetInterpolation(DA,DA,Mat*,Vec*);

196: EXTERN int DAGetAdicArray(DA,PetscTruth,void**,void**,int*);
197: EXTERN int DARestoreAdicArray(DA,PetscTruth,void**,void**,int*);
198: EXTERN int DAGetAdicMFArray(DA,PetscTruth,void**,void**,int*);
199: EXTERN int DARestoreAdicMFArray(DA,PetscTruth,void**,void**,int*);
200: EXTERN int DAGetArray(DA,PetscTruth,void**);
201: EXTERN int DARestoreArray(DA,PetscTruth,void**);
202: EXTERN int ad_DAGetArray(DA,PetscTruth,void**);
203: EXTERN int ad_DARestoreArray(DA,PetscTruth,void**);
204: EXTERN int admf_DAGetArray(DA,PetscTruth,void**);
205: EXTERN int admf_DARestoreArray(DA,PetscTruth,void**);

207:  #include petscpf.h
208: EXTERN int DACreatePF(DA,PF*);

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

216:    Level: beginner

218:   Concepts: multi-component, LNKS solvers

220: .seealso:  VecPackCreate(), VecPackDestroy()
221: S*/
222: typedef struct _p_VecPack *VecPack;

224: EXTERN int VecPackCreate(MPI_Comm,VecPack*);
225: EXTERN int VecPackDestroy(VecPack);
226: EXTERN int VecPackAddArray(VecPack,int);
227: EXTERN int VecPackAddDA(VecPack,DA);
228: EXTERN int VecPackAddVecScatter(VecPack,VecScatter);
229: EXTERN int VecPackScatter(VecPack,Vec,...);
230: EXTERN int VecPackGather(VecPack,Vec,...);
231: EXTERN int VecPackGetAccess(VecPack,Vec,...);
232: EXTERN int VecPackRestoreAccess(VecPack,Vec,...);
233: EXTERN int VecPackGetLocalVectors(VecPack,...);
234: EXTERN int VecPackGetEntries(VecPack,...);
235: EXTERN int VecPackRestoreLocalVectors(VecPack,...);
236: EXTERN int VecPackCreateGlobalVector(VecPack,Vec*);
237: EXTERN int VecPackGetGlobalIndices(VecPack,...);
238: EXTERN int VecPackRefine(VecPack,MPI_Comm,VecPack*);
239: EXTERN int VecPackGetInterpolation(VecPack,VecPack,Mat*,Vec*);

241:  #include petscsnes.h

243: /*S
244:      DM - Abstract PETSc object that manages an abstract grid object
245:           
246:    Level: intermediate

248:   Concepts: grids, grid refinement

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

252: .seealso:  VecPackCreate(), DA, VecPack
253: S*/
254: typedef struct _p_DM* DM;

256: EXTERN int DMView(DM,PetscViewer);
257: EXTERN int DMDestroy(DM);
258: EXTERN int DMCreateGlobalVector(DM,Vec*);
259: EXTERN int DMGetColoring(DM,ISColoringType,ISColoring*);
260: EXTERN int DMGetMatrix(DM,MatType,Mat*);
261: EXTERN int DMGetInterpolation(DM,DM,Mat*,Vec*);
262: EXTERN int DMRefine(DM,MPI_Comm,DM*);
263: EXTERN int DMGetInterpolationScale(DM,DM,Mat,Vec*);

265: /*S
266:      DMMG -  Data structure to easily manage multi-level non-linear solvers on grids managed by DM
267:           
268:    Level: intermediate

270:   Concepts: multigrid, Newton-multigrid

272: .seealso:  VecPackCreate(), DA, VecPack, DM, DMMGCreate()
273: S*/
274: typedef struct _p_DMMG *DMMG;
275: struct _p_DMMG {
276:   DM         dm;                   /* grid information for this level */
277:   Vec        x,b,r;                /* global vectors used in multigrid preconditioner for this level*/
278:   Mat        J;                    /* matrix on this level */
279:   Mat        R;                    /* restriction to next coarser level (not defined on level 0) */
280:   int        nlevels;              /* number of levels above this one (total number of levels on level 0)*/
281:   MPI_Comm   comm;
282:   int        (*solve)(DMMG*,int);
283:   void       *user;
284:   PetscTruth galerkin;             /* for A_c = R*A*R^T */

286:   /* SLES only */
287:   SLES       sles;
288:   int        (*rhs)(DMMG,Vec);
289:   PetscTruth matricesset;              /* User had called DMMGSetSLES() and the matrices have been computed */

291:   /* SNES only */
292:   Mat           B;
293:   Vec           Rscale;                /* scaling to restriction before computing Jacobian */
294:   int           (*computejacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
295:   int           (*computefunction)(SNES,Vec,Vec,void*);

297:   MatFDColoring    fdcoloring;            /* only used with FD coloring for Jacobian */
298:   SNES             snes;
299:   int              (*initialguess)(SNES,Vec,void*);
300:   Vec              w,work1,work2;         /* global vectors */
301:   Vec              lwork1;
302: };

304: EXTERN int DMMGCreate(MPI_Comm,int,void*,DMMG**);
305: EXTERN int DMMGDestroy(DMMG*);
306: EXTERN int DMMGSetUp(DMMG*);
307: EXTERN int DMMGSetSLES(DMMG*,int (*)(DMMG,Vec),int (*)(DMMG,Mat));
308: EXTERN int DMMGSetSNES(DMMG*,int (*)(SNES,Vec,Vec,void*),int (*)(SNES,Vec,Mat*,Mat*,MatStructure*,void*));
309: EXTERN int DMMGSetInitialGuess(DMMG*,int (*)(SNES,Vec,void*));
310: EXTERN int DMMGView(DMMG*,PetscViewer);
311: EXTERN int DMMGSolve(DMMG*);
312: EXTERN int DMMGSetUseMatrixFree(DMMG*);
313: EXTERN int DMMGSetDM(DMMG*,DM);
314: EXTERN int DMMGSetUpLevel(DMMG*,SLES,int);
315: EXTERN int DMMGSetUseGalerkinCoarse(DMMG*);

317: EXTERN int DMMGSetSNESLocal_Private(DMMG*,DALocalFunction1,DALocalFunction1,DALocalFunction1,DALocalFunction1);
318: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
319: #  define DMMGSetSNESLocal(dmmg,function,jacobian,ad_function,admf_function) 
320:   DMMGSetSNESLocal_Private(dmmg,(DALocalFunction1)function,(DALocalFunction1)jacobian,(DALocalFunction1)(ad_function),(DALocalFunction1)(admf_function))
321: #else
322: #  define DMMGSetSNESLocal(dmmg,function,jacobian,ad_function,admf_function) DMMGSetSNESLocal_Private(dmmg,(DALocalFunction1)function,(DALocalFunction1)jacobian,(DALocalFunction1)0,(DALocalFunction1)0)
323: #endif

325: EXTERN int DMMGSetSNESLocali_Private(DMMG*,int (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*),int (*)(DALocalInfo*,MatStencil*,void*,void*,void*),int (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
326: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
327: #  define DMMGSetSNESLocali(dmmg,function,ad_function,admf_function) DMMGSetSNESLocali_Private(dmmg,(int (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*))function,(int (*)(DALocalInfo*,MatStencil*,void*,void*,void*))(ad_function),(int (*)(DALocalInfo*,MatStencil*,void*,void*,void*))(admf_function))
328: #else
329: #  define DMMGSetSNESLocali(dmmg,function,ad_function,admf_function) DMMGSetSNESLocali_Private(dmmg,(int (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*))function,0,0)
330: #endif

332: #define DMMGGetb(ctx)              (ctx)[(ctx)[0]->nlevels-1]->b
333: #define DMMGGetr(ctx)              (ctx)[(ctx)[0]->nlevels-1]->r

335: /*MC
336:    DMMGGetx - Returns the solution vector from a DMMG solve on the finest grid

338:    Synopsis:
339:    Vec DMMGGetx(DMMG *dmmg)

341:    Not Collective, but resulting vector is parallel

343:    Input Parameters:
344: .   dmmg - DMMG solve context

346:    Level: intermediate

348:    Fortran Usage:
349: .     DMMGGetx(DMMG dmmg,Vec x,int ierr)

351: .seealso: DMMGCreate(), DMMGSetSNES(), DMMGSetSLES(), DMMGSetSNESLocal()

353: M*/
354: #define DMMGGetx(ctx)              (ctx)[(ctx)[0]->nlevels-1]->x

356: #define DMMGGetJ(ctx)              (ctx)[(ctx)[0]->nlevels-1]->J
357: #define DMMGGetB(ctx)              (ctx)[(ctx)[0]->nlevels-1]->B
358: #define DMMGGetFine(ctx)           (ctx)[(ctx)[0]->nlevels-1]
359: #define DMMGGetSLES(ctx)           (ctx)[(ctx)[0]->nlevels-1]->sles
360: #define DMMGGetSNES(ctx)           (ctx)[(ctx)[0]->nlevels-1]->snes
361: #define DMMGGetDA(ctx)             (DA)((ctx)[(ctx)[0]->nlevels-1]->dm)
362: #define DMMGGetVecPack(ctx)        (VecPack)((ctx)[(ctx)[0]->nlevels-1]->dm)
363: #define DMMGGetUser(ctx,level)     ((ctx)[levels]->user)
364: #define DMMGSetUser(ctx,level,usr) 0,(ctx)[levels]->user = usr
365: #define DMMGGetLevels(ctx)         (ctx)[0]->nlevels

367: #endif