Actual source code: dmcoordinates.c
1: #include <petsc/private/dmimpl.h>
3: #include <petscdmplex.h>
4: #include <petscsf.h>
6: PetscErrorCode DMRestrictHook_Coordinates(DM dm, DM dmc, void *ctx)
7: {
8: DM dm_coord, dmc_coord;
9: Vec coords, ccoords;
10: Mat inject;
11: PetscFunctionBegin;
12: PetscCall(DMGetCoordinateDM(dm, &dm_coord));
13: PetscCall(DMGetCoordinateDM(dmc, &dmc_coord));
14: PetscCall(DMGetCoordinates(dm, &coords));
15: PetscCall(DMGetCoordinates(dmc, &ccoords));
16: if (coords && !ccoords) {
17: PetscCall(DMCreateGlobalVector(dmc_coord, &ccoords));
18: PetscCall(PetscObjectSetName((PetscObject)ccoords, "coordinates"));
19: PetscCall(DMCreateInjection(dmc_coord, dm_coord, &inject));
20: PetscCall(MatRestrict(inject, coords, ccoords));
21: PetscCall(MatDestroy(&inject));
22: PetscCall(DMSetCoordinates(dmc, ccoords));
23: PetscCall(VecDestroy(&ccoords));
24: }
25: PetscFunctionReturn(PETSC_SUCCESS);
26: }
28: static PetscErrorCode DMSubDomainHook_Coordinates(DM dm, DM subdm, void *ctx)
29: {
30: DM dm_coord, subdm_coord;
31: Vec coords, ccoords, clcoords;
32: VecScatter *scat_i, *scat_g;
33: PetscFunctionBegin;
34: PetscCall(DMGetCoordinateDM(dm, &dm_coord));
35: PetscCall(DMGetCoordinateDM(subdm, &subdm_coord));
36: PetscCall(DMGetCoordinates(dm, &coords));
37: PetscCall(DMGetCoordinates(subdm, &ccoords));
38: if (coords && !ccoords) {
39: PetscCall(DMCreateGlobalVector(subdm_coord, &ccoords));
40: PetscCall(PetscObjectSetName((PetscObject)ccoords, "coordinates"));
41: PetscCall(DMCreateLocalVector(subdm_coord, &clcoords));
42: PetscCall(PetscObjectSetName((PetscObject)clcoords, "coordinates"));
43: PetscCall(DMCreateDomainDecompositionScatters(dm_coord, 1, &subdm_coord, NULL, &scat_i, &scat_g));
44: PetscCall(VecScatterBegin(scat_i[0], coords, ccoords, INSERT_VALUES, SCATTER_FORWARD));
45: PetscCall(VecScatterEnd(scat_i[0], coords, ccoords, INSERT_VALUES, SCATTER_FORWARD));
46: PetscCall(VecScatterBegin(scat_g[0], coords, clcoords, INSERT_VALUES, SCATTER_FORWARD));
47: PetscCall(VecScatterEnd(scat_g[0], coords, clcoords, INSERT_VALUES, SCATTER_FORWARD));
48: PetscCall(DMSetCoordinates(subdm, ccoords));
49: PetscCall(DMSetCoordinatesLocal(subdm, clcoords));
50: PetscCall(VecScatterDestroy(&scat_i[0]));
51: PetscCall(VecScatterDestroy(&scat_g[0]));
52: PetscCall(VecDestroy(&ccoords));
53: PetscCall(VecDestroy(&clcoords));
54: PetscCall(PetscFree(scat_i));
55: PetscCall(PetscFree(scat_g));
56: }
57: PetscFunctionReturn(PETSC_SUCCESS);
58: }
60: /*@
61: DMGetCoordinateDM - Gets the `DM` that prescribes coordinate layout and scatters between global and local coordinates
63: Collective
65: Input Parameter:
66: . dm - the `DM`
68: Output Parameter:
69: . cdm - coordinate `DM`
71: Level: intermediate
73: .seealso: `DM`, `DMSetCoordinateDM()`, `DMSetCoordinates()`, `DMSetCoordinatesLocal()`, `DMGetCoordinates()`, `DMGetCoordinatesLocal()`, `DMGSetCellCoordinateDM()`,
75: @*/
76: PetscErrorCode DMGetCoordinateDM(DM dm, DM *cdm)
77: {
78: PetscFunctionBegin;
80: PetscAssertPointer(cdm, 2);
81: if (!dm->coordinates[0].dm) {
82: DM cdm;
84: PetscUseTypeMethod(dm, createcoordinatedm, &cdm);
85: PetscCall(PetscObjectSetName((PetscObject)cdm, "coordinateDM"));
86: /* Just in case the DM sets the coordinate DM when creating it (DMP4est can do this, because it may not setup
87: * until the call to CreateCoordinateDM) */
88: PetscCall(DMDestroy(&dm->coordinates[0].dm));
89: dm->coordinates[0].dm = cdm;
90: }
91: *cdm = dm->coordinates[0].dm;
92: PetscFunctionReturn(PETSC_SUCCESS);
93: }
95: /*@
96: DMSetCoordinateDM - Sets the `DM` that prescribes coordinate layout and scatters between global and local coordinates
98: Logically Collective
100: Input Parameters:
101: + dm - the `DM`
102: - cdm - coordinate `DM`
104: Level: intermediate
106: .seealso: `DM`, `DMGetCoordinateDM()`, `DMSetCoordinates()`, `DMGetCellCoordinateDM()`, `DMSetCoordinatesLocal()`, `DMGetCoordinates()`, `DMGetCoordinatesLocal()`,
107: `DMGSetCellCoordinateDM()`
108: @*/
109: PetscErrorCode DMSetCoordinateDM(DM dm, DM cdm)
110: {
111: PetscFunctionBegin;
114: PetscCall(PetscObjectReference((PetscObject)cdm));
115: PetscCall(DMDestroy(&dm->coordinates[0].dm));
116: dm->coordinates[0].dm = cdm;
117: PetscFunctionReturn(PETSC_SUCCESS);
118: }
120: /*@
121: DMGetCellCoordinateDM - Gets the `DM` that prescribes cellwise coordinate layout and scatters between global and local cellwise coordinates
123: Collective
125: Input Parameter:
126: . dm - the `DM`
128: Output Parameter:
129: . cdm - cellwise coordinate `DM`, or `NULL` if they are not defined
131: Level: intermediate
133: Note:
134: Call `DMLocalizeCoordinates()` to automatically create cellwise coordinates for periodic geometries.
136: .seealso: `DM`, `DMSetCellCoordinateDM()`, `DMSetCellCoordinates()`, `DMSetCellCoordinatesLocal()`, `DMGetCellCoordinates()`, `DMGetCellCoordinatesLocal()`,
137: `DMLocalizeCoordinates()`, `DMSetCoordinateDM()`, `DMGetCoordinateDM()`
138: @*/
139: PetscErrorCode DMGetCellCoordinateDM(DM dm, DM *cdm)
140: {
141: PetscFunctionBegin;
143: PetscAssertPointer(cdm, 2);
144: *cdm = dm->coordinates[1].dm;
145: PetscFunctionReturn(PETSC_SUCCESS);
146: }
148: /*@
149: DMSetCellCoordinateDM - Sets the `DM` that prescribes cellwise coordinate layout and scatters between global and local cellwise coordinates
151: Logically Collective
153: Input Parameters:
154: + dm - the `DM`
155: - cdm - cellwise coordinate `DM`
157: Level: intermediate
159: Note:
160: As opposed to `DMSetCoordinateDM()` these coordinates are useful for discontinuous Galerkin methods since they support coordinate fields that are discontinuous at cell boundaries.
162: .seealso: `DMGetCellCoordinateDM()`, `DMSetCellCoordinates()`, `DMSetCellCoordinatesLocal()`, `DMGetCellCoordinates()`, `DMGetCellCoordinatesLocal()`,
163: `DMSetCoordinateDM()`, `DMGetCoordinateDM()`
164: @*/
165: PetscErrorCode DMSetCellCoordinateDM(DM dm, DM cdm)
166: {
167: PetscInt dim;
169: PetscFunctionBegin;
171: if (cdm) {
173: PetscCall(DMGetCoordinateDim(dm, &dim));
174: dm->coordinates[1].dim = dim;
175: }
176: PetscCall(PetscObjectReference((PetscObject)cdm));
177: PetscCall(DMDestroy(&dm->coordinates[1].dm));
178: dm->coordinates[1].dm = cdm;
179: PetscFunctionReturn(PETSC_SUCCESS);
180: }
182: /*@
183: DMGetCoordinateDim - Retrieve the dimension of the embedding space for coordinate values. For example a mesh on the surface of a sphere would have a 3 dimensional embedding space
185: Not Collective
187: Input Parameter:
188: . dm - The `DM` object
190: Output Parameter:
191: . dim - The embedding dimension
193: Level: intermediate
195: .seealso: `DM`, `DMSetCoordinateDim()`, `DMGetCoordinateSection()`, `DMGetCoordinateDM()`, `DMGetLocalSection()`, `DMSetLocalSection()`
196: @*/
197: PetscErrorCode DMGetCoordinateDim(DM dm, PetscInt *dim)
198: {
199: PetscFunctionBegin;
201: PetscAssertPointer(dim, 2);
202: if (dm->coordinates[0].dim == PETSC_DEFAULT) dm->coordinates[0].dim = dm->dim;
203: *dim = dm->coordinates[0].dim;
204: PetscFunctionReturn(PETSC_SUCCESS);
205: }
207: /*@
208: DMSetCoordinateDim - Set the dimension of the embedding space for coordinate values.
210: Not Collective
212: Input Parameters:
213: + dm - The `DM` object
214: - dim - The embedding dimension
216: Level: intermediate
218: .seealso: `DM`, `DMGetCoordinateDim()`, `DMSetCoordinateSection()`, `DMGetCoordinateSection()`, `DMGetLocalSection()`, `DMSetLocalSection()`
219: @*/
220: PetscErrorCode DMSetCoordinateDim(DM dm, PetscInt dim)
221: {
222: PetscDS ds;
223: PetscInt Nds, n;
225: PetscFunctionBegin;
227: dm->coordinates[0].dim = dim;
228: if (dm->dim >= 0) {
229: PetscCall(DMGetNumDS(dm, &Nds));
230: for (n = 0; n < Nds; ++n) {
231: PetscCall(DMGetRegionNumDS(dm, n, NULL, NULL, &ds, NULL));
232: PetscCall(PetscDSSetCoordinateDimension(ds, dim));
233: }
234: }
235: PetscFunctionReturn(PETSC_SUCCESS);
236: }
238: /*@
239: DMGetCoordinateSection - Retrieve the `PetscSection` of coordinate values over the mesh.
241: Collective
243: Input Parameter:
244: . dm - The `DM` object
246: Output Parameter:
247: . section - The `PetscSection` object
249: Level: intermediate
251: Note:
252: This just retrieves the local section from the coordinate `DM`. In other words,
253: .vb
254: DMGetCoordinateDM(dm, &cdm);
255: DMGetLocalSection(cdm, §ion);
256: .ve
258: .seealso: `DM`, `DMGetCoordinateDM()`, `DMGetLocalSection()`, `DMSetLocalSection()`
259: @*/
260: PetscErrorCode DMGetCoordinateSection(DM dm, PetscSection *section)
261: {
262: DM cdm;
264: PetscFunctionBegin;
266: PetscAssertPointer(section, 2);
267: PetscCall(DMGetCoordinateDM(dm, &cdm));
268: PetscCall(DMGetLocalSection(cdm, section));
269: PetscFunctionReturn(PETSC_SUCCESS);
270: }
272: /*@
273: DMSetCoordinateSection - Set the `PetscSection` of coordinate values over the mesh.
275: Not Collective
277: Input Parameters:
278: + dm - The `DM` object
279: . dim - The embedding dimension, or `PETSC_DETERMINE`
280: - section - The `PetscSection` object
282: Level: intermediate
284: .seealso: `DM`, `DMGetCoordinateDim()`, `DMGetCoordinateSection()`, `DMGetLocalSection()`, `DMSetLocalSection()`
285: @*/
286: PetscErrorCode DMSetCoordinateSection(DM dm, PetscInt dim, PetscSection section)
287: {
288: DM cdm;
290: PetscFunctionBegin;
293: PetscCall(DMGetCoordinateDM(dm, &cdm));
294: PetscCall(DMSetLocalSection(cdm, section));
295: if (dim == PETSC_DETERMINE) {
296: PetscInt d = PETSC_DEFAULT;
297: PetscInt pStart, pEnd, vStart, vEnd, v, dd;
299: PetscCall(PetscSectionGetChart(section, &pStart, &pEnd));
300: PetscCall(DMGetDimPoints(dm, 0, &vStart, &vEnd));
301: pStart = PetscMax(vStart, pStart);
302: pEnd = PetscMin(vEnd, pEnd);
303: for (v = pStart; v < pEnd; ++v) {
304: PetscCall(PetscSectionGetDof(section, v, &dd));
305: if (dd) {
306: d = dd;
307: break;
308: }
309: }
310: if (d >= 0) PetscCall(DMSetCoordinateDim(dm, d));
311: }
312: PetscFunctionReturn(PETSC_SUCCESS);
313: }
315: /*@
316: DMGetCellCoordinateSection - Retrieve the `PetscSection` of cellwise coordinate values over the mesh.
318: Collective
320: Input Parameter:
321: . dm - The `DM` object
323: Output Parameter:
324: . section - The `PetscSection` object, or `NULL` if no cellwise coordinates are defined
326: Level: intermediate
328: Note:
329: This just retrieves the local section from the cell coordinate `DM`. In other words,
330: .vb
331: DMGetCellCoordinateDM(dm, &cdm);
332: DMGetLocalSection(cdm, §ion);
333: .ve
335: .seealso: `DM`, `DMGetCoordinateSection()`, `DMSetCellCoordinateSection()`, `DMGetCellCoordinateDM()`, `DMGetCoordinateDM()`, `DMGetLocalSection()`, `DMSetLocalSection()`
336: @*/
337: PetscErrorCode DMGetCellCoordinateSection(DM dm, PetscSection *section)
338: {
339: DM cdm;
341: PetscFunctionBegin;
343: PetscAssertPointer(section, 2);
344: *section = NULL;
345: PetscCall(DMGetCellCoordinateDM(dm, &cdm));
346: if (cdm) PetscCall(DMGetLocalSection(cdm, section));
347: PetscFunctionReturn(PETSC_SUCCESS);
348: }
350: /*@
351: DMSetCellCoordinateSection - Set the `PetscSection` of cellwise coordinate values over the mesh.
353: Not Collective
355: Input Parameters:
356: + dm - The `DM` object
357: . dim - The embedding dimension, or `PETSC_DETERMINE`
358: - section - The `PetscSection` object for a cellwise layout
360: Level: intermediate
362: .seealso: `DM`, `DMGetCoordinateDim()`, `DMSetCoordinateSection()`, `DMGetCellCoordinateSection()`, `DMGetCoordinateSection()`, `DMGetCellCoordinateDM()`, `DMGetLocalSection()`, `DMSetLocalSection()`
363: @*/
364: PetscErrorCode DMSetCellCoordinateSection(DM dm, PetscInt dim, PetscSection section)
365: {
366: DM cdm;
368: PetscFunctionBegin;
371: PetscCall(DMGetCellCoordinateDM(dm, &cdm));
372: PetscCheck(cdm, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "No DM defined for cellwise coordinates");
373: PetscCall(DMSetLocalSection(cdm, section));
374: if (dim == PETSC_DETERMINE) {
375: PetscInt d = PETSC_DEFAULT;
376: PetscInt pStart, pEnd, vStart, vEnd, v, dd;
378: PetscCall(PetscSectionGetChart(section, &pStart, &pEnd));
379: PetscCall(DMGetDimPoints(dm, 0, &vStart, &vEnd));
380: pStart = PetscMax(vStart, pStart);
381: pEnd = PetscMin(vEnd, pEnd);
382: for (v = pStart; v < pEnd; ++v) {
383: PetscCall(PetscSectionGetDof(section, v, &dd));
384: if (dd) {
385: d = dd;
386: break;
387: }
388: }
389: if (d >= 0) PetscCall(DMSetCoordinateDim(dm, d));
390: }
391: PetscFunctionReturn(PETSC_SUCCESS);
392: }
394: /*@
395: DMGetCoordinates - Gets a global vector with the coordinates associated with the `DM`.
397: Collective
399: Input Parameter:
400: . dm - the `DM`
402: Output Parameter:
403: . c - global coordinate vector
405: Level: intermediate
407: Notes:
408: This is a borrowed reference, so the user should NOT destroy this vector. When the `DM` is
409: destroyed `c` will no longer be valid.
411: Each process has only the locally-owned portion of the global coordinates (does NOT have the ghost coordinates).
413: For `DMDA`, in two and three dimensions coordinates are interlaced (x_0,y_0,x_1,y_1,...)
414: and (x_0,y_0,z_0,x_1,y_1,z_1...)
416: .seealso: `DM`, `DMDA`, `DMSetCoordinates()`, `DMGetCoordinatesLocal()`, `DMGetCoordinateDM()`, `DMDASetUniformCoordinates()`
417: @*/
418: PetscErrorCode DMGetCoordinates(DM dm, Vec *c)
419: {
420: PetscFunctionBegin;
422: PetscAssertPointer(c, 2);
423: if (!dm->coordinates[0].x && dm->coordinates[0].xl) {
424: DM cdm = NULL;
426: PetscCall(DMGetCoordinateDM(dm, &cdm));
427: PetscCall(DMCreateGlobalVector(cdm, &dm->coordinates[0].x));
428: PetscCall(PetscObjectSetName((PetscObject)dm->coordinates[0].x, "coordinates"));
429: PetscCall(DMLocalToGlobalBegin(cdm, dm->coordinates[0].xl, INSERT_VALUES, dm->coordinates[0].x));
430: PetscCall(DMLocalToGlobalEnd(cdm, dm->coordinates[0].xl, INSERT_VALUES, dm->coordinates[0].x));
431: }
432: *c = dm->coordinates[0].x;
433: PetscFunctionReturn(PETSC_SUCCESS);
434: }
436: /*@
437: DMSetCoordinates - Sets into the `DM` a global vector that holds the coordinates
439: Collective
441: Input Parameters:
442: + dm - the `DM`
443: - c - coordinate vector
445: Level: intermediate
447: Notes:
448: The coordinates do not include those for ghost points, which are in the local vector.
450: The vector `c` can be destroyed after the call
452: .seealso: `DM`, `DMSetCoordinatesLocal()`, `DMGetCoordinates()`, `DMGetCoordinatesLocal()`, `DMGetCoordinateDM()`, `DMDASetUniformCoordinates()`
453: @*/
454: PetscErrorCode DMSetCoordinates(DM dm, Vec c)
455: {
456: PetscFunctionBegin;
459: PetscCall(PetscObjectReference((PetscObject)c));
460: PetscCall(VecDestroy(&dm->coordinates[0].x));
461: dm->coordinates[0].x = c;
462: PetscCall(VecDestroy(&dm->coordinates[0].xl));
463: PetscCall(DMCoarsenHookAdd(dm, DMRestrictHook_Coordinates, NULL, NULL));
464: PetscCall(DMSubDomainHookAdd(dm, DMSubDomainHook_Coordinates, NULL, NULL));
465: PetscFunctionReturn(PETSC_SUCCESS);
466: }
468: /*@
469: DMGetCellCoordinates - Gets a global vector with the cellwise coordinates associated with the `DM`.
471: Collective
473: Input Parameter:
474: . dm - the `DM`
476: Output Parameter:
477: . c - global coordinate vector
479: Level: intermediate
481: Notes:
482: This is a borrowed reference, so the user should NOT destroy this vector. When the `DM` is
483: destroyed `c` will no longer be valid.
485: Each process has only the locally-owned portion of the global coordinates (does NOT have the ghost coordinates).
487: .seealso: `DM`, `DMGetCoordinates()`, `DMSetCellCoordinates()`, `DMGetCellCoordinatesLocal()`, `DMGetCellCoordinateDM()`
488: @*/
489: PetscErrorCode DMGetCellCoordinates(DM dm, Vec *c)
490: {
491: PetscFunctionBegin;
493: PetscAssertPointer(c, 2);
494: if (!dm->coordinates[1].x && dm->coordinates[1].xl) {
495: DM cdm = NULL;
497: PetscCall(DMGetCellCoordinateDM(dm, &cdm));
498: PetscCall(DMCreateGlobalVector(cdm, &dm->coordinates[1].x));
499: PetscCall(PetscObjectSetName((PetscObject)dm->coordinates[1].x, "DG coordinates"));
500: PetscCall(DMLocalToGlobalBegin(cdm, dm->coordinates[1].xl, INSERT_VALUES, dm->coordinates[1].x));
501: PetscCall(DMLocalToGlobalEnd(cdm, dm->coordinates[1].xl, INSERT_VALUES, dm->coordinates[1].x));
502: }
503: *c = dm->coordinates[1].x;
504: PetscFunctionReturn(PETSC_SUCCESS);
505: }
507: /*@
508: DMSetCellCoordinates - Sets into the `DM` a global vector that holds the cellwise coordinates
510: Collective
512: Input Parameters:
513: + dm - the `DM`
514: - c - cellwise coordinate vector
516: Level: intermediate
518: Notes:
519: The coordinates do not include those for ghost points, which are in the local vector.
521: The vector `c` should be destroyed by the caller.
523: .seealso: `DM`, `DMGetCoordinates()`, `DMSetCellCoordinatesLocal()`, `DMGetCellCoordinates()`, `DMGetCellCoordinatesLocal()`, `DMGetCellCoordinateDM()`
524: @*/
525: PetscErrorCode DMSetCellCoordinates(DM dm, Vec c)
526: {
527: PetscFunctionBegin;
530: PetscCall(PetscObjectReference((PetscObject)c));
531: PetscCall(VecDestroy(&dm->coordinates[1].x));
532: dm->coordinates[1].x = c;
533: PetscCall(VecDestroy(&dm->coordinates[1].xl));
534: PetscFunctionReturn(PETSC_SUCCESS);
535: }
537: /*@
538: DMGetCoordinatesLocalSetUp - Prepares a local vector of coordinates, so that `DMGetCoordinatesLocalNoncollective()` can be used as non-collective afterwards.
540: Collective
542: Input Parameter:
543: . dm - the `DM`
545: Level: advanced
547: .seealso: `DM`, `DMSetCoordinates()`, `DMGetCoordinatesLocalNoncollective()`
548: @*/
549: PetscErrorCode DMGetCoordinatesLocalSetUp(DM dm)
550: {
551: PetscFunctionBegin;
553: if (!dm->coordinates[0].xl && dm->coordinates[0].x) {
554: DM cdm = NULL;
555: PetscInt bs;
557: PetscCall(DMGetCoordinateDM(dm, &cdm));
558: PetscCall(DMCreateLocalVector(cdm, &dm->coordinates[0].xl));
559: // If the size of the vector is 0, it will not get the right block size
560: PetscCall(VecGetBlockSize(dm->coordinates[0].x, &bs));
561: PetscCall(VecSetBlockSize(dm->coordinates[0].xl, bs));
562: PetscCall(PetscObjectSetName((PetscObject)dm->coordinates[0].xl, "coordinates"));
563: PetscCall(DMGlobalToLocalBegin(cdm, dm->coordinates[0].x, INSERT_VALUES, dm->coordinates[0].xl));
564: PetscCall(DMGlobalToLocalEnd(cdm, dm->coordinates[0].x, INSERT_VALUES, dm->coordinates[0].xl));
565: }
566: PetscFunctionReturn(PETSC_SUCCESS);
567: }
569: /*@
570: DMGetCoordinatesLocal - Gets a local vector with the coordinates associated with the `DM`.
572: Collective the first time it is called
574: Input Parameter:
575: . dm - the `DM`
577: Output Parameter:
578: . c - coordinate vector
580: Level: intermediate
582: Notes:
583: This is a borrowed reference, so the user should NOT destroy `c`
585: Each process has the local and ghost coordinates
587: For `DMDA`, in two and three dimensions coordinates are interlaced (x_0,y_0,x_1,y_1,...)
588: and (x_0,y_0,z_0,x_1,y_1,z_1...)
590: .seealso: `DM`, `DMSetCoordinatesLocal()`, `DMGetCoordinates()`, `DMSetCoordinates()`, `DMGetCoordinateDM()`, `DMGetCoordinatesLocalNoncollective()`
591: @*/
592: PetscErrorCode DMGetCoordinatesLocal(DM dm, Vec *c)
593: {
594: PetscFunctionBegin;
596: PetscAssertPointer(c, 2);
597: PetscCall(DMGetCoordinatesLocalSetUp(dm));
598: *c = dm->coordinates[0].xl;
599: PetscFunctionReturn(PETSC_SUCCESS);
600: }
602: /*@
603: DMGetCoordinatesLocalNoncollective - Non-collective version of `DMGetCoordinatesLocal()`. Fails if global coordinates have been set and `DMGetCoordinatesLocalSetUp()` not called.
605: Not Collective
607: Input Parameter:
608: . dm - the `DM`
610: Output Parameter:
611: . c - coordinate vector
613: Level: advanced
615: Note:
616: A previous call to `DMGetCoordinatesLocal()` or `DMGetCoordinatesLocalSetUp()` ensures that a call to this function will not error.
618: .seealso: `DM`, `DMGetCoordinatesLocalSetUp()`, `DMGetCoordinatesLocal()`, `DMSetCoordinatesLocal()`, `DMGetCoordinates()`, `DMSetCoordinates()`, `DMGetCoordinateDM()`
619: @*/
620: PetscErrorCode DMGetCoordinatesLocalNoncollective(DM dm, Vec *c)
621: {
622: PetscFunctionBegin;
624: PetscAssertPointer(c, 2);
625: PetscCheck(dm->coordinates[0].xl || !dm->coordinates[0].x, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "DMGetCoordinatesLocalSetUp() has not been called");
626: *c = dm->coordinates[0].xl;
627: PetscFunctionReturn(PETSC_SUCCESS);
628: }
630: /*@
631: DMGetCoordinatesLocalTuple - Gets a local vector with the coordinates of specified points and the section describing its layout.
633: Not Collective
635: Input Parameters:
636: + dm - the `DM`
637: - p - the `IS` of points whose coordinates will be returned
639: Output Parameters:
640: + pCoordSection - the `PetscSection` describing the layout of pCoord, i.e. each point corresponds to one point in `p`, and DOFs correspond to coordinates
641: - pCoord - the `Vec` with coordinates of points in `p`
643: Level: advanced
645: Notes:
646: `DMGetCoordinatesLocalSetUp()` must be called first. This function employs `DMGetCoordinatesLocalNoncollective()` so it is not collective.
648: This creates a new vector, so the user SHOULD destroy this vector
650: Each process has the local and ghost coordinates
652: For `DMDA`, in two and three dimensions coordinates are interlaced (x_0,y_0,x_1,y_1,...)
653: and (x_0,y_0,z_0,x_1,y_1,z_1...)
655: .seealso: `DM`, `DMDA`, `DMSetCoordinatesLocal()`, `DMGetCoordinatesLocal()`, `DMGetCoordinatesLocalNoncollective()`, `DMGetCoordinatesLocalSetUp()`, `DMGetCoordinates()`, `DMSetCoordinates()`, `DMGetCoordinateDM()`
656: @*/
657: PetscErrorCode DMGetCoordinatesLocalTuple(DM dm, IS p, PetscSection *pCoordSection, Vec *pCoord)
658: {
659: DM cdm;
660: PetscSection cs, newcs;
661: Vec coords;
662: const PetscScalar *arr;
663: PetscScalar *newarr = NULL;
664: PetscInt n;
666: PetscFunctionBegin;
669: if (pCoordSection) PetscAssertPointer(pCoordSection, 3);
670: if (pCoord) PetscAssertPointer(pCoord, 4);
671: PetscCall(DMGetCoordinateDM(dm, &cdm));
672: PetscCall(DMGetLocalSection(cdm, &cs));
673: PetscCall(DMGetCoordinatesLocal(dm, &coords));
674: PetscCheck(coords, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "DMGetCoordinatesLocalSetUp() has not been called or coordinates not set");
675: PetscCheck(cdm && cs, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "DM not supported");
676: PetscCall(VecGetArrayRead(coords, &arr));
677: PetscCall(PetscSectionExtractDofsFromArray(cs, MPIU_SCALAR, arr, p, &newcs, pCoord ? ((void **)&newarr) : NULL));
678: PetscCall(VecRestoreArrayRead(coords, &arr));
679: if (pCoord) {
680: PetscCall(PetscSectionGetStorageSize(newcs, &n));
681: /* set array in two steps to mimic PETSC_OWN_POINTER */
682: PetscCall(VecCreateSeqWithArray(PetscObjectComm((PetscObject)p), 1, n, NULL, pCoord));
683: PetscCall(VecReplaceArray(*pCoord, newarr));
684: } else {
685: PetscCall(PetscFree(newarr));
686: }
687: if (pCoordSection) {
688: *pCoordSection = newcs;
689: } else PetscCall(PetscSectionDestroy(&newcs));
690: PetscFunctionReturn(PETSC_SUCCESS);
691: }
693: /*@
694: DMSetCoordinatesLocal - Sets into the `DM` a local vector, including ghost points, that holds the coordinates
696: Not Collective
698: Input Parameters:
699: + dm - the `DM`
700: - c - coordinate vector
702: Level: intermediate
704: Notes:
705: The coordinates of ghost points can be set using `DMSetCoordinates()`
706: followed by `DMGetCoordinatesLocal()`. This is intended to enable the
707: setting of ghost coordinates outside of the domain.
709: The vector `c` should be destroyed by the caller.
711: .seealso: `DM`, `DMGetCoordinatesLocal()`, `DMSetCoordinates()`, `DMGetCoordinates()`, `DMGetCoordinateDM()`
712: @*/
713: PetscErrorCode DMSetCoordinatesLocal(DM dm, Vec c)
714: {
715: PetscFunctionBegin;
718: PetscCall(PetscObjectReference((PetscObject)c));
719: PetscCall(VecDestroy(&dm->coordinates[0].xl));
720: dm->coordinates[0].xl = c;
721: PetscCall(VecDestroy(&dm->coordinates[0].x));
722: PetscFunctionReturn(PETSC_SUCCESS);
723: }
725: /*@
726: DMGetCellCoordinatesLocalSetUp - Prepares a local vector of cellwise coordinates, so that `DMGetCellCoordinatesLocalNoncollective()` can be used as non-collective afterwards.
728: Collective
730: Input Parameter:
731: . dm - the `DM`
733: Level: advanced
735: .seealso: `DM`, `DMGetCellCoordinatesLocalNoncollective()`
736: @*/
737: PetscErrorCode DMGetCellCoordinatesLocalSetUp(DM dm)
738: {
739: PetscFunctionBegin;
741: if (!dm->coordinates[1].xl && dm->coordinates[1].x) {
742: DM cdm = NULL;
744: PetscCall(DMGetCellCoordinateDM(dm, &cdm));
745: PetscCall(DMCreateLocalVector(cdm, &dm->coordinates[1].xl));
746: PetscCall(PetscObjectSetName((PetscObject)dm->coordinates[1].xl, "DG coordinates"));
747: PetscCall(DMGlobalToLocalBegin(cdm, dm->coordinates[1].x, INSERT_VALUES, dm->coordinates[1].xl));
748: PetscCall(DMGlobalToLocalEnd(cdm, dm->coordinates[1].x, INSERT_VALUES, dm->coordinates[1].xl));
749: }
750: PetscFunctionReturn(PETSC_SUCCESS);
751: }
753: /*@
754: DMGetCellCoordinatesLocal - Gets a local vector with the cellwise coordinates associated with the `DM`.
756: Collective
758: Input Parameter:
759: . dm - the `DM`
761: Output Parameter:
762: . c - coordinate vector
764: Level: intermediate
766: Notes:
767: This is a borrowed reference, so the user should NOT destroy this vector
769: Each process has the local and ghost coordinates
771: .seealso: `DM`, `DMSetCellCoordinatesLocal()`, `DMGetCellCoordinates()`, `DMSetCellCoordinates()`, `DMGetCellCoordinateDM()`, `DMGetCellCoordinatesLocalNoncollective()`
772: @*/
773: PetscErrorCode DMGetCellCoordinatesLocal(DM dm, Vec *c)
774: {
775: PetscFunctionBegin;
777: PetscAssertPointer(c, 2);
778: PetscCall(DMGetCellCoordinatesLocalSetUp(dm));
779: *c = dm->coordinates[1].xl;
780: PetscFunctionReturn(PETSC_SUCCESS);
781: }
783: /*@
784: DMGetCellCoordinatesLocalNoncollective - Non-collective version of `DMGetCellCoordinatesLocal()`. Fails if global cellwise coordinates have been set and `DMGetCellCoordinatesLocalSetUp()` not called.
786: Not Collective
788: Input Parameter:
789: . dm - the `DM`
791: Output Parameter:
792: . c - cellwise coordinate vector
794: Level: advanced
796: .seealso: `DM`, `DMGetCellCoordinatesLocalSetUp()`, `DMGetCellCoordinatesLocal()`, `DMSetCellCoordinatesLocal()`, `DMGetCellCoordinates()`, `DMSetCellCoordinates()`, `DMGetCellCoordinateDM()`
797: @*/
798: PetscErrorCode DMGetCellCoordinatesLocalNoncollective(DM dm, Vec *c)
799: {
800: PetscFunctionBegin;
802: PetscAssertPointer(c, 2);
803: PetscCheck(dm->coordinates[1].xl || !dm->coordinates[1].x, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "DMGetCellCoordinatesLocalSetUp() has not been called");
804: *c = dm->coordinates[1].xl;
805: PetscFunctionReturn(PETSC_SUCCESS);
806: }
808: /*@
809: DMSetCellCoordinatesLocal - Sets into the `DM` a local vector including ghost points that holds the cellwise coordinates
811: Not Collective
813: Input Parameters:
814: + dm - the `DM`
815: - c - cellwise coordinate vector
817: Level: intermediate
819: Notes:
820: The coordinates of ghost points can be set using `DMSetCoordinates()`
821: followed by `DMGetCoordinatesLocal()`. This is intended to enable the
822: setting of ghost coordinates outside of the domain.
824: The vector `c` should be destroyed by the caller.
826: .seealso: `DM`, `DMGetCellCoordinatesLocal()`, `DMSetCellCoordinates()`, `DMGetCellCoordinates()`, `DMGetCellCoordinateDM()`
827: @*/
828: PetscErrorCode DMSetCellCoordinatesLocal(DM dm, Vec c)
829: {
830: PetscFunctionBegin;
833: PetscCall(PetscObjectReference((PetscObject)c));
834: PetscCall(VecDestroy(&dm->coordinates[1].xl));
835: dm->coordinates[1].xl = c;
836: PetscCall(VecDestroy(&dm->coordinates[1].x));
837: PetscFunctionReturn(PETSC_SUCCESS);
838: }
840: PetscErrorCode DMGetCoordinateField(DM dm, DMField *field)
841: {
842: PetscFunctionBegin;
844: PetscAssertPointer(field, 2);
845: if (!dm->coordinates[0].field) {
846: if (dm->ops->createcoordinatefield) PetscCall((*dm->ops->createcoordinatefield)(dm, &dm->coordinates[0].field));
847: }
848: *field = dm->coordinates[0].field;
849: PetscFunctionReturn(PETSC_SUCCESS);
850: }
852: PetscErrorCode DMSetCoordinateField(DM dm, DMField field)
853: {
854: PetscFunctionBegin;
857: PetscCall(PetscObjectReference((PetscObject)field));
858: PetscCall(DMFieldDestroy(&dm->coordinates[0].field));
859: dm->coordinates[0].field = field;
860: PetscFunctionReturn(PETSC_SUCCESS);
861: }
863: /*@
864: DMGetLocalBoundingBox - Returns the bounding box for the piece of the `DM` on this process.
866: Not Collective
868: Input Parameter:
869: . dm - the `DM`
871: Output Parameters:
872: + lmin - local minimum coordinates (length coord dim, optional)
873: - lmax - local maximum coordinates (length coord dim, optional)
875: Level: beginner
877: Note:
878: If the `DM` is a `DMDA` and has no coordinates, the index bounds are returned instead.
880: .seealso: `DM`, `DMGetCoordinates()`, `DMGetCoordinatesLocal()`, `DMGetBoundingBox()`
881: @*/
882: PetscErrorCode DMGetLocalBoundingBox(DM dm, PetscReal lmin[], PetscReal lmax[])
883: {
884: Vec coords = NULL;
885: PetscReal min[3] = {PETSC_MAX_REAL, PETSC_MAX_REAL, PETSC_MAX_REAL};
886: PetscReal max[3] = {PETSC_MIN_REAL, PETSC_MIN_REAL, PETSC_MIN_REAL};
887: PetscInt cdim, i, j;
889: PetscFunctionBegin;
891: PetscCall(DMGetCoordinateDim(dm, &cdim));
892: PetscCall(DMGetCoordinatesLocal(dm, &coords));
893: if (coords) {
894: const PetscScalar *local_coords;
895: PetscInt N, Ni;
897: for (j = cdim; j < 3; ++j) {
898: min[j] = 0;
899: max[j] = 0;
900: }
901: PetscCall(VecGetArrayRead(coords, &local_coords));
902: PetscCall(VecGetLocalSize(coords, &N));
903: Ni = N / cdim;
904: for (i = 0; i < Ni; ++i) {
905: for (j = 0; j < cdim; ++j) {
906: min[j] = PetscMin(min[j], PetscRealPart(local_coords[i * cdim + j]));
907: max[j] = PetscMax(max[j], PetscRealPart(local_coords[i * cdim + j]));
908: }
909: }
910: PetscCall(VecRestoreArrayRead(coords, &local_coords));
911: PetscCall(DMGetCellCoordinatesLocal(dm, &coords));
912: if (coords) {
913: PetscCall(VecGetArrayRead(coords, &local_coords));
914: PetscCall(VecGetLocalSize(coords, &N));
915: Ni = N / cdim;
916: for (i = 0; i < Ni; ++i) {
917: for (j = 0; j < cdim; ++j) {
918: min[j] = PetscMin(min[j], PetscRealPart(local_coords[i * cdim + j]));
919: max[j] = PetscMax(max[j], PetscRealPart(local_coords[i * cdim + j]));
920: }
921: }
922: PetscCall(VecRestoreArrayRead(coords, &local_coords));
923: }
924: } else {
925: PetscBool isda;
927: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMDA, &isda));
928: if (isda) PetscCall(DMGetLocalBoundingIndices_DMDA(dm, min, max));
929: }
930: if (lmin) PetscCall(PetscArraycpy(lmin, min, cdim));
931: if (lmax) PetscCall(PetscArraycpy(lmax, max, cdim));
932: PetscFunctionReturn(PETSC_SUCCESS);
933: }
935: /*@
936: DMGetBoundingBox - Returns the global bounding box for the `DM`.
938: Collective
940: Input Parameter:
941: . dm - the `DM`
943: Output Parameters:
944: + gmin - global minimum coordinates (length coord dim, optional)
945: - gmax - global maximum coordinates (length coord dim, optional)
947: Level: beginner
949: .seealso: `DM`, `DMGetLocalBoundingBox()`, `DMGetCoordinates()`, `DMGetCoordinatesLocal()`
950: @*/
951: PetscErrorCode DMGetBoundingBox(DM dm, PetscReal gmin[], PetscReal gmax[])
952: {
953: PetscReal lmin[3], lmax[3];
954: PetscInt cdim;
955: PetscMPIInt count;
957: PetscFunctionBegin;
959: PetscCall(DMGetCoordinateDim(dm, &cdim));
960: PetscCall(PetscMPIIntCast(cdim, &count));
961: PetscCall(DMGetLocalBoundingBox(dm, lmin, lmax));
962: if (gmin) PetscCall(MPIU_Allreduce(lmin, gmin, count, MPIU_REAL, MPIU_MIN, PetscObjectComm((PetscObject)dm)));
963: if (gmax) PetscCall(MPIU_Allreduce(lmax, gmax, count, MPIU_REAL, MPIU_MAX, PetscObjectComm((PetscObject)dm)));
964: PetscFunctionReturn(PETSC_SUCCESS);
965: }
967: static void evaluate_coordinates(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar xnew[])
968: {
969: for (PetscInt i = 0; i < dim; i++) xnew[i] = x[i];
970: }
972: /*@
973: DMProjectCoordinates - Project coordinates to a different space
975: Input Parameters:
976: + dm - The `DM` object
977: - disc - The new coordinate discretization or `NULL` to ensure a coordinate discretization exists
979: Level: intermediate
981: Notes:
982: A `PetscFE` defines an approximation space using a `PetscSpace`, which represents the basis functions, and a `PetscDualSpace`, which defines the interpolation operation
983: in the space.
985: This function takes the current mesh coordinates, which are discretized using some `PetscFE` space, and projects this function into a new `PetscFE` space.
986: The coordinate projection is done on the continuous coordinates, and if possible, the discontinuous coordinates are also updated.
988: Developer Note:
989: With more effort, we could directly project the discontinuous coordinates also.
991: .seealso: `DM`, `PetscFE`, `DMGetCoordinateField()`
992: @*/
993: PetscErrorCode DMProjectCoordinates(DM dm, PetscFE disc)
994: {
995: PetscFE discOld;
996: PetscClassId classid;
997: DM cdmOld, cdmNew;
998: Vec coordsOld, coordsNew;
999: PetscBool same_space = PETSC_TRUE;
1000: const char *prefix;
1002: PetscFunctionBegin;
1006: PetscCall(DMGetCoordinateDM(dm, &cdmOld));
1007: /* Check current discretization is compatible */
1008: PetscCall(DMGetField(cdmOld, 0, NULL, (PetscObject *)&discOld));
1009: PetscCall(PetscObjectGetClassId((PetscObject)discOld, &classid));
1010: if (classid != PETSCFE_CLASSID) {
1011: if (classid == PETSC_CONTAINER_CLASSID) {
1012: PetscFE feLinear;
1013: DMPolytopeType ct;
1014: PetscInt dim, dE, cStart, cEnd, ctTmp;
1016: /* Assume linear vertex coordinates */
1017: PetscCall(DMGetDimension(dm, &dim));
1018: PetscCall(DMGetCoordinateDim(dm, &dE));
1019: PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd));
1020: if (cEnd > cStart) PetscCall(DMPlexGetCellType(dm, cStart, &ct));
1021: else ct = DM_POLYTOPE_UNKNOWN;
1022: ctTmp = (PetscInt)ct;
1023: PetscCall(MPIU_Allreduce(MPI_IN_PLACE, &ctTmp, 1, MPIU_INT, MPI_MIN, PetscObjectComm((PetscObject)dm)));
1024: ct = (DMPolytopeType)ctTmp;
1025: PetscCall(PetscFECreateLagrangeByCell(PETSC_COMM_SELF, dim, dE, ct, 1, -1, &feLinear));
1026: PetscCall(DMSetField(cdmOld, 0, NULL, (PetscObject)feLinear));
1027: PetscCall(PetscFEDestroy(&feLinear));
1028: PetscCall(DMCreateDS(cdmOld));
1029: PetscCall(DMGetField(cdmOld, 0, NULL, (PetscObject *)&discOld));
1030: } else {
1031: const char *discname;
1033: PetscCall(PetscObjectGetType((PetscObject)discOld, &discname));
1034: SETERRQ(PetscObjectComm((PetscObject)discOld), PETSC_ERR_SUP, "Discretization type %s not supported", discname);
1035: }
1036: }
1037: if (!disc) PetscFunctionReturn(PETSC_SUCCESS);
1038: { // Check if the new space is the same as the old modulo quadrature
1039: PetscDualSpace dsOld, ds;
1040: PetscCall(PetscFEGetDualSpace(discOld, &dsOld));
1041: PetscCall(PetscFEGetDualSpace(disc, &ds));
1042: PetscCall(PetscDualSpaceEqual(dsOld, ds, &same_space));
1043: }
1044: /* Make a fresh clone of the coordinate DM */
1045: PetscCall(DMClone(cdmOld, &cdmNew));
1046: cdmNew->cloneOpts = PETSC_TRUE;
1047: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)cdmOld, &prefix));
1048: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)cdmNew, prefix));
1049: PetscCall(DMSetField(cdmNew, 0, NULL, (PetscObject)disc));
1050: PetscCall(DMCreateDS(cdmNew));
1051: if (cdmOld->periodic.setup) {
1052: cdmNew->periodic.setup = cdmOld->periodic.setup;
1053: PetscCall(cdmNew->periodic.setup(cdmNew));
1054: }
1055: if (dm->setfromoptionscalled) PetscCall(DMSetFromOptions(cdmNew));
1056: PetscCall(DMGetCoordinates(dm, &coordsOld));
1057: PetscCall(DMCreateGlobalVector(cdmNew, &coordsNew));
1058: if (same_space) {
1059: // Need to copy so that the new vector has the right dm
1060: PetscCall(VecCopy(coordsOld, coordsNew));
1061: } else { // Project the coordinate vector from old to new space
1062: void (*funcs[])(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[]) = {evaluate_coordinates};
1063: // We can't call DMProjectField directly because it depends on KSP for DMGlobalToLocalSolve(), but we can use the core strategy
1064: Vec X_new_loc;
1065: PetscCall(DMCreateLocalVector(cdmNew, &X_new_loc));
1066: PetscCall(DMSetCoordinateDM(cdmNew, cdmOld));
1067: // See DMPlexRemapGeometry() for a similar pattern handling the coordinate field
1068: DMField cf;
1069: PetscCall(DMGetCoordinateField(dm, &cf));
1070: cdmNew->coordinates[0].field = cf;
1071: PetscCall(DMProjectFieldLocal(cdmNew, 0.0, NULL, funcs, INSERT_VALUES, X_new_loc));
1072: cdmNew->coordinates[0].field = NULL;
1073: PetscCall(DMSetCoordinateDM(cdmNew, NULL));
1074: PetscCall(DMLocalToGlobal(cdmNew, X_new_loc, INSERT_VALUES, coordsNew));
1075: PetscCall(VecDestroy(&X_new_loc));
1076: }
1077: /* Set new coordinate structures */
1078: PetscCall(DMSetCoordinateField(dm, NULL));
1079: PetscCall(DMSetCoordinateDM(dm, cdmNew));
1080: PetscCall(DMSetCoordinates(dm, coordsNew));
1081: PetscCall(VecDestroy(&coordsNew));
1082: PetscCall(DMDestroy(&cdmNew));
1083: PetscFunctionReturn(PETSC_SUCCESS);
1084: }
1086: /*@
1087: DMLocatePoints - Locate the points in `v` in the mesh and return a `PetscSF` of the containing cells
1089: Collective
1091: Input Parameters:
1092: + dm - The `DM`
1093: - ltype - The type of point location, e.g. `DM_POINTLOCATION_NONE` or `DM_POINTLOCATION_NEAREST`
1095: Input/Output Parameters:
1096: + v - The `Vec` of points, on output contains the nearest mesh points to the given points if `DM_POINTLOCATION_NEAREST` is used
1097: - cellSF - Points to either `NULL`, or a `PetscSF` with guesses for which cells contain each point;
1098: on output, the `PetscSF` containing the MPI ranks and local indices of the containing points
1100: Level: developer
1102: Notes:
1103: To do a search of the local cells of the mesh, `v` should have `PETSC_COMM_SELF` as its communicator.
1104: To do a search of all the cells in the distributed mesh, `v` should have the same MPI communicator as `dm`.
1106: Points will only be located in owned cells, not overlap cells arising from `DMPlexDistribute()` or other overlapping distributions.
1108: If *cellSF is `NULL` on input, a `PetscSF` will be created.
1109: If *cellSF is not `NULL` on input, it should point to an existing `PetscSF`, whose graph will be used as initial guesses.
1111: An array that maps each point to its containing cell can be obtained with
1112: .vb
1113: const PetscSFNode *cells;
1114: PetscInt nFound;
1115: const PetscInt *found;
1117: PetscSFGetGraph(cellSF,NULL,&nFound,&found,&cells);
1118: .ve
1120: Where cells[i].rank is the MPI rank of the process owning the cell containing point found[i] (or i if found == NULL), and cells[i].index is
1121: the index of the cell in its MPI process' local numbering. This rank is in the communicator for `v`, so if `v` is on `PETSC_COMM_SELF` then the rank will always be 0.
1123: .seealso: `DM`, `DMSetCoordinates()`, `DMSetCoordinatesLocal()`, `DMGetCoordinates()`, `DMGetCoordinatesLocal()`, `DMPointLocationType`
1124: @*/
1125: PetscErrorCode DMLocatePoints(DM dm, Vec v, DMPointLocationType ltype, PetscSF *cellSF)
1126: {
1127: PetscFunctionBegin;
1130: PetscAssertPointer(cellSF, 4);
1131: if (*cellSF) {
1132: PetscMPIInt result;
1135: PetscCallMPI(MPI_Comm_compare(PetscObjectComm((PetscObject)v), PetscObjectComm((PetscObject)*cellSF), &result));
1136: PetscCheck(result == MPI_IDENT || result == MPI_CONGRUENT, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "cellSF must have a communicator congruent to v's");
1137: } else {
1138: PetscCall(PetscSFCreate(PetscObjectComm((PetscObject)v), cellSF));
1139: }
1140: PetscCall(PetscLogEventBegin(DM_LocatePoints, dm, 0, 0, 0));
1141: PetscUseTypeMethod(dm, locatepoints, v, ltype, *cellSF);
1142: PetscCall(PetscLogEventEnd(DM_LocatePoints, dm, 0, 0, 0));
1143: PetscFunctionReturn(PETSC_SUCCESS);
1144: }