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 PETSCDM_DLLEXPORT 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 PETSCDM_DLLEXPORT DASetElementType(DA,DAElementType);
85: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetElements(DA,PetscInt *,const PetscInt*[]);
86: EXTERN PetscErrorCode PETSCDM_DLLEXPORT 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 PETSCDM_DLLEXPORT DACreate1d(MPI_Comm,DAPeriodicType,PetscInt,PetscInt,PetscInt,PetscInt*,DA *);
98: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DACreate2d(MPI_Comm,DAPeriodicType,DAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,PetscInt*,DA *);
99: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DACreate3d(MPI_Comm,DAPeriodicType,DAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,PetscInt*,PetscInt*,DA*);
100: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DACreate(MPI_Comm,PetscInt,DAPeriodicType,DAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,PetscInt*,PetscInt*,DA*);
101: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DADestroy(DA);
102: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAView(DA,PetscViewer);
103: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAPrintHelp(DA);
105: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGlobalToLocalBegin(DA,Vec,InsertMode,Vec);
106: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGlobalToLocalEnd(DA,Vec,InsertMode,Vec);
107: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGlobalToNaturalBegin(DA,Vec,InsertMode,Vec);
108: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGlobalToNaturalEnd(DA,Vec,InsertMode,Vec);
109: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DANaturalToGlobalBegin(DA,Vec,InsertMode,Vec);
110: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DANaturalToGlobalEnd(DA,Vec,InsertMode,Vec);
111: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DALocalToLocalBegin(DA,Vec,InsertMode,Vec);
112: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DALocalToLocalEnd(DA,Vec,InsertMode,Vec);
113: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DALocalToGlobal(DA,Vec,InsertMode,Vec);
114: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DALocalToGlobalBegin(DA,Vec,Vec);
115: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DALocalToGlobalEnd(DA,Vec,Vec);
116: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetOwnershipRange(DA,PetscInt **,PetscInt **,PetscInt **);
117: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DACreateGlobalVector(DA,Vec *);
118: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DACreateNaturalVector(DA,Vec *);
119: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DACreateLocalVector(DA,Vec *);
120: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetLocalVector(DA,Vec *);
121: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DARestoreLocalVector(DA,Vec *);
122: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetGlobalVector(DA,Vec *);
123: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DARestoreGlobalVector(DA,Vec *);
124: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DALoad(PetscViewer,PetscInt,PetscInt,PetscInt,DA *);
125: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetCorners(DA,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
126: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetGhostCorners(DA,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
127: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetInfo(DA,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,DAPeriodicType*,DAStencilType*);
128: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetProcessorSubset(DA,DADirection,PetscInt,MPI_Comm*);
129: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DARefine(DA,MPI_Comm,DA*);
131: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGlobalToNaturalAllCreate(DA,VecScatter*);
132: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DANaturalAllToGlobalCreate(DA,VecScatter*);
134: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetGlobalIndices(DA,PetscInt*,PetscInt**);
135: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetISLocalToGlobalMapping(DA,ISLocalToGlobalMapping*);
136: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetISLocalToGlobalMappingBlck(DA,ISLocalToGlobalMapping*);
138: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetScatter(DA,VecScatter*,VecScatter*,VecScatter*);
140: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetAO(DA,AO*);
141: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetCoordinates(DA,Vec);
142: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetCoordinates(DA,Vec *);
143: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetGhostedCoordinates(DA,Vec *);
144: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetCoordinateDA(DA,DA *);
145: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetUniformCoordinates(DA,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal);
146: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetFieldName(DA,PetscInt,const char[]);
147: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetFieldName(DA,PetscInt,char **);
149: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAVecGetArray(DA,Vec,void *);
150: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAVecRestoreArray(DA,Vec,void *);
152: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAVecGetArrayDOF(DA,Vec,void *);
153: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAVecRestoreArrayDOF(DA,Vec,void *);
155: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASplitComm2d(MPI_Comm,PetscInt,PetscInt,PetscInt,MPI_Comm*);
157: /*S
158: SDA - This provides a simplified interface to the DA distributed
159: array object in PETSc. This is intended for people who are
160: NOT using PETSc vectors or objects but just want to distribute
161: simple rectangular arrays amoung a number of procesors and have
162: PETSc handle moving the ghost-values when needed.
164: In certain applications this can serve as a replacement for
165: BlockComm (which is apparently being phased out?).
168: Level: beginner
170: Concepts: simplified distributed array
172: .seealso: SDACreate1d(), SDACreate2d(), SDACreate3d(), SDADestroy(), DA, SDALocalToLocalBegin(),
173: SDALocalToLocalEnd(), SDAGetCorners(), SDAGetGhostCorners()
174: S*/
175: typedef struct _n_SDA* SDA;
177: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SDACreate3d(MPI_Comm,DAPeriodicType,DAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,PetscInt*,PetscInt*,SDA*);
178: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SDACreate2d(MPI_Comm,DAPeriodicType,DAStencilType,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,PetscInt*,SDA*);
179: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SDACreate1d(MPI_Comm,DAPeriodicType,PetscInt,PetscInt,PetscInt,PetscInt*,SDA*);
180: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SDADestroy(SDA);
181: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SDALocalToLocalBegin(SDA,PetscScalar*,InsertMode,PetscScalar*);
182: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SDALocalToLocalEnd(SDA,PetscScalar*,InsertMode,PetscScalar*);
183: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SDAGetCorners(SDA,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
184: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SDAGetGhostCorners(SDA,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
185: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SDAArrayView(SDA,PetscScalar*,PetscViewer);
187: EXTERN PetscErrorCode PETSCDM_DLLEXPORT MatRegisterDAAD(void);
188: EXTERN PetscErrorCode PETSCDM_DLLEXPORT MatCreateDAAD(DA,Mat*);
190: /*S
191: DALocalInfo - C struct that contains information about a structured grid and a processors logical
192: location in it.
194: Level: beginner
196: Concepts: distributed array
198: Developer note: Then entries in this struct are int instead of PetscInt so that the elements may
199: be extracted in Fortran as if from an integer array
201: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DADestroy(), DA, DAGetLocalInfo(), DAGetInfo()
202: S*/
203: typedef struct {
204: PetscInt dim,dof,sw;
205: PetscInt mx,my,mz; /* global number of grid points in each direction */
206: PetscInt xs,ys,zs; /* starting pointd of this processor, excluding ghosts */
207: PetscInt xm,ym,zm; /* number of grid points on this processor, excluding ghosts */
208: PetscInt gxs,gys,gzs; /* starting point of this processor including ghosts */
209: PetscInt gxm,gym,gzm; /* number of grid points on this processor including ghosts */
210: DAPeriodicType pt;
211: DAStencilType st;
212: DA da;
213: } DALocalInfo;
215: /*MC
216: DAForEachPointBegin2d - Starts a loop over the local part of a two dimensional DA
218: Synopsis:
219: void DAForEachPointBegin2d(DALocalInfo *info,PetscInt i,PetscInt j);
220:
221: Level: intermediate
223: .seealso: DAForEachPointEnd2d(), DAVecGetArray()
224: M*/
225: #define DAForEachPointBegin2d(info,i,j) {\
226: PetscInt _xints = info->xs,_xinte = info->xs+info->xm,_yints = info->ys,_yinte = info->ys+info->ym;\
227: for (j=_yints; j<_yinte; j++) {\
228: for (i=_xints; i<_xinte; i++) {\
230: /*MC
231: DAForEachPointEnd2d - Ends a loop over the local part of a two dimensional DA
233: Synopsis:
234: void DAForEachPointEnd2d;
235:
236: Level: intermediate
238: .seealso: DAForEachPointBegin2d(), DAVecGetArray()
239: M*/
240: #define DAForEachPointEnd2d }}}
242: /*MC
243: DACoor2d - Structure for holding 2d (x and y) coordinates.
245: Level: intermediate
247: Sample Usage:
248: DACoor2d **coors;
249: Vec vcoors;
250: DA cda;
252: DAGetCoordinates(da,&vcoors);
253: DAGetCoordinateDA(da,&cda);
254: DAVecGetArray(dac,vcoors,&coors);
255: DAGetCorners(dac,&mstart,&nstart,0,&m,&n,0)
256: for (i=mstart; i<mstart+m; i++) {
257: for (j=nstart; j<nstart+n; j++) {
258: x = coors[j][i].x;
259: y = coors[j][i].y;
260: ......
261: }
262: }
263: DAVecRestoreArray(dac,vcoors,&coors);
265: .seealso: DACoor3d, DAForEachPointBegin(), DAGetCoordinateDA(), DAGetCoordinates(), DAGetGhostCoordinates()
266: M*/
267: typedef struct {PetscScalar x,y;} DACoor2d;
269: /*MC
270: DACoor3d - Structure for holding 3d (x, y and z) coordinates.
272: Level: intermediate
274: Sample Usage:
275: DACoor3d **coors;
276: Vec vcoors;
277: DA cda;
279: DAGetCoordinates(da,&vcoors);
280: DAGetCoordinateDA(da,&cda);
281: DAVecGetArray(dac,vcoors,&coors);
282: DAGetCorners(dac,&mstart,&nstart,&pstart,&m,&n,&p)
283: for (i=mstart; i<mstart+m; i++) {
284: for (j=nstart; j<nstart+n; j++) {
285: for (k=pstart; k<pstart+p; k++) {
286: x = coors[k][j][i].x;
287: y = coors[k][j][i].y;
288: z = coors[k][j][i].z;
289: ......
290: }
291: }
292: DAVecRestoreArray(dac,vcoors,&coors);
294: .seealso: DACoor2d, DAForEachPointBegin(), DAGetCoordinateDA(), DAGetCoordinates(), DAGetGhostCoordinates()
295: M*/
296: typedef struct {PetscScalar x,y,z;} DACoor3d;
297:
298: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetLocalInfo(DA,DALocalInfo*);
299: typedef PetscErrorCode (*DALocalFunction1)(DALocalInfo*,void*,void*,void*);
300: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAFormFunction1(DA,Vec,Vec,void*);
301: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAFormFunctioni1(DA,PetscInt,Vec,PetscScalar*,void*);
302: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAFormFunctionib1(DA,PetscInt,Vec,PetscScalar*,void*);
303: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAComputeJacobian1WithAdic(DA,Vec,Mat,void*);
304: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAComputeJacobian1WithAdifor(DA,Vec,Mat,void*);
305: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAMultiplyByJacobian1WithAdic(DA,Vec,Vec,Vec,void*);
306: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAMultiplyByJacobian1WithAdifor(DA,Vec,Vec,Vec,void*);
307: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAMultiplyByJacobian1WithAD(DA,Vec,Vec,Vec,void*);
308: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAComputeJacobian1(DA,Vec,Mat,void*);
309: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetLocalFunction(DA,DALocalFunction1*);
310: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetLocalFunction(DA,DALocalFunction1);
311: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetLocalFunctioni(DA,PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*));
312: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetLocalFunctionib(DA,PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*));
313: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetLocalJacobian(DA,DALocalFunction1);
314: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetLocalAdicFunction_Private(DA,DALocalFunction1);
316: /*MC
317: DASetLocalAdicFunction - Caches in a DA a local function computed by ADIC/ADIFOR
319: Collective on DA
321: Synopsis:
322: PetscErrorCode DASetLocalAdicFunction(DA da,DALocalFunction1 ad_lf)
323:
324: Input Parameter:
325: + da - initial distributed array
326: - ad_lf - the local function as computed by ADIC/ADIFOR
328: Level: intermediate
330: .keywords: distributed array, refine
332: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DADestroy(), DAGetLocalFunction(), DASetLocalFunction(),
333: DASetLocalJacobian()
334: M*/
335: #if defined(PETSC_HAVE_ADIC)
336: # define DASetLocalAdicFunction(a,d) DASetLocalAdicFunction_Private(a,(DALocalFunction1)d)
337: #else
338: # define DASetLocalAdicFunction(a,d) DASetLocalAdicFunction_Private(a,0)
339: #endif
341: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetLocalAdicMFFunction_Private(DA,DALocalFunction1);
342: #if defined(PETSC_HAVE_ADIC)
343: # define DASetLocalAdicMFFunction(a,d) DASetLocalAdicMFFunction_Private(a,(DALocalFunction1)d)
344: #else
345: # define DASetLocalAdicMFFunction(a,d) DASetLocalAdicMFFunction_Private(a,0)
346: #endif
347: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetLocalAdicFunctioni_Private(DA,PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
348: #if defined(PETSC_HAVE_ADIC)
349: # define DASetLocalAdicFunctioni(a,d) DASetLocalAdicFunctioni_Private(a,(PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
350: #else
351: # define DASetLocalAdicFunctioni(a,d) DASetLocalAdicFunctioni_Private(a,0)
352: #endif
353: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetLocalAdicMFFunctioni_Private(DA,PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
354: #if defined(PETSC_HAVE_ADIC)
355: # define DASetLocalAdicMFFunctioni(a,d) DASetLocalAdicMFFunctioni_Private(a,(PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
356: #else
357: # define DASetLocalAdicMFFunctioni(a,d) DASetLocalAdicMFFunctioni_Private(a,0)
358: #endif
360: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetLocalAdicFunctionib_Private(DA,PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
361: #if defined(PETSC_HAVE_ADIC)
362: # define DASetLocalAdicFunctionib(a,d) DASetLocalAdicFunctionib_Private(a,(PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
363: #else
364: # define DASetLocalAdicFunctionib(a,d) DASetLocalAdicFunctionib_Private(a,0)
365: #endif
366: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetLocalAdicMFFunctionib_Private(DA,PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*));
367: #if defined(PETSC_HAVE_ADIC)
368: # define DASetLocalAdicMFFunctionib(a,d) DASetLocalAdicMFFunctionib_Private(a,(PetscErrorCode (*)(DALocalInfo*,MatStencil*,void*,void*,void*))d)
369: #else
370: # define DASetLocalAdicMFFunctionib(a,d) DASetLocalAdicMFFunctionib_Private(a,0)
371: #endif
373: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAFormFunctioniTest1(DA,void*);
375: #include petscmat.h
376: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetColoring(DA,ISColoringType,ISColoring *);
377: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetMatrix(DA, MatType,Mat *);
378: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetGetMatrix(DA,PetscErrorCode (*)(DA, MatType,Mat *));
379: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetInterpolation(DA,DA,Mat*,Vec*);
380: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetInjection(DA,DA,VecScatter*);
381: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetBlockFills(DA,PetscInt*,PetscInt*);
382: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetMatPreallocateOnly(DA,PetscTruth);
383: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DASetRefinementFactor(DA,PetscInt,PetscInt,PetscInt);
384: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetRefinementFactor(DA,PetscInt*,PetscInt*,PetscInt*);
386: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetAdicArray(DA,PetscTruth,void**,void**,PetscInt*);
387: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DARestoreAdicArray(DA,PetscTruth,void**,void**,PetscInt*);
388: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetAdicMFArray(DA,PetscTruth,void**,void**,PetscInt*);
389: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetAdicMFArray4(DA,PetscTruth,void**,void**,PetscInt*);
390: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetAdicMFArray9(DA,PetscTruth,void**,void**,PetscInt*);
391: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetAdicMFArrayb(DA,PetscTruth,void**,void**,PetscInt*);
392: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DARestoreAdicMFArray(DA,PetscTruth,void**,void**,PetscInt*);
393: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DAGetArray(DA,PetscTruth,void**);
394: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DARestoreArray(DA,PetscTruth,void**);
395: EXTERN PetscErrorCode PETSCDM_DLLEXPORT ad_DAGetArray(DA,PetscTruth,void**);
396: EXTERN PetscErrorCode PETSCDM_DLLEXPORT ad_DARestoreArray(DA,PetscTruth,void**);
397: EXTERN PetscErrorCode PETSCDM_DLLEXPORT admf_DAGetArray(DA,PetscTruth,void**);
398: EXTERN PetscErrorCode PETSCDM_DLLEXPORT admf_DARestoreArray(DA,PetscTruth,void**);
400: #include petscpf.h
401: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DACreatePF(DA,PF*);
403: /*S
404: VecPack - Abstract PETSc object that manages treating several distinct vectors as if they
405: were one. The VecPack routines allow one to manage a nonlinear solver that works on a
406: vector that consists of several distinct parts. This is mostly used for LNKS solvers,
407: that is design optimization problems that are written as a nonlinear system
409: Level: beginner
411: Concepts: multi-component, LNKS solvers
413: .seealso: VecPackCreate(), VecPackDestroy(), DM
414: S*/
415: typedef struct _p_VecPack* VecPack;
417: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackCreate(MPI_Comm,VecPack*);
418: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackDestroy(VecPack);
419: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackAddArray(VecPack,PetscInt);
420: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackAddDA(VecPack,DA);
421: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackAddVecScatter(VecPack,VecScatter);
422: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackScatter(VecPack,Vec,...);
423: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackGather(VecPack,Vec,...);
424: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackGetAccess(VecPack,Vec,...);
425: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackRestoreAccess(VecPack,Vec,...);
426: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackGetLocalVectors(VecPack,...);
427: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackGetEntries(VecPack,...);
428: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackRestoreLocalVectors(VecPack,...);
429: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackCreateGlobalVector(VecPack,Vec*);
430: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackGetGlobalIndices(VecPack,...);
431: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackRefine(VecPack,MPI_Comm,VecPack*);
432: EXTERN PetscErrorCode PETSCDM_DLLEXPORT VecPackGetInterpolation(VecPack,VecPack,Mat*,Vec*);
434: /*S
435: Slice - Abstract PETSc object that manages distributed field data for a simple unstructured matrix
437: Level: beginner
439: Concepts: distributed array
441: .seealso: DACreate1d(), DACreate2d(), DACreate3d(), DADestroy(), VecScatter, DACreate(), VecPackCreate(), VecPack
442: S*/
443: typedef struct _p_Sliced* Sliced;
445: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SlicedView(Sliced,PetscViewer);
446: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SlicedCreate(MPI_Comm,Sliced*);
447: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SlicedDestroy(Sliced);
448: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SlicedCreateGlobalVector(Sliced,Vec*);
449: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SlicedGetMatrix(Sliced, MatType,Mat*);
450: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SlicedGetGlobalIndices(Sliced,PetscInt*[]);
451: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SlicedSetPreallocation(Sliced,PetscInt,const PetscInt[],PetscInt,const PetscInt[]);
452: EXTERN PetscErrorCode PETSCDM_DLLEXPORT SlicedSetGhosts(Sliced,PetscInt,PetscInt,PetscInt,const PetscInt[]);
454: /*S
455: DM - Abstract PETSc object that manages an abstract grid object
456:
457: Level: intermediate
459: Concepts: grids, grid refinement
461: Notes: The DA object and the VecPack object are examples of DMs
463: Though the DA objects require the petscsnes.h include files the DM library is
464: NOT dependent on the SNES or KSP library. In fact, the KSP and SNES libraries depend on
465: DM. (This is not great design, but not trivial to fix).
467: .seealso: VecPackCreate(), DA, VecPack
468: S*/
469: typedef struct _p_DM* DM;
471: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DMView(DM,PetscViewer);
472: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DMDestroy(DM);
473: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DMCreateGlobalVector(DM,Vec*);
474: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DMGetColoring(DM,ISColoringType,ISColoring*);
475: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DMGetMatrix(DM, MatType,Mat*);
476: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DMGetInterpolation(DM,DM,Mat*,Vec*);
477: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DMGetInjection(DM,DM,VecScatter*);
478: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DMRefine(DM,MPI_Comm,DM*);
479: EXTERN PetscErrorCode PETSCDM_DLLEXPORT DMGetInterpolationScale(DM,DM,Mat,Vec*);
481: typedef struct NLF_DAAD* NLF;
485: #endif