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,int,int,int,int,int,int,int,int,int *,int *,int *,DA *);
72: EXTERN int DADestroy(DA);
73: EXTERN int DAView(DA,PetscViewer);
75: EXTERN int DAPrintHelp(DA);
77: EXTERN int DAGlobalToLocalBegin(DA,Vec,InsertMode,Vec);
78: EXTERN int DAGlobalToLocalEnd(DA,Vec,InsertMode,Vec);
79: EXTERN int DAGlobalToNaturalBegin(DA,Vec,InsertMode,Vec);
80: EXTERN int DAGlobalToNaturalEnd(DA,Vec,InsertMode,Vec);
81: EXTERN int DANaturalToGlobalBegin(DA,Vec,InsertMode,Vec);
82: EXTERN int DANaturalToGlobalEnd(DA,Vec,InsertMode,Vec);
83: EXTERN int DALocalToLocalBegin(DA,Vec,InsertMode,Vec);
84: EXTERN int DALocalToLocalEnd(DA,Vec,InsertMode,Vec);
85: EXTERN int DALocalToGlobal(DA,Vec,InsertMode,Vec);
86: EXTERN int DALocalToGlobalBegin(DA,Vec,Vec);
87: EXTERN int DALocalToGlobalEnd(DA,Vec,Vec);
88: EXTERN int DAGetOwnershipRange(DA,int **,int **,int **);
89: EXTERN int DACreateGlobalVector(DA,Vec *);
90: EXTERN int DACreateNaturalVector(DA,Vec *);
91: EXTERN int DACreateLocalVector(DA,Vec *);
92: EXTERN int DAGetLocalVector(DA,Vec *);
93: EXTERN int DARestoreLocalVector(DA,Vec *);
94: EXTERN int DAGetGlobalVector(DA,Vec *);
95: EXTERN int DARestoreGlobalVector(DA,Vec *);
96: EXTERN int DALoad(PetscViewer,int,int,int,DA *);
97: EXTERN int DAGetCorners(DA,int*,int*,int*,int*,int*,int*);
98: EXTERN int DAGetGhostCorners(DA,int*,int*,int*,int*,int*,int*);
99: EXTERN int DAGetInfo(DA,int*,int*,int*,int*,int*,int*,int*,int*,int*,DAPeriodicType*,DAStencilType*);
100: EXTERN int DAGetProcessorSubset(DA,DADirection,int,MPI_Comm*);
101: EXTERN int DARefine(DA,MPI_Comm,DA*);
103: EXTERN int DAGlobalToNaturalAllCreate(DA,VecScatter*);
104: EXTERN int DANaturalAllToGlobalCreate(DA,VecScatter*);
106: EXTERN int DAGetGlobalIndices(DA,int*,int**);
107: EXTERN int DAGetISLocalToGlobalMapping(DA,ISLocalToGlobalMapping*);
108: EXTERN int DAGetISLocalToGlobalMappingBlck(DA,ISLocalToGlobalMapping*);
110: EXTERN int DAGetScatter(DA,VecScatter*,VecScatter*,VecScatter*);
112: EXTERN int DAGetAO(DA,AO*);
113: EXTERN int DASetCoordinates(DA,Vec);
114: EXTERN int DAGetCoordinates(DA,Vec *);
115: EXTERN int DASetUniformCoordinates(DA,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal);
116: EXTERN int DASetFieldName(DA,int,const char[]);
117: EXTERN int DAGetFieldName(DA,int,char **);
119: EXTERN int DAVecGetArray(DA,Vec,void **);
120: EXTERN int DAVecRestoreArray(DA,Vec,void **);
122: EXTERN int DASplitComm2d(MPI_Comm,int,int,int,MPI_Comm*);
124: EXTERN int MatRegisterDAAD(void);
125: EXTERN int MatCreateDAAD(DA,Mat*);
127: /*S
128: DALocalInfo - C struct that contains information about a structured grid and a processors logical
129: location in it.
131: Level: beginner
133: Concepts: distributed array
135: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DADestroy(), DA, DAGetLocalInfo(), DAGetInfo()
136: S*/
137: typedef struct {
138: int dim,dof,sw;
139: DAPeriodicType pt;
140: DAStencilType st;
141: int mx,my,mz; /* global number of grid points in each direction */
142: int xs,ys,zs; /* starting point of this processor, excluding ghosts */
143: int xm,ym,zm; /* number of grid points on this processor, excluding ghosts */
144: int gxs,gys,gzs; /* starting point of this processor including ghosts */
145: int gxm,gym,gzm; /* number of grid points on this processor including ghosts */
146: DA da;
147: } DALocalInfo;
149: EXTERN int DAGetLocalInfo(DA,DALocalInfo*);
150: typedef int (*DALocalFunction1)(DALocalInfo*,void*,void*,void*);
151: EXTERN int DAFormFunction1(DA,Vec,Vec,void*);
152: EXTERN int DAFormFunctioni1(DA,int,Vec,PetscScalar*,void*);
153: EXTERN int DAComputeJacobian1WithAdic(DA,Vec,Mat,void*);
154: EXTERN int DAComputeJacobian1WithAdifor(DA,Vec,Mat,void*);
155: EXTERN int DAMultiplyByJacobian1WithAdic(DA,Vec,Vec,Vec,void*);
156: EXTERN int DAMultiplyByJacobian1WithAdifor(DA,Vec,Vec,Vec,void*);
157: EXTERN int DAMultiplyByJacobian1WithAD(DA,Vec,Vec,Vec,void*);
158: EXTERN int DAComputeJacobian1(DA,Vec,Mat,void*);
159: EXTERN int DAGetLocalFunction(DA,DALocalFunction1*);
160: EXTERN int DASetLocalFunction(DA,DALocalFunction1);
161: EXTERN int DASetLocalFunctioni(DA,int (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*));
162: EXTERN int DASetLocalJacobian(DA,DALocalFunction1);
163: EXTERN int DASetLocalAdicFunction_Private(DA,DALocalFunction1);
164: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
165: # define DASetLocalAdicFunction(a,d) DASetLocalAdicFunction_Private(a,(DALocalFunction1)d)
166: #else
167: # define DASetLocalAdicFunction(a,d) DASetLocalAdicFunction_Private(a,0)
168: #endif
169: EXTERN int DASetLocalAdicMFFunction_Private(DA,DALocalFunction1);
170: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
171: # define DASetLocalAdicMFFunction(a,d) DASetLocalAdicMFFunction_Private(a,(DALocalFunction1)d)
172: #else
173: # define DASetLocalAdicMFFunction(a,d) DASetLocalAdicMFFunction_Private(a,0)
174: #endif
175: EXTERN int DASetLocalAdicFunctioni_Private(DA,int (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
176: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
177: # define DASetLocalAdicFunctioni(a,d) DASetLocalAdicFunctioni_Private(a,(int (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
178: #else
179: # define DASetLocalAdicFunctioni(a,d) DASetLocalAdicFunctioni_Private(a,0)
180: #endif
181: EXTERN int DASetLocalAdicMFFunctioni_Private(DA,int (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
182: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
183: # define DASetLocalAdicMFFunctioni(a,d) DASetLocalAdicMFFunctioni_Private(a,(int (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
184: #else
185: # define DASetLocalAdicMFFunctioni(a,d) DASetLocalAdicMFFunctioni_Private(a,0)
186: #endif
187: EXTERN int DAFormFunctioniTest1(DA,void*);
190: #include petscmat.h
191: EXTERN int DAGetColoring(DA,ISColoringType,ISColoring *);
192: EXTERN int DAGetMatrix(DA,MatType,Mat *);
193: EXTERN int DAGetInterpolation(DA,DA,Mat*,Vec*);
195: EXTERN int DAGetAdicArray(DA,PetscTruth,void**,void**,int*);
196: EXTERN int DARestoreAdicArray(DA,PetscTruth,void**,void**,int*);
197: EXTERN int DAGetAdicMFArray(DA,PetscTruth,void**,void**,int*);
198: EXTERN int DARestoreAdicMFArray(DA,PetscTruth,void**,void**,int*);
199: EXTERN int DAGetArray(DA,PetscTruth,void**);
200: EXTERN int DARestoreArray(DA,PetscTruth,void**);
201: EXTERN int ad_DAGetArray(DA,PetscTruth,void**);
202: EXTERN int ad_DARestoreArray(DA,PetscTruth,void**);
203: EXTERN int admf_DAGetArray(DA,PetscTruth,void**);
204: EXTERN int admf_DARestoreArray(DA,PetscTruth,void**);
206: #include petscpf.h
207: EXTERN int DACreatePF(DA,PF*);
209: /*S
210: VecPack - Abstract PETSc object that manages treating several distinct vectors as if they
211: were one. The VecPack routines allow one to manage a nonlinear solver that works on a
212: vector that consists of several distinct parts. This is mostly used for LNKS solvers,
213: that is design optimization problems that are written as a nonlinear system
215: Level: beginner
217: Concepts: multi-component, LNKS solvers
219: .seealso: VecPackCreate(), VecPackDestroy()
220: S*/
221: typedef struct _p_VecPack *VecPack;
223: EXTERN int VecPackCreate(MPI_Comm,VecPack*);
224: EXTERN int VecPackDestroy(VecPack);
225: EXTERN int VecPackAddArray(VecPack,int);
226: EXTERN int VecPackAddDA(VecPack,DA);
227: EXTERN int VecPackAddVecScatter(VecPack,VecScatter);
228: EXTERN int VecPackScatter(VecPack,Vec,...);
229: EXTERN int VecPackGather(VecPack,Vec,...);
230: EXTERN int VecPackGetAccess(VecPack,Vec,...);
231: EXTERN int VecPackRestoreAccess(VecPack,Vec,...);
232: EXTERN int VecPackGetLocalVectors(VecPack,...);
233: EXTERN int VecPackGetEntries(VecPack,...);
234: EXTERN int VecPackRestoreLocalVectors(VecPack,...);
235: EXTERN int VecPackCreateGlobalVector(VecPack,Vec*);
236: EXTERN int VecPackGetGlobalIndices(VecPack,...);
237: EXTERN int VecPackRefine(VecPack,MPI_Comm,VecPack*);
238: EXTERN int VecPackGetInterpolation(VecPack,VecPack,Mat*,Vec*);
240: #include petscsnes.h
242: /*S
243: DM - Abstract PETSc object that manages an abstract grid object
244:
245: Level: intermediate
247: Concepts: grids, grid refinement
249: Notes: The DA object and the VecPack object are examples of DMs
251: .seealso: VecPackCreate(), DA, VecPack
252: S*/
253: typedef struct _p_DM* DM;
255: EXTERN int DMView(DM,PetscViewer);
256: EXTERN int DMDestroy(DM);
257: EXTERN int DMCreateGlobalVector(DM,Vec*);
258: EXTERN int DMGetColoring(DM,ISColoringType,ISColoring*);
259: EXTERN int DMGetMatrix(DM,MatType,Mat*);
260: EXTERN int DMGetInterpolation(DM,DM,Mat*,Vec*);
261: EXTERN int DMRefine(DM,MPI_Comm,DM*);
262: EXTERN int DMGetInterpolationScale(DM,DM,Mat,Vec*);
264: /*S
265: DMMG - Data structure to easily manage multi-level non-linear solvers on grids managed by DM
266:
267: Level: intermediate
269: Concepts: multigrid, Newton-multigrid
271: .seealso: VecPackCreate(), DA, VecPack, DM, DMMGCreate()
272: S*/
273: typedef struct _p_DMMG *DMMG;
274: struct _p_DMMG {
275: DM dm; /* grid information for this level */
276: Vec x,b,r; /* global vectors used in multigrid preconditioner for this level*/
277: Mat J; /* matrix on this level */
278: Mat R; /* restriction to next coarser level (not defined on level 0) */
279: int nlevels; /* number of levels above this one (total number of levels on level 0)*/
280: MPI_Comm comm;
281: int (*solve)(DMMG*,int);
282: void *user;
283: PetscTruth galerkin; /* for A_c = R*A*R^T */
285: /* SLES only */
286: SLES sles;
287: int (*rhs)(DMMG,Vec);
288: PetscTruth matricesset; /* User had called DMMGSetSLES() and the matrices have been computed */
290: /* SNES only */
291: Mat B;
292: Vec Rscale; /* scaling to restriction before computing Jacobian */
293: int (*computejacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
294: int (*computefunction)(SNES,Vec,Vec,void*);
296: MatFDColoring fdcoloring; /* only used with FD coloring for Jacobian */
297: SNES snes;
298: int (*initialguess)(SNES,Vec,void*);
299: Vec w,work1,work2; /* global vectors */
300: Vec lwork1;
301: };
303: EXTERN int DMMGCreate(MPI_Comm,int,void*,DMMG**);
304: EXTERN int DMMGDestroy(DMMG*);
305: EXTERN int DMMGSetUp(DMMG*);
306: EXTERN int DMMGSetSLES(DMMG*,int (*)(DMMG,Vec),int (*)(DMMG,Mat));
307: EXTERN int DMMGSetSNES(DMMG*,int (*)(SNES,Vec,Vec,void*),int (*)(SNES,Vec,Mat*,Mat*,MatStructure*,void*));
308: EXTERN int DMMGSetInitialGuess(DMMG*,int (*)(SNES,Vec,void*));
309: EXTERN int DMMGView(DMMG*,PetscViewer);
310: EXTERN int DMMGSolve(DMMG*);
311: EXTERN int DMMGSetUseMatrixFree(DMMG*);
312: EXTERN int DMMGSetDM(DMMG*,DM);
313: EXTERN int DMMGSetUpLevel(DMMG*,SLES,int);
314: EXTERN int DMMGSetUseGalerkinCoarse(DMMG*);
316: EXTERN int DMMGSetSNESLocal_Private(DMMG*,DALocalFunction1,DALocalFunction1,DALocalFunction1,DALocalFunction1);
317: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
318: # define DMMGSetSNESLocal(dmmg,function,jacobian,ad_function,admf_function)
319: DMMGSetSNESLocal_Private(dmmg,(DALocalFunction1)function,(DALocalFunction1)jacobian,(DALocalFunction1)(ad_function),(DALocalFunction1)(admf_function))
320: #else
321: # define DMMGSetSNESLocal(dmmg,function,jacobian,ad_function,admf_function) DMMGSetSNESLocal_Private(dmmg,(DALocalFunction1)function,(DALocalFunction1)jacobian,(DALocalFunction1)0,(DALocalFunction1)0)
322: #endif
324: EXTERN int DMMGSetSNESLocali_Private(DMMG*,int (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*),int (*)(DALocalInfo*,MatStencil*,void*,void*,void*),int (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
325: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
326: # 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))
327: #else
328: # define DMMGSetSNESLocali(dmmg,function,ad_function,admf_function) DMMGSetSNESLocali_Private(dmmg,(int (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*))function,0,0)
329: #endif
331: #define DMMGGetb(ctx) (ctx)[(ctx)[0]->nlevels-1]->b
332: #define DMMGGetr(ctx) (ctx)[(ctx)[0]->nlevels-1]->r
334: /*MC
335: DMMGGetx - Returns the solution vector from a DMMG solve on the finest grid
337: Synopsis:
338: Vec DMMGGetx(DMMG *dmmg)
340: Not Collective, but resulting vector is parallel
342: Input Parameters:
343: . dmmg - DMMG solve context
345: Level: intermediate
347: Fortran Usage:
348: . DMMGGetx(DMMG dmmg,Vec x,int ierr)
350: .seealso: DMMGCreate(), DMMGSetSNES(), DMMGSetSLES(), DMMGSetSNESLocal()
352: M*/
353: #define DMMGGetx(ctx) (ctx)[(ctx)[0]->nlevels-1]->x
355: #define DMMGGetJ(ctx) (ctx)[(ctx)[0]->nlevels-1]->J
356: #define DMMGGetB(ctx) (ctx)[(ctx)[0]->nlevels-1]->B
357: #define DMMGGetFine(ctx) (ctx)[(ctx)[0]->nlevels-1]
358: #define DMMGGetSLES(ctx) (ctx)[(ctx)[0]->nlevels-1]->sles
359: #define DMMGGetSNES(ctx) (ctx)[(ctx)[0]->nlevels-1]->snes
360: #define DMMGGetDA(ctx) (DA)((ctx)[(ctx)[0]->nlevels-1]->dm)
361: #define DMMGGetVecPack(ctx) (VecPack)((ctx)[(ctx)[0]->nlevels-1]->dm)
362: #define DMMGGetUser(ctx,level) ((ctx)[levels]->user)
363: #define DMMGSetUser(ctx,level,usr) 0,(ctx)[levels]->user = usr
364: #define DMMGGetLevels(ctx) (ctx)[0]->nlevels
366: #endif