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 northeast, 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: extern int DA_COOKIE;
64: /* Logging support; why is this done this way? Should be changed to the pattern used by other objects */
65: enum {DA_GlobalToLocal, DA_LocalToGlobal, DA_LocalADFunction,DA_MAX_EVENTS};
66: extern int DAEvents[DA_MAX_EVENTS];
67: #define DALogEventBegin(e,o1,o2,o3,o4) PetscLogEventBegin(DAEvents[e],o1,o2,o3,o4)
68: #define DALogEventEnd(e,o1,o2,o3,o4) PetscLogEventEnd(DAEvents[e],o1,o2,o3,o4)
70: EXTERN int DACreate1d(MPI_Comm,DAPeriodicType,int,int,int,int*,DA *);
71: EXTERN int DACreate2d(MPI_Comm,DAPeriodicType,DAStencilType,int,int,int,int,int,int,int*,int*,DA *);
72: EXTERN int DACreate3d(MPI_Comm,DAPeriodicType,DAStencilType,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: #define DAForEachPointBegin2d(info,i,j) {\
151: int _xints = info->xs,_xinte = info->xs+info->xm,_yints = info->ys,_yinte = info->ys+info->ym;\
152: for (j=_yints; j<_yinte; j++) {\
153: for (i=_xints; i<_xinte; i++) {\
155: #define DAForEachPointEnd2d }}}
157:
159: EXTERN int DAGetLocalInfo(DA,DALocalInfo*);
160: typedef int (*DALocalFunction1)(DALocalInfo*,void*,void*,void*);
161: EXTERN int DAFormFunction1(DA,Vec,Vec,void*);
162: EXTERN int DAFormFunctioni1(DA,int,Vec,PetscScalar*,void*);
163: EXTERN int DAComputeJacobian1WithAdic(DA,Vec,Mat,void*);
164: EXTERN int DAComputeJacobian1WithAdifor(DA,Vec,Mat,void*);
165: EXTERN int DAMultiplyByJacobian1WithAdic(DA,Vec,Vec,Vec,void*);
166: EXTERN int DAMultiplyByJacobian1WithAdifor(DA,Vec,Vec,Vec,void*);
167: EXTERN int DAMultiplyByJacobian1WithAD(DA,Vec,Vec,Vec,void*);
168: EXTERN int DAComputeJacobian1(DA,Vec,Mat,void*);
169: EXTERN int DAGetLocalFunction(DA,DALocalFunction1*);
170: EXTERN int DASetLocalFunction(DA,DALocalFunction1);
171: EXTERN int DASetLocalFunctioni(DA,int (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*));
172: EXTERN int DASetLocalJacobian(DA,DALocalFunction1);
173: EXTERN int DASetLocalAdicFunction_Private(DA,DALocalFunction1);
175: /*MC
176: DASetLocalAdicFunction - Caches in a DA a local function computed by ADIC/ADIFOR
178: Collective on DA
180: Synopsis:
181: int int DASetLocalAdicFunction(DA da,DALocalFunction1 ad_lf)
182:
183: Input Parameter:
184: + da - initial distributed array
185: - ad_lf - the local function as computed by ADIC/ADIFOR
187: Level: intermediate
189: .keywords: distributed array, refine
191: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DADestroy(), DAGetLocalFunction(), DASetLocalFunction(),
192: DASetLocalJacobian()
193: M*/
194: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
195: # define DASetLocalAdicFunction(a,d) DASetLocalAdicFunction_Private(a,(DALocalFunction1)d)
196: #else
197: # define DASetLocalAdicFunction(a,d) DASetLocalAdicFunction_Private(a,0)
198: #endif
200: EXTERN int DASetLocalAdicMFFunction_Private(DA,DALocalFunction1);
201: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
202: # define DASetLocalAdicMFFunction(a,d) DASetLocalAdicMFFunction_Private(a,(DALocalFunction1)d)
203: #else
204: # define DASetLocalAdicMFFunction(a,d) DASetLocalAdicMFFunction_Private(a,0)
205: #endif
206: EXTERN int DASetLocalAdicFunctioni_Private(DA,int (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
207: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
208: # define DASetLocalAdicFunctioni(a,d) DASetLocalAdicFunctioni_Private(a,(int (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
209: #else
210: # define DASetLocalAdicFunctioni(a,d) DASetLocalAdicFunctioni_Private(a,0)
211: #endif
212: EXTERN int DASetLocalAdicMFFunctioni_Private(DA,int (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
213: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
214: # define DASetLocalAdicMFFunctioni(a,d) DASetLocalAdicMFFunctioni_Private(a,(int (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
215: #else
216: # define DASetLocalAdicMFFunctioni(a,d) DASetLocalAdicMFFunctioni_Private(a,0)
217: #endif
218: EXTERN int DAFormFunctioniTest1(DA,void*);
221: #include petscmat.h
222: EXTERN int DAGetColoring(DA,ISColoringType,ISColoring *);
223: EXTERN int DAGetMatrix(DA,MatType,Mat *);
224: EXTERN int DASetGetMatrix(DA,int (*)(DA,MatType,Mat *));
225: EXTERN int DAGetInterpolation(DA,DA,Mat*,Vec*);
226: EXTERN int DAGetInjection(DA,DA,VecScatter*);
227: EXTERN int DASetBlockFills(DA,int*,int*);
229: EXTERN int DAGetAdicArray(DA,PetscTruth,void**,void**,int*);
230: EXTERN int DARestoreAdicArray(DA,PetscTruth,void**,void**,int*);
231: EXTERN int DAGetAdicMFArray(DA,PetscTruth,void**,void**,int*);
232: EXTERN int DARestoreAdicMFArray(DA,PetscTruth,void**,void**,int*);
233: EXTERN int DAGetArray(DA,PetscTruth,void**);
234: EXTERN int DARestoreArray(DA,PetscTruth,void**);
235: EXTERN int ad_DAGetArray(DA,PetscTruth,void**);
236: EXTERN int ad_DARestoreArray(DA,PetscTruth,void**);
237: EXTERN int admf_DAGetArray(DA,PetscTruth,void**);
238: EXTERN int admf_DARestoreArray(DA,PetscTruth,void**);
240: #include petscpf.h
241: EXTERN int DACreatePF(DA,PF*);
243: /*S
244: VecPack - Abstract PETSc object that manages treating several distinct vectors as if they
245: were one. The VecPack routines allow one to manage a nonlinear solver that works on a
246: vector that consists of several distinct parts. This is mostly used for LNKS solvers,
247: that is design optimization problems that are written as a nonlinear system
249: Level: beginner
251: Concepts: multi-component, LNKS solvers
253: .seealso: VecPackCreate(), VecPackDestroy()
254: S*/
255: typedef struct _p_VecPack *VecPack;
257: EXTERN int VecPackCreate(MPI_Comm,VecPack*);
258: EXTERN int VecPackDestroy(VecPack);
259: EXTERN int VecPackAddArray(VecPack,int);
260: EXTERN int VecPackAddDA(VecPack,DA);
261: EXTERN int VecPackAddVecScatter(VecPack,VecScatter);
262: EXTERN int VecPackScatter(VecPack,Vec,...);
263: EXTERN int VecPackGather(VecPack,Vec,...);
264: EXTERN int VecPackGetAccess(VecPack,Vec,...);
265: EXTERN int VecPackRestoreAccess(VecPack,Vec,...);
266: EXTERN int VecPackGetLocalVectors(VecPack,...);
267: EXTERN int VecPackGetEntries(VecPack,...);
268: EXTERN int VecPackRestoreLocalVectors(VecPack,...);
269: EXTERN int VecPackCreateGlobalVector(VecPack,Vec*);
270: EXTERN int VecPackGetGlobalIndices(VecPack,...);
271: EXTERN int VecPackRefine(VecPack,MPI_Comm,VecPack*);
272: EXTERN int VecPackGetInterpolation(VecPack,VecPack,Mat*,Vec*);
274: #include petscsnes.h
276: /*S
277: DM - Abstract PETSc object that manages an abstract grid object
278:
279: Level: intermediate
281: Concepts: grids, grid refinement
283: Notes: The DA object and the VecPack object are examples of DMs
285: .seealso: VecPackCreate(), DA, VecPack
286: S*/
287: typedef struct _p_DM* DM;
289: EXTERN int DMView(DM,PetscViewer);
290: EXTERN int DMDestroy(DM);
291: EXTERN int DMCreateGlobalVector(DM,Vec*);
292: EXTERN int DMGetColoring(DM,ISColoringType,ISColoring*);
293: EXTERN int DMGetMatrix(DM,MatType,Mat*);
294: EXTERN int DMGetInterpolation(DM,DM,Mat*,Vec*);
295: EXTERN int DMGetInjection(DM,DM,VecScatter*);
296: EXTERN int DMRefine(DM,MPI_Comm,DM*);
297: EXTERN int DMGetInterpolationScale(DM,DM,Mat,Vec*);
299: typedef struct NLF_DAAD* NLF;
301: /*S
302: DMMG - Data structure to easily manage multi-level non-linear solvers on grids managed by DM
303:
304: Level: intermediate
306: Concepts: multigrid, Newton-multigrid
308: .seealso: VecPackCreate(), DA, VecPack, DM, DMMGCreate()
309: S*/
310: typedef struct _p_DMMG *DMMG;
311: struct _p_DMMG {
312: DM dm; /* grid information for this level */
313: Vec x,b,r; /* global vectors used in multigrid preconditioner for this level*/
314: Mat J; /* matrix on this level */
315: Mat R; /* restriction to next coarser level (not defined on level 0) */
316: int nlevels; /* number of levels above this one (total number of levels on level 0)*/
317: MPI_Comm comm;
318: int (*solve)(DMMG*,int);
319: void *user;
320: PetscTruth galerkin; /* for A_c = R*A*R^T */
322: /* SLES only */
323: SLES sles;
324: int (*rhs)(DMMG,Vec);
325: PetscTruth matricesset; /* User had called DMMGSetSLES() and the matrices have been computed */
327: /* SNES only */
328: Mat B;
329: Vec Rscale; /* scaling to restriction before computing Jacobian */
330: int (*computejacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
331: int (*computefunction)(SNES,Vec,Vec,void*);
333: PetscTruth updatejacobian; /* compute new Jacobian when DMMGComputeJacobian_Multigrid() is called */
334: int updatejacobianperiod; /* how often, inside a SNES, the Jacobian is recomputed */
336: MatFDColoring fdcoloring; /* only used with FD coloring for Jacobian */
337: SNES snes;
338: int (*initialguess)(SNES,Vec,void*);
339: Vec w,work1,work2; /* global vectors */
340: Vec lwork1;
342: /* FAS only */
343: NLF nlf; /* FAS smoother object */
344: VecScatter inject; /* inject from this level to the next coarsest */
345: PetscTruth monitor,monitorall;
346: int presmooth,postsmooth,coarsesmooth;
347: PetscReal rtol,atol,rrtol; /* convergence tolerance */
348:
349: };
351: EXTERN int DMMGCreate(MPI_Comm,int,void*,DMMG**);
352: EXTERN int DMMGDestroy(DMMG*);
353: EXTERN int DMMGSetUp(DMMG*);
354: EXTERN int DMMGSetSLES(DMMG*,int (*)(DMMG,Vec),int (*)(DMMG,Mat));
355: EXTERN int DMMGSetSNES(DMMG*,int (*)(SNES,Vec,Vec,void*),int (*)(SNES,Vec,Mat*,Mat*,MatStructure*,void*));
356: EXTERN int DMMGSetInitialGuess(DMMG*,int (*)(SNES,Vec,void*));
357: EXTERN int DMMGView(DMMG*,PetscViewer);
358: EXTERN int DMMGSolve(DMMG*);
359: EXTERN int DMMGSetUseMatrixFree(DMMG*);
360: EXTERN int DMMGSetDM(DMMG*,DM);
361: EXTERN int DMMGSetUpLevel(DMMG*,SLES,int);
362: EXTERN int DMMGSetUseGalerkinCoarse(DMMG*);
364: EXTERN int DMMGSetSNESLocal_Private(DMMG*,DALocalFunction1,DALocalFunction1,DALocalFunction1,DALocalFunction1);
365: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
366: # define DMMGSetSNESLocal(dmmg,function,jacobian,ad_function,admf_function) \
367: DMMGSetSNESLocal_Private(dmmg,(DALocalFunction1)function,(DALocalFunction1)jacobian,(DALocalFunction1)(ad_function),(DALocalFunction1)(admf_function))
368: #else
369: # define DMMGSetSNESLocal(dmmg,function,jacobian,ad_function,admf_function) DMMGSetSNESLocal_Private(dmmg,(DALocalFunction1)function,(DALocalFunction1)jacobian,(DALocalFunction1)0,(DALocalFunction1)0)
370: #endif
372: EXTERN int DMMGSetSNESLocali_Private(DMMG*,int (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*),int (*)(DALocalInfo*,MatStencil*,void*,void*,void*),int (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
373: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
374: # 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))
375: #else
376: # define DMMGSetSNESLocali(dmmg,function,ad_function,admf_function) DMMGSetSNESLocali_Private(dmmg,(int (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*))function,0,0)
377: #endif
379: #define DMMGGetb(ctx) (ctx)[(ctx)[0]->nlevels-1]->b
380: #define DMMGGetr(ctx) (ctx)[(ctx)[0]->nlevels-1]->r
382: /*MC
383: DMMGGetx - Returns the solution vector from a DMMG solve on the finest grid
385: Synopsis:
386: Vec DMMGGetx(DMMG *dmmg)
388: Not Collective, but resulting vector is parallel
390: Input Parameters:
391: . dmmg - DMMG solve context
393: Level: intermediate
395: Fortran Usage:
396: . DMMGGetx(DMMG dmmg,Vec x,int ierr)
398: .seealso: DMMGCreate(), DMMGSetSNES(), DMMGSetSLES(), DMMGSetSNESLocal()
400: M*/
401: #define DMMGGetx(ctx) (ctx)[(ctx)[0]->nlevels-1]->x
403: #define DMMGGetJ(ctx) (ctx)[(ctx)[0]->nlevels-1]->J
404: #define DMMGGetComm(ctx) (ctx)[(ctx)[0]->nlevels-1]->comm
405: #define DMMGGetB(ctx) (ctx)[(ctx)[0]->nlevels-1]->B
406: #define DMMGGetFine(ctx) (ctx)[(ctx)[0]->nlevels-1]
407: #define DMMGGetSLES(ctx) (ctx)[(ctx)[0]->nlevels-1]->sles
408: #define DMMGGetSNES(ctx) (ctx)[(ctx)[0]->nlevels-1]->snes
409: #define DMMGGetDA(ctx) (DA)((ctx)[(ctx)[0]->nlevels-1]->dm)
410: #define DMMGGetVecPack(ctx) (VecPack)((ctx)[(ctx)[0]->nlevels-1]->dm)
411: #define DMMGGetUser(ctx,level) ((ctx)[levels]->user)
412: #define DMMGSetUser(ctx,level,usr) 0,(ctx)[levels]->user = usr
413: #define DMMGGetLevels(ctx) (ctx)[0]->nlevels
415: #endif