Actual source code: isltog.c
1: #include <petsc/private/isimpl.h>
2: #include <petsc/private/hashmapi.h>
3: #include <petscsf.h>
4: #include <petscviewer.h>
5: #include <petscbt.h>
7: PetscClassId IS_LTOGM_CLASSID;
8: static PetscErrorCode ISLocalToGlobalMappingSetUpBlockInfo_Private(ISLocalToGlobalMapping);
10: typedef struct {
11: PetscInt *globals;
12: } ISLocalToGlobalMapping_Basic;
14: typedef struct {
15: PetscHMapI globalht;
16: } ISLocalToGlobalMapping_Hash;
18: /*@
19: ISGetPointRange - Returns a description of the points in an `IS` suitable for traversal
21: Not Collective
23: Input Parameter:
24: . pointIS - The `IS` object
26: Output Parameters:
27: + pStart - The first index, see notes
28: . pEnd - One past the last index, see notes
29: - points - The indices, see notes
31: Level: intermediate
33: Notes:
34: If the `IS` contains contiguous indices in an `ISSTRIDE`, then the indices are contained in [pStart, pEnd) and points = `NULL`.
35: Otherwise, `pStart = 0`, `pEnd = numIndices`, and points is an array of the indices. This supports the following pattern
36: .vb
37: ISGetPointRange(is, &pStart, &pEnd, &points);
38: for (p = pStart; p < pEnd; ++p) {
39: const PetscInt point = points ? points[p] : p;
40: // use point
41: }
42: ISRestorePointRange(is, &pstart, &pEnd, &points);
43: .ve
44: Hence the same code can be written for `pointIS` being a `ISSTRIDE` or `ISGENERAL`
46: .seealso: [](sec_scatter), `IS`, `ISRestorePointRange()`, `ISGetPointSubrange()`, `ISGetIndices()`, `ISCreateStride()`
47: @*/
48: PetscErrorCode ISGetPointRange(IS pointIS, PetscInt *pStart, PetscInt *pEnd, const PetscInt *points[])
49: {
50: PetscInt numCells, step = 1;
51: PetscBool isStride;
53: PetscFunctionBeginHot;
54: *pStart = 0;
55: *points = NULL;
56: PetscCall(ISGetLocalSize(pointIS, &numCells));
57: PetscCall(PetscObjectTypeCompare((PetscObject)pointIS, ISSTRIDE, &isStride));
58: if (isStride) PetscCall(ISStrideGetInfo(pointIS, pStart, &step));
59: *pEnd = *pStart + numCells;
60: if (!isStride || step != 1) PetscCall(ISGetIndices(pointIS, points));
61: PetscFunctionReturn(PETSC_SUCCESS);
62: }
64: /*@
65: ISRestorePointRange - Destroys the traversal description created with `ISGetPointRange()`
67: Not Collective
69: Input Parameters:
70: + pointIS - The `IS` object
71: . pStart - The first index, from `ISGetPointRange()`
72: . pEnd - One past the last index, from `ISGetPointRange()`
73: - points - The indices, from `ISGetPointRange()`
75: Level: intermediate
77: .seealso: [](sec_scatter), `IS`, `ISGetPointRange()`, `ISGetPointSubrange()`, `ISGetIndices()`, `ISCreateStride()`
78: @*/
79: PetscErrorCode ISRestorePointRange(IS pointIS, PetscInt *pStart, PetscInt *pEnd, const PetscInt *points[])
80: {
81: PetscInt step = 1;
82: PetscBool isStride;
84: PetscFunctionBeginHot;
85: PetscCall(PetscObjectTypeCompare((PetscObject)pointIS, ISSTRIDE, &isStride));
86: if (isStride) PetscCall(ISStrideGetInfo(pointIS, pStart, &step));
87: if (!isStride || step != 1) PetscCall(ISGetIndices(pointIS, points));
88: PetscFunctionReturn(PETSC_SUCCESS);
89: }
91: /*@
92: ISGetPointSubrange - Configures the input `IS` to be a subrange for the traversal information given
94: Not Collective
96: Input Parameters:
97: + subpointIS - The `IS` object to be configured
98: . pStart - The first index of the subrange
99: . pEnd - One past the last index for the subrange
100: - points - The indices for the entire range, from `ISGetPointRange()`
102: Output Parameters:
103: . subpointIS - The `IS` object now configured to be a subrange
105: Level: intermediate
107: Note:
108: The input `IS` will now respond properly to calls to `ISGetPointRange()` and return the subrange.
110: .seealso: [](sec_scatter), `IS`, `ISGetPointRange()`, `ISRestorePointRange()`, `ISGetIndices()`, `ISCreateStride()`
111: @*/
112: PetscErrorCode ISGetPointSubrange(IS subpointIS, PetscInt pStart, PetscInt pEnd, const PetscInt points[])
113: {
114: PetscFunctionBeginHot;
115: if (points) {
116: PetscCall(ISSetType(subpointIS, ISGENERAL));
117: PetscCall(ISGeneralSetIndices(subpointIS, pEnd - pStart, &points[pStart], PETSC_USE_POINTER));
118: } else {
119: PetscCall(ISSetType(subpointIS, ISSTRIDE));
120: PetscCall(ISStrideSetStride(subpointIS, pEnd - pStart, pStart, 1));
121: }
122: PetscFunctionReturn(PETSC_SUCCESS);
123: }
125: /*
126: Creates the global mapping information in the ISLocalToGlobalMapping structure
128: If the user has not selected how to handle the global to local mapping then use HASH for "large" problems
129: */
130: static PetscErrorCode ISGlobalToLocalMappingSetUp(ISLocalToGlobalMapping mapping)
131: {
132: PetscInt i, *idx = mapping->indices, n = mapping->n, end, start;
134: PetscFunctionBegin;
135: if (mapping->data) PetscFunctionReturn(PETSC_SUCCESS);
136: end = 0;
137: start = PETSC_INT_MAX;
139: for (i = 0; i < n; i++) {
140: if (idx[i] < 0) continue;
141: if (idx[i] < start) start = idx[i];
142: if (idx[i] > end) end = idx[i];
143: }
144: if (start > end) {
145: start = 0;
146: end = -1;
147: }
148: mapping->globalstart = start;
149: mapping->globalend = end;
150: if (!((PetscObject)mapping)->type_name) {
151: if ((end - start) > PetscMax(4 * n, 1000000)) {
152: PetscCall(ISLocalToGlobalMappingSetType(mapping, ISLOCALTOGLOBALMAPPINGHASH));
153: } else {
154: PetscCall(ISLocalToGlobalMappingSetType(mapping, ISLOCALTOGLOBALMAPPINGBASIC));
155: }
156: }
157: PetscTryTypeMethod(mapping, globaltolocalmappingsetup);
158: PetscFunctionReturn(PETSC_SUCCESS);
159: }
161: static PetscErrorCode ISGlobalToLocalMappingSetUp_Basic(ISLocalToGlobalMapping mapping)
162: {
163: PetscInt i, *idx = mapping->indices, n = mapping->n, end, start, *globals;
164: ISLocalToGlobalMapping_Basic *map;
166: PetscFunctionBegin;
167: start = mapping->globalstart;
168: end = mapping->globalend;
169: PetscCall(PetscNew(&map));
170: PetscCall(PetscMalloc1(end - start + 2, &globals));
171: map->globals = globals;
172: for (i = 0; i < end - start + 1; i++) globals[i] = -1;
173: for (i = 0; i < n; i++) {
174: if (idx[i] < 0) continue;
175: globals[idx[i] - start] = i;
176: }
177: mapping->data = (void *)map;
178: PetscFunctionReturn(PETSC_SUCCESS);
179: }
181: static PetscErrorCode ISGlobalToLocalMappingSetUp_Hash(ISLocalToGlobalMapping mapping)
182: {
183: PetscInt i, *idx = mapping->indices, n = mapping->n;
184: ISLocalToGlobalMapping_Hash *map;
186: PetscFunctionBegin;
187: PetscCall(PetscNew(&map));
188: PetscCall(PetscHMapICreate(&map->globalht));
189: for (i = 0; i < n; i++) {
190: if (idx[i] < 0) continue;
191: PetscCall(PetscHMapISet(map->globalht, idx[i], i));
192: }
193: mapping->data = (void *)map;
194: PetscFunctionReturn(PETSC_SUCCESS);
195: }
197: static PetscErrorCode ISLocalToGlobalMappingDestroy_Basic(ISLocalToGlobalMapping mapping)
198: {
199: ISLocalToGlobalMapping_Basic *map = (ISLocalToGlobalMapping_Basic *)mapping->data;
201: PetscFunctionBegin;
202: if (!map) PetscFunctionReturn(PETSC_SUCCESS);
203: PetscCall(PetscFree(map->globals));
204: PetscCall(PetscFree(mapping->data));
205: PetscFunctionReturn(PETSC_SUCCESS);
206: }
208: static PetscErrorCode ISLocalToGlobalMappingDestroy_Hash(ISLocalToGlobalMapping mapping)
209: {
210: ISLocalToGlobalMapping_Hash *map = (ISLocalToGlobalMapping_Hash *)mapping->data;
212: PetscFunctionBegin;
213: if (!map) PetscFunctionReturn(PETSC_SUCCESS);
214: PetscCall(PetscHMapIDestroy(&map->globalht));
215: PetscCall(PetscFree(mapping->data));
216: PetscFunctionReturn(PETSC_SUCCESS);
217: }
219: static PetscErrorCode ISLocalToGlobalMappingResetBlockInfo_Private(ISLocalToGlobalMapping mapping)
220: {
221: PetscFunctionBegin;
222: PetscCall(PetscFree(mapping->info_procs));
223: PetscCall(PetscFree(mapping->info_numprocs));
224: if (mapping->info_indices) {
225: for (PetscInt i = 0; i < mapping->info_nproc; i++) PetscCall(PetscFree(mapping->info_indices[i]));
226: PetscCall(PetscFree(mapping->info_indices));
227: }
228: if (mapping->info_nodei) PetscCall(PetscFree(mapping->info_nodei[0]));
229: PetscCall(PetscFree2(mapping->info_nodec, mapping->info_nodei));
230: PetscCall(PetscSFDestroy(&mapping->multileaves_sf));
231: PetscFunctionReturn(PETSC_SUCCESS);
232: }
234: #define GTOLTYPE _Basic
235: #define GTOLNAME _Basic
236: #define GTOLBS mapping->bs
237: #define GTOL(g, local) \
238: do { \
239: local = map->globals[g / bs - start]; \
240: if (local >= 0) local = bs * local + (g % bs); \
241: } while (0)
243: #include <../src/vec/is/utils/isltog.h>
245: #define GTOLTYPE _Basic
246: #define GTOLNAME Block_Basic
247: #define GTOLBS 1
248: #define GTOL(g, local) \
249: do { \
250: local = map->globals[g - start]; \
251: } while (0)
252: #include <../src/vec/is/utils/isltog.h>
254: #define GTOLTYPE _Hash
255: #define GTOLNAME _Hash
256: #define GTOLBS mapping->bs
257: #define GTOL(g, local) \
258: do { \
259: (void)PetscHMapIGet(map->globalht, g / bs, &local); \
260: if (local >= 0) local = bs * local + (g % bs); \
261: } while (0)
262: #include <../src/vec/is/utils/isltog.h>
264: #define GTOLTYPE _Hash
265: #define GTOLNAME Block_Hash
266: #define GTOLBS 1
267: #define GTOL(g, local) \
268: do { \
269: (void)PetscHMapIGet(map->globalht, g, &local); \
270: } while (0)
271: #include <../src/vec/is/utils/isltog.h>
273: /*@
274: ISLocalToGlobalMappingDuplicate - Duplicates the local to global mapping object
276: Not Collective
278: Input Parameter:
279: . ltog - local to global mapping
281: Output Parameter:
282: . nltog - the duplicated local to global mapping
284: Level: advanced
286: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`
287: @*/
288: PetscErrorCode ISLocalToGlobalMappingDuplicate(ISLocalToGlobalMapping ltog, ISLocalToGlobalMapping *nltog)
289: {
290: ISLocalToGlobalMappingType l2gtype;
292: PetscFunctionBegin;
294: PetscCall(ISLocalToGlobalMappingCreate(PetscObjectComm((PetscObject)ltog), ltog->bs, ltog->n, ltog->indices, PETSC_COPY_VALUES, nltog));
295: PetscCall(ISLocalToGlobalMappingGetType(ltog, &l2gtype));
296: PetscCall(ISLocalToGlobalMappingSetType(*nltog, l2gtype));
297: PetscFunctionReturn(PETSC_SUCCESS);
298: }
300: /*@
301: ISLocalToGlobalMappingGetSize - Gets the local size of a local to global mapping
303: Not Collective
305: Input Parameter:
306: . mapping - local to global mapping
308: Output Parameter:
309: . n - the number of entries in the local mapping, `ISLocalToGlobalMappingGetIndices()` returns an array of this length
311: Level: advanced
313: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`
314: @*/
315: PetscErrorCode ISLocalToGlobalMappingGetSize(ISLocalToGlobalMapping mapping, PetscInt *n)
316: {
317: PetscFunctionBegin;
319: PetscAssertPointer(n, 2);
320: *n = mapping->bs * mapping->n;
321: PetscFunctionReturn(PETSC_SUCCESS);
322: }
324: /*@
325: ISLocalToGlobalMappingViewFromOptions - View an `ISLocalToGlobalMapping` based on values in the options database
327: Collective
329: Input Parameters:
330: + A - the local to global mapping object
331: . obj - optional object that provides the options prefix used for the options database query, pass `NULL` to use the options prefix of `A`
332: - name - command line option
334: Options Database Key:
335: . -name viewer_specification - See `PetscOptionsCreateViewer()` for the values of `viewer_specification`
337: Level: intermediate
339: Note:
340: This checks the options database, creates the viewer on-the-fly, uses it and then destroys it. Hence it should not be called in heavily used routines,
341: rather `PetscOptionsCreateViewer()` should be used to construct the viewer once which can then be utilized in the heavily used routine.
343: .seealso: [](sec_scatter), `PetscViewer`, `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingView()`, `PetscObjectViewFromOptions()`, `ISLocalToGlobalMappingCreate()`, `PetscOptionsCreateViewer()`
344: @*/
345: PetscErrorCode ISLocalToGlobalMappingViewFromOptions(ISLocalToGlobalMapping A, PetscObject obj, const char name[])
346: {
347: PetscFunctionBegin;
349: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
350: PetscFunctionReturn(PETSC_SUCCESS);
351: }
353: /*@
354: ISLocalToGlobalMappingView - View a local to global mapping
356: Collective on viewer
358: Input Parameters:
359: + mapping - local to global mapping
360: - viewer - viewer
362: Level: intermediate
364: .seealso: [](sec_scatter), `PetscViewer`, `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`
365: @*/
366: PetscErrorCode ISLocalToGlobalMappingView(ISLocalToGlobalMapping mapping, PetscViewer viewer)
367: {
368: PetscBool isascii, isbinary;
369: PetscViewerFormat format;
371: PetscFunctionBegin;
373: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)mapping), &viewer));
376: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
377: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
378: PetscCall(PetscViewerGetFormat(viewer, &format));
379: if (isascii) {
380: if (format == PETSC_VIEWER_ASCII_MATLAB) {
381: const PetscInt *idxs;
382: IS is;
383: const char *name = ((PetscObject)mapping)->name;
384: char iname[PETSC_MAX_PATH_LEN];
386: PetscCall(PetscSNPrintf(iname, sizeof(iname), "%sl2g", name ? name : ""));
387: PetscCall(ISLocalToGlobalMappingGetIndices(mapping, &idxs));
388: PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)viewer), mapping->n * mapping->bs, idxs, PETSC_USE_POINTER, &is));
389: PetscCall(PetscObjectSetName((PetscObject)is, iname));
390: PetscCall(ISView(is, viewer));
391: PetscCall(ISLocalToGlobalMappingRestoreIndices(mapping, &idxs));
392: PetscCall(ISDestroy(&is));
393: } else {
394: PetscMPIInt rank;
396: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)mapping), &rank));
397: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)mapping, viewer));
398: PetscCall(PetscViewerASCIIPushSynchronized(viewer));
399: for (PetscInt i = 0; i < mapping->n; i++) {
400: PetscInt bs = mapping->bs, g = mapping->indices[i];
401: if (bs == 1) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d] %" PetscInt_FMT " %" PetscInt_FMT "\n", rank, i, g));
402: else PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d] %" PetscInt_FMT ":%" PetscInt_FMT " %" PetscInt_FMT ":%" PetscInt_FMT "\n", rank, i * bs, (i + 1) * bs, g * bs, (g + 1) * bs));
403: }
404: PetscCall(PetscViewerFlush(viewer));
405: PetscCall(PetscViewerASCIIPopSynchronized(viewer));
406: }
407: } else if (isbinary) {
408: PetscBool skipHeader;
410: PetscCall(PetscViewerSetUp(viewer));
411: PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &skipHeader));
412: if (!skipHeader) {
413: PetscMPIInt size;
414: PetscInt tr[3];
416: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)viewer), &size));
417: tr[0] = IS_LTOGM_FILE_CLASSID;
418: tr[1] = mapping->bs;
419: tr[2] = size;
420: PetscCall(PetscViewerBinaryWrite(viewer, tr, 3, PETSC_INT));
421: PetscCall(PetscViewerBinaryWriteAll(viewer, &mapping->n, 1, PETSC_DETERMINE, PETSC_DETERMINE, PETSC_INT));
422: }
423: /* write block indices */
424: PetscCall(PetscViewerBinaryWriteAll(viewer, mapping->indices, mapping->n, PETSC_DETERMINE, PETSC_DETERMINE, PETSC_INT));
425: }
426: PetscFunctionReturn(PETSC_SUCCESS);
427: }
429: /*@
430: ISLocalToGlobalMappingLoad - Loads a local-to-global mapping that has been stored in binary format.
432: Collective on viewer
434: Input Parameters:
435: + mapping - the newly loaded map, this needs to have been created with `ISLocalToGlobalMappingCreate()` or some related function before a call to `ISLocalToGlobalMappingLoad()`
436: - viewer - binary file viewer, obtained from `PetscViewerBinaryOpen()`
438: Level: intermediate
440: .seealso: [](sec_scatter), `PetscViewer`, `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingView()`, `ISLocalToGlobalMappingCreate()`
441: @*/
442: PetscErrorCode ISLocalToGlobalMappingLoad(ISLocalToGlobalMapping mapping, PetscViewer viewer)
443: {
444: PetscBool isbinary, skipHeader;
446: PetscFunctionBegin;
449: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
450: PetscCheck(isbinary, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Invalid viewer of type %s", ((PetscObject)viewer)->type_name);
452: /* reset previous data */
453: PetscCall(ISLocalToGlobalMappingResetBlockInfo_Private(mapping));
455: PetscCall(PetscViewerSetUp(viewer));
456: PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &skipHeader));
458: /* When skipping header, it assumes bs and n have been already set */
459: if (!skipHeader) {
460: MPI_Comm comm = PetscObjectComm((PetscObject)viewer);
461: PetscInt tr[3], nold = mapping->n, *sizes, nmaps = PETSC_DECIDE, st = 0;
463: PetscCall(PetscViewerBinaryRead(viewer, tr, 3, NULL, PETSC_INT));
464: PetscCheck(tr[0] == IS_LTOGM_FILE_CLASSID, PETSC_COMM_SELF, PETSC_ERR_FILE_UNEXPECTED, "Not a local-to-global map next in file");
466: mapping->bs = tr[1];
467: PetscCall(PetscMalloc1(tr[2], &sizes));
468: PetscCall(PetscViewerBinaryRead(viewer, sizes, tr[2], NULL, PETSC_INT));
470: /* consume the input, read multiple maps per process if needed */
471: PetscCall(PetscSplitOwnership(comm, &nmaps, &tr[2]));
472: PetscCallMPI(MPI_Exscan(&nmaps, &st, 1, MPIU_INT, MPI_SUM, comm));
473: mapping->n = 0;
474: for (PetscInt i = st; i < st + nmaps; i++) mapping->n += sizes[i];
475: PetscCall(PetscFree(sizes));
477: if (nold != mapping->n) {
478: if (mapping->dealloc_indices) PetscCall(PetscFree(mapping->indices));
479: mapping->indices = NULL;
480: }
481: }
483: /* read indices */
484: if (mapping->n && !mapping->indices) {
485: PetscCall(PetscMalloc1(mapping->n, &mapping->indices));
486: mapping->dealloc_indices = PETSC_TRUE;
487: }
488: PetscCall(PetscViewerBinaryReadAll(viewer, mapping->indices, mapping->n, PETSC_DETERMINE, PETSC_DETERMINE, PETSC_INT));
489: PetscFunctionReturn(PETSC_SUCCESS);
490: }
492: /*@
493: ISLocalToGlobalMappingCreateIS - Creates a mapping between a local (0 to n)
494: ordering and a global parallel ordering.
496: Not Collective
498: Input Parameter:
499: . is - index set containing the global numbers for each local number
501: Output Parameter:
502: . mapping - new mapping data structure
504: Level: advanced
506: Note:
507: the block size of the `IS` determines the block size of the mapping
509: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingSetFromOptions()`
510: @*/
511: PetscErrorCode ISLocalToGlobalMappingCreateIS(IS is, ISLocalToGlobalMapping *mapping)
512: {
513: PetscInt n, bs;
514: const PetscInt *indices;
515: MPI_Comm comm;
516: PetscBool isblock;
518: PetscFunctionBegin;
520: PetscAssertPointer(mapping, 2);
522: PetscCall(PetscObjectGetComm((PetscObject)is, &comm));
523: PetscCall(ISGetLocalSize(is, &n));
524: PetscCall(PetscObjectTypeCompare((PetscObject)is, ISBLOCK, &isblock));
525: if (!isblock) {
526: PetscCall(ISGetIndices(is, &indices));
527: PetscCall(ISLocalToGlobalMappingCreate(comm, 1, n, indices, PETSC_COPY_VALUES, mapping));
528: PetscCall(ISRestoreIndices(is, &indices));
529: } else {
530: PetscCall(ISGetBlockSize(is, &bs));
531: PetscCall(ISBlockGetIndices(is, &indices));
532: PetscCall(ISLocalToGlobalMappingCreate(comm, bs, n / bs, indices, PETSC_COPY_VALUES, mapping));
533: PetscCall(ISBlockRestoreIndices(is, &indices));
534: }
535: PetscFunctionReturn(PETSC_SUCCESS);
536: }
538: /*@
539: ISLocalToGlobalMappingCreateSF - Creates a mapping between a local (0 to n) ordering and a global parallel ordering induced by a star forest.
541: Collective
543: Input Parameters:
544: + sf - star forest mapping contiguous local indices to (rank, offset)
545: - start - first global index on this process, or `PETSC_DECIDE` to compute contiguous global numbering automatically
547: Output Parameter:
548: . mapping - new mapping data structure
550: Level: advanced
552: Note:
553: If a process calls this function with `start` = `PETSC_DECIDE` then all processes must, otherwise the program will hang.
555: .seealso: [](sec_scatter), `PetscSF`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingSetFromOptions()`
556: @*/
557: PetscErrorCode ISLocalToGlobalMappingCreateSF(PetscSF sf, PetscInt start, ISLocalToGlobalMapping *mapping)
558: {
559: PetscInt i, maxlocal, nroots, nleaves, *globals, *ltog;
560: MPI_Comm comm;
562: PetscFunctionBegin;
564: PetscAssertPointer(mapping, 3);
565: PetscCall(PetscObjectGetComm((PetscObject)sf, &comm));
566: PetscCall(PetscSFGetGraph(sf, &nroots, &nleaves, NULL, NULL));
567: if (start == PETSC_DECIDE) {
568: start = 0;
569: PetscCallMPI(MPI_Exscan(&nroots, &start, 1, MPIU_INT, MPI_SUM, comm));
570: } else PetscCheck(start >= 0, comm, PETSC_ERR_ARG_OUTOFRANGE, "start must be nonnegative or PETSC_DECIDE");
571: PetscCall(PetscSFGetLeafRange(sf, NULL, &maxlocal));
572: ++maxlocal;
573: PetscCall(PetscMalloc1(nroots, &globals));
574: PetscCall(PetscMalloc1(maxlocal, <og));
575: for (i = 0; i < nroots; i++) globals[i] = start + i;
576: for (i = 0; i < maxlocal; i++) ltog[i] = -1;
577: PetscCall(PetscSFBcastBegin(sf, MPIU_INT, globals, ltog, MPI_REPLACE));
578: PetscCall(PetscSFBcastEnd(sf, MPIU_INT, globals, ltog, MPI_REPLACE));
579: PetscCall(ISLocalToGlobalMappingCreate(comm, 1, maxlocal, ltog, PETSC_OWN_POINTER, mapping));
580: PetscCall(PetscFree(globals));
581: PetscFunctionReturn(PETSC_SUCCESS);
582: }
584: /*@
585: ISLocalToGlobalMappingSetBlockSize - Sets the blocksize of the mapping
587: Not Collective
589: Input Parameters:
590: + mapping - mapping data structure
591: - bs - the blocksize
593: Level: advanced
595: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`
596: @*/
597: PetscErrorCode ISLocalToGlobalMappingSetBlockSize(ISLocalToGlobalMapping mapping, PetscInt bs)
598: {
599: PetscInt *nid;
600: const PetscInt *oid;
601: PetscInt i, cn, on, obs, nn;
603: PetscFunctionBegin;
605: PetscCheck(bs >= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid block size %" PetscInt_FMT, bs);
606: if (bs == mapping->bs) PetscFunctionReturn(PETSC_SUCCESS);
607: on = mapping->n;
608: obs = mapping->bs;
609: oid = mapping->indices;
610: nn = (on * obs) / bs;
611: PetscCheck((on * obs) % bs == 0, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Block size %" PetscInt_FMT " is inconsistent with block size %" PetscInt_FMT " and number of block indices %" PetscInt_FMT, bs, obs, on);
613: PetscCall(PetscMalloc1(nn, &nid));
614: PetscCall(ISLocalToGlobalMappingGetIndices(mapping, &oid));
615: for (i = 0; i < nn; i++) {
616: PetscInt j;
617: for (j = 0, cn = 0; j < bs - 1; j++) {
618: if (oid[i * bs + j] < 0) {
619: cn++;
620: continue;
621: }
622: PetscCheck(oid[i * bs + j] == oid[i * bs + j + 1] - 1, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Block sizes %" PetscInt_FMT " and %" PetscInt_FMT " are incompatible with the block indices: non consecutive indices %" PetscInt_FMT " %" PetscInt_FMT, bs, obs, oid[i * bs + j], oid[i * bs + j + 1]);
623: }
624: if (oid[i * bs + j] < 0) cn++;
625: if (cn) {
626: PetscCheck(cn == bs, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Block sizes %" PetscInt_FMT " and %" PetscInt_FMT " are incompatible with the block indices: invalid number of negative entries in block %" PetscInt_FMT, bs, obs, cn);
627: nid[i] = -1;
628: } else {
629: nid[i] = oid[i * bs] / bs;
630: }
631: }
632: PetscCall(ISLocalToGlobalMappingRestoreIndices(mapping, &oid));
634: mapping->n = nn;
635: mapping->bs = bs;
636: PetscCall(PetscFree(mapping->indices));
637: mapping->indices = nid;
638: mapping->globalstart = 0;
639: mapping->globalend = 0;
641: /* reset the cached information */
642: PetscCall(ISLocalToGlobalMappingResetBlockInfo_Private(mapping));
643: PetscTryTypeMethod(mapping, destroy);
644: PetscFunctionReturn(PETSC_SUCCESS);
645: }
647: /*@
648: ISLocalToGlobalMappingGetBlockSize - Gets the blocksize of the mapping
649: ordering and a global parallel ordering.
651: Not Collective
653: Input Parameter:
654: . mapping - mapping data structure
656: Output Parameter:
657: . bs - the blocksize
659: Level: advanced
661: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`
662: @*/
663: PetscErrorCode ISLocalToGlobalMappingGetBlockSize(ISLocalToGlobalMapping mapping, PetscInt *bs)
664: {
665: PetscFunctionBegin;
667: *bs = mapping->bs;
668: PetscFunctionReturn(PETSC_SUCCESS);
669: }
671: /*@
672: ISLocalToGlobalMappingCreate - Creates a mapping between a local (0 to n)
673: ordering and a global parallel ordering.
675: Not Collective, but communicator may have more than one process
677: Input Parameters:
678: + comm - MPI communicator
679: . bs - the block size
680: . n - the number of local elements divided by the block size, or equivalently the number of block indices
681: . indices - the global index for each local element, these do not need to be in increasing order (sorted), these values should not be scaled (i.e. multiplied) by the blocksize bs
682: - mode - see PetscCopyMode
684: Output Parameter:
685: . mapping - new mapping data structure
687: Level: advanced
689: Notes:
690: There is one integer value in indices per block and it represents the actual indices bs*idx + j, where j=0,..,bs-1
692: For "small" problems when using `ISGlobalToLocalMappingApply()` and `ISGlobalToLocalMappingApplyBlock()`, the `ISLocalToGlobalMappingType`
693: of `ISLOCALTOGLOBALMAPPINGBASIC` will be used; this uses more memory but is faster; this approach is not scalable for extremely large mappings.
695: For large problems `ISLOCALTOGLOBALMAPPINGHASH` is used, this is scalable.
696: Use `ISLocalToGlobalMappingSetType()` or call `ISLocalToGlobalMappingSetFromOptions()` with the option
697: `-islocaltoglobalmapping_type` <`basic`,`hash`> to control which is used.
699: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingSetFromOptions()`,
700: `ISLOCALTOGLOBALMAPPINGBASIC`, `ISLOCALTOGLOBALMAPPINGHASH`,
701: `ISLocalToGlobalMappingSetType()`, `ISLocalToGlobalMappingType`
702: @*/
703: PetscErrorCode ISLocalToGlobalMappingCreate(MPI_Comm comm, PetscInt bs, PetscInt n, const PetscInt indices[], PetscCopyMode mode, ISLocalToGlobalMapping *mapping)
704: {
705: PetscInt *in;
707: PetscFunctionBegin;
708: if (n) PetscAssertPointer(indices, 4);
709: PetscAssertPointer(mapping, 6);
711: *mapping = NULL;
712: PetscCall(ISInitializePackage());
714: PetscCall(PetscHeaderCreate(*mapping, IS_LTOGM_CLASSID, "ISLocalToGlobalMapping", "Local to global mapping", "IS", comm, ISLocalToGlobalMappingDestroy, ISLocalToGlobalMappingView));
715: (*mapping)->n = n;
716: (*mapping)->bs = bs;
717: if (mode == PETSC_COPY_VALUES) {
718: PetscCall(PetscMalloc1(n, &in));
719: PetscCall(PetscArraycpy(in, indices, n));
720: (*mapping)->indices = in;
721: (*mapping)->dealloc_indices = PETSC_TRUE;
722: } else if (mode == PETSC_OWN_POINTER) {
723: (*mapping)->indices = (PetscInt *)indices;
724: (*mapping)->dealloc_indices = PETSC_TRUE;
725: } else if (mode == PETSC_USE_POINTER) {
726: (*mapping)->indices = (PetscInt *)indices;
727: } else SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "Invalid mode %d", mode);
728: PetscFunctionReturn(PETSC_SUCCESS);
729: }
731: PetscFunctionList ISLocalToGlobalMappingList = NULL;
733: /*@
734: ISLocalToGlobalMappingSetFromOptions - Set mapping options from the options database.
736: Not Collective
738: Input Parameter:
739: . mapping - mapping data structure
741: Options Database Key:
742: . -islocaltoglobalmapping_type (basic|hash) - nonscalable and scalable versions
744: Level: advanced
746: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`,
747: `ISLocalToGlobalMappingCreateIS()`, `ISLOCALTOGLOBALMAPPINGBASIC`,
748: `ISLOCALTOGLOBALMAPPINGHASH`, `ISLocalToGlobalMappingSetType()`, `ISLocalToGlobalMappingType`
749: @*/
750: PetscErrorCode ISLocalToGlobalMappingSetFromOptions(ISLocalToGlobalMapping mapping)
751: {
752: char type[256];
753: ISLocalToGlobalMappingType defaulttype = "Not set";
754: PetscBool flg;
756: PetscFunctionBegin;
758: PetscCall(ISLocalToGlobalMappingRegisterAll());
759: PetscObjectOptionsBegin((PetscObject)mapping);
760: PetscCall(PetscOptionsFList("-islocaltoglobalmapping_type", "ISLocalToGlobalMapping method", "ISLocalToGlobalMappingSetType", ISLocalToGlobalMappingList, ((PetscObject)mapping)->type_name ? ((PetscObject)mapping)->type_name : defaulttype, type, sizeof(type), &flg));
761: if (flg) PetscCall(ISLocalToGlobalMappingSetType(mapping, type));
762: PetscOptionsEnd();
763: PetscFunctionReturn(PETSC_SUCCESS);
764: }
766: /*@
767: ISLocalToGlobalMappingDestroy - Destroys a mapping between a local (0 to n)
768: ordering and a global parallel ordering.
770: Not Collective
772: Input Parameter:
773: . mapping - mapping data structure
775: Level: advanced
777: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingCreate()`
778: @*/
779: PetscErrorCode ISLocalToGlobalMappingDestroy(ISLocalToGlobalMapping *mapping)
780: {
781: PetscFunctionBegin;
782: if (!*mapping) PetscFunctionReturn(PETSC_SUCCESS);
784: if (--((PetscObject)*mapping)->refct > 0) {
785: *mapping = NULL;
786: PetscFunctionReturn(PETSC_SUCCESS);
787: }
788: if ((*mapping)->dealloc_indices) PetscCall(PetscFree((*mapping)->indices));
789: PetscCall(ISLocalToGlobalMappingResetBlockInfo_Private(*mapping));
790: PetscTryTypeMethod(*mapping, destroy);
791: PetscCall(PetscHeaderDestroy(mapping));
792: PetscFunctionReturn(PETSC_SUCCESS);
793: }
795: /*@
796: ISLocalToGlobalMappingApplyIS - Creates from an `IS` in the local numbering
797: a new index set using the global numbering defined in an `ISLocalToGlobalMapping`
798: context.
800: Collective
802: Input Parameters:
803: + mapping - mapping between local and global numbering
804: - is - index set in local numbering
806: Output Parameter:
807: . newis - index set in global numbering
809: Level: advanced
811: Note:
812: The output `IS` will have the same communicator as the input `IS` as well as the same block size.
814: .seealso: [](sec_scatter), `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingCreate()`,
815: `ISLocalToGlobalMappingDestroy()`, `ISGlobalToLocalMappingApply()`
816: @*/
817: PetscErrorCode ISLocalToGlobalMappingApplyIS(ISLocalToGlobalMapping mapping, IS is, IS *newis)
818: {
819: PetscInt n, *idxout, bs;
820: const PetscInt *idxin;
822: PetscFunctionBegin;
825: PetscAssertPointer(newis, 3);
827: PetscCall(ISGetLocalSize(is, &n));
828: PetscCall(ISGetBlockSize(is, &bs));
829: PetscCall(ISGetIndices(is, &idxin));
830: PetscCall(PetscMalloc1(n, &idxout));
831: PetscCall(ISLocalToGlobalMappingApply(mapping, n, idxin, idxout));
832: PetscCall(ISRestoreIndices(is, &idxin));
833: PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)is), n, idxout, PETSC_OWN_POINTER, newis));
834: PetscCall(ISSetBlockSize(*newis, bs));
835: PetscFunctionReturn(PETSC_SUCCESS);
836: }
838: /*@
839: ISLocalToGlobalMappingApply - Takes a list of integers in a local numbering
840: and converts them to the global numbering.
842: Not Collective
844: Input Parameters:
845: + mapping - the local to global mapping context
846: . N - number of integers
847: - in - input indices in local numbering
849: Output Parameter:
850: . out - indices in global numbering
852: Level: advanced
854: Note:
855: The `in` and `out` array parameters may be identical.
857: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingApplyBlock()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingDestroy()`,
858: `ISLocalToGlobalMappingApplyIS()`, `AOCreateBasic()`, `AOApplicationToPetsc()`,
859: `AOPetscToApplication()`, `ISGlobalToLocalMappingApply()`
860: @*/
861: PetscErrorCode ISLocalToGlobalMappingApply(ISLocalToGlobalMapping mapping, PetscInt N, const PetscInt in[], PetscInt out[])
862: {
863: PetscInt i, bs, Nmax;
865: PetscFunctionBegin;
867: bs = mapping->bs;
868: Nmax = bs * mapping->n;
869: if (bs == 1) {
870: const PetscInt *idx = mapping->indices;
871: for (i = 0; i < N; i++) {
872: if (in[i] < 0) {
873: out[i] = in[i];
874: continue;
875: }
876: PetscCheck(in[i] < Nmax, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Local index %" PetscInt_FMT " too large %" PetscInt_FMT " (max) at %" PetscInt_FMT, in[i], Nmax - 1, i);
877: out[i] = idx[in[i]];
878: }
879: } else {
880: const PetscInt *idx = mapping->indices;
881: for (i = 0; i < N; i++) {
882: if (in[i] < 0) {
883: out[i] = in[i];
884: continue;
885: }
886: PetscCheck(in[i] < Nmax, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Local index %" PetscInt_FMT " too large %" PetscInt_FMT " (max) at %" PetscInt_FMT, in[i], Nmax - 1, i);
887: out[i] = idx[in[i] / bs] * bs + (in[i] % bs);
888: }
889: }
890: PetscFunctionReturn(PETSC_SUCCESS);
891: }
893: /*@
894: ISLocalToGlobalMappingApplyBlock - Takes a list of integers in a local block numbering and converts them to the global block numbering
896: Not Collective
898: Input Parameters:
899: + mapping - the local to global mapping context
900: . N - number of integers
901: - in - input indices in local block numbering
903: Output Parameter:
904: . out - indices in global block numbering
906: Example:
907: If the index values are {0,1,6,7} set with a call to `ISLocalToGlobalMappingCreate`(`PETSC_COMM_SELF`,2,2,{0,3}) then the mapping applied to 0
908: (the first block) would produce 0 and the mapping applied to 1 (the second block) would produce 3.
910: Level: advanced
912: Note:
913: The `in` and `out` array parameters may be identical.
915: .seealso: [](sec_scatter), `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingDestroy()`,
916: `ISLocalToGlobalMappingApplyIS()`, `AOCreateBasic()`, `AOApplicationToPetsc()`,
917: `AOPetscToApplication()`, `ISGlobalToLocalMappingApply()`
918: @*/
919: PetscErrorCode ISLocalToGlobalMappingApplyBlock(ISLocalToGlobalMapping mapping, PetscInt N, const PetscInt in[], PetscInt out[])
920: {
921: PetscInt i, Nmax;
922: const PetscInt *idx;
924: PetscFunctionBegin;
926: Nmax = mapping->n;
927: idx = mapping->indices;
928: for (i = 0; i < N; i++) {
929: if (in[i] < 0) {
930: out[i] = in[i];
931: continue;
932: }
933: PetscCheck(in[i] < Nmax, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Local block index %" PetscInt_FMT " too large %" PetscInt_FMT " (max) at %" PetscInt_FMT, in[i], Nmax - 1, i);
934: out[i] = idx[in[i]];
935: }
936: PetscFunctionReturn(PETSC_SUCCESS);
937: }
939: /*@
940: ISGlobalToLocalMappingApply - Provides the local numbering for a list of integers
941: specified with a global numbering.
943: Not Collective
945: Input Parameters:
946: + mapping - mapping between local and global numbering
947: . type - `IS_GTOLM_MASK` - maps global indices with no local value to -1 in the output list (i.e., mask them)
948: `IS_GTOLM_DROP` - drops the indices with no local value from the output list
949: . n - number of global indices to map
950: - idx - global indices to map
952: Output Parameters:
953: + nout - number of indices in output array (if type == `IS_GTOLM_MASK` then nout = n)
954: - idxout - local index of each global index, one must pass in an array long enough
955: to hold all the indices. You can call `ISGlobalToLocalMappingApply()` with
956: idxout == NULL to determine the required length (returned in nout)
957: and then allocate the required space and call `ISGlobalToLocalMappingApply()`
958: a second time to set the values.
960: Level: advanced
962: Notes:
963: Either `nout` or `idxout` may be `NULL`. `idx` and `idxout` may be identical.
965: For "small" problems when using `ISGlobalToLocalMappingApply()` and `ISGlobalToLocalMappingApplyBlock()`, the `ISLocalToGlobalMappingType` of
966: `ISLOCALTOGLOBALMAPPINGBASIC` will be used;
967: this uses more memory but is faster; this approach is not scalable for extremely large mappings. For large problems `ISLOCALTOGLOBALMAPPINGHASH` is used, this is scalable.
968: Use `ISLocalToGlobalMappingSetType()` or call `ISLocalToGlobalMappingSetFromOptions()` with the option -islocaltoglobalmapping_type <basic,hash> to control which is used.
970: Developer Notes:
971: The manual page states that `idx` and `idxout` may be identical but the calling
972: sequence declares `idx` as const so it cannot be the same as `idxout`.
974: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingApply()`, `ISGlobalToLocalMappingApplyBlock()`, `ISLocalToGlobalMappingCreate()`,
975: `ISLocalToGlobalMappingDestroy()`
976: @*/
977: PetscErrorCode ISGlobalToLocalMappingApply(ISLocalToGlobalMapping mapping, ISGlobalToLocalMappingMode type, PetscInt n, const PetscInt idx[], PetscInt *nout, PetscInt idxout[])
978: {
979: PetscFunctionBegin;
981: if (!mapping->data) PetscCall(ISGlobalToLocalMappingSetUp(mapping));
982: PetscUseTypeMethod(mapping, globaltolocalmappingapply, type, n, idx, nout, idxout);
983: PetscFunctionReturn(PETSC_SUCCESS);
984: }
986: /*@
987: ISGlobalToLocalMappingApplyIS - Creates from an `IS` in the global numbering
988: a new index set using the local numbering defined in an `ISLocalToGlobalMapping`
989: context.
991: Not Collective
993: Input Parameters:
994: + mapping - mapping between local and global numbering
995: . type - `IS_GTOLM_MASK` - maps global indices with no local value to -1 in the output list (i.e., mask them)
996: `IS_GTOLM_DROP` - drops the indices with no local value from the output list
997: - is - index set in global numbering
999: Output Parameter:
1000: . newis - index set in local numbering
1002: Level: advanced
1004: Notes:
1005: The output `IS` will be sequential, as it encodes a purely local operation
1007: If `type` is `IS_GTOLM_MASK`, `newis` will have the same block size as `is`
1009: .seealso: [](sec_scatter), `ISGlobalToLocalMapping`, `ISGlobalToLocalMappingApply()`, `ISLocalToGlobalMappingCreate()`,
1010: `ISLocalToGlobalMappingDestroy()`
1011: @*/
1012: PetscErrorCode ISGlobalToLocalMappingApplyIS(ISLocalToGlobalMapping mapping, ISGlobalToLocalMappingMode type, IS is, IS *newis)
1013: {
1014: PetscInt n, nout, *idxout, bs;
1015: const PetscInt *idxin;
1017: PetscFunctionBegin;
1020: PetscAssertPointer(newis, 4);
1022: PetscCall(ISGetLocalSize(is, &n));
1023: PetscCall(ISGetIndices(is, &idxin));
1024: if (type == IS_GTOLM_MASK) PetscCall(PetscMalloc1(n, &idxout));
1025: else {
1026: PetscCall(ISGlobalToLocalMappingApply(mapping, type, n, idxin, &nout, NULL));
1027: PetscCall(PetscMalloc1(nout, &idxout));
1028: }
1029: PetscCall(ISGlobalToLocalMappingApply(mapping, type, n, idxin, &nout, idxout));
1030: PetscCall(ISRestoreIndices(is, &idxin));
1031: PetscCall(ISCreateGeneral(PETSC_COMM_SELF, nout, idxout, PETSC_OWN_POINTER, newis));
1032: if (type == IS_GTOLM_MASK) {
1033: PetscCall(ISGetBlockSize(is, &bs));
1034: PetscCall(ISSetBlockSize(*newis, bs));
1035: }
1036: PetscFunctionReturn(PETSC_SUCCESS);
1037: }
1039: /*@
1040: ISGlobalToLocalMappingApplyBlock - Provides the local block numbering for a list of integers
1041: specified with a block global numbering.
1043: Not Collective
1045: Input Parameters:
1046: + mapping - mapping between local and global numbering
1047: . type - `IS_GTOLM_MASK` - maps global indices with no local value to -1 in the output list (i.e., mask them)
1048: `IS_GTOLM_DROP` - drops the indices with no local value from the output list
1049: . n - number of global indices to map
1050: - idx - global indices to map
1052: Output Parameters:
1053: + nout - number of indices in output array (if type == `IS_GTOLM_MASK` then nout = n)
1054: - idxout - local index of each global index, one must pass in an array long enough
1055: to hold all the indices. You can call `ISGlobalToLocalMappingApplyBlock()` with
1056: idxout == NULL to determine the required length (returned in nout)
1057: and then allocate the required space and call `ISGlobalToLocalMappingApplyBlock()`
1058: a second time to set the values.
1060: Level: advanced
1062: Notes:
1063: Either `nout` or `idxout` may be `NULL`. `idx` and `idxout` may be identical.
1065: For "small" problems when using `ISGlobalToLocalMappingApply()` and `ISGlobalToLocalMappingApplyBlock()`, the `ISLocalToGlobalMappingType` of
1066: `ISLOCALTOGLOBALMAPPINGBASIC` will be used;
1067: this uses more memory but is faster; this approach is not scalable for extremely large mappings. For large problems `ISLOCALTOGLOBALMAPPINGHASH` is used, this is scalable.
1068: Use `ISLocalToGlobalMappingSetType()` or call `ISLocalToGlobalMappingSetFromOptions()` with the option -islocaltoglobalmapping_type <basic,hash> to control which is used.
1070: Developer Notes:
1071: The manual page states that `idx` and `idxout` may be identical but the calling
1072: sequence declares `idx` as const so it cannot be the same as `idxout`.
1074: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingApply()`, `ISGlobalToLocalMappingApply()`, `ISLocalToGlobalMappingCreate()`,
1075: `ISLocalToGlobalMappingDestroy()`
1076: @*/
1077: PetscErrorCode ISGlobalToLocalMappingApplyBlock(ISLocalToGlobalMapping mapping, ISGlobalToLocalMappingMode type, PetscInt n, const PetscInt idx[], PetscInt *nout, PetscInt idxout[])
1078: {
1079: PetscFunctionBegin;
1081: if (!mapping->data) PetscCall(ISGlobalToLocalMappingSetUp(mapping));
1082: PetscUseTypeMethod(mapping, globaltolocalmappingapplyblock, type, n, idx, nout, idxout);
1083: PetscFunctionReturn(PETSC_SUCCESS);
1084: }
1086: /*@
1087: ISLocalToGlobalMappingGetBlockInfo - Gets the neighbor information
1089: Collective the first time it is called
1091: Input Parameter:
1092: . mapping - the mapping from local to global indexing
1094: Output Parameters:
1095: + nproc - number of processes that are connected to the calling process
1096: . procs - neighboring processes
1097: . numprocs - number of block indices for each process
1098: - indices - block indices (in local numbering) shared with neighbors (sorted by global numbering)
1100: Level: advanced
1102: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1103: `ISLocalToGlobalMappingRestoreBlockInfo()`, `ISLocalToGlobalMappingGetBlockMultiLeavesSF()`
1104: @*/
1105: PetscErrorCode ISLocalToGlobalMappingGetBlockInfo(ISLocalToGlobalMapping mapping, PetscInt *nproc, PetscInt *procs[], PetscInt *numprocs[], PetscInt **indices[])
1106: {
1107: PetscFunctionBegin;
1109: PetscCall(ISLocalToGlobalMappingSetUpBlockInfo_Private(mapping));
1110: if (nproc) *nproc = mapping->info_nproc;
1111: if (procs) *procs = mapping->info_procs;
1112: if (numprocs) *numprocs = mapping->info_numprocs;
1113: if (indices) *indices = mapping->info_indices;
1114: PetscFunctionReturn(PETSC_SUCCESS);
1115: }
1117: /*@
1118: ISLocalToGlobalMappingGetBlockNodeInfo - Gets the neighbor information for each local block index
1120: Collective the first time it is called
1122: Input Parameter:
1123: . mapping - the mapping from local to global indexing
1125: Output Parameters:
1126: + n - number of local block nodes
1127: . n_procs - an array storing the number of processes for each local block node (including self)
1128: - procs - the processes' rank for each local block node (sorted, self is first)
1130: Level: advanced
1132: Notes:
1133: The user needs to call `ISLocalToGlobalMappingRestoreBlockNodeInfo()` when the data is no longer needed.
1134: The information returned by this function complements that of `ISLocalToGlobalMappingGetBlockInfo()`.
1135: The latter only provides local information, and the neighboring information
1136: cannot be inferred in the general case, unless the mapping is locally one-to-one on each process.
1138: .seealso: `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1139: `ISLocalToGlobalMappingGetBlockInfo()`, `ISLocalToGlobalMappingRestoreBlockNodeInfo()`, `ISLocalToGlobalMappingGetNodeInfo()`
1140: @*/
1141: PetscErrorCode ISLocalToGlobalMappingGetBlockNodeInfo(ISLocalToGlobalMapping mapping, PetscInt *n, PetscInt *n_procs[], PetscInt **procs[])
1142: {
1143: PetscFunctionBegin;
1145: PetscCall(ISLocalToGlobalMappingSetUpBlockInfo_Private(mapping));
1146: if (n) *n = mapping->n;
1147: if (n_procs) *n_procs = mapping->info_nodec;
1148: if (procs) *procs = mapping->info_nodei;
1149: PetscFunctionReturn(PETSC_SUCCESS);
1150: }
1152: /*@
1153: ISLocalToGlobalMappingRestoreBlockNodeInfo - Frees the memory allocated by `ISLocalToGlobalMappingGetBlockNodeInfo()`
1155: Not Collective
1157: Input Parameters:
1158: + mapping - the mapping from local to global indexing
1159: . n - number of local block nodes
1160: . n_procs - an array storing the number of processes for each local block nodes (including self)
1161: - procs - the processes' rank for each local block node (sorted, self is first)
1163: Level: advanced
1165: .seealso: `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1166: `ISLocalToGlobalMappingGetBlockNodeInfo()`
1167: @*/
1168: PetscErrorCode ISLocalToGlobalMappingRestoreBlockNodeInfo(ISLocalToGlobalMapping mapping, PetscInt *n, PetscInt *n_procs[], PetscInt **procs[])
1169: {
1170: PetscFunctionBegin;
1172: if (n) *n = 0;
1173: if (n_procs) *n_procs = NULL;
1174: if (procs) *procs = NULL;
1175: PetscFunctionReturn(PETSC_SUCCESS);
1176: }
1178: /*@
1179: ISLocalToGlobalMappingGetBlockMultiLeavesSF - Get the star-forest to communicate multi-leaf block data
1181: Collective the first time it is called
1183: Input Parameter:
1184: . mapping - the mapping from local to global indexing
1186: Output Parameter:
1187: . mlsf - the `PetscSF`
1189: Level: advanced
1191: Notes:
1192: The returned star forest is suitable to exchange local information with other processes sharing the same global block index.
1193: For example, suppose a mapping with two processes has been created with
1194: .vb
1195: rank 0 global block indices: [0, 1, 2]
1196: rank 1 global block indices: [2, 3, 4]
1197: .ve
1198: and we want to share the local information
1199: .vb
1200: rank 0 data: [-1, -2, -3]
1201: rank 1 data: [1, 2, 3]
1202: .ve
1203: then, the broadcasting action of `mlsf` will allow to collect
1204: .vb
1205: rank 0 mlleafdata: [-1, -2, -3, 3]
1206: rank 1 mlleafdata: [-3, 3, 1, 2]
1207: .ve
1208: Use ``ISLocalToGlobalMappingGetBlockNodeInfo()`` to index into the multi-leaf data.
1210: .seealso: [](sec_scatter), `ISLocalToGlobalMappingGetBlockNodeInfo()`, `PetscSF`
1211: @*/
1212: PetscErrorCode ISLocalToGlobalMappingGetBlockMultiLeavesSF(ISLocalToGlobalMapping mapping, PetscSF *mlsf)
1213: {
1214: PetscFunctionBegin;
1216: PetscAssertPointer(mlsf, 2);
1217: PetscCall(ISLocalToGlobalMappingSetUpBlockInfo_Private(mapping));
1218: *mlsf = mapping->multileaves_sf;
1219: PetscFunctionReturn(PETSC_SUCCESS);
1220: }
1222: static PetscErrorCode ISLocalToGlobalMappingSetUpBlockInfo_Private(ISLocalToGlobalMapping mapping)
1223: {
1224: PetscSF sf, sf2, imsf, msf;
1225: MPI_Comm comm;
1226: const PetscSFNode *sfnode;
1227: PetscSFNode *newsfnode;
1228: PetscLayout layout;
1229: PetscHMapI neighs;
1230: PetscHashIter iter;
1231: PetscBool missing;
1232: const PetscInt *gidxs, *rootdegree;
1233: PetscInt *mask, *mrootdata, *leafdata, *newleafdata, *leafrd, *tmpg;
1234: PetscInt nroots, nleaves, newnleaves, bs, i, j, m, mnroots, p;
1235: PetscMPIInt rank, size;
1237: PetscFunctionBegin;
1238: if (mapping->multileaves_sf) PetscFunctionReturn(PETSC_SUCCESS);
1239: PetscCall(PetscObjectGetComm((PetscObject)mapping, &comm));
1240: PetscCallMPI(MPI_Comm_size(comm, &size));
1241: PetscCallMPI(MPI_Comm_rank(comm, &rank));
1243: /* Get mapping indices */
1244: PetscCall(ISLocalToGlobalMappingGetBlockSize(mapping, &bs));
1245: PetscCall(ISLocalToGlobalMappingGetBlockIndices(mapping, &gidxs));
1246: PetscCall(ISLocalToGlobalMappingGetSize(mapping, &nleaves));
1247: nleaves /= bs;
1249: /* Create layout for global indices */
1250: for (i = 0, m = 0; i < nleaves; i++) m = PetscMax(m, gidxs[i]);
1251: PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &m, 1, MPIU_INT, MPI_MAX, comm));
1252: PetscCall(PetscLayoutCreate(comm, &layout));
1253: PetscCall(PetscLayoutSetSize(layout, m + 1));
1254: PetscCall(PetscLayoutSetUp(layout));
1256: /* Create SF to share global indices */
1257: PetscCall(PetscSFCreate(comm, &sf));
1258: PetscCall(PetscSFSetGraphLayout(sf, layout, nleaves, NULL, PETSC_OWN_POINTER, gidxs));
1259: PetscCall(PetscSFSetUp(sf));
1260: PetscCall(PetscLayoutDestroy(&layout));
1262: /* communicate root degree to leaves */
1263: PetscCall(PetscSFGetGraph(sf, &nroots, NULL, NULL, &sfnode));
1264: PetscCall(PetscSFComputeDegreeBegin(sf, &rootdegree));
1265: PetscCall(PetscSFComputeDegreeEnd(sf, &rootdegree));
1266: for (i = 0, mnroots = 0; i < nroots; i++) mnroots += rootdegree[i];
1267: PetscCall(PetscMalloc3(2 * PetscMax(mnroots, nroots), &mrootdata, 2 * nleaves, &leafdata, nleaves, &leafrd));
1268: for (i = 0, m = 0; i < nroots; i++) {
1269: mrootdata[2 * i + 0] = rootdegree[i];
1270: mrootdata[2 * i + 1] = m;
1271: m += rootdegree[i];
1272: }
1273: PetscCall(PetscSFBcastBegin(sf, MPIU_2INT, mrootdata, leafdata, MPI_REPLACE));
1274: PetscCall(PetscSFBcastEnd(sf, MPIU_2INT, mrootdata, leafdata, MPI_REPLACE));
1276: /* allocate enough space to store ranks */
1277: for (i = 0, newnleaves = 0; i < nleaves; i++) {
1278: newnleaves += leafdata[2 * i];
1279: leafrd[i] = leafdata[2 * i];
1280: }
1282: /* create new SF nodes to collect multi-root data at leaves */
1283: PetscCall(PetscMalloc1(newnleaves, &newsfnode));
1284: for (i = 0, m = 0; i < nleaves; i++) {
1285: for (j = 0; j < leafrd[i]; j++) {
1286: newsfnode[m].rank = sfnode[i].rank;
1287: newsfnode[m].index = leafdata[2 * i + 1] + j;
1288: m++;
1289: }
1290: }
1292: /* gather ranks at multi roots */
1293: for (i = 0; i < mnroots; i++) mrootdata[i] = -1;
1294: for (i = 0; i < nleaves; i++) leafdata[i] = rank;
1296: PetscCall(PetscSFGatherBegin(sf, MPIU_INT, leafdata, mrootdata));
1297: PetscCall(PetscSFGatherEnd(sf, MPIU_INT, leafdata, mrootdata));
1299: /* from multi-roots to multi-leaves */
1300: PetscCall(PetscSFCreate(comm, &sf2));
1301: PetscCall(PetscSFSetGraph(sf2, mnroots, newnleaves, NULL, PETSC_OWN_POINTER, newsfnode, PETSC_OWN_POINTER));
1302: PetscCall(PetscSFSetUp(sf2));
1304: /* broadcast multi-root data to multi-leaves */
1305: PetscCall(PetscMalloc1(newnleaves, &newleafdata));
1306: PetscCall(PetscSFBcastBegin(sf2, MPIU_INT, mrootdata, newleafdata, MPI_REPLACE));
1307: PetscCall(PetscSFBcastEnd(sf2, MPIU_INT, mrootdata, newleafdata, MPI_REPLACE));
1309: /* sort sharing ranks */
1310: for (i = 0, m = 0; i < nleaves; i++) {
1311: PetscCall(PetscSortInt(leafrd[i], newleafdata + m));
1312: m += leafrd[i];
1313: }
1315: /* Number of neighbors and their ranks */
1316: PetscCall(PetscHMapICreate(&neighs));
1317: for (i = 0; i < newnleaves; i++) PetscCall(PetscHMapIPut(neighs, newleafdata[i], &iter, &missing));
1318: PetscCall(PetscHMapIGetSize(neighs, &mapping->info_nproc));
1319: PetscCall(PetscMalloc1(mapping->info_nproc, &mapping->info_procs));
1320: PetscCall(PetscHMapIGetKeys(neighs, (i = 0, &i), mapping->info_procs));
1321: for (i = 0; i < mapping->info_nproc; i++) { /* put info for self first */
1322: if (mapping->info_procs[i] == rank) {
1323: PetscInt newr = mapping->info_procs[0];
1325: mapping->info_procs[0] = rank;
1326: mapping->info_procs[i] = newr;
1327: break;
1328: }
1329: }
1330: if (mapping->info_nproc) PetscCall(PetscSortInt(mapping->info_nproc - 1, mapping->info_procs + 1));
1331: PetscCall(PetscHMapIDestroy(&neighs));
1333: /* collect info data */
1334: PetscCall(PetscMalloc1(mapping->info_nproc, &mapping->info_numprocs));
1335: PetscCall(PetscMalloc1(mapping->info_nproc, &mapping->info_indices));
1336: for (i = 0; i < mapping->info_nproc; i++) mapping->info_indices[i] = NULL;
1338: PetscCall(PetscMalloc1(nleaves, &mask));
1339: PetscCall(PetscMalloc1(nleaves, &tmpg));
1340: for (p = 0; p < mapping->info_nproc; p++) {
1341: PetscInt *tmp, trank = mapping->info_procs[p];
1343: PetscCall(PetscMemzero(mask, nleaves * sizeof(*mask)));
1344: for (i = 0, m = 0; i < nleaves; i++) {
1345: for (j = 0; j < leafrd[i]; j++) {
1346: if (newleafdata[m] == trank) mask[i]++;
1347: if (!p && newleafdata[m] != rank) mask[i]++;
1348: m++;
1349: }
1350: }
1351: for (i = 0, m = 0; i < nleaves; i++)
1352: if (mask[i] > (!p ? 1 : 0)) m++;
1354: PetscCall(PetscMalloc1(m, &tmp));
1355: for (i = 0, m = 0; i < nleaves; i++)
1356: if (mask[i] > (!p ? 1 : 0)) {
1357: tmp[m] = i;
1358: tmpg[m] = gidxs[i];
1359: m++;
1360: }
1361: PetscCall(PetscSortIntWithArray(m, tmpg, tmp));
1362: mapping->info_indices[p] = tmp;
1363: mapping->info_numprocs[p] = m;
1364: }
1366: /* Node info */
1367: PetscCall(PetscMalloc2(nleaves, &mapping->info_nodec, nleaves + 1, &mapping->info_nodei));
1368: PetscCall(PetscArraycpy(mapping->info_nodec, leafrd, nleaves));
1369: PetscCall(PetscMalloc1(newnleaves, &mapping->info_nodei[0]));
1370: for (i = 0; i < nleaves - 1; i++) mapping->info_nodei[i + 1] = mapping->info_nodei[i] + mapping->info_nodec[i];
1371: PetscCall(PetscArraycpy(mapping->info_nodei[0], newleafdata, newnleaves));
1373: /* Create SF from leaves to multi-leaves */
1374: PetscCall(PetscSFGetMultiSF(sf, &msf));
1375: PetscCall(PetscSFCreateInverseSF(msf, &imsf));
1376: PetscCall(PetscSFCompose(imsf, sf2, &mapping->multileaves_sf));
1377: PetscCall(PetscSFDestroy(&imsf));
1378: PetscCall(PetscSFDestroy(&sf));
1379: PetscCall(PetscSFDestroy(&sf2));
1381: PetscCall(ISLocalToGlobalMappingRestoreBlockIndices(mapping, &gidxs));
1382: PetscCall(PetscFree(tmpg));
1383: PetscCall(PetscFree(mask));
1384: PetscCall(PetscFree3(mrootdata, leafdata, leafrd));
1385: PetscCall(PetscFree(newleafdata));
1386: PetscFunctionReturn(PETSC_SUCCESS);
1387: }
1389: /*@
1390: ISLocalToGlobalMappingRestoreBlockInfo - Frees the memory allocated by `ISLocalToGlobalMappingGetBlockInfo()`
1392: Not Collective
1394: Input Parameters:
1395: + mapping - the mapping from local to global indexing
1396: . nproc - number of processes that are connected to the calling process
1397: . procs - neighboring processes
1398: . numprocs - number of block indices for each process
1399: - indices - block indices (in local numbering) shared with neighbors (sorted by global numbering)
1401: Level: advanced
1403: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1404: `ISLocalToGlobalMappingGetInfo()`
1405: @*/
1406: PetscErrorCode ISLocalToGlobalMappingRestoreBlockInfo(ISLocalToGlobalMapping mapping, PetscInt *nproc, PetscInt *procs[], PetscInt *numprocs[], PetscInt **indices[])
1407: {
1408: PetscFunctionBegin;
1410: if (nproc) *nproc = 0;
1411: if (procs) *procs = NULL;
1412: if (numprocs) *numprocs = NULL;
1413: if (indices) *indices = NULL;
1414: PetscFunctionReturn(PETSC_SUCCESS);
1415: }
1417: /*@
1418: ISLocalToGlobalMappingGetInfo - Gets the neighbor information for each process
1420: Collective the first time it is called
1422: Input Parameter:
1423: . mapping - the mapping from local to global indexing
1425: Output Parameters:
1426: + nproc - number of processes that are connected to the calling process
1427: . procs - neighboring processes
1428: . numprocs - number of indices for each process
1429: - indices - indices (in local numbering) shared with neighbors (sorted by global numbering)
1431: Level: advanced
1433: Note:
1434: The user needs to call `ISLocalToGlobalMappingRestoreInfo()` when the data is no longer needed.
1436: Fortran Notes:
1437: There is no `ISLocalToGlobalMappingRestoreInfo()` in Fortran. You must make sure that
1438: `procs`[], `numprocs`[] and `indices`[][] are large enough arrays, either by allocating them
1439: dynamically or defining static ones large enough.
1441: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1442: `ISLocalToGlobalMappingRestoreInfo()`, `ISLocalToGlobalMappingGetNodeInfo()`
1443: @*/
1444: PetscErrorCode ISLocalToGlobalMappingGetInfo(ISLocalToGlobalMapping mapping, PetscInt *nproc, PetscInt *procs[], PetscInt *numprocs[], PetscInt **indices[])
1445: {
1446: PetscInt **bindices = NULL, *bnumprocs = NULL, bs, i, j, k, n, *bprocs;
1448: PetscFunctionBegin;
1450: bs = mapping->bs;
1451: PetscCall(ISLocalToGlobalMappingGetBlockInfo(mapping, &n, &bprocs, &bnumprocs, &bindices));
1452: if (bs > 1) { /* we need to expand the cached info */
1453: if (indices) PetscCall(PetscCalloc1(n, indices));
1454: if (numprocs) PetscCall(PetscCalloc1(n, numprocs));
1455: if (indices || numprocs) {
1456: for (i = 0; i < n; i++) {
1457: if (indices) {
1458: PetscCall(PetscMalloc1(bs * bnumprocs[i], &(*indices)[i]));
1459: for (j = 0; j < bnumprocs[i]; j++) {
1460: for (k = 0; k < bs; k++) (*indices)[i][j * bs + k] = bs * bindices[i][j] + k;
1461: }
1462: }
1463: if (numprocs) (*numprocs)[i] = bnumprocs[i] * bs;
1464: }
1465: }
1466: } else {
1467: if (numprocs) *numprocs = bnumprocs;
1468: if (indices) *indices = bindices;
1469: }
1470: if (nproc) *nproc = n;
1471: if (procs) *procs = bprocs;
1472: PetscFunctionReturn(PETSC_SUCCESS);
1473: }
1475: /*@
1476: ISLocalToGlobalMappingRestoreInfo - Frees the memory allocated by `ISLocalToGlobalMappingGetInfo()`
1478: Not Collective
1480: Input Parameters:
1481: + mapping - the mapping from local to global indexing
1482: . nproc - number of processes that are connected to the calling process
1483: . procs - neighboring processes
1484: . numprocs - number of indices for each process
1485: - indices - indices (in local numbering) shared with neighbors (sorted by global numbering)
1487: Level: advanced
1489: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1490: `ISLocalToGlobalMappingGetInfo()`
1491: @*/
1492: PetscErrorCode ISLocalToGlobalMappingRestoreInfo(ISLocalToGlobalMapping mapping, PetscInt *nproc, PetscInt *procs[], PetscInt *numprocs[], PetscInt **indices[])
1493: {
1494: PetscFunctionBegin;
1496: if (mapping->bs > 1) {
1497: if (numprocs) PetscCall(PetscFree(*numprocs));
1498: if (indices) {
1499: if (*indices)
1500: for (PetscInt i = 0; i < *nproc; i++) PetscCall(PetscFree((*indices)[i]));
1501: PetscCall(PetscFree(*indices));
1502: }
1503: }
1504: PetscFunctionReturn(PETSC_SUCCESS);
1505: }
1507: /*@
1508: ISLocalToGlobalMappingGetNodeInfo - Gets the neighbor information of local nodes
1510: Collective the first time it is called
1512: Input Parameter:
1513: . mapping - the mapping from local to global indexing
1515: Output Parameters:
1516: + n - number of local nodes
1517: . n_procs - an array storing the number of processes for each local node (including self)
1518: - procs - the processes' rank for each local node (sorted, self is first)
1520: Level: advanced
1522: Note:
1523: The user needs to call `ISLocalToGlobalMappingRestoreNodeInfo()` when the data is no longer needed.
1525: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1526: `ISLocalToGlobalMappingGetInfo()`, `ISLocalToGlobalMappingRestoreNodeInfo()`, `ISLocalToGlobalMappingGetBlockNodeInfo()`
1527: @*/
1528: PetscErrorCode ISLocalToGlobalMappingGetNodeInfo(ISLocalToGlobalMapping mapping, PetscInt *n, PetscInt *n_procs[], PetscInt **procs[])
1529: {
1530: PetscInt **bprocs = NULL, *bn_procs = NULL, bs, i, j, k, bn;
1532: PetscFunctionBegin;
1534: bs = mapping->bs;
1535: PetscCall(ISLocalToGlobalMappingGetBlockNodeInfo(mapping, &bn, &bn_procs, &bprocs));
1536: if (bs > 1) { /* we need to expand the cached info */
1537: PetscInt *tn_procs;
1538: PetscInt c;
1540: PetscCall(PetscMalloc1(bn * bs, &tn_procs));
1541: for (i = 0, c = 0; i < bn; i++) {
1542: for (k = 0; k < bs; k++) tn_procs[i * bs + k] = bn_procs[i];
1543: c += bs * bn_procs[i];
1544: }
1545: if (n) *n = bn * bs;
1546: if (procs) {
1547: PetscInt **tprocs;
1548: PetscInt tn = bn * bs;
1550: PetscCall(PetscMalloc1(tn, &tprocs));
1551: if (tn) PetscCall(PetscMalloc1(c, &tprocs[0]));
1552: for (i = 0; i < tn - 1; i++) tprocs[i + 1] = tprocs[i] + tn_procs[i];
1553: for (i = 0; i < bn; i++) {
1554: for (k = 0; k < bs; k++) {
1555: for (j = 0; j < bn_procs[i]; j++) tprocs[i * bs + k][j] = bprocs[i][j];
1556: }
1557: }
1558: *procs = tprocs;
1559: }
1560: if (n_procs) *n_procs = tn_procs;
1561: else PetscCall(PetscFree(tn_procs));
1562: } else {
1563: if (n) *n = bn;
1564: if (n_procs) *n_procs = bn_procs;
1565: if (procs) *procs = bprocs;
1566: }
1567: PetscFunctionReturn(PETSC_SUCCESS);
1568: }
1570: /*@
1571: ISLocalToGlobalMappingRestoreNodeInfo - Frees the memory allocated by `ISLocalToGlobalMappingGetNodeInfo()`
1573: Not Collective
1575: Input Parameters:
1576: + mapping - the mapping from local to global indexing
1577: . n - number of local nodes
1578: . n_procs - an array storing the number of processes for each local node (including self)
1579: - procs - the processes' rank for each local node (sorted, self is first)
1581: Level: advanced
1583: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1584: `ISLocalToGlobalMappingGetInfo()`
1585: @*/
1586: PetscErrorCode ISLocalToGlobalMappingRestoreNodeInfo(ISLocalToGlobalMapping mapping, PetscInt *n, PetscInt *n_procs[], PetscInt **procs[])
1587: {
1588: PetscFunctionBegin;
1590: if (mapping->bs > 1) {
1591: if (n_procs) PetscCall(PetscFree(*n_procs));
1592: if (procs) {
1593: if (*procs) PetscCall(PetscFree((*procs)[0]));
1594: PetscCall(PetscFree(*procs));
1595: }
1596: }
1597: PetscCall(ISLocalToGlobalMappingRestoreBlockNodeInfo(mapping, n, n_procs, procs));
1598: PetscFunctionReturn(PETSC_SUCCESS);
1599: }
1601: /*@
1602: ISLocalToGlobalMappingGetIndices - Get global indices for every local point that is mapped
1604: Not Collective
1606: Input Parameter:
1607: . ltog - local to global mapping
1609: Output Parameter:
1610: . array - array of indices, the length of this array may be obtained with `ISLocalToGlobalMappingGetSize()`
1612: Level: advanced
1614: Note:
1615: `ISLocalToGlobalMappingGetSize()` returns the length the this array
1617: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingRestoreIndices()`,
1618: `ISLocalToGlobalMappingGetBlockIndices()`, `ISLocalToGlobalMappingRestoreBlockIndices()`
1619: @*/
1620: PetscErrorCode ISLocalToGlobalMappingGetIndices(ISLocalToGlobalMapping ltog, const PetscInt *array[])
1621: {
1622: PetscFunctionBegin;
1624: PetscAssertPointer(array, 2);
1625: if (ltog->bs == 1) {
1626: *array = ltog->indices;
1627: } else {
1628: PetscInt *jj, k, i, j, n = ltog->n, bs = ltog->bs;
1629: const PetscInt *ii;
1631: PetscCall(PetscMalloc1(bs * n, &jj));
1632: *array = jj;
1633: k = 0;
1634: ii = ltog->indices;
1635: for (i = 0; i < n; i++)
1636: for (j = 0; j < bs; j++) jj[k++] = bs * ii[i] + j;
1637: }
1638: PetscFunctionReturn(PETSC_SUCCESS);
1639: }
1641: /*@
1642: ISLocalToGlobalMappingRestoreIndices - Restore indices obtained with `ISLocalToGlobalMappingGetIndices()`
1644: Not Collective
1646: Input Parameters:
1647: + ltog - local to global mapping
1648: - array - array of indices
1650: Level: advanced
1652: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingGetIndices()`
1653: @*/
1654: PetscErrorCode ISLocalToGlobalMappingRestoreIndices(ISLocalToGlobalMapping ltog, const PetscInt *array[])
1655: {
1656: PetscFunctionBegin;
1658: PetscAssertPointer(array, 2);
1659: PetscCheck(ltog->bs != 1 || *array == ltog->indices, PETSC_COMM_SELF, PETSC_ERR_ARG_BADPTR, "Trying to return mismatched pointer");
1660: if (ltog->bs > 1) PetscCall(PetscFree(*(void **)array));
1661: PetscFunctionReturn(PETSC_SUCCESS);
1662: }
1664: /*@
1665: ISLocalToGlobalMappingGetBlockIndices - Get global indices for every local block in a `ISLocalToGlobalMapping`
1667: Not Collective
1669: Input Parameter:
1670: . ltog - local to global mapping
1672: Output Parameter:
1673: . array - array of indices
1675: Level: advanced
1677: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingApply()`,
1678: `ISLocalToGlobalMappingRestoreBlockIndices()`
1679: @*/
1680: PetscErrorCode ISLocalToGlobalMappingGetBlockIndices(ISLocalToGlobalMapping ltog, const PetscInt *array[])
1681: {
1682: PetscFunctionBegin;
1684: PetscAssertPointer(array, 2);
1685: *array = ltog->indices;
1686: PetscFunctionReturn(PETSC_SUCCESS);
1687: }
1689: /*@
1690: ISLocalToGlobalMappingRestoreBlockIndices - Restore indices obtained with `ISLocalToGlobalMappingGetBlockIndices()`
1692: Not Collective
1694: Input Parameters:
1695: + ltog - local to global mapping
1696: - array - array of indices
1698: Level: advanced
1700: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingGetIndices()`
1701: @*/
1702: PetscErrorCode ISLocalToGlobalMappingRestoreBlockIndices(ISLocalToGlobalMapping ltog, const PetscInt *array[])
1703: {
1704: PetscFunctionBegin;
1706: PetscAssertPointer(array, 2);
1707: PetscCheck(*array == ltog->indices, PETSC_COMM_SELF, PETSC_ERR_ARG_BADPTR, "Trying to return mismatched pointer");
1708: *array = NULL;
1709: PetscFunctionReturn(PETSC_SUCCESS);
1710: }
1712: /*@
1713: ISLocalToGlobalMappingConcatenate - Create a new mapping that concatenates a list of mappings
1715: Not Collective
1717: Input Parameters:
1718: + comm - communicator for the new mapping, must contain the communicator of every mapping to concatenate
1719: . n - number of mappings to concatenate
1720: - ltogs - local to global mappings
1722: Output Parameter:
1723: . ltogcat - new mapping
1725: Level: advanced
1727: Note:
1728: This currently always returns a mapping with block size of 1
1730: Developer Notes:
1731: If all the input mapping have the same block size we could easily handle that as a special case
1733: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingCreate()`
1734: @*/
1735: PetscErrorCode ISLocalToGlobalMappingConcatenate(MPI_Comm comm, PetscInt n, const ISLocalToGlobalMapping ltogs[], ISLocalToGlobalMapping *ltogcat)
1736: {
1737: PetscInt i, cnt, m, *idx;
1739: PetscFunctionBegin;
1740: PetscCheck(n >= 0, comm, PETSC_ERR_ARG_OUTOFRANGE, "Must have a non-negative number of mappings, given %" PetscInt_FMT, n);
1741: if (n > 0) PetscAssertPointer(ltogs, 3);
1743: PetscAssertPointer(ltogcat, 4);
1744: for (cnt = 0, i = 0; i < n; i++) {
1745: PetscCall(ISLocalToGlobalMappingGetSize(ltogs[i], &m));
1746: cnt += m;
1747: }
1748: PetscCall(PetscMalloc1(cnt, &idx));
1749: for (cnt = 0, i = 0; i < n; i++) {
1750: const PetscInt *subidx;
1751: PetscCall(ISLocalToGlobalMappingGetSize(ltogs[i], &m));
1752: PetscCall(ISLocalToGlobalMappingGetIndices(ltogs[i], &subidx));
1753: PetscCall(PetscArraycpy(&idx[cnt], subidx, m));
1754: PetscCall(ISLocalToGlobalMappingRestoreIndices(ltogs[i], &subidx));
1755: cnt += m;
1756: }
1757: PetscCall(ISLocalToGlobalMappingCreate(comm, 1, cnt, idx, PETSC_OWN_POINTER, ltogcat));
1758: PetscFunctionReturn(PETSC_SUCCESS);
1759: }
1761: /*MC
1762: ISLOCALTOGLOBALMAPPINGBASIC - basic implementation of the `ISLocalToGlobalMapping` object. When `ISGlobalToLocalMappingApply()` is
1763: used this is good for only small and moderate size problems.
1765: Options Database Key:
1766: . -islocaltoglobalmapping_type basic - select this method
1768: Level: beginner
1770: Developer Note:
1771: This stores all the mapping information on each MPI rank.
1773: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingSetType()`, `ISLOCALTOGLOBALMAPPINGHASH`
1774: M*/
1775: PETSC_EXTERN PetscErrorCode ISLocalToGlobalMappingCreate_Basic(ISLocalToGlobalMapping ltog)
1776: {
1777: PetscFunctionBegin;
1778: ltog->ops->globaltolocalmappingapply = ISGlobalToLocalMappingApply_Basic;
1779: ltog->ops->globaltolocalmappingsetup = ISGlobalToLocalMappingSetUp_Basic;
1780: ltog->ops->globaltolocalmappingapplyblock = ISGlobalToLocalMappingApplyBlock_Basic;
1781: ltog->ops->destroy = ISLocalToGlobalMappingDestroy_Basic;
1782: PetscFunctionReturn(PETSC_SUCCESS);
1783: }
1785: /*MC
1786: ISLOCALTOGLOBALMAPPINGHASH - hash implementation of the `ISLocalToGlobalMapping` object. When `ISGlobalToLocalMappingApply()` is
1787: used this is good for large memory problems.
1789: Options Database Key:
1790: . -islocaltoglobalmapping_type hash - select this method
1792: Level: beginner
1794: Note:
1795: This is selected automatically for large problems if the user does not set the type.
1797: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingSetType()`, `ISLOCALTOGLOBALMAPPINGBASIC`
1798: M*/
1799: PETSC_EXTERN PetscErrorCode ISLocalToGlobalMappingCreate_Hash(ISLocalToGlobalMapping ltog)
1800: {
1801: PetscFunctionBegin;
1802: ltog->ops->globaltolocalmappingapply = ISGlobalToLocalMappingApply_Hash;
1803: ltog->ops->globaltolocalmappingsetup = ISGlobalToLocalMappingSetUp_Hash;
1804: ltog->ops->globaltolocalmappingapplyblock = ISGlobalToLocalMappingApplyBlock_Hash;
1805: ltog->ops->destroy = ISLocalToGlobalMappingDestroy_Hash;
1806: PetscFunctionReturn(PETSC_SUCCESS);
1807: }
1809: /*@
1810: ISLocalToGlobalMappingRegister - Registers a method for applying a global to local mapping with an `ISLocalToGlobalMapping`
1812: Not Collective, No Fortran Support
1814: Input Parameters:
1815: + sname - name of a new method
1816: - function - routine to create method context
1818: Example Usage:
1819: .vb
1820: ISLocalToGlobalMappingRegister("my_mapper", MyCreate);
1821: .ve
1823: Then, your mapping can be chosen with the procedural interface via
1824: .vb
1825: ISLocalToGlobalMappingSetType(ltog, "my_mapper")
1826: .ve
1827: or at runtime via the option
1828: .vb
1829: -islocaltoglobalmapping_type my_mapper
1830: .ve
1832: Level: advanced
1834: Note:
1835: `ISLocalToGlobalMappingRegister()` may be called multiple times to add several user-defined mappings.
1837: .seealso: [](sec_scatter), `ISLocalToGlobalMappingRegisterAll()`, `ISLocalToGlobalMappingRegisterDestroy()`, `ISLOCALTOGLOBALMAPPINGBASIC`,
1838: `ISLOCALTOGLOBALMAPPINGHASH`, `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingApply()`
1839: @*/
1840: PetscErrorCode ISLocalToGlobalMappingRegister(const char sname[], PetscErrorCode (*function)(ISLocalToGlobalMapping))
1841: {
1842: PetscFunctionBegin;
1843: PetscCall(ISInitializePackage());
1844: PetscCall(PetscFunctionListAdd(&ISLocalToGlobalMappingList, sname, function));
1845: PetscFunctionReturn(PETSC_SUCCESS);
1846: }
1848: /*@
1849: ISLocalToGlobalMappingSetType - Sets the implementation type `ISLocalToGlobalMapping` will use
1851: Logically Collective
1853: Input Parameters:
1854: + ltog - the `ISLocalToGlobalMapping` object
1855: - type - a known method
1857: Options Database Key:
1858: . -islocaltoglobalmapping_type (basic|hash) - Sets the method used for applying the mapping, see `ISLocalToGlobalMappingType`
1860: Level: intermediate
1862: Notes:
1863: See `ISLocalToGlobalMappingType` for available methods
1865: Normally, it is best to use the `ISLocalToGlobalMappingSetFromOptions()` command and
1866: then set the `ISLocalToGlobalMappingType` from the options database rather than by using
1867: this routine.
1869: Developer Notes:
1870: `ISLocalToGlobalMappingRegister()` is used to add new types to `ISLocalToGlobalMappingList` from which they
1871: are accessed by `ISLocalToGlobalMappingSetType()`.
1873: .seealso: [](sec_scatter), `ISLocalToGlobalMappingType`, `ISLocalToGlobalMappingRegister()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingGetType()`
1874: @*/
1875: PetscErrorCode ISLocalToGlobalMappingSetType(ISLocalToGlobalMapping ltog, ISLocalToGlobalMappingType type)
1876: {
1877: PetscBool match;
1878: PetscErrorCode (*r)(ISLocalToGlobalMapping) = NULL;
1880: PetscFunctionBegin;
1882: if (type) PetscAssertPointer(type, 2);
1884: PetscCall(PetscObjectTypeCompare((PetscObject)ltog, type, &match));
1885: if (match) PetscFunctionReturn(PETSC_SUCCESS);
1887: /* L2G maps defer type setup at globaltolocal calls, allow passing NULL here */
1888: if (type) {
1889: PetscCall(PetscFunctionListFind(ISLocalToGlobalMappingList, type, &r));
1890: PetscCheck(r, PetscObjectComm((PetscObject)ltog), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unable to find requested ISLocalToGlobalMapping type %s", type);
1891: }
1892: /* Destroy the previous private LTOG context */
1893: PetscTryTypeMethod(ltog, destroy);
1894: ltog->ops->destroy = NULL;
1896: PetscCall(PetscObjectChangeTypeName((PetscObject)ltog, type));
1897: if (r) PetscCall((*r)(ltog));
1898: PetscFunctionReturn(PETSC_SUCCESS);
1899: }
1901: /*@
1902: ISLocalToGlobalMappingGetType - Get the type of the `ISLocalToGlobalMapping`
1904: Not Collective
1906: Input Parameter:
1907: . ltog - the `ISLocalToGlobalMapping` object
1909: Output Parameter:
1910: . type - the type
1912: Level: intermediate
1914: Note:
1915: `type` should not be retained for later use as it will be an invalid pointer if the `ISLocalToGlobalMappingType` of `ltog` is changed.
1917: .seealso: [](sec_scatter), `ISLocalToGlobalMappingType`, `ISLocalToGlobalMappingRegister()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingSetType()`,
1918: `PetscObjectTypeCompare()`, `PetscObjectTypeCompareAny()`
1919: @*/
1920: PetscErrorCode ISLocalToGlobalMappingGetType(ISLocalToGlobalMapping ltog, ISLocalToGlobalMappingType *type)
1921: {
1922: PetscFunctionBegin;
1924: PetscAssertPointer(type, 2);
1925: *type = ((PetscObject)ltog)->type_name;
1926: PetscFunctionReturn(PETSC_SUCCESS);
1927: }
1929: PetscBool ISLocalToGlobalMappingRegisterAllCalled = PETSC_FALSE;
1931: /*@
1932: ISLocalToGlobalMappingRegisterAll - Registers all of the local to global mapping components in the `IS` package.
1934: Not Collective
1936: Level: advanced
1938: .seealso: [](sec_scatter), `ISRegister()`, `ISLocalToGlobalRegister()`
1939: @*/
1940: PetscErrorCode ISLocalToGlobalMappingRegisterAll(void)
1941: {
1942: PetscFunctionBegin;
1943: if (ISLocalToGlobalMappingRegisterAllCalled) PetscFunctionReturn(PETSC_SUCCESS);
1944: ISLocalToGlobalMappingRegisterAllCalled = PETSC_TRUE;
1945: PetscCall(ISLocalToGlobalMappingRegister(ISLOCALTOGLOBALMAPPINGBASIC, ISLocalToGlobalMappingCreate_Basic));
1946: PetscCall(ISLocalToGlobalMappingRegister(ISLOCALTOGLOBALMAPPINGHASH, ISLocalToGlobalMappingCreate_Hash));
1947: PetscFunctionReturn(PETSC_SUCCESS);
1948: }