Actual source code: petscda.h
1: /*
2: Regular array object, for easy parallelism of simple grid
3: problems on regular distributed arrays.
4: */
7: #include petscvec.h
8: #include petscao.h
11: /*S
12: DA - Abstract PETSc object that manages distributed field data for a single structured grid
14: Level: beginner
16: Concepts: distributed array
18: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DADestroy(), VecScatter, DACreate(), DM, VecPack
19: S*/
20: typedef struct _p_DA* DA;
22: /*E
23: DAStencilType - Determines if the stencil extends only along the coordinate directions, or also
24: to the northeast, northwest etc
26: Level: beginner
28: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DA, DACreate()
29: E*/
30: typedef enum { DA_STENCIL_STAR,DA_STENCIL_BOX } DAStencilType;
32: /*MC
33: DA_STENCIL_STAR - "Star"-type stencil. In logical grid coordinates, only (i,j,k), (i+s,j,k), (i,j+s,k),
34: (i,j,k+s) are in the stencil NOT, for example, (i+s,j+s,k)
36: Level: beginner
38: .seealso: DA_STENCIL_BOX, DAStencilType
39: M*/
41: /*MC
42: DA_STENCIL_Box - "Box"-type stencil. In logical grid coordinates, any of (i,j,k), (i+s,j+r,k+t) may
43: be in the stencil.
45: Level: beginner
47: .seealso: DA_STENCIL_STAR, DAStencilType
48: M*/
50: /*E
51: DAPeriodicType - Is the domain periodic in one or more directions
53: Level: beginner
55: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DA, DACreate()
56: E*/
57: typedef enum { DA_NONPERIODIC,DA_XPERIODIC,DA_YPERIODIC,DA_XYPERIODIC,
58: DA_XYZPERIODIC,DA_XZPERIODIC,DA_YZPERIODIC,DA_ZPERIODIC}
59: DAPeriodicType;
61: /*E
62: DAInterpolationType - Defines the type of interpolation that will be returned by
63: DAGetInterpolation.
65: Level: beginner
67: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DA, DAGetInterpolation(), DASetInterpolationType(), DACreate()
68: E*/
69: typedef enum { DA_Q0, DA_Q1 } DAInterpolationType;
71: EXTERN PetscErrorCode DASetInterpolationType(DA,DAInterpolationType);
73: /*E
74: DAElementType - Defines the type of elements that will be returned by
75: DAGetElements.
77: Level: beginner
79: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DA, DAGetInterpolation(), DASetInterpolationType(),
80: DASetElementType(), DAGetElements(), DARestoreElements(), DACreate()
81: E*/
82: typedef enum { DA_ELEMENT_P1, DA_ELEMENT_Q1 } DAElementType;
84: EXTERN PetscErrorCode DASetElementType(DA,DAElementType);
85: EXTERN PetscErrorCode DAGetElements(DA,PetscInt *,const PetscInt*[]);
86: EXTERN PetscErrorCode DARestoreElements(DA,PetscInt *,const PetscInt*[]);
89: #define DAXPeriodic(pt) ((pt)==DA_XPERIODIC||(pt)==DA_XYPERIODIC||(pt)==DA_XZPERIODIC||(pt)==DA_XYZPERIODIC)
90: #define DAYPeriodic(pt) ((pt)==DA_YPERIODIC||(pt)==DA_XYPERIODIC||(pt)==DA_YZPERIODIC||(pt)==DA_XYZPERIODIC)
91: #define DAZPeriodic(pt) ((pt)==DA_ZPERIODIC||(pt)==DA_XZPERIODIC||(pt)==DA_YZPERIODIC||(pt)==DA_XYZPERIODIC)
93: typedef enum { DA_X,DA_Y,DA_Z } DADirection;
97: EXTERN PetscErrorCode DACreate1d(MPI_Comm,DAPeriodicType,PetscInt,PetscInt,PetscInt,PetscInt*,DA *);
98: EXTERN PetscErrorCode DACreate2d(MPI_Comm,DAPeriodicType,DAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,PetscInt*,DA *);
99: EXTERN PetscErrorCode DACreate3d(MPI_Comm,DAPeriodicType,DAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,PetscInt*,PetscInt*,DA*);
100: EXTERN PetscErrorCode DACreate(MPI_Comm,PetscInt,DAPeriodicType,DAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,PetscInt*,PetscInt*,DA*);
101: EXTERN PetscErrorCode DADestroy(DA);
102: EXTERN PetscErrorCode DAView(DA,PetscViewer);
104: EXTERN PetscErrorCode DAGlobalToLocalBegin(DA,Vec,InsertMode,Vec);
105: EXTERN PetscErrorCode DAGlobalToLocalEnd(DA,Vec,InsertMode,Vec);
106: EXTERN PetscErrorCode DAGlobalToNaturalBegin(DA,Vec,InsertMode,Vec);
107: EXTERN PetscErrorCode DAGlobalToNaturalEnd(DA,Vec,InsertMode,Vec);
108: EXTERN PetscErrorCode DANaturalToGlobalBegin(DA,Vec,InsertMode,Vec);
109: EXTERN PetscErrorCode DANaturalToGlobalEnd(DA,Vec,InsertMode,Vec);
110: EXTERN PetscErrorCode DALocalToLocalBegin(DA,Vec,InsertMode,Vec);
111: EXTERN PetscErrorCode DALocalToLocalEnd(DA,Vec,InsertMode,Vec);
112: EXTERN PetscErrorCode DALocalToGlobal(DA,Vec,InsertMode,Vec);
113: EXTERN PetscErrorCode DALocalToGlobalBegin(DA,Vec,Vec);
114: EXTERN PetscErrorCode DALocalToGlobalEnd(DA,Vec,Vec);
115: EXTERN PetscErrorCode DAGetOwnershipRange(DA,PetscInt **,PetscInt **,PetscInt **);
116: EXTERN PetscErrorCode DACreateGlobalVector(DA,Vec *);
117: EXTERN PetscErrorCode DACreateNaturalVector(DA,Vec *);
118: EXTERN PetscErrorCode DACreateLocalVector(DA,Vec *);
119: EXTERN PetscErrorCode DAGetLocalVector(DA,Vec *);
120: EXTERN PetscErrorCode DARestoreLocalVector(DA,Vec *);
121: EXTERN PetscErrorCode DAGetGlobalVector(DA,Vec *);
122: EXTERN PetscErrorCode DARestoreGlobalVector(DA,Vec *);
123: EXTERN PetscErrorCode DALoad(PetscViewer,PetscInt,PetscInt,PetscInt,DA *);
124: EXTERN PetscErrorCode DAGetCorners(DA,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
125: EXTERN PetscErrorCode DAGetGhostCorners(DA,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
126: EXTERN PetscErrorCode DAGetInfo(DA,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,DAPeriodicType*,DAStencilType*);
127: EXTERN PetscErrorCode DAGetProcessorSubset(DA,DADirection,PetscInt,MPI_Comm*);
128: EXTERN PetscErrorCode DARefine(DA,MPI_Comm,DA*);
130: EXTERN PetscErrorCode DAGlobalToNaturalAllCreate(DA,VecScatter*);
131: EXTERN PetscErrorCode DANaturalAllToGlobalCreate(DA,VecScatter*);
133: EXTERN PetscErrorCode DAGetGlobalIndices(DA,PetscInt*,PetscInt**);
134: EXTERN PetscErrorCode DAGetISLocalToGlobalMapping(DA,ISLocalToGlobalMapping*);
135: EXTERN PetscErrorCode DAGetISLocalToGlobalMappingBlck(DA,ISLocalToGlobalMapping*);
137: EXTERN PetscErrorCode DAGetScatter(DA,VecScatter*,VecScatter*,VecScatter*);
139: EXTERN PetscErrorCode DAGetAO(DA,AO*);
140: EXTERN PetscErrorCode DASetCoordinates(DA,Vec);
141: EXTERN PetscErrorCode DAGetCoordinates(DA,Vec *);
142: EXTERN PetscErrorCode DAGetGhostedCoordinates(DA,Vec *);
143: EXTERN PetscErrorCode DAGetCoordinateDA(DA,DA *);
144: EXTERN PetscErrorCode DASetUniformCoordinates(DA,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal);
145: EXTERN PetscErrorCode DASetFieldName(DA,PetscInt,const char[]);
146: EXTERN PetscErrorCode DAGetFieldName(DA,PetscInt,char **);
148: EXTERN PetscErrorCode DAVecGetArray(DA,Vec,void *);
149: EXTERN PetscErrorCode DAVecRestoreArray(DA,Vec,void *);
151: EXTERN PetscErrorCode DAVecGetArrayDOF(DA,Vec,void *);
152: EXTERN PetscErrorCode DAVecRestoreArrayDOF(DA,Vec,void *);
154: EXTERN PetscErrorCode DASplitComm2d(MPI_Comm,PetscInt,PetscInt,PetscInt,MPI_Comm*);
156: /*S
157: SDA - This provides a simplified interface to the DA distributed
158: array object in PETSc. This is intended for people who are
159: NOT using PETSc vectors or objects but just want to distribute
160: simple rectangular arrays amoung a number of procesors and have
161: PETSc handle moving the ghost-values when needed.
163: In certain applications this can serve as a replacement for
164: BlockComm (which is apparently being phased out?).
167: Level: beginner
169: Concepts: simplified distributed array
171: .seealso: SDACreate1d(), SDACreate2d(), SDACreate3d(), SDADestroy(), DA, SDALocalToLocalBegin(),
172: SDALocalToLocalEnd(), SDAGetCorners(), SDAGetGhostCorners()
173: S*/
174: typedef struct _n_SDA* SDA;
176: EXTERN PetscErrorCode SDACreate3d(MPI_Comm,DAPeriodicType,DAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,PetscInt*,PetscInt*,SDA*);
177: EXTERN PetscErrorCode SDACreate2d(MPI_Comm,DAPeriodicType,DAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,PetscInt*,SDA*);
178: EXTERN PetscErrorCode SDACreate1d(MPI_Comm,DAPeriodicType,PetscInt,PetscInt,PetscInt,PetscInt*,SDA*);
179: EXTERN PetscErrorCode SDADestroy(SDA);
180: EXTERN PetscErrorCode SDALocalToLocalBegin(SDA,PetscScalar*,InsertMode,PetscScalar*);
181: EXTERN PetscErrorCode SDALocalToLocalEnd(SDA,PetscScalar*,InsertMode,PetscScalar*);
182: EXTERN PetscErrorCode SDAGetCorners(SDA,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
183: EXTERN PetscErrorCode SDAGetGhostCorners(SDA,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
184: EXTERN PetscErrorCode SDAArrayView(SDA,PetscScalar*,PetscViewer);
186: EXTERN PetscErrorCode MatRegisterDAAD(void);
187: EXTERN PetscErrorCode MatCreateDAAD(DA,Mat*);
189: /*S
190: DALocalInfo - C struct that contains information about a structured grid and a processors logical
191: location in it.
193: Level: beginner
195: Concepts: distributed array
197: Developer note: Then entries in this struct are int instead of PetscInt so that the elements may
198: be extracted in Fortran as if from an integer array
200: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DADestroy(), DA, DAGetLocalInfo(), DAGetInfo()
201: S*/
202: typedef struct {
203: PetscInt dim,dof,sw;
204: PetscInt mx,my,mz; /* global number of grid points in each direction */
205: PetscInt xs,ys,zs; /* starting pointd of this processor, excluding ghosts */
206: PetscInt xm,ym,zm; /* number of grid points on this processor, excluding ghosts */
207: PetscInt gxs,gys,gzs; /* starting point of this processor including ghosts */
208: PetscInt gxm,gym,gzm; /* number of grid points on this processor including ghosts */
209: DAPeriodicType pt;
210: DAStencilType st;
211: DA da;
212: } DALocalInfo;
214: /*MC
215: DAForEachPointBegin2d - Starts a loop over the local part of a two dimensional DA
217: Synopsis:
218: void DAForEachPointBegin2d(DALocalInfo *info,PetscInt i,PetscInt j);
219:
220: Level: intermediate
222: .seealso: DAForEachPointEnd2d(), DAVecGetArray()
223: M*/
224: #define DAForEachPointBegin2d(info,i,j) {\
225: PetscInt _xints = info->xs,_xinte = info->xs+info->xm,_yints = info->ys,_yinte = info->ys+info->ym;\
226: for (j=_yints; j<_yinte; j++) {\
227: for (i=_xints; i<_xinte; i++) {\
229: /*MC
230: DAForEachPointEnd2d - Ends a loop over the local part of a two dimensional DA
232: Synopsis:
233: void DAForEachPointEnd2d;
234:
235: Level: intermediate
237: .seealso: DAForEachPointBegin2d(), DAVecGetArray()
238: M*/
239: #define DAForEachPointEnd2d }}}
241: /*MC
242: DACoor2d - Structure for holding 2d (x and y) coordinates.
244: Level: intermediate
246: Sample Usage:
247: DACoor2d **coors;
248: Vec vcoors;
249: DA cda;
251: DAGetCoordinates(da,&vcoors);
252: DAGetCoordinateDA(da,&cda);
253: DAVecGetArray(dac,vcoors,&coors);
254: DAGetCorners(dac,&mstart,&nstart,0,&m,&n,0)
255: for (i=mstart; i<mstart+m; i++) {
256: for (j=nstart; j<nstart+n; j++) {
257: x = coors[j][i].x;
258: y = coors[j][i].y;
259: ......
260: }
261: }
262: DAVecRestoreArray(dac,vcoors,&coors);
264: .seealso: DACoor3d, DAForEachPointBegin(), DAGetCoordinateDA(), DAGetCoordinates(), DAGetGhostCoordinates()
265: M*/
266: typedef struct {PetscScalar x,y;} DACoor2d;
268: /*MC
269: DACoor3d - Structure for holding 3d (x, y and z) coordinates.
271: Level: intermediate
273: Sample Usage:
274: DACoor3d **coors;
275: Vec vcoors;
276: DA cda;
278: DAGetCoordinates(da,&vcoors);
279: DAGetCoordinateDA(da,&cda);
280: DAVecGetArray(dac,vcoors,&coors);
281: DAGetCorners(dac,&mstart,&nstart,&pstart,&m,&n,&p)
282: for (i=mstart; i<mstart+m; i++) {
283: for (j=nstart; j<nstart+n; j++) {
284: for (k=pstart; k<pstart+p; k++) {
285: x = coors[k][j][i].x;
286: y = coors[k][j][i].y;
287: z = coors[k][j][i].z;
288: ......
289: }
290: }
291: DAVecRestoreArray(dac,vcoors,&coors);
293: .seealso: DACoor2d, DAForEachPointBegin(), DAGetCoordinateDA(), DAGetCoordinates(), DAGetGhostCoordinates()
294: M*/
295: typedef struct {PetscScalar x,y,z;} DACoor3d;
296:
297: EXTERN PetscErrorCode DAGetLocalInfo(DA,DALocalInfo*);
298: typedef PetscErrorCode (*DALocalFunction1)(DALocalInfo*,void*,void*,void*);
299: EXTERN PetscErrorCode DAFormFunctionLocal(DA, DALocalFunction1, Vec, Vec, void *);
300: EXTERN PetscErrorCode DAFormFunctionLocalGhost(DA, DALocalFunction1, Vec, Vec, void *);
301: EXTERN PetscErrorCode DAFormJacobianLocal(DA, DALocalFunction1, Vec, Mat, void *);
302: EXTERN PetscErrorCode DAFormFunction1(DA,Vec,Vec,void*);
303: EXTERN PetscErrorCode DAFormFunction(DA,PetscErrorCode (*)(void),Vec,Vec,void*);
304: EXTERN PetscErrorCode DAFormFunctioni1(DA,PetscInt,Vec,PetscScalar*,void*);
305: EXTERN PetscErrorCode DAFormFunctionib1(DA,PetscInt,Vec,PetscScalar*,void*);
306: EXTERN PetscErrorCode DAComputeJacobian1WithAdic(DA,Vec,Mat,void*);
307: EXTERN PetscErrorCode DAComputeJacobian1WithAdifor(DA,Vec,Mat,void*);
308: EXTERN PetscErrorCode DAMultiplyByJacobian1WithAdic(DA,Vec,Vec,Vec,void*);
309: EXTERN PetscErrorCode DAMultiplyByJacobian1WithAdifor(DA,Vec,Vec,Vec,void*);
310: EXTERN PetscErrorCode DAMultiplyByJacobian1WithAD(DA,Vec,Vec,Vec,void*);
311: EXTERN PetscErrorCode DAComputeJacobian1(DA,Vec,Mat,void*);
312: EXTERN PetscErrorCode DAGetLocalFunction(DA,DALocalFunction1*);
313: EXTERN PetscErrorCode DASetLocalFunction(DA,DALocalFunction1);
314: EXTERN PetscErrorCode DASetLocalFunctioni(DA,PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*));
315: EXTERN PetscErrorCode DASetLocalFunctionib(DA,PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*));
316: EXTERN PetscErrorCode DASetLocalJacobian(DA,DALocalFunction1);
317: EXTERN PetscErrorCode DASetLocalAdicFunction_Private(DA,DALocalFunction1);
319: /*MC
320: DASetLocalAdicFunction - Caches in a DA a local function computed by ADIC/ADIFOR
322: Collective on DA
324: Synopsis:
325: PetscErrorCode DASetLocalAdicFunction(DA da,DALocalFunction1 ad_lf)
326:
327: Input Parameter:
328: + da - initial distributed array
329: - ad_lf - the local function as computed by ADIC/ADIFOR
331: Level: intermediate
333: .keywords: distributed array, refine
335: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DADestroy(), DAGetLocalFunction(), DASetLocalFunction(),
336: DASetLocalJacobian()
337: M*/
338: #if defined(PETSC_HAVE_ADIC)
339: # define DASetLocalAdicFunction(a,d) DASetLocalAdicFunction_Private(a,(DALocalFunction1)d)
340: #else
341: # define DASetLocalAdicFunction(a,d) DASetLocalAdicFunction_Private(a,0)
342: #endif
344: EXTERN PetscErrorCode DASetLocalAdicMFFunction_Private(DA,DALocalFunction1);
345: #if defined(PETSC_HAVE_ADIC)
346: # define DASetLocalAdicMFFunction(a,d) DASetLocalAdicMFFunction_Private(a,(DALocalFunction1)d)
347: #else
348: # define DASetLocalAdicMFFunction(a,d) DASetLocalAdicMFFunction_Private(a,0)
349: #endif
350: EXTERN PetscErrorCode DASetLocalAdicFunctioni_Private(DA,PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
351: #if defined(PETSC_HAVE_ADIC)
352: # define DASetLocalAdicFunctioni(a,d) DASetLocalAdicFunctioni_Private(a,(PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
353: #else
354: # define DASetLocalAdicFunctioni(a,d) DASetLocalAdicFunctioni_Private(a,0)
355: #endif
356: EXTERN PetscErrorCode DASetLocalAdicMFFunctioni_Private(DA,PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
357: #if defined(PETSC_HAVE_ADIC)
358: # define DASetLocalAdicMFFunctioni(a,d) DASetLocalAdicMFFunctioni_Private(a,(PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
359: #else
360: # define DASetLocalAdicMFFunctioni(a,d) DASetLocalAdicMFFunctioni_Private(a,0)
361: #endif
363: EXTERN PetscErrorCode DASetLocalAdicFunctionib_Private(DA,PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
364: #if defined(PETSC_HAVE_ADIC)
365: # define DASetLocalAdicFunctionib(a,d) DASetLocalAdicFunctionib_Private(a,(PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
366: #else
367: # define DASetLocalAdicFunctionib(a,d) DASetLocalAdicFunctionib_Private(a,0)
368: #endif
369: EXTERN PetscErrorCode DASetLocalAdicMFFunctionib_Private(DA,PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
370: #if defined(PETSC_HAVE_ADIC)
371: # define DASetLocalAdicMFFunctionib(a,d) DASetLocalAdicMFFunctionib_Private(a,(PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
372: #else
373: # define DASetLocalAdicMFFunctionib(a,d) DASetLocalAdicMFFunctionib_Private(a,0)
374: #endif
376: EXTERN PetscErrorCode DAFormFunctioniTest1(DA,void*);
378: #include petscmat.h
379: EXTERN PetscErrorCode DAGetColoring(DA,ISColoringType,ISColoring *);
380: EXTERN PetscErrorCode DAGetMatrix(DA, MatType,Mat *);
381: EXTERN PetscErrorCode DASetGetMatrix(DA,PetscErrorCode (*)(DA, MatType,Mat *));
382: EXTERN PetscErrorCode DAGetInterpolation(DA,DA,Mat*,Vec*);
383: EXTERN PetscErrorCode DAGetInjection(DA,DA,VecScatter*);
384: EXTERN PetscErrorCode DASetBlockFills(DA,PetscInt*,PetscInt*);
385: EXTERN PetscErrorCode DASetMatPreallocateOnly(DA,PetscTruth);
386: EXTERN PetscErrorCode DASetRefinementFactor(DA,PetscInt,PetscInt,PetscInt);
387: EXTERN PetscErrorCode DAGetRefinementFactor(DA,PetscInt*,PetscInt*,PetscInt*);
389: EXTERN PetscErrorCode DAGetAdicArray(DA,PetscTruth,void**,void**,PetscInt*);
390: EXTERN PetscErrorCode DARestoreAdicArray(DA,PetscTruth,void**,void**,PetscInt*);
391: EXTERN PetscErrorCode DAGetAdicMFArray(DA,PetscTruth,void**,void**,PetscInt*);
392: EXTERN PetscErrorCode DAGetAdicMFArray4(DA,PetscTruth,void**,void**,PetscInt*);
393: EXTERN PetscErrorCode DAGetAdicMFArray9(DA,PetscTruth,void**,void**,PetscInt*);
394: EXTERN PetscErrorCode DAGetAdicMFArrayb(DA,PetscTruth,void**,void**,PetscInt*);
395: EXTERN PetscErrorCode DARestoreAdicMFArray(DA,PetscTruth,void**,void**,PetscInt*);
396: EXTERN PetscErrorCode DAGetArray(DA,PetscTruth,void**);
397: EXTERN PetscErrorCode DARestoreArray(DA,PetscTruth,void**);
398: EXTERN PetscErrorCode ad_DAGetArray(DA,PetscTruth,void**);
399: EXTERN PetscErrorCode ad_DARestoreArray(DA,PetscTruth,void**);
400: EXTERN PetscErrorCode admf_DAGetArray(DA,PetscTruth,void**);
401: EXTERN PetscErrorCode admf_DARestoreArray(DA,PetscTruth,void**);
403: #include petscpf.h
404: EXTERN PetscErrorCode DACreatePF(DA,PF*);
406: /*S
407: VecPack - Abstract PETSc object that manages treating several distinct vectors as if they
408: were one. The VecPack routines allow one to manage a nonlinear solver that works on a
409: vector that consists of several distinct parts. This is mostly used for LNKS solvers,
410: that is design optimization problems that are written as a nonlinear system
412: Level: beginner
414: Concepts: multi-component, LNKS solvers
416: .seealso: VecPackCreate(), VecPackDestroy(), DM
417: S*/
418: typedef struct _p_VecPack* VecPack;
420: EXTERN PetscErrorCode VecPackCreate(MPI_Comm,VecPack*);
421: EXTERN PetscErrorCode VecPackDestroy(VecPack);
422: EXTERN PetscErrorCode VecPackAddArray(VecPack,PetscMPIInt,PetscInt);
423: EXTERN PetscErrorCode VecPackAddDA(VecPack,DA);
424: EXTERN PetscErrorCode VecPackAddVecScatter(VecPack,VecScatter);
425: EXTERN PetscErrorCode VecPackScatter(VecPack,Vec,...);
426: EXTERN PetscErrorCode VecPackGather(VecPack,Vec,...);
427: EXTERN PetscErrorCode VecPackGetAccess(VecPack,Vec,...);
428: EXTERN PetscErrorCode VecPackRestoreAccess(VecPack,Vec,...);
429: EXTERN PetscErrorCode VecPackGetLocalVectors(VecPack,...);
430: EXTERN PetscErrorCode VecPackGetEntries(VecPack,...);
431: EXTERN PetscErrorCode VecPackRestoreLocalVectors(VecPack,...);
432: EXTERN PetscErrorCode VecPackCreateGlobalVector(VecPack,Vec*);
433: EXTERN PetscErrorCode VecPackGetGlobalIndices(VecPack,...);
434: EXTERN PetscErrorCode VecPackRefine(VecPack,MPI_Comm,VecPack*);
435: EXTERN PetscErrorCode VecPackGetInterpolation(VecPack,VecPack,Mat*,Vec*);
436: EXTERN PetscErrorCode VecPackGetMatrix(VecPack,MatType,Mat*);
437: EXTERN PetscErrorCode VecPackGetColoring(VecPack,ISColoringType,ISColoring*);
439: /*S
440: Slice - Abstract PETSc object that manages distributed field data for a simple unstructured matrix
442: Level: beginner
444: Concepts: distributed array
446: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DADestroy(), VecScatter, DACreate(), VecPackCreate(), VecPack
447: S*/
448: typedef struct _p_Sliced* Sliced;
450: EXTERN PetscErrorCode SlicedView(Sliced,PetscViewer);
451: EXTERN PetscErrorCode SlicedCreate(MPI_Comm,Sliced*);
452: EXTERN PetscErrorCode SlicedDestroy(Sliced);
453: EXTERN PetscErrorCode SlicedCreateGlobalVector(Sliced,Vec*);
454: EXTERN PetscErrorCode SlicedGetMatrix(Sliced, MatType,Mat*);
455: EXTERN PetscErrorCode SlicedGetGlobalIndices(Sliced,PetscInt*[]);
456: EXTERN PetscErrorCode SlicedSetPreallocation(Sliced,PetscInt,const PetscInt[],PetscInt,const PetscInt[]);
457: EXTERN PetscErrorCode SlicedSetGhosts(Sliced,PetscInt,PetscInt,PetscInt,const PetscInt[]);
459: /*S
460: DM - Abstract PETSc object that manages an abstract grid object
461:
462: Level: intermediate
464: Concepts: grids, grid refinement
466: Notes: The DA object and the VecPack object are examples of DMs
468: Though the DA objects require the petscsnes.h include files the DM library is
469: NOT dependent on the SNES or KSP library. In fact, the KSP and SNES libraries depend on
470: DM. (This is not great design, but not trivial to fix).
472: .seealso: VecPackCreate(), DA, VecPack
473: S*/
474: typedef struct _p_DM* DM;
476: EXTERN PetscErrorCode DMView(DM,PetscViewer);
477: EXTERN PetscErrorCode DMDestroy(DM);
478: EXTERN PetscErrorCode DMCreateGlobalVector(DM,Vec*);
479: EXTERN PetscErrorCode DMGetColoring(DM,ISColoringType,ISColoring*);
480: EXTERN PetscErrorCode DMGetMatrix(DM, MatType,Mat*);
481: EXTERN PetscErrorCode DMGetInterpolation(DM,DM,Mat*,Vec*);
482: EXTERN PetscErrorCode DMGetInjection(DM,DM,VecScatter*);
483: EXTERN PetscErrorCode DMRefine(DM,MPI_Comm,DM*);
484: EXTERN PetscErrorCode DMCoarsen(DM,MPI_Comm,DM*);
485: EXTERN PetscErrorCode DMRefineHierarchy(DM,PetscInt,DM**);
486: EXTERN PetscErrorCode DMCoarsenHierarchy(DM,PetscInt,DM**);
487: EXTERN PetscErrorCode DMGetInterpolationScale(DM,DM,Mat,Vec*);
489: typedef struct NLF_DAAD* NLF;
491: #include <petscbag.h>
493: EXTERN PetscErrorCode PetscViewerBinaryMatlabOpen(MPI_Comm, const char [], PetscViewer*);
494: EXTERN PetscErrorCode PetscViewerBinaryMatlabDestroy(PetscViewer);
495: EXTERN PetscErrorCode PetscViewerBinaryMatlabOutputBag(PetscViewer, const char [], PetscBag);
496: EXTERN PetscErrorCode PetscViewerBinaryMatlabOutputVec(PetscViewer, const char [], Vec);
497: EXTERN PetscErrorCode PetscViewerBinaryMatlabOutputVecDA(PetscViewer, const char [], Vec, DA);
500: #endif