Actual source code: dmimpl.h
petsc-3.10.5 2019-03-28
3: #if !defined(_DMIMPL_H)
4: #define _DMIMPL_H
6: #include <petscdm.h>
7: #include <petsc/private/petscimpl.h>
8: #include <petsc/private/petscdsimpl.h>
9: #include <petsc/private/isimpl.h>
11: PETSC_EXTERN PetscBool DMRegisterAllCalled;
12: PETSC_EXTERN PetscErrorCode DMRegisterAll(void);
13: typedef PetscErrorCode (*NullSpaceFunc)(DM dm, PetscInt field, MatNullSpace *nullSpace);
15: typedef struct _DMOps *DMOps;
16: struct _DMOps {
17: PetscErrorCode (*view)(DM,PetscViewer);
18: PetscErrorCode (*load)(DM,PetscViewer);
19: PetscErrorCode (*clone)(DM,DM*);
20: PetscErrorCode (*setfromoptions)(PetscOptionItems*,DM);
21: PetscErrorCode (*setup)(DM);
22: PetscErrorCode (*createdefaultsection)(DM);
23: PetscErrorCode (*createdefaultconstraints)(DM);
24: PetscErrorCode (*createglobalvector)(DM,Vec*);
25: PetscErrorCode (*createlocalvector)(DM,Vec*);
26: PetscErrorCode (*getlocaltoglobalmapping)(DM);
27: PetscErrorCode (*createfieldis)(DM,PetscInt*,char***,IS**);
28: PetscErrorCode (*createcoordinatedm)(DM,DM*);
29: PetscErrorCode (*createcoordinatefield)(DM,DMField*);
31: PetscErrorCode (*getcoloring)(DM,ISColoringType,ISColoring*);
32: PetscErrorCode (*creatematrix)(DM, Mat*);
33: PetscErrorCode (*createinterpolation)(DM,DM,Mat*,Vec*);
34: PetscErrorCode (*createrestriction)(DM,DM,Mat*);
35: PetscErrorCode (*createmassmatrix)(DM,DM,Mat*);
36: PetscErrorCode (*getaggregates)(DM,DM,Mat*);
37: PetscErrorCode (*hascreateinjection)(DM,PetscBool*);
38: PetscErrorCode (*getinjection)(DM,DM,Mat*);
40: PetscErrorCode (*refine)(DM,MPI_Comm,DM*);
41: PetscErrorCode (*coarsen)(DM,MPI_Comm,DM*);
42: PetscErrorCode (*refinehierarchy)(DM,PetscInt,DM*);
43: PetscErrorCode (*coarsenhierarchy)(DM,PetscInt,DM*);
44: PetscErrorCode (*adaptlabel)(DM,DMLabel,DM*);
45: PetscErrorCode (*adaptmetric)(DM,Vec,DMLabel,DM*);
47: PetscErrorCode (*globaltolocalbegin)(DM,Vec,InsertMode,Vec);
48: PetscErrorCode (*globaltolocalend)(DM,Vec,InsertMode,Vec);
49: PetscErrorCode (*localtoglobalbegin)(DM,Vec,InsertMode,Vec);
50: PetscErrorCode (*localtoglobalend)(DM,Vec,InsertMode,Vec);
51: PetscErrorCode (*localtolocalbegin)(DM,Vec,InsertMode,Vec);
52: PetscErrorCode (*localtolocalend)(DM,Vec,InsertMode,Vec);
54: PetscErrorCode (*destroy)(DM);
56: PetscErrorCode (*computevariablebounds)(DM,Vec,Vec);
58: PetscErrorCode (*createsubdm)(DM,PetscInt,const PetscInt*,IS*,DM*);
59: PetscErrorCode (*createsuperdm)(DM*,PetscInt,IS**,DM*);
60: PetscErrorCode (*createfielddecomposition)(DM,PetscInt*,char***,IS**,DM**);
61: PetscErrorCode (*createdomaindecomposition)(DM,PetscInt*,char***,IS**,IS**,DM**);
62: PetscErrorCode (*createddscatters)(DM,PetscInt,DM*,VecScatter**,VecScatter**,VecScatter**);
64: PetscErrorCode (*getdimpoints)(DM,PetscInt,PetscInt*,PetscInt*);
65: PetscErrorCode (*locatepoints)(DM,Vec,DMPointLocationType,PetscSF);
66: PetscErrorCode (*getneighbors)(DM,PetscInt*,const PetscMPIInt**);
67: PetscErrorCode (*getboundingbox)(DM,PetscReal*,PetscReal*);
68: PetscErrorCode (*getlocalboundingbox)(DM,PetscReal*,PetscReal*);
69: PetscErrorCode (*locatepointssubdomain)(DM,Vec,PetscMPIInt**);
71: PetscErrorCode (*projectfunctionlocal)(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void **,InsertMode,Vec);
72: PetscErrorCode (*projectfunctionlabellocal)(DM,PetscReal,DMLabel,PetscInt,const PetscInt[],PetscInt,const PetscInt[],PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void **,InsertMode,Vec);
73: PetscErrorCode (*projectfieldlocal)(DM,PetscReal,Vec,void(**)(PetscInt,PetscInt,PetscInt,const PetscInt[],const PetscInt[],const PetscScalar[],const PetscScalar[],const PetscScalar[],const PetscInt[],const PetscInt[],const PetscScalar[],const PetscScalar[],const PetscScalar[],PetscReal,const PetscReal[],PetscInt,const PetscScalar[],PetscScalar[]),InsertMode,Vec);
74: PetscErrorCode (*projectfieldlabellocal)(DM,PetscReal,DMLabel,PetscInt,const PetscInt[],PetscInt,const PetscInt[],Vec,void(**)(PetscInt,PetscInt,PetscInt,const PetscInt[],const PetscInt[],const PetscScalar[],const PetscScalar[],const PetscScalar[],const PetscInt[],const PetscInt[],const PetscScalar[],const PetscScalar[],const PetscScalar[],PetscReal,const PetscReal[],PetscInt,const PetscScalar[],PetscScalar[]),InsertMode,Vec);
75: PetscErrorCode (*computel2diff)(DM,PetscReal,PetscErrorCode(**)(PetscInt, PetscReal,const PetscReal [], PetscInt, PetscScalar *, void *), void **, Vec, PetscReal *);
76: PetscErrorCode (*computel2gradientdiff)(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal [],const PetscReal[],PetscInt, PetscScalar *,void *),void **,Vec,const PetscReal[],PetscReal *);
77: PetscErrorCode (*computel2fielddiff)(DM,PetscReal,PetscErrorCode(**)(PetscInt, PetscReal,const PetscReal [], PetscInt, PetscScalar *, void *), void **, Vec, PetscReal *);
79: PetscErrorCode (*getcompatibility)(DM,DM,PetscBool*,PetscBool*);
80: };
82: PETSC_EXTERN PetscErrorCode DMLocalizeCoordinate_Internal(DM, PetscInt, const PetscScalar[], const PetscScalar[], PetscScalar[]);
83: PETSC_EXTERN PetscErrorCode DMLocalizeCoordinateReal_Internal(DM, PetscInt, const PetscReal[], const PetscReal[], PetscReal[]);
84: PETSC_EXTERN PetscErrorCode DMLocalizeAddCoordinate_Internal(DM, PetscInt, const PetscScalar[], const PetscScalar[], PetscScalar[]);
86: typedef struct _DMCoarsenHookLink *DMCoarsenHookLink;
87: struct _DMCoarsenHookLink {
88: PetscErrorCode (*coarsenhook)(DM,DM,void*); /* Run once, when coarse DM is created */
89: PetscErrorCode (*restricthook)(DM,Mat,Vec,Mat,DM,void*); /* Run each time a new problem is restricted to a coarse grid */
90: void *ctx;
91: DMCoarsenHookLink next;
92: };
94: typedef struct _DMRefineHookLink *DMRefineHookLink;
95: struct _DMRefineHookLink {
96: PetscErrorCode (*refinehook)(DM,DM,void*); /* Run once, when a fine DM is created */
97: PetscErrorCode (*interphook)(DM,Mat,DM,void*); /* Run each time a new problem is interpolated to a fine grid */
98: void *ctx;
99: DMRefineHookLink next;
100: };
102: typedef struct _DMSubDomainHookLink *DMSubDomainHookLink;
103: struct _DMSubDomainHookLink {
104: PetscErrorCode (*ddhook)(DM,DM,void*);
105: PetscErrorCode (*restricthook)(DM,VecScatter,VecScatter,DM,void*);
106: void *ctx;
107: DMSubDomainHookLink next;
108: };
110: typedef struct _DMGlobalToLocalHookLink *DMGlobalToLocalHookLink;
111: struct _DMGlobalToLocalHookLink {
112: PetscErrorCode (*beginhook)(DM,Vec,InsertMode,Vec,void*);
113: PetscErrorCode (*endhook)(DM,Vec,InsertMode,Vec,void*);
114: void *ctx;
115: DMGlobalToLocalHookLink next;
116: };
118: typedef struct _DMLocalToGlobalHookLink *DMLocalToGlobalHookLink;
119: struct _DMLocalToGlobalHookLink {
120: PetscErrorCode (*beginhook)(DM,Vec,InsertMode,Vec,void*);
121: PetscErrorCode (*endhook)(DM,Vec,InsertMode,Vec,void*);
122: void *ctx;
123: DMLocalToGlobalHookLink next;
124: };
126: typedef enum {DMVEC_STATUS_IN,DMVEC_STATUS_OUT} DMVecStatus;
127: typedef struct _DMNamedVecLink *DMNamedVecLink;
128: struct _DMNamedVecLink {
129: Vec X;
130: char *name;
131: DMVecStatus status;
132: DMNamedVecLink next;
133: };
135: typedef struct _DMWorkLink *DMWorkLink;
136: struct _DMWorkLink {
137: size_t bytes;
138: void *mem;
139: DMWorkLink next;
140: };
142: #define DM_MAX_WORK_VECTORS 100 /* work vectors available to users via DMGetGlobalVector(), DMGetLocalVector() */
144: struct _n_DMLabelLink {
145: DMLabel label;
146: PetscBool output;
147: struct _n_DMLabelLink *next;
148: };
149: typedef struct _n_DMLabelLink *DMLabelLink;
151: struct _n_DMLabelLinkList {
152: PetscInt refct;
153: DMLabelLink next;
154: };
155: typedef struct _n_DMLabelLinkList *DMLabelLinkList;
157: typedef struct _n_Boundary *DMBoundary;
159: struct _n_Boundary {
160: DSBoundary dsboundary;
161: DMLabel label;
162: DMBoundary next;
163: };
165: PETSC_EXTERN PetscErrorCode DMDestroyLabelLinkList(DM);
167: struct _p_DM {
168: PETSCHEADER(struct _DMOps);
169: Vec localin[DM_MAX_WORK_VECTORS],localout[DM_MAX_WORK_VECTORS];
170: Vec globalin[DM_MAX_WORK_VECTORS],globalout[DM_MAX_WORK_VECTORS];
171: DMNamedVecLink namedglobal;
172: DMNamedVecLink namedlocal;
173: DMWorkLink workin,workout;
174: DMLabelLinkList labels; /* Linked list of labels */
175: DMLabel depthLabel; /* Optimized access to depth label */
176: void *ctx; /* a user context */
177: PetscErrorCode (*ctxdestroy)(void**);
178: Vec x; /* location at which the functions/Jacobian are computed */
179: ISColoringType coloringtype;
180: MatFDColoring fd;
181: VecType vectype; /* type of vector created with DMCreateLocalVector() and DMCreateGlobalVector() */
182: MatType mattype; /* type of matrix created with DMCreateMatrix() */
183: PetscInt bs;
184: ISLocalToGlobalMapping ltogmap;
185: PetscBool prealloc_only; /* Flag indicating the DMCreateMatrix() should only preallocate, not fill the matrix */
186: PetscBool structure_only; /* Flag indicating the DMCreateMatrix() create matrix structure without values */
187: PetscInt levelup,leveldown; /* if the DM has been obtained by refining (or coarsening) this indicates how many times that process has been used to generate this DM */
188: PetscBool setupcalled; /* Indicates that the DM has been set up, methods that modify a DM such that a fresh setup is required should reset this flag */
189: void *data;
190: /* Hierarchy / Submeshes */
191: DM coarseMesh;
192: DM fineMesh;
193: DMCoarsenHookLink coarsenhook; /* For transfering auxiliary problem data to coarser grids */
194: DMRefineHookLink refinehook;
195: DMSubDomainHookLink subdomainhook;
196: DMGlobalToLocalHookLink gtolhook;
197: DMLocalToGlobalHookLink ltoghook;
198: /* Topology */
199: PetscInt dim; /* The topological dimension */
200: /* Flexible communication */
201: PetscSF sfMigration; /* SF for point distribution created during distribution */
202: PetscSF sf; /* SF for parallel point overlap */
203: PetscSF defaultSF; /* SF for parallel dof overlap using default section */
204: PetscSF sfNatural; /* SF mapping to the "natural" ordering */
205: PetscBool useNatural; /* Create the natural SF */
206: /* Allows a non-standard data layout */
207: PetscSection defaultSection; /* Layout for local vectors */
208: PetscSection defaultGlobalSection; /* Layout for global vectors */
209: PetscLayout map;
210: /* Constraints */
211: PetscSection defaultConstraintSection;
212: Mat defaultConstraintMat;
213: /* Coordinates */
214: PetscInt dimEmbed; /* The dimension of the embedding space */
215: DM coordinateDM; /* Layout for coordinates (default section) */
216: Vec coordinates; /* Coordinate values in global vector */
217: Vec coordinatesLocal; /* Coordinate values in local vector */
218: PetscBool periodic; /* Is the DM periodic? */
219: DMField coordinateField; /* Coordinates as an abstract field */
220: PetscReal *L, *maxCell; /* Size of periodic box and max cell size for determining periodicity */
221: DMBoundaryType *bdtype; /* Indicates type of topological boundary */
222: /* Null spaces -- of course I should make this have a variable number of fields */
223: /* I now believe this might not be the right way: see below */
224: NullSpaceFunc nullspaceConstructors[10];
225: /* Fields are represented by objects */
226: PetscDS prob;
227: DMBoundary boundary; /* List of boundary conditions */
228: /* Output structures */
229: DM dmBC; /* The DM with boundary conditions in the global DM */
230: PetscInt outputSequenceNum; /* The current sequence number for output */
231: PetscReal outputSequenceVal; /* The current sequence value for output */
233: PetscObject dmksp,dmsnes,dmts;
234: };
236: PETSC_EXTERN PetscLogEvent DM_Convert;
237: PETSC_EXTERN PetscLogEvent DM_GlobalToLocal;
238: PETSC_EXTERN PetscLogEvent DM_LocalToGlobal;
239: PETSC_EXTERN PetscLogEvent DM_LocatePoints;
240: PETSC_EXTERN PetscLogEvent DM_Coarsen;
241: PETSC_EXTERN PetscLogEvent DM_Refine;
242: PETSC_EXTERN PetscLogEvent DM_CreateInterpolation;
243: PETSC_EXTERN PetscLogEvent DM_CreateRestriction;
245: PETSC_EXTERN PetscErrorCode DMCreateGlobalVector_Section_Private(DM,Vec*);
246: PETSC_EXTERN PetscErrorCode DMCreateLocalVector_Section_Private(DM,Vec*);
247: PETSC_EXTERN PetscErrorCode DMCreateSubDM_Section_Private(DM,PetscInt,const PetscInt[],IS*,DM*);
248: PETSC_EXTERN PetscErrorCode DMCreateSuperDM_Section_Private(DM[],PetscInt,IS**,DM*);
250: /*
252: Composite Vectors
254: Single global representation
255: Individual global representations
256: Single local representation
257: Individual local representations
259: Subsets of individual as a single????? Do we handle this by having DMComposite inside composite??????
261: DM da_u, da_v, da_p
263: DM dm_velocities
265: DM dm
267: DMDACreate(,&da_u);
268: DMDACreate(,&da_v);
269: DMCompositeCreate(,&dm_velocities);
270: DMCompositeAddDM(dm_velocities,(DM)du);
271: DMCompositeAddDM(dm_velocities,(DM)dv);
273: DMDACreate(,&da_p);
274: DMCompositeCreate(,&dm_velocities);
275: DMCompositeAddDM(dm,(DM)dm_velocities);
276: DMCompositeAddDM(dm,(DM)dm_p);
279: Access parts of composite vectors (Composite only)
280: ---------------------------------
281: DMCompositeGetAccess - access the global vector as subvectors and array (for redundant arrays)
282: ADD for local vector -
284: Element access
285: --------------
286: From global vectors
287: -DAVecGetArray - for DMDA
288: -VecGetArray - for DMSliced
289: ADD for DMComposite??? maybe
291: From individual vector
292: -DAVecGetArray - for DMDA
293: -VecGetArray -for sliced
294: ADD for DMComposite??? maybe
296: From single local vector
297: ADD * single local vector as arrays?
299: Communication
300: -------------
301: DMGlobalToLocal - global vector to single local vector
303: DMCompositeScatter/Gather - direct to individual local vectors and arrays CHANGE name closer to GlobalToLocal?
305: Obtaining vectors
306: -----------------
307: DMCreateGlobal/Local
308: DMGetGlobal/Local
309: DMCompositeGetLocalVectors - gives individual local work vectors and arrays
312: ????? individual global vectors ????
314: */
316: #if defined(PETSC_HAVE_HDF5)
317: PETSC_EXTERN PetscErrorCode DMSequenceLoad_HDF5_Internal(DM, const char *, PetscInt, PetscScalar *, PetscViewer);
318: #endif
320: PETSC_STATIC_INLINE PetscErrorCode DMGetLocalOffset_Private(DM dm, PetscInt point, PetscInt *start, PetscInt *end)
321: {
323: #if defined(PETSC_USE_DEBUG)
324: {
325: PetscInt dof;
327: *start = *end = 0; /* Silence overzealous compiler warning */
328: if (!dm->defaultSection) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "DM must have a default Section, see DMSetDefaultSection()");
329: PetscSectionGetOffset(dm->defaultSection, point, start);
330: PetscSectionGetDof(dm->defaultSection, point, &dof);
331: *end = *start + dof;
332: }
333: #else
334: {
335: const PetscSection s = dm->defaultSection;
336: *start = s->atlasOff[point - s->pStart];
337: *end = *start + s->atlasDof[point - s->pStart];
338: }
339: #endif
340: return(0);
341: }
343: PETSC_STATIC_INLINE PetscErrorCode DMGetLocalFieldOffset_Private(DM dm, PetscInt point, PetscInt field, PetscInt *start, PetscInt *end)
344: {
346: #if defined(PETSC_USE_DEBUG)
347: {
348: PetscInt dof;
350: *start = *end = 0; /* Silence overzealous compiler warning */
351: if (!dm->defaultSection) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "DM must have a default Section, see DMSetDefaultSection()");
352: PetscSectionGetFieldOffset(dm->defaultSection, point, field, start);
353: PetscSectionGetFieldDof(dm->defaultSection, point, field, &dof);
354: *end = *start + dof;
355: }
356: #else
357: {
358: const PetscSection s = dm->defaultSection->field[field];
359: *start = s->atlasOff[point - s->pStart];
360: *end = *start + s->atlasDof[point - s->pStart];
361: }
362: #endif
363: return(0);
364: }
366: PETSC_STATIC_INLINE PetscErrorCode DMGetGlobalOffset_Private(DM dm, PetscInt point, PetscInt *start, PetscInt *end)
367: {
369: #if defined(PETSC_USE_DEBUG)
370: {
372: PetscInt dof,cdof;
373: *start = *end = 0; /* Silence overzealous compiler warning */
374: if (!dm->defaultSection) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "DM must have a default Section, see DMSetDefaultSection()");
375: if (!dm->defaultGlobalSection) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "DM must have a default global Section. It will be created automatically by DMGetDefaultGlobalSection()");
376: PetscSectionGetOffset(dm->defaultGlobalSection, point, start);
377: PetscSectionGetDof(dm->defaultGlobalSection, point, &dof);
378: PetscSectionGetConstraintDof(dm->defaultGlobalSection, point, &cdof);
379: *end = *start + dof - cdof + (dof < 0 ? 1 : 0);
380: }
381: #else
382: {
383: const PetscSection s = dm->defaultGlobalSection;
384: const PetscInt dof = s->atlasDof[point - s->pStart];
385: const PetscInt cdof = s->bc ? s->bc->atlasDof[point - s->bc->pStart] : 0;
386: *start = s->atlasOff[point - s->pStart];
387: *end = *start + dof - cdof + (dof < 0 ? 1 : 0);
388: }
389: #endif
390: return(0);
391: }
393: PETSC_STATIC_INLINE PetscErrorCode DMGetGlobalFieldOffset_Private(DM dm, PetscInt point, PetscInt field, PetscInt *start, PetscInt *end)
394: {
396: #if defined(PETSC_USE_DEBUG)
397: {
398: PetscInt loff, lfoff, fdof, fcdof, ffcdof, f;
400: *start = *end = 0; /* Silence overzealous compiler warning */
401: if (!dm->defaultSection) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "DM must have a default Section, see DMSetDefaultSection()");
402: if (!dm->defaultGlobalSection) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "DM must have a default global Section. It will be crated automatically by DMGetDefaultGlobalSection()");
403: PetscSectionGetOffset(dm->defaultGlobalSection, point, start);
404: PetscSectionGetOffset(dm->defaultSection, point, &loff);
405: PetscSectionGetFieldOffset(dm->defaultSection, point, field, &lfoff);
406: PetscSectionGetFieldDof(dm->defaultSection, point, field, &fdof);
407: PetscSectionGetFieldConstraintDof(dm->defaultSection, point, field, &fcdof);
408: *start = *start < 0 ? *start - (lfoff-loff) : *start + lfoff-loff;
409: for (f = 0; f < field; ++f) {
410: PetscSectionGetFieldConstraintDof(dm->defaultSection, point, f, &ffcdof);
411: *start = *start < 0 ? *start + ffcdof : *start - ffcdof;
412: }
413: *end = *start < 0 ? *start - (fdof-fcdof) : *start + fdof-fcdof;
414: }
415: #else
416: {
417: const PetscSection s = dm->defaultSection;
418: const PetscSection fs = dm->defaultSection->field[field];
419: const PetscSection gs = dm->defaultGlobalSection;
420: const PetscInt loff = s->atlasOff[point - s->pStart];
421: const PetscInt goff = gs->atlasOff[point - s->pStart];
422: const PetscInt lfoff = fs->atlasOff[point - s->pStart];
423: const PetscInt fdof = fs->atlasDof[point - s->pStart];
424: const PetscInt fcdof = fs->bc ? fs->bc->atlasDof[point - fs->bc->pStart] : 0;
425: PetscInt ffcdof = 0, f;
427: for (f = 0; f < field; ++f) {
428: const PetscSection ffs = dm->defaultSection->field[f];
429: ffcdof += ffs->bc ? ffs->bc->atlasDof[point - ffs->bc->pStart] : 0;
430: }
431: *start = goff + (goff < 0 ? loff-lfoff + ffcdof : lfoff-loff - ffcdof);
432: *end = *start < 0 ? *start - (fdof-fcdof) : *start + fdof-fcdof;
433: }
434: #endif
435: return(0);
436: }
438: #endif