Actual source code: isltog.c
1: #include <petsc/private/isimpl.h>
2: #include <petsc/private/hashmapi.h>
3: #include <petscsf.h>
5: PetscClassId IS_LTOGM_CLASSID;
6: static PetscErrorCode ISLocalToGlobalMappingSetUpBlockInfo_Private(ISLocalToGlobalMapping);
8: typedef struct {
9: PetscInt *globals;
10: } ISLocalToGlobalMapping_Basic;
12: typedef struct {
13: PetscHMapI globalht;
14: } ISLocalToGlobalMapping_Hash;
16: /*@
17: ISGetPointRange - Returns a description of the points in an `IS` suitable for traversal
19: Not Collective
21: Input Parameter:
22: . pointIS - The `IS` object
24: Output Parameters:
25: + pStart - The first index, see notes
26: . pEnd - One past the last index, see notes
27: - points - The indices, see notes
29: Level: intermediate
31: Notes:
32: If the `IS` contains contiguous indices in an `ISSTRIDE`, then the indices are contained in [pStart, pEnd) and points = `NULL`.
33: Otherwise, `pStart = 0`, `pEnd = numIndices`, and points is an array of the indices. This supports the following pattern
34: .vb
35: ISGetPointRange(is, &pStart, &pEnd, &points);
36: for (p = pStart; p < pEnd; ++p) {
37: const PetscInt point = points ? points[p] : p;
38: // use point
39: }
40: ISRestorePointRange(is, &pstart, &pEnd, &points);
41: .ve
42: Hence the same code can be written for `pointIS` being a `ISSTRIDE` or `ISGENERAL`
44: .seealso: [](sec_scatter), `IS`, `ISRestorePointRange()`, `ISGetPointSubrange()`, `ISGetIndices()`, `ISCreateStride()`
45: @*/
46: PetscErrorCode ISGetPointRange(IS pointIS, PetscInt *pStart, PetscInt *pEnd, const PetscInt *points[])
47: {
48: PetscInt numCells, step = 1;
49: PetscBool isStride;
51: PetscFunctionBeginHot;
52: *pStart = 0;
53: *points = NULL;
54: PetscCall(ISGetLocalSize(pointIS, &numCells));
55: PetscCall(PetscObjectTypeCompare((PetscObject)pointIS, ISSTRIDE, &isStride));
56: if (isStride) PetscCall(ISStrideGetInfo(pointIS, pStart, &step));
57: *pEnd = *pStart + numCells;
58: if (!isStride || step != 1) PetscCall(ISGetIndices(pointIS, points));
59: PetscFunctionReturn(PETSC_SUCCESS);
60: }
62: /*@
63: ISRestorePointRange - Destroys the traversal description created with `ISGetPointRange()`
65: Not Collective
67: Input Parameters:
68: + pointIS - The `IS` object
69: . pStart - The first index, from `ISGetPointRange()`
70: . pEnd - One past the last index, from `ISGetPointRange()`
71: - points - The indices, from `ISGetPointRange()`
73: Level: intermediate
75: .seealso: [](sec_scatter), `IS`, `ISGetPointRange()`, `ISGetPointSubrange()`, `ISGetIndices()`, `ISCreateStride()`
76: @*/
77: PetscErrorCode ISRestorePointRange(IS pointIS, PetscInt *pStart, PetscInt *pEnd, const PetscInt *points[])
78: {
79: PetscInt step = 1;
80: PetscBool isStride;
82: PetscFunctionBeginHot;
83: PetscCall(PetscObjectTypeCompare((PetscObject)pointIS, ISSTRIDE, &isStride));
84: if (isStride) PetscCall(ISStrideGetInfo(pointIS, pStart, &step));
85: if (!isStride || step != 1) PetscCall(ISGetIndices(pointIS, points));
86: PetscFunctionReturn(PETSC_SUCCESS);
87: }
89: /*@
90: ISGetPointSubrange - Configures the input `IS` to be a subrange for the traversal information given
92: Not Collective
94: Input Parameters:
95: + subpointIS - The `IS` object to be configured
96: . pStart - The first index of the subrange
97: . pEnd - One past the last index for the subrange
98: - points - The indices for the entire range, from `ISGetPointRange()`
100: Output Parameters:
101: . subpointIS - The `IS` object now configured to be a subrange
103: Level: intermediate
105: Note:
106: The input `IS` will now respond properly to calls to `ISGetPointRange()` and return the subrange.
108: .seealso: [](sec_scatter), `IS`, `ISGetPointRange()`, `ISRestorePointRange()`, `ISGetIndices()`, `ISCreateStride()`
109: @*/
110: PetscErrorCode ISGetPointSubrange(IS subpointIS, PetscInt pStart, PetscInt pEnd, const PetscInt points[])
111: {
112: PetscFunctionBeginHot;
113: if (points) {
114: PetscCall(ISSetType(subpointIS, ISGENERAL));
115: PetscCall(ISGeneralSetIndices(subpointIS, pEnd - pStart, &points[pStart], PETSC_USE_POINTER));
116: } else {
117: PetscCall(ISSetType(subpointIS, ISSTRIDE));
118: PetscCall(ISStrideSetStride(subpointIS, pEnd - pStart, pStart, 1));
119: }
120: PetscFunctionReturn(PETSC_SUCCESS);
121: }
123: /*
124: Creates the global mapping information in the ISLocalToGlobalMapping structure
126: If the user has not selected how to handle the global to local mapping then use HASH for "large" problems
127: */
128: static PetscErrorCode ISGlobalToLocalMappingSetUp(ISLocalToGlobalMapping mapping)
129: {
130: PetscInt i, *idx = mapping->indices, n = mapping->n, end, start;
132: PetscFunctionBegin;
133: if (mapping->data) PetscFunctionReturn(PETSC_SUCCESS);
134: end = 0;
135: start = PETSC_INT_MAX;
137: for (i = 0; i < n; i++) {
138: if (idx[i] < 0) continue;
139: if (idx[i] < start) start = idx[i];
140: if (idx[i] > end) end = idx[i];
141: }
142: if (start > end) {
143: start = 0;
144: end = -1;
145: }
146: mapping->globalstart = start;
147: mapping->globalend = end;
148: if (!((PetscObject)mapping)->type_name) {
149: if ((end - start) > PetscMax(4 * n, 1000000)) {
150: PetscCall(ISLocalToGlobalMappingSetType(mapping, ISLOCALTOGLOBALMAPPINGHASH));
151: } else {
152: PetscCall(ISLocalToGlobalMappingSetType(mapping, ISLOCALTOGLOBALMAPPINGBASIC));
153: }
154: }
155: PetscTryTypeMethod(mapping, globaltolocalmappingsetup);
156: PetscFunctionReturn(PETSC_SUCCESS);
157: }
159: static PetscErrorCode ISGlobalToLocalMappingSetUp_Basic(ISLocalToGlobalMapping mapping)
160: {
161: PetscInt i, *idx = mapping->indices, n = mapping->n, end, start, *globals;
162: ISLocalToGlobalMapping_Basic *map;
164: PetscFunctionBegin;
165: start = mapping->globalstart;
166: end = mapping->globalend;
167: PetscCall(PetscNew(&map));
168: PetscCall(PetscMalloc1(end - start + 2, &globals));
169: map->globals = globals;
170: for (i = 0; i < end - start + 1; i++) globals[i] = -1;
171: for (i = 0; i < n; i++) {
172: if (idx[i] < 0) continue;
173: globals[idx[i] - start] = i;
174: }
175: mapping->data = (void *)map;
176: PetscFunctionReturn(PETSC_SUCCESS);
177: }
179: static PetscErrorCode ISGlobalToLocalMappingSetUp_Hash(ISLocalToGlobalMapping mapping)
180: {
181: PetscInt i, *idx = mapping->indices, n = mapping->n;
182: ISLocalToGlobalMapping_Hash *map;
184: PetscFunctionBegin;
185: PetscCall(PetscNew(&map));
186: PetscCall(PetscHMapICreate(&map->globalht));
187: for (i = 0; i < n; i++) {
188: if (idx[i] < 0) continue;
189: PetscCall(PetscHMapISet(map->globalht, idx[i], i));
190: }
191: mapping->data = (void *)map;
192: PetscFunctionReturn(PETSC_SUCCESS);
193: }
195: static PetscErrorCode ISLocalToGlobalMappingDestroy_Basic(ISLocalToGlobalMapping mapping)
196: {
197: ISLocalToGlobalMapping_Basic *map = (ISLocalToGlobalMapping_Basic *)mapping->data;
199: PetscFunctionBegin;
200: if (!map) PetscFunctionReturn(PETSC_SUCCESS);
201: PetscCall(PetscFree(map->globals));
202: PetscCall(PetscFree(mapping->data));
203: PetscFunctionReturn(PETSC_SUCCESS);
204: }
206: static PetscErrorCode ISLocalToGlobalMappingDestroy_Hash(ISLocalToGlobalMapping mapping)
207: {
208: ISLocalToGlobalMapping_Hash *map = (ISLocalToGlobalMapping_Hash *)mapping->data;
210: PetscFunctionBegin;
211: if (!map) PetscFunctionReturn(PETSC_SUCCESS);
212: PetscCall(PetscHMapIDestroy(&map->globalht));
213: PetscCall(PetscFree(mapping->data));
214: PetscFunctionReturn(PETSC_SUCCESS);
215: }
217: static PetscErrorCode ISLocalToGlobalMappingResetBlockInfo_Private(ISLocalToGlobalMapping mapping)
218: {
219: PetscFunctionBegin;
220: PetscCall(PetscFree(mapping->info_procs));
221: PetscCall(PetscFree(mapping->info_numprocs));
222: if (mapping->info_indices) {
223: for (PetscInt i = 0; i < mapping->info_nproc; i++) PetscCall(PetscFree(mapping->info_indices[i]));
224: PetscCall(PetscFree(mapping->info_indices));
225: }
226: if (mapping->info_nodei) PetscCall(PetscFree(mapping->info_nodei[0]));
227: PetscCall(PetscFree2(mapping->info_nodec, mapping->info_nodei));
228: PetscCall(PetscSFDestroy(&mapping->multileaves_sf));
229: PetscFunctionReturn(PETSC_SUCCESS);
230: }
232: #define GTOLTYPE _Basic
233: #define GTOLNAME _Basic
234: #define GTOLBS mapping->bs
235: #define GTOL(g, local) \
236: do { \
237: local = map->globals[(g) / bs - start]; \
238: if ((local) >= 0) local = bs * (local) + ((g) % bs); \
239: } while (0)
241: #include <../src/vec/is/utils/isltog.h>
243: #define GTOLTYPE _Basic
244: #define GTOLNAME Block_Basic
245: #define GTOLBS 1
246: #define GTOL(g, local) \
247: do { \
248: local = map->globals[(g) - start]; \
249: } while (0)
250: #include <../src/vec/is/utils/isltog.h>
252: #define GTOLTYPE _Hash
253: #define GTOLNAME _Hash
254: #define GTOLBS mapping->bs
255: #define GTOL(g, local) \
256: do { \
257: (void)PetscHMapIGet(map->globalht, (g) / bs, &(local)); \
258: if ((local) >= 0) local = bs * (local) + ((g) % bs); \
259: } while (0)
260: #include <../src/vec/is/utils/isltog.h>
262: #define GTOLTYPE _Hash
263: #define GTOLNAME Block_Hash
264: #define GTOLBS 1
265: #define GTOL(g, local) \
266: do { \
267: (void)PetscHMapIGet(map->globalht, g, &(local)); \
268: } while (0)
269: #include <../src/vec/is/utils/isltog.h>
271: /*@
272: ISLocalToGlobalMappingDuplicate - Duplicates the local to global mapping object
274: Not Collective
276: Input Parameter:
277: . ltog - local to global mapping
279: Output Parameter:
280: . nltog - the duplicated local to global mapping
282: Level: advanced
284: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`
285: @*/
286: PetscErrorCode ISLocalToGlobalMappingDuplicate(ISLocalToGlobalMapping ltog, ISLocalToGlobalMapping *nltog)
287: {
288: ISLocalToGlobalMappingType l2gtype;
290: PetscFunctionBegin;
292: PetscCall(ISLocalToGlobalMappingCreate(PetscObjectComm((PetscObject)ltog), ltog->bs, ltog->n, ltog->indices, PETSC_COPY_VALUES, nltog));
293: PetscCall(ISLocalToGlobalMappingGetType(ltog, &l2gtype));
294: PetscCall(ISLocalToGlobalMappingSetType(*nltog, l2gtype));
295: PetscFunctionReturn(PETSC_SUCCESS);
296: }
298: /*@
299: ISLocalToGlobalMappingGetSize - Gets the local size of a local to global mapping
301: Not Collective
303: Input Parameter:
304: . mapping - local to global mapping
306: Output Parameter:
307: . n - the number of entries in the local mapping, `ISLocalToGlobalMappingGetIndices()` returns an array of this length
309: Level: advanced
311: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`
312: @*/
313: PetscErrorCode ISLocalToGlobalMappingGetSize(ISLocalToGlobalMapping mapping, PetscInt *n)
314: {
315: PetscFunctionBegin;
317: PetscAssertPointer(n, 2);
318: *n = mapping->bs * mapping->n;
319: PetscFunctionReturn(PETSC_SUCCESS);
320: }
322: /*@
323: ISLocalToGlobalMappingViewFromOptions - View an `ISLocalToGlobalMapping` based on values in the options database
325: Collective
327: Input Parameters:
328: + A - the local to global mapping object
329: . obj - optional object that provides the options prefix used for the options database query, pass `NULL` to use the options prefix of `A`
330: - name - command line option
332: Options Database Key:
333: . -name viewer_specification - See `PetscOptionsCreateViewer()` for the values of `viewer_specification`
335: Level: intermediate
337: Note:
338: 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,
339: rather `PetscOptionsCreateViewer()` should be used to construct the viewer once which can then be utilized in the heavily used routine.
341: .seealso: [](sec_scatter), `PetscViewer`, `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingView()`, `PetscObjectViewFromOptions()`, `ISLocalToGlobalMappingCreate()`, `PetscOptionsCreateViewer()`
342: @*/
343: PetscErrorCode ISLocalToGlobalMappingViewFromOptions(ISLocalToGlobalMapping A, PetscObject obj, const char name[])
344: {
345: PetscFunctionBegin;
347: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
348: PetscFunctionReturn(PETSC_SUCCESS);
349: }
351: /*@
352: ISLocalToGlobalMappingView - View a local to global mapping
354: Collective on viewer
356: Input Parameters:
357: + mapping - local to global mapping
358: - viewer - viewer
360: Level: intermediate
362: .seealso: [](sec_scatter), `PetscViewer`, `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`
363: @*/
364: PetscErrorCode ISLocalToGlobalMappingView(ISLocalToGlobalMapping mapping, PetscViewer viewer)
365: {
366: PetscBool isascii, isbinary;
367: PetscViewerFormat format;
369: PetscFunctionBegin;
371: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)mapping), &viewer));
374: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
375: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
376: PetscCall(PetscViewerGetFormat(viewer, &format));
377: if (isascii) {
378: if (format == PETSC_VIEWER_ASCII_MATLAB) {
379: const PetscInt *idxs;
380: IS is;
381: const char *name = ((PetscObject)mapping)->name;
382: char iname[PETSC_MAX_PATH_LEN];
384: PetscCall(PetscSNPrintf(iname, sizeof(iname), "%sl2g", name ? name : ""));
385: PetscCall(ISLocalToGlobalMappingGetIndices(mapping, &idxs));
386: PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)viewer), mapping->n * mapping->bs, idxs, PETSC_USE_POINTER, &is));
387: PetscCall(PetscObjectSetName((PetscObject)is, iname));
388: PetscCall(ISView(is, viewer));
389: PetscCall(ISLocalToGlobalMappingRestoreIndices(mapping, &idxs));
390: PetscCall(ISDestroy(&is));
391: } else {
392: PetscMPIInt rank;
394: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)mapping), &rank));
395: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)mapping, viewer));
396: PetscCall(PetscViewerASCIIPushSynchronized(viewer));
397: for (PetscInt i = 0; i < mapping->n; i++) {
398: PetscInt bs = mapping->bs, g = mapping->indices[i];
399: if (bs == 1) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d] %" PetscInt_FMT " %" PetscInt_FMT "\n", rank, i, g));
400: 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));
401: }
402: PetscCall(PetscViewerFlush(viewer));
403: PetscCall(PetscViewerASCIIPopSynchronized(viewer));
404: }
405: } else if (isbinary) {
406: PetscBool skipHeader;
408: PetscCall(PetscViewerSetUp(viewer));
409: PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &skipHeader));
410: if (!skipHeader) {
411: PetscMPIInt size;
412: PetscInt tr[3];
414: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)viewer), &size));
415: tr[0] = IS_LTOGM_FILE_CLASSID;
416: tr[1] = mapping->bs;
417: tr[2] = size;
418: PetscCall(PetscViewerBinaryWrite(viewer, tr, 3, PETSC_INT));
419: PetscCall(PetscViewerBinaryWriteAll(viewer, &mapping->n, 1, PETSC_DETERMINE, PETSC_DETERMINE, PETSC_INT));
420: }
421: /* write block indices */
422: PetscCall(PetscViewerBinaryWriteAll(viewer, mapping->indices, mapping->n, PETSC_DETERMINE, PETSC_DETERMINE, PETSC_INT));
423: }
424: PetscFunctionReturn(PETSC_SUCCESS);
425: }
427: /*@
428: ISLocalToGlobalMappingLoad - Loads a local-to-global mapping that has been stored in binary format.
430: Collective on viewer
432: Input Parameters:
433: + mapping - the newly loaded map, this needs to have been created with `ISLocalToGlobalMappingCreate()` or some related function before a call to `ISLocalToGlobalMappingLoad()`
434: - viewer - binary file viewer, obtained from `PetscViewerBinaryOpen()`
436: Level: intermediate
438: .seealso: [](sec_scatter), `PetscViewer`, `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingView()`, `ISLocalToGlobalMappingCreate()`
439: @*/
440: PetscErrorCode ISLocalToGlobalMappingLoad(ISLocalToGlobalMapping mapping, PetscViewer viewer)
441: {
442: PetscBool isbinary, skipHeader;
444: PetscFunctionBegin;
447: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
448: PetscCheck(isbinary, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Invalid viewer of type %s", ((PetscObject)viewer)->type_name);
450: /* reset previous data */
451: PetscCall(ISLocalToGlobalMappingResetBlockInfo_Private(mapping));
453: PetscCall(PetscViewerSetUp(viewer));
454: PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &skipHeader));
456: /* When skipping header, it assumes bs and n have been already set */
457: if (!skipHeader) {
458: MPI_Comm comm = PetscObjectComm((PetscObject)viewer);
459: PetscInt tr[3], nold = mapping->n, *sizes, nmaps = PETSC_DECIDE, st = 0;
461: PetscCall(PetscViewerBinaryRead(viewer, tr, 3, NULL, PETSC_INT));
462: PetscCheck(tr[0] == IS_LTOGM_FILE_CLASSID, PETSC_COMM_SELF, PETSC_ERR_FILE_UNEXPECTED, "Not a local-to-global map next in file");
464: mapping->bs = tr[1];
465: PetscCall(PetscMalloc1(tr[2], &sizes));
466: PetscCall(PetscViewerBinaryRead(viewer, sizes, tr[2], NULL, PETSC_INT));
468: /* consume the input, read multiple maps per process if needed */
469: PetscCall(PetscSplitOwnership(comm, &nmaps, &tr[2]));
470: PetscCallMPI(MPI_Exscan(&nmaps, &st, 1, MPIU_INT, MPI_SUM, comm));
471: mapping->n = 0;
472: for (PetscInt i = st; i < st + nmaps; i++) mapping->n += sizes[i];
473: PetscCall(PetscFree(sizes));
475: if (nold != mapping->n) {
476: if (mapping->dealloc_indices) PetscCall(PetscFree(mapping->indices));
477: mapping->indices = NULL;
478: }
479: }
481: /* read indices */
482: if (mapping->n && !mapping->indices) {
483: PetscCall(PetscMalloc1(mapping->n, &mapping->indices));
484: mapping->dealloc_indices = PETSC_TRUE;
485: }
486: PetscCall(PetscViewerBinaryReadAll(viewer, mapping->indices, mapping->n, PETSC_DETERMINE, PETSC_DETERMINE, PETSC_INT));
487: PetscFunctionReturn(PETSC_SUCCESS);
488: }
490: /*@
491: ISLocalToGlobalMappingCreateIS - Creates a mapping between a local (0 to n)
492: ordering and a global parallel ordering.
494: Not Collective
496: Input Parameter:
497: . is - index set containing the global numbers for each local number
499: Output Parameter:
500: . mapping - new mapping data structure
502: Level: advanced
504: Note:
505: the block size of the `IS` determines the block size of the mapping
507: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingSetFromOptions()`
508: @*/
509: PetscErrorCode ISLocalToGlobalMappingCreateIS(IS is, ISLocalToGlobalMapping *mapping)
510: {
511: PetscInt n, bs;
512: const PetscInt *indices;
513: MPI_Comm comm;
514: PetscBool isblock;
516: PetscFunctionBegin;
518: PetscAssertPointer(mapping, 2);
520: PetscCall(PetscObjectGetComm((PetscObject)is, &comm));
521: PetscCall(ISGetLocalSize(is, &n));
522: PetscCall(PetscObjectTypeCompare((PetscObject)is, ISBLOCK, &isblock));
523: if (!isblock) {
524: PetscCall(ISGetIndices(is, &indices));
525: PetscCall(ISLocalToGlobalMappingCreate(comm, 1, n, indices, PETSC_COPY_VALUES, mapping));
526: PetscCall(ISRestoreIndices(is, &indices));
527: } else {
528: PetscCall(ISGetBlockSize(is, &bs));
529: PetscCall(ISBlockGetIndices(is, &indices));
530: PetscCall(ISLocalToGlobalMappingCreate(comm, bs, n / bs, indices, PETSC_COPY_VALUES, mapping));
531: PetscCall(ISBlockRestoreIndices(is, &indices));
532: }
533: PetscFunctionReturn(PETSC_SUCCESS);
534: }
536: /*@
537: ISLocalToGlobalMappingCreateSF - Creates a mapping between a local (0 to n) ordering and a global parallel ordering induced by a star forest.
539: Collective
541: Input Parameters:
542: + sf - star forest mapping contiguous local indices to (rank, offset)
543: - start - first global index on this process, or `PETSC_DECIDE` to compute contiguous global numbering automatically
545: Output Parameter:
546: . mapping - new mapping data structure
548: Level: advanced
550: Note:
551: If a process calls this function with `start` = `PETSC_DECIDE` then all processes must, otherwise the program will hang.
553: .seealso: [](sec_scatter), `PetscSF`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingSetFromOptions()`
554: @*/
555: PetscErrorCode ISLocalToGlobalMappingCreateSF(PetscSF sf, PetscInt start, ISLocalToGlobalMapping *mapping)
556: {
557: PetscInt i, maxlocal, nroots, nleaves, *globals, *ltog;
558: MPI_Comm comm;
560: PetscFunctionBegin;
562: PetscAssertPointer(mapping, 3);
563: PetscCall(PetscObjectGetComm((PetscObject)sf, &comm));
564: PetscCall(PetscSFGetGraph(sf, &nroots, &nleaves, NULL, NULL));
565: if (start == PETSC_DECIDE) {
566: start = 0;
567: PetscCallMPI(MPI_Exscan(&nroots, &start, 1, MPIU_INT, MPI_SUM, comm));
568: } else PetscCheck(start >= 0, comm, PETSC_ERR_ARG_OUTOFRANGE, "start must be nonnegative or PETSC_DECIDE");
569: PetscCall(PetscSFGetLeafRange(sf, NULL, &maxlocal));
570: ++maxlocal;
571: PetscCall(PetscMalloc1(nroots, &globals));
572: PetscCall(PetscMalloc1(maxlocal, <og));
573: for (i = 0; i < nroots; i++) globals[i] = start + i;
574: for (i = 0; i < maxlocal; i++) ltog[i] = -1;
575: PetscCall(PetscSFBcastBegin(sf, MPIU_INT, globals, ltog, MPI_REPLACE));
576: PetscCall(PetscSFBcastEnd(sf, MPIU_INT, globals, ltog, MPI_REPLACE));
577: PetscCall(ISLocalToGlobalMappingCreate(comm, 1, maxlocal, ltog, PETSC_OWN_POINTER, mapping));
578: PetscCall(PetscFree(globals));
579: PetscFunctionReturn(PETSC_SUCCESS);
580: }
582: /*@
583: ISLocalToGlobalMappingSetBlockSize - Sets the blocksize of the mapping
585: Not Collective
587: Input Parameters:
588: + mapping - mapping data structure
589: - bs - the blocksize
591: Level: advanced
593: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`
594: @*/
595: PetscErrorCode ISLocalToGlobalMappingSetBlockSize(ISLocalToGlobalMapping mapping, PetscInt bs)
596: {
597: PetscInt *nid;
598: const PetscInt *oid;
599: PetscInt i, cn, on, obs, nn;
601: PetscFunctionBegin;
603: PetscCheck(bs >= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid block size %" PetscInt_FMT, bs);
604: if (bs == mapping->bs) PetscFunctionReturn(PETSC_SUCCESS);
605: on = mapping->n;
606: obs = mapping->bs;
607: oid = mapping->indices;
608: nn = (on * obs) / bs;
609: 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);
611: PetscCall(PetscMalloc1(nn, &nid));
612: PetscCall(ISLocalToGlobalMappingGetIndices(mapping, &oid));
613: for (i = 0; i < nn; i++) {
614: PetscInt j;
615: for (j = 0, cn = 0; j < bs - 1; j++) {
616: if (oid[i * bs + j] < 0) {
617: cn++;
618: continue;
619: }
620: 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]);
621: }
622: if (oid[i * bs + j] < 0) cn++;
623: if (cn) {
624: 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);
625: nid[i] = -1;
626: } else {
627: nid[i] = oid[i * bs] / bs;
628: }
629: }
630: PetscCall(ISLocalToGlobalMappingRestoreIndices(mapping, &oid));
632: mapping->n = nn;
633: mapping->bs = bs;
634: PetscCall(PetscFree(mapping->indices));
635: mapping->indices = nid;
636: mapping->globalstart = 0;
637: mapping->globalend = 0;
639: /* reset the cached information */
640: PetscCall(ISLocalToGlobalMappingResetBlockInfo_Private(mapping));
641: PetscTryTypeMethod(mapping, destroy);
642: PetscFunctionReturn(PETSC_SUCCESS);
643: }
645: /*@
646: ISLocalToGlobalMappingGetBlockSize - Gets the blocksize of the mapping
647: ordering and a global parallel ordering.
649: Not Collective
651: Input Parameter:
652: . mapping - mapping data structure
654: Output Parameter:
655: . bs - the blocksize
657: Level: advanced
659: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`
660: @*/
661: PetscErrorCode ISLocalToGlobalMappingGetBlockSize(ISLocalToGlobalMapping mapping, PetscInt *bs)
662: {
663: PetscFunctionBegin;
665: *bs = mapping->bs;
666: PetscFunctionReturn(PETSC_SUCCESS);
667: }
669: /*@
670: ISLocalToGlobalMappingCreate - Creates a mapping between a local (0 to n)
671: ordering and a global parallel ordering.
673: Not Collective, but communicator may have more than one process
675: Input Parameters:
676: + comm - MPI communicator
677: . bs - the block size
678: . n - the number of local elements divided by the block size, or equivalently the number of block indices
679: . 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
680: - mode - see PetscCopyMode
682: Output Parameter:
683: . mapping - new mapping data structure
685: Level: advanced
687: Notes:
688: There is one integer value in indices per block and it represents the actual indices bs*idx + j, where j=0,..,bs-1
690: For "small" problems when using `ISGlobalToLocalMappingApply()` and `ISGlobalToLocalMappingApplyBlock()`, the `ISLocalToGlobalMappingType`
691: of `ISLOCALTOGLOBALMAPPINGBASIC` will be used; this uses more memory but is faster; this approach is not scalable for extremely large mappings.
693: For large problems `ISLOCALTOGLOBALMAPPINGHASH` is used, this is scalable.
694: Use `ISLocalToGlobalMappingSetType()` or call `ISLocalToGlobalMappingSetFromOptions()` with the option
695: `-islocaltoglobalmapping_type` <`basic`,`hash`> to control which is used.
697: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingSetFromOptions()`,
698: `ISLOCALTOGLOBALMAPPINGBASIC`, `ISLOCALTOGLOBALMAPPINGHASH`,
699: `ISLocalToGlobalMappingSetType()`, `ISLocalToGlobalMappingType`
700: @*/
701: PetscErrorCode ISLocalToGlobalMappingCreate(MPI_Comm comm, PetscInt bs, PetscInt n, const PetscInt indices[], PetscCopyMode mode, ISLocalToGlobalMapping *mapping)
702: {
703: PetscInt *in;
705: PetscFunctionBegin;
706: if (n) PetscAssertPointer(indices, 4);
707: PetscAssertPointer(mapping, 6);
709: *mapping = NULL;
710: PetscCall(ISInitializePackage());
712: PetscCall(PetscHeaderCreate(*mapping, IS_LTOGM_CLASSID, "ISLocalToGlobalMapping", "Local to global mapping", "IS", comm, ISLocalToGlobalMappingDestroy, ISLocalToGlobalMappingView));
713: (*mapping)->n = n;
714: (*mapping)->bs = bs;
715: if (mode == PETSC_COPY_VALUES) {
716: PetscCall(PetscMalloc1(n, &in));
717: PetscCall(PetscArraycpy(in, indices, n));
718: (*mapping)->indices = in;
719: (*mapping)->dealloc_indices = PETSC_TRUE;
720: } else if (mode == PETSC_OWN_POINTER) {
721: (*mapping)->indices = (PetscInt *)indices;
722: (*mapping)->dealloc_indices = PETSC_TRUE;
723: } else if (mode == PETSC_USE_POINTER) {
724: (*mapping)->indices = (PetscInt *)indices;
725: } else SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "Invalid mode %d", mode);
726: PetscFunctionReturn(PETSC_SUCCESS);
727: }
729: PetscFunctionList ISLocalToGlobalMappingList = NULL;
731: /*@
732: ISLocalToGlobalMappingSetFromOptions - Set mapping options from the options database.
734: Not Collective
736: Input Parameter:
737: . mapping - mapping data structure
739: Options Database Key:
740: . -islocaltoglobalmapping_type (basic|hash) - nonscalable and scalable versions
742: Level: advanced
744: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`,
745: `ISLocalToGlobalMappingCreateIS()`, `ISLOCALTOGLOBALMAPPINGBASIC`,
746: `ISLOCALTOGLOBALMAPPINGHASH`, `ISLocalToGlobalMappingSetType()`, `ISLocalToGlobalMappingType`
747: @*/
748: PetscErrorCode ISLocalToGlobalMappingSetFromOptions(ISLocalToGlobalMapping mapping)
749: {
750: char type[256];
751: ISLocalToGlobalMappingType defaulttype = "Not set";
752: PetscBool flg;
754: PetscFunctionBegin;
756: PetscCall(ISLocalToGlobalMappingRegisterAll());
757: PetscObjectOptionsBegin((PetscObject)mapping);
758: PetscCall(PetscOptionsFList("-islocaltoglobalmapping_type", "ISLocalToGlobalMapping method", "ISLocalToGlobalMappingSetType", ISLocalToGlobalMappingList, ((PetscObject)mapping)->type_name ? ((PetscObject)mapping)->type_name : defaulttype, type, sizeof(type), &flg));
759: if (flg) PetscCall(ISLocalToGlobalMappingSetType(mapping, type));
760: PetscOptionsEnd();
761: PetscFunctionReturn(PETSC_SUCCESS);
762: }
764: /*@
765: ISLocalToGlobalMappingDestroy - Destroys a mapping between a local (0 to n)
766: ordering and a global parallel ordering.
768: Not Collective
770: Input Parameter:
771: . mapping - mapping data structure
773: Level: advanced
775: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingCreate()`
776: @*/
777: PetscErrorCode ISLocalToGlobalMappingDestroy(ISLocalToGlobalMapping *mapping)
778: {
779: PetscFunctionBegin;
780: if (!*mapping) PetscFunctionReturn(PETSC_SUCCESS);
782: if (--((PetscObject)*mapping)->refct > 0) {
783: *mapping = NULL;
784: PetscFunctionReturn(PETSC_SUCCESS);
785: }
786: if ((*mapping)->dealloc_indices) PetscCall(PetscFree((*mapping)->indices));
787: PetscCall(ISLocalToGlobalMappingResetBlockInfo_Private(*mapping));
788: PetscTryTypeMethod(*mapping, destroy);
789: PetscCall(PetscHeaderDestroy(mapping));
790: PetscFunctionReturn(PETSC_SUCCESS);
791: }
793: /*@
794: ISLocalToGlobalMappingApplyIS - Creates from an `IS` in the local numbering
795: a new index set using the global numbering defined in an `ISLocalToGlobalMapping`
796: context.
798: Collective
800: Input Parameters:
801: + mapping - mapping between local and global numbering
802: - is - index set in local numbering
804: Output Parameter:
805: . newis - index set in global numbering
807: Level: advanced
809: Note:
810: The output `IS` will have the same communicator as the input `IS` as well as the same block size.
812: .seealso: [](sec_scatter), `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingCreate()`,
813: `ISLocalToGlobalMappingDestroy()`, `ISGlobalToLocalMappingApply()`
814: @*/
815: PetscErrorCode ISLocalToGlobalMappingApplyIS(ISLocalToGlobalMapping mapping, IS is, IS *newis)
816: {
817: PetscInt n, *idxout, bs;
818: const PetscInt *idxin;
820: PetscFunctionBegin;
823: PetscAssertPointer(newis, 3);
825: PetscCall(ISGetLocalSize(is, &n));
826: PetscCall(ISGetBlockSize(is, &bs));
827: PetscCall(ISGetIndices(is, &idxin));
828: PetscCall(PetscMalloc1(n, &idxout));
829: PetscCall(ISLocalToGlobalMappingApply(mapping, n, idxin, idxout));
830: PetscCall(ISRestoreIndices(is, &idxin));
831: PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)is), n, idxout, PETSC_OWN_POINTER, newis));
832: PetscCall(ISSetBlockSize(*newis, bs));
833: PetscFunctionReturn(PETSC_SUCCESS);
834: }
836: /*@
837: ISLocalToGlobalMappingApply - Takes a list of integers in a local numbering
838: and converts them to the global numbering.
840: Not Collective
842: Input Parameters:
843: + mapping - the local to global mapping context
844: . N - number of integers
845: - in - input indices in local numbering
847: Output Parameter:
848: . out - indices in global numbering
850: Level: advanced
852: Note:
853: The `in` and `out` array parameters may be identical.
855: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingApplyBlock()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingDestroy()`,
856: `ISLocalToGlobalMappingApplyIS()`, `AOCreateBasic()`, `AOApplicationToPetsc()`,
857: `AOPetscToApplication()`, `ISGlobalToLocalMappingApply()`
858: @*/
859: PetscErrorCode ISLocalToGlobalMappingApply(ISLocalToGlobalMapping mapping, PetscInt N, const PetscInt in[], PetscInt out[])
860: {
861: PetscInt i, bs, Nmax;
863: PetscFunctionBegin;
865: bs = mapping->bs;
866: Nmax = bs * mapping->n;
867: if (bs == 1) {
868: const PetscInt *idx = mapping->indices;
869: for (i = 0; i < N; i++) {
870: if (in[i] < 0) {
871: out[i] = in[i];
872: continue;
873: }
874: 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);
875: out[i] = idx[in[i]];
876: }
877: } else {
878: const PetscInt *idx = mapping->indices;
879: for (i = 0; i < N; i++) {
880: if (in[i] < 0) {
881: out[i] = in[i];
882: continue;
883: }
884: 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);
885: out[i] = idx[in[i] / bs] * bs + (in[i] % bs);
886: }
887: }
888: PetscFunctionReturn(PETSC_SUCCESS);
889: }
891: /*@
892: ISLocalToGlobalMappingApplyBlock - Takes a list of integers in a local block numbering and converts them to the global block numbering
894: Not Collective
896: Input Parameters:
897: + mapping - the local to global mapping context
898: . N - number of integers
899: - in - input indices in local block numbering
901: Output Parameter:
902: . out - indices in global block numbering
904: Example:
905: 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
906: (the first block) would produce 0 and the mapping applied to 1 (the second block) would produce 3.
908: Level: advanced
910: Note:
911: The `in` and `out` array parameters may be identical.
913: .seealso: [](sec_scatter), `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingDestroy()`,
914: `ISLocalToGlobalMappingApplyIS()`, `AOCreateBasic()`, `AOApplicationToPetsc()`,
915: `AOPetscToApplication()`, `ISGlobalToLocalMappingApply()`
916: @*/
917: PetscErrorCode ISLocalToGlobalMappingApplyBlock(ISLocalToGlobalMapping mapping, PetscInt N, const PetscInt in[], PetscInt out[])
918: {
919: PetscInt i, Nmax;
920: const PetscInt *idx;
922: PetscFunctionBegin;
924: Nmax = mapping->n;
925: idx = mapping->indices;
926: for (i = 0; i < N; i++) {
927: if (in[i] < 0) {
928: out[i] = in[i];
929: continue;
930: }
931: 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);
932: out[i] = idx[in[i]];
933: }
934: PetscFunctionReturn(PETSC_SUCCESS);
935: }
937: /*@
938: ISGlobalToLocalMappingApply - Provides the local numbering for a list of integers
939: specified with a global numbering.
941: Not Collective
943: Input Parameters:
944: + mapping - mapping between local and global numbering
945: . type - `IS_GTOLM_MASK` - maps global indices with no local value to -1 in the output list (i.e., mask them)
946: `IS_GTOLM_DROP` - drops the indices with no local value from the output list
947: . n - number of global indices to map
948: - idx - global indices to map
950: Output Parameters:
951: + nout - number of indices in output array (if type == `IS_GTOLM_MASK` then nout = n)
952: - idxout - local index of each global index, one must pass in an array long enough
953: to hold all the indices. You can call `ISGlobalToLocalMappingApply()` with
954: idxout == NULL to determine the required length (returned in nout)
955: and then allocate the required space and call `ISGlobalToLocalMappingApply()`
956: a second time to set the values.
958: Level: advanced
960: Notes:
961: Either `nout` or `idxout` may be `NULL`. `idx` and `idxout` may be identical.
963: For "small" problems when using `ISGlobalToLocalMappingApply()` and `ISGlobalToLocalMappingApplyBlock()`, the `ISLocalToGlobalMappingType` of
964: `ISLOCALTOGLOBALMAPPINGBASIC` will be used;
965: 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.
966: Use `ISLocalToGlobalMappingSetType()` or call `ISLocalToGlobalMappingSetFromOptions()` with the option -islocaltoglobalmapping_type <basic,hash> to control which is used.
968: Developer Notes:
969: The manual page states that `idx` and `idxout` may be identical but the calling
970: sequence declares `idx` as const so it cannot be the same as `idxout`.
972: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingApply()`, `ISGlobalToLocalMappingApplyBlock()`, `ISLocalToGlobalMappingCreate()`,
973: `ISLocalToGlobalMappingDestroy()`
974: @*/
975: PetscErrorCode ISGlobalToLocalMappingApply(ISLocalToGlobalMapping mapping, ISGlobalToLocalMappingMode type, PetscInt n, const PetscInt idx[], PetscInt *nout, PetscInt idxout[])
976: {
977: PetscFunctionBegin;
979: if (!mapping->data) PetscCall(ISGlobalToLocalMappingSetUp(mapping));
980: PetscUseTypeMethod(mapping, globaltolocalmappingapply, type, n, idx, nout, idxout);
981: PetscFunctionReturn(PETSC_SUCCESS);
982: }
984: /*@
985: ISGlobalToLocalMappingApplyIS - Creates from an `IS` in the global numbering
986: a new index set using the local numbering defined in an `ISLocalToGlobalMapping`
987: context.
989: Not Collective
991: Input Parameters:
992: + mapping - mapping between local and global numbering
993: . type - `IS_GTOLM_MASK` - maps global indices with no local value to -1 in the output list (i.e., mask them)
994: `IS_GTOLM_DROP` - drops the indices with no local value from the output list
995: - is - index set in global numbering
997: Output Parameter:
998: . newis - index set in local numbering
1000: Level: advanced
1002: Notes:
1003: The output `IS` will be sequential, as it encodes a purely local operation
1005: If `type` is `IS_GTOLM_MASK`, `newis` will have the same block size as `is`
1007: .seealso: [](sec_scatter), `ISGlobalToLocalMapping`, `ISGlobalToLocalMappingApply()`, `ISLocalToGlobalMappingCreate()`,
1008: `ISLocalToGlobalMappingDestroy()`
1009: @*/
1010: PetscErrorCode ISGlobalToLocalMappingApplyIS(ISLocalToGlobalMapping mapping, ISGlobalToLocalMappingMode type, IS is, IS *newis)
1011: {
1012: PetscInt n, nout, *idxout, bs;
1013: const PetscInt *idxin;
1015: PetscFunctionBegin;
1018: PetscAssertPointer(newis, 4);
1020: PetscCall(ISGetLocalSize(is, &n));
1021: PetscCall(ISGetIndices(is, &idxin));
1022: if (type == IS_GTOLM_MASK) PetscCall(PetscMalloc1(n, &idxout));
1023: else {
1024: PetscCall(ISGlobalToLocalMappingApply(mapping, type, n, idxin, &nout, NULL));
1025: PetscCall(PetscMalloc1(nout, &idxout));
1026: }
1027: PetscCall(ISGlobalToLocalMappingApply(mapping, type, n, idxin, &nout, idxout));
1028: PetscCall(ISRestoreIndices(is, &idxin));
1029: PetscCall(ISCreateGeneral(PETSC_COMM_SELF, nout, idxout, PETSC_OWN_POINTER, newis));
1030: if (type == IS_GTOLM_MASK) {
1031: PetscCall(ISGetBlockSize(is, &bs));
1032: PetscCall(ISSetBlockSize(*newis, bs));
1033: }
1034: PetscFunctionReturn(PETSC_SUCCESS);
1035: }
1037: /*@
1038: ISGlobalToLocalMappingApplyBlock - Provides the local block numbering for a list of integers
1039: specified with a block global numbering.
1041: Not Collective
1043: Input Parameters:
1044: + mapping - mapping between local and global numbering
1045: . type - `IS_GTOLM_MASK` - maps global indices with no local value to -1 in the output list (i.e., mask them)
1046: `IS_GTOLM_DROP` - drops the indices with no local value from the output list
1047: . n - number of global indices to map
1048: - idx - global indices to map
1050: Output Parameters:
1051: + nout - number of indices in output array (if type == `IS_GTOLM_MASK` then nout = n)
1052: - idxout - local index of each global index, one must pass in an array long enough
1053: to hold all the indices. You can call `ISGlobalToLocalMappingApplyBlock()` with
1054: idxout == NULL to determine the required length (returned in nout)
1055: and then allocate the required space and call `ISGlobalToLocalMappingApplyBlock()`
1056: a second time to set the values.
1058: Level: advanced
1060: Notes:
1061: Either `nout` or `idxout` may be `NULL`. `idx` and `idxout` may be identical.
1063: For "small" problems when using `ISGlobalToLocalMappingApply()` and `ISGlobalToLocalMappingApplyBlock()`, the `ISLocalToGlobalMappingType` of
1064: `ISLOCALTOGLOBALMAPPINGBASIC` will be used;
1065: 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.
1066: Use `ISLocalToGlobalMappingSetType()` or call `ISLocalToGlobalMappingSetFromOptions()` with the option -islocaltoglobalmapping_type <basic,hash> to control which is used.
1068: Developer Notes:
1069: The manual page states that `idx` and `idxout` may be identical but the calling
1070: sequence declares `idx` as const so it cannot be the same as `idxout`.
1072: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingApply()`, `ISGlobalToLocalMappingApply()`, `ISLocalToGlobalMappingCreate()`,
1073: `ISLocalToGlobalMappingDestroy()`
1074: @*/
1075: PetscErrorCode ISGlobalToLocalMappingApplyBlock(ISLocalToGlobalMapping mapping, ISGlobalToLocalMappingMode type, PetscInt n, const PetscInt idx[], PetscInt *nout, PetscInt idxout[])
1076: {
1077: PetscFunctionBegin;
1079: if (!mapping->data) PetscCall(ISGlobalToLocalMappingSetUp(mapping));
1080: PetscUseTypeMethod(mapping, globaltolocalmappingapplyblock, type, n, idx, nout, idxout);
1081: PetscFunctionReturn(PETSC_SUCCESS);
1082: }
1084: /*@
1085: ISLocalToGlobalMappingGetBlockInfo - Gets the neighbor information
1087: Collective the first time it is called
1089: Input Parameter:
1090: . mapping - the mapping from local to global indexing
1092: Output Parameters:
1093: + nproc - number of processes that are connected to the calling process
1094: . procs - neighboring processes
1095: . numprocs - number of block indices for each process
1096: - indices - block indices (in local numbering) shared with neighbors (sorted by global numbering)
1098: Level: advanced
1100: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1101: `ISLocalToGlobalMappingRestoreBlockInfo()`, `ISLocalToGlobalMappingGetBlockMultiLeavesSF()`
1102: @*/
1103: PetscErrorCode ISLocalToGlobalMappingGetBlockInfo(ISLocalToGlobalMapping mapping, PetscInt *nproc, PetscInt *procs[], PetscInt *numprocs[], PetscInt **indices[])
1104: {
1105: PetscFunctionBegin;
1107: PetscCall(ISLocalToGlobalMappingSetUpBlockInfo_Private(mapping));
1108: if (nproc) *nproc = mapping->info_nproc;
1109: if (procs) *procs = mapping->info_procs;
1110: if (numprocs) *numprocs = mapping->info_numprocs;
1111: if (indices) *indices = mapping->info_indices;
1112: PetscFunctionReturn(PETSC_SUCCESS);
1113: }
1115: /*@
1116: ISLocalToGlobalMappingGetBlockNodeInfo - Gets the neighbor information for each local block index
1118: Collective the first time it is called
1120: Input Parameter:
1121: . mapping - the mapping from local to global indexing
1123: Output Parameters:
1124: + n - number of local block nodes
1125: . n_procs - an array storing the number of processes for each local block node (including self)
1126: - procs - the processes' rank for each local block node (sorted, self is first)
1128: Level: advanced
1130: Notes:
1131: The user needs to call `ISLocalToGlobalMappingRestoreBlockNodeInfo()` when the data is no longer needed.
1132: The information returned by this function complements that of `ISLocalToGlobalMappingGetBlockInfo()`.
1133: The latter only provides local information, and the neighboring information
1134: cannot be inferred in the general case, unless the mapping is locally one-to-one on each process.
1136: .seealso: `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1137: `ISLocalToGlobalMappingGetBlockInfo()`, `ISLocalToGlobalMappingRestoreBlockNodeInfo()`, `ISLocalToGlobalMappingGetNodeInfo()`
1138: @*/
1139: PetscErrorCode ISLocalToGlobalMappingGetBlockNodeInfo(ISLocalToGlobalMapping mapping, PetscInt *n, PetscInt *n_procs[], PetscInt **procs[])
1140: {
1141: PetscFunctionBegin;
1143: PetscCall(ISLocalToGlobalMappingSetUpBlockInfo_Private(mapping));
1144: if (n) *n = mapping->n;
1145: if (n_procs) *n_procs = mapping->info_nodec;
1146: if (procs) *procs = mapping->info_nodei;
1147: PetscFunctionReturn(PETSC_SUCCESS);
1148: }
1150: /*@
1151: ISLocalToGlobalMappingRestoreBlockNodeInfo - Frees the memory allocated by `ISLocalToGlobalMappingGetBlockNodeInfo()`
1153: Not Collective
1155: Input Parameters:
1156: + mapping - the mapping from local to global indexing
1157: . n - number of local block nodes
1158: . n_procs - an array storing the number of processes for each local block nodes (including self)
1159: - procs - the processes' rank for each local block node (sorted, self is first)
1161: Level: advanced
1163: .seealso: `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1164: `ISLocalToGlobalMappingGetBlockNodeInfo()`
1165: @*/
1166: PetscErrorCode ISLocalToGlobalMappingRestoreBlockNodeInfo(ISLocalToGlobalMapping mapping, PetscInt *n, PetscInt *n_procs[], PetscInt **procs[])
1167: {
1168: PetscFunctionBegin;
1170: if (n) *n = 0;
1171: if (n_procs) *n_procs = NULL;
1172: if (procs) *procs = NULL;
1173: PetscFunctionReturn(PETSC_SUCCESS);
1174: }
1176: /*@
1177: ISLocalToGlobalMappingGetBlockMultiLeavesSF - Get the star-forest to communicate multi-leaf block data
1179: Collective the first time it is called
1181: Input Parameter:
1182: . mapping - the mapping from local to global indexing
1184: Output Parameter:
1185: . mlsf - the `PetscSF`
1187: Level: advanced
1189: Notes:
1190: The returned star forest is suitable to exchange local information with other processes sharing the same global block index.
1191: For example, suppose a mapping with two processes has been created with
1192: .vb
1193: rank 0 global block indices: [0, 1, 2]
1194: rank 1 global block indices: [2, 3, 4]
1195: .ve
1196: and we want to share the local information
1197: .vb
1198: rank 0 data: [-1, -2, -3]
1199: rank 1 data: [1, 2, 3]
1200: .ve
1201: then, the broadcasting action of `mlsf` will allow to collect
1202: .vb
1203: rank 0 mlleafdata: [-1, -2, -3, 3]
1204: rank 1 mlleafdata: [-3, 3, 1, 2]
1205: .ve
1206: Use ``ISLocalToGlobalMappingGetBlockNodeInfo()`` to index into the multi-leaf data.
1208: .seealso: [](sec_scatter), `ISLocalToGlobalMappingGetBlockNodeInfo()`, `PetscSF`
1209: @*/
1210: PetscErrorCode ISLocalToGlobalMappingGetBlockMultiLeavesSF(ISLocalToGlobalMapping mapping, PetscSF *mlsf)
1211: {
1212: PetscFunctionBegin;
1214: PetscAssertPointer(mlsf, 2);
1215: PetscCall(ISLocalToGlobalMappingSetUpBlockInfo_Private(mapping));
1216: *mlsf = mapping->multileaves_sf;
1217: PetscFunctionReturn(PETSC_SUCCESS);
1218: }
1220: static PetscErrorCode ISLocalToGlobalMappingSetUpBlockInfo_Private(ISLocalToGlobalMapping mapping)
1221: {
1222: PetscSF sf, sf2, imsf, msf;
1223: MPI_Comm comm;
1224: const PetscSFNode *sfnode;
1225: PetscSFNode *newsfnode;
1226: PetscLayout layout;
1227: PetscHMapI neighs;
1228: PetscHashIter iter;
1229: PetscBool missing;
1230: const PetscInt *gidxs, *rootdegree;
1231: PetscInt *mask, *mrootdata, *leafdata, *newleafdata, *leafrd, *tmpg;
1232: PetscInt nroots, nleaves, newnleaves, bs, i, j, m, mnroots, p;
1233: PetscMPIInt rank, size;
1235: PetscFunctionBegin;
1236: if (mapping->multileaves_sf) PetscFunctionReturn(PETSC_SUCCESS);
1237: PetscCall(PetscObjectGetComm((PetscObject)mapping, &comm));
1238: PetscCallMPI(MPI_Comm_size(comm, &size));
1239: PetscCallMPI(MPI_Comm_rank(comm, &rank));
1241: /* Get mapping indices */
1242: PetscCall(ISLocalToGlobalMappingGetBlockSize(mapping, &bs));
1243: PetscCall(ISLocalToGlobalMappingGetBlockIndices(mapping, &gidxs));
1244: PetscCall(ISLocalToGlobalMappingGetSize(mapping, &nleaves));
1245: nleaves /= bs;
1247: /* Create layout for global indices */
1248: for (i = 0, m = 0; i < nleaves; i++) m = PetscMax(m, gidxs[i]);
1249: PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &m, 1, MPIU_INT, MPI_MAX, comm));
1250: PetscCall(PetscLayoutCreate(comm, &layout));
1251: PetscCall(PetscLayoutSetSize(layout, m + 1));
1252: PetscCall(PetscLayoutSetUp(layout));
1254: /* Create SF to share global indices */
1255: PetscCall(PetscSFCreate(comm, &sf));
1256: PetscCall(PetscSFSetGraphLayout(sf, layout, nleaves, NULL, PETSC_OWN_POINTER, gidxs));
1257: PetscCall(PetscSFSetUp(sf));
1258: PetscCall(PetscLayoutDestroy(&layout));
1260: /* communicate root degree to leaves */
1261: PetscCall(PetscSFGetGraph(sf, &nroots, NULL, NULL, &sfnode));
1262: PetscCall(PetscSFComputeDegreeBegin(sf, &rootdegree));
1263: PetscCall(PetscSFComputeDegreeEnd(sf, &rootdegree));
1264: for (i = 0, mnroots = 0; i < nroots; i++) mnroots += rootdegree[i];
1265: PetscCall(PetscMalloc3(2 * PetscMax(mnroots, nroots), &mrootdata, 2 * nleaves, &leafdata, nleaves, &leafrd));
1266: for (i = 0, m = 0; i < nroots; i++) {
1267: mrootdata[2 * i + 0] = rootdegree[i];
1268: mrootdata[2 * i + 1] = m;
1269: m += rootdegree[i];
1270: }
1271: PetscCall(PetscSFBcastBegin(sf, MPIU_2INT, mrootdata, leafdata, MPI_REPLACE));
1272: PetscCall(PetscSFBcastEnd(sf, MPIU_2INT, mrootdata, leafdata, MPI_REPLACE));
1274: /* allocate enough space to store ranks */
1275: for (i = 0, newnleaves = 0; i < nleaves; i++) {
1276: newnleaves += leafdata[2 * i];
1277: leafrd[i] = leafdata[2 * i];
1278: }
1280: /* create new SF nodes to collect multi-root data at leaves */
1281: PetscCall(PetscMalloc1(newnleaves, &newsfnode));
1282: for (i = 0, m = 0; i < nleaves; i++) {
1283: for (j = 0; j < leafrd[i]; j++) {
1284: newsfnode[m].rank = sfnode[i].rank;
1285: newsfnode[m].index = leafdata[2 * i + 1] + j;
1286: m++;
1287: }
1288: }
1290: /* gather ranks at multi roots */
1291: for (i = 0; i < mnroots; i++) mrootdata[i] = -1;
1292: for (i = 0; i < nleaves; i++) leafdata[i] = rank;
1294: PetscCall(PetscSFGatherBegin(sf, MPIU_INT, leafdata, mrootdata));
1295: PetscCall(PetscSFGatherEnd(sf, MPIU_INT, leafdata, mrootdata));
1297: /* from multi-roots to multi-leaves */
1298: PetscCall(PetscSFCreate(comm, &sf2));
1299: PetscCall(PetscSFSetGraph(sf2, mnroots, newnleaves, NULL, PETSC_OWN_POINTER, newsfnode, PETSC_OWN_POINTER));
1300: PetscCall(PetscSFSetUp(sf2));
1302: /* broadcast multi-root data to multi-leaves */
1303: PetscCall(PetscMalloc1(newnleaves, &newleafdata));
1304: PetscCall(PetscSFBcastBegin(sf2, MPIU_INT, mrootdata, newleafdata, MPI_REPLACE));
1305: PetscCall(PetscSFBcastEnd(sf2, MPIU_INT, mrootdata, newleafdata, MPI_REPLACE));
1307: /* sort sharing ranks */
1308: for (i = 0, m = 0; i < nleaves; i++) {
1309: PetscCall(PetscSortInt(leafrd[i], newleafdata + m));
1310: m += leafrd[i];
1311: }
1313: /* Number of neighbors and their ranks */
1314: PetscCall(PetscHMapICreate(&neighs));
1315: for (i = 0; i < newnleaves; i++) PetscCall(PetscHMapIPut(neighs, newleafdata[i], &iter, &missing));
1316: PetscCall(PetscHMapIGetSize(neighs, &mapping->info_nproc));
1317: PetscCall(PetscMalloc1(mapping->info_nproc, &mapping->info_procs));
1318: PetscCall(PetscHMapIGetKeys(neighs, (i = 0, &i), mapping->info_procs));
1319: for (i = 0; i < mapping->info_nproc; i++) { /* put info for self first */
1320: if (mapping->info_procs[i] == rank) {
1321: PetscInt newr = mapping->info_procs[0];
1323: mapping->info_procs[0] = rank;
1324: mapping->info_procs[i] = newr;
1325: break;
1326: }
1327: }
1328: if (mapping->info_nproc) PetscCall(PetscSortInt(mapping->info_nproc - 1, mapping->info_procs + 1));
1329: PetscCall(PetscHMapIDestroy(&neighs));
1331: /* collect info data */
1332: PetscCall(PetscMalloc1(mapping->info_nproc, &mapping->info_numprocs));
1333: PetscCall(PetscMalloc1(mapping->info_nproc, &mapping->info_indices));
1334: for (i = 0; i < mapping->info_nproc; i++) mapping->info_indices[i] = NULL;
1336: PetscCall(PetscMalloc1(nleaves, &mask));
1337: PetscCall(PetscMalloc1(nleaves, &tmpg));
1338: for (p = 0; p < mapping->info_nproc; p++) {
1339: PetscInt *tmp, trank = mapping->info_procs[p];
1341: PetscCall(PetscMemzero(mask, nleaves * sizeof(*mask)));
1342: for (i = 0, m = 0; i < nleaves; i++) {
1343: for (j = 0; j < leafrd[i]; j++) {
1344: if (newleafdata[m] == trank) mask[i]++;
1345: if (!p && newleafdata[m] != rank) mask[i]++;
1346: m++;
1347: }
1348: }
1349: for (i = 0, m = 0; i < nleaves; i++)
1350: if (mask[i] > (!p ? 1 : 0)) m++;
1352: PetscCall(PetscMalloc1(m, &tmp));
1353: for (i = 0, m = 0; i < nleaves; i++)
1354: if (mask[i] > (!p ? 1 : 0)) {
1355: tmp[m] = i;
1356: tmpg[m] = gidxs[i];
1357: m++;
1358: }
1359: PetscCall(PetscSortIntWithArray(m, tmpg, tmp));
1360: mapping->info_indices[p] = tmp;
1361: mapping->info_numprocs[p] = m;
1362: }
1364: /* Node info */
1365: PetscCall(PetscMalloc2(nleaves, &mapping->info_nodec, nleaves + 1, &mapping->info_nodei));
1366: PetscCall(PetscArraycpy(mapping->info_nodec, leafrd, nleaves));
1367: PetscCall(PetscMalloc1(newnleaves, &mapping->info_nodei[0]));
1368: for (i = 0; i < nleaves - 1; i++) mapping->info_nodei[i + 1] = mapping->info_nodei[i] + mapping->info_nodec[i];
1369: PetscCall(PetscArraycpy(mapping->info_nodei[0], newleafdata, newnleaves));
1371: /* Create SF from leaves to multi-leaves */
1372: PetscCall(PetscSFGetMultiSF(sf, &msf));
1373: PetscCall(PetscSFCreateInverseSF(msf, &imsf));
1374: PetscCall(PetscSFCompose(imsf, sf2, &mapping->multileaves_sf));
1375: PetscCall(PetscSFDestroy(&imsf));
1376: PetscCall(PetscSFDestroy(&sf));
1377: PetscCall(PetscSFDestroy(&sf2));
1379: PetscCall(ISLocalToGlobalMappingRestoreBlockIndices(mapping, &gidxs));
1380: PetscCall(PetscFree(tmpg));
1381: PetscCall(PetscFree(mask));
1382: PetscCall(PetscFree3(mrootdata, leafdata, leafrd));
1383: PetscCall(PetscFree(newleafdata));
1384: PetscFunctionReturn(PETSC_SUCCESS);
1385: }
1387: /*@
1388: ISLocalToGlobalMappingRestoreBlockInfo - Frees the memory allocated by `ISLocalToGlobalMappingGetBlockInfo()`
1390: Not Collective
1392: Input Parameters:
1393: + mapping - the mapping from local to global indexing
1394: . nproc - number of processes that are connected to the calling process
1395: . procs - neighboring processes
1396: . numprocs - number of block indices for each process
1397: - indices - block indices (in local numbering) shared with neighbors (sorted by global numbering)
1399: Level: advanced
1401: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1402: `ISLocalToGlobalMappingGetInfo()`
1403: @*/
1404: PetscErrorCode ISLocalToGlobalMappingRestoreBlockInfo(ISLocalToGlobalMapping mapping, PetscInt *nproc, PetscInt *procs[], PetscInt *numprocs[], PetscInt **indices[])
1405: {
1406: PetscFunctionBegin;
1408: if (nproc) *nproc = 0;
1409: if (procs) *procs = NULL;
1410: if (numprocs) *numprocs = NULL;
1411: if (indices) *indices = NULL;
1412: PetscFunctionReturn(PETSC_SUCCESS);
1413: }
1415: /*@
1416: ISLocalToGlobalMappingGetInfo - Gets the neighbor information for each process
1418: Collective the first time it is called
1420: Input Parameter:
1421: . mapping - the mapping from local to global indexing
1423: Output Parameters:
1424: + nproc - number of processes that are connected to the calling process
1425: . procs - neighboring processes
1426: . numprocs - number of indices for each process
1427: - indices - indices (in local numbering) shared with neighbors (sorted by global numbering)
1429: Level: advanced
1431: Note:
1432: The user needs to call `ISLocalToGlobalMappingRestoreInfo()` when the data is no longer needed.
1434: Fortran Notes:
1435: There is no `ISLocalToGlobalMappingRestoreInfo()` in Fortran. You must make sure that
1436: `procs`[], `numprocs`[] and `indices`[][] are large enough arrays, either by allocating them
1437: dynamically or defining static ones large enough.
1439: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1440: `ISLocalToGlobalMappingRestoreInfo()`, `ISLocalToGlobalMappingGetNodeInfo()`
1441: @*/
1442: PetscErrorCode ISLocalToGlobalMappingGetInfo(ISLocalToGlobalMapping mapping, PetscInt *nproc, PetscInt *procs[], PetscInt *numprocs[], PetscInt **indices[])
1443: {
1444: PetscInt **bindices = NULL, *bnumprocs = NULL, bs, i, j, k, n, *bprocs;
1446: PetscFunctionBegin;
1448: bs = mapping->bs;
1449: PetscCall(ISLocalToGlobalMappingGetBlockInfo(mapping, &n, &bprocs, &bnumprocs, &bindices));
1450: if (bs > 1) { /* we need to expand the cached info */
1451: if (indices) PetscCall(PetscCalloc1(n, indices));
1452: if (numprocs) PetscCall(PetscCalloc1(n, numprocs));
1453: if (indices || numprocs) {
1454: for (i = 0; i < n; i++) {
1455: if (indices) {
1456: PetscCall(PetscMalloc1(bs * bnumprocs[i], &(*indices)[i]));
1457: for (j = 0; j < bnumprocs[i]; j++) {
1458: for (k = 0; k < bs; k++) (*indices)[i][j * bs + k] = bs * bindices[i][j] + k;
1459: }
1460: }
1461: if (numprocs) (*numprocs)[i] = bnumprocs[i] * bs;
1462: }
1463: }
1464: } else {
1465: if (numprocs) *numprocs = bnumprocs;
1466: if (indices) *indices = bindices;
1467: }
1468: if (nproc) *nproc = n;
1469: if (procs) *procs = bprocs;
1470: PetscFunctionReturn(PETSC_SUCCESS);
1471: }
1473: /*@
1474: ISLocalToGlobalMappingRestoreInfo - Frees the memory allocated by `ISLocalToGlobalMappingGetInfo()`
1476: Not Collective
1478: Input Parameters:
1479: + mapping - the mapping from local to global indexing
1480: . nproc - number of processes that are connected to the calling process
1481: . procs - neighboring processes
1482: . numprocs - number of indices for each process
1483: - indices - indices (in local numbering) shared with neighbors (sorted by global numbering)
1485: Level: advanced
1487: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1488: `ISLocalToGlobalMappingGetInfo()`
1489: @*/
1490: PetscErrorCode ISLocalToGlobalMappingRestoreInfo(ISLocalToGlobalMapping mapping, PetscInt *nproc, PetscInt *procs[], PetscInt *numprocs[], PetscInt **indices[])
1491: {
1492: PetscFunctionBegin;
1494: if (mapping->bs > 1) {
1495: if (numprocs) PetscCall(PetscFree(*numprocs));
1496: if (indices) {
1497: if (*indices)
1498: for (PetscInt i = 0; i < *nproc; i++) PetscCall(PetscFree((*indices)[i]));
1499: PetscCall(PetscFree(*indices));
1500: }
1501: }
1502: PetscFunctionReturn(PETSC_SUCCESS);
1503: }
1505: /*@
1506: ISLocalToGlobalMappingGetNodeInfo - Gets the neighbor information of local nodes
1508: Collective the first time it is called
1510: Input Parameter:
1511: . mapping - the mapping from local to global indexing
1513: Output Parameters:
1514: + n - number of local nodes
1515: . n_procs - an array storing the number of processes for each local node (including self)
1516: - procs - the processes' rank for each local node (sorted, self is first)
1518: Level: advanced
1520: Note:
1521: The user needs to call `ISLocalToGlobalMappingRestoreNodeInfo()` when the data is no longer needed.
1523: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1524: `ISLocalToGlobalMappingGetInfo()`, `ISLocalToGlobalMappingRestoreNodeInfo()`, `ISLocalToGlobalMappingGetBlockNodeInfo()`
1525: @*/
1526: PetscErrorCode ISLocalToGlobalMappingGetNodeInfo(ISLocalToGlobalMapping mapping, PetscInt *n, PetscInt *n_procs[], PetscInt **procs[])
1527: {
1528: PetscInt **bprocs = NULL, *bn_procs = NULL, bs, i, j, k, bn;
1530: PetscFunctionBegin;
1532: bs = mapping->bs;
1533: PetscCall(ISLocalToGlobalMappingGetBlockNodeInfo(mapping, &bn, &bn_procs, &bprocs));
1534: if (bs > 1) { /* we need to expand the cached info */
1535: PetscInt *tn_procs;
1536: PetscInt c;
1538: PetscCall(PetscMalloc1(bn * bs, &tn_procs));
1539: for (i = 0, c = 0; i < bn; i++) {
1540: for (k = 0; k < bs; k++) tn_procs[i * bs + k] = bn_procs[i];
1541: c += bs * bn_procs[i];
1542: }
1543: if (n) *n = bn * bs;
1544: if (procs) {
1545: PetscInt **tprocs;
1546: PetscInt tn = bn * bs;
1548: PetscCall(PetscMalloc1(tn, &tprocs));
1549: if (tn) PetscCall(PetscMalloc1(c, &tprocs[0]));
1550: for (i = 0; i < tn - 1; i++) tprocs[i + 1] = tprocs[i] + tn_procs[i];
1551: for (i = 0; i < bn; i++) {
1552: for (k = 0; k < bs; k++) {
1553: for (j = 0; j < bn_procs[i]; j++) tprocs[i * bs + k][j] = bprocs[i][j];
1554: }
1555: }
1556: *procs = tprocs;
1557: }
1558: if (n_procs) *n_procs = tn_procs;
1559: else PetscCall(PetscFree(tn_procs));
1560: } else {
1561: if (n) *n = bn;
1562: if (n_procs) *n_procs = bn_procs;
1563: if (procs) *procs = bprocs;
1564: }
1565: PetscFunctionReturn(PETSC_SUCCESS);
1566: }
1568: /*@
1569: ISLocalToGlobalMappingRestoreNodeInfo - Frees the memory allocated by `ISLocalToGlobalMappingGetNodeInfo()`
1571: Not Collective
1573: Input Parameters:
1574: + mapping - the mapping from local to global indexing
1575: . n - number of local nodes
1576: . n_procs - an array storing the number of processes for each local node (including self)
1577: - procs - the processes' rank for each local node (sorted, self is first)
1579: Level: advanced
1581: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1582: `ISLocalToGlobalMappingGetInfo()`
1583: @*/
1584: PetscErrorCode ISLocalToGlobalMappingRestoreNodeInfo(ISLocalToGlobalMapping mapping, PetscInt *n, PetscInt *n_procs[], PetscInt **procs[])
1585: {
1586: PetscFunctionBegin;
1588: if (mapping->bs > 1) {
1589: if (n_procs) PetscCall(PetscFree(*n_procs));
1590: if (procs) {
1591: if (*procs) PetscCall(PetscFree((*procs)[0]));
1592: PetscCall(PetscFree(*procs));
1593: }
1594: }
1595: PetscCall(ISLocalToGlobalMappingRestoreBlockNodeInfo(mapping, n, n_procs, procs));
1596: PetscFunctionReturn(PETSC_SUCCESS);
1597: }
1599: /*@
1600: ISLocalToGlobalMappingGetIndices - Get global indices for every local point that is mapped
1602: Not Collective
1604: Input Parameter:
1605: . ltog - local to global mapping
1607: Output Parameter:
1608: . array - array of indices, the length of this array may be obtained with `ISLocalToGlobalMappingGetSize()`
1610: Level: advanced
1612: Note:
1613: `ISLocalToGlobalMappingGetSize()` returns the length the this array
1615: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingRestoreIndices()`,
1616: `ISLocalToGlobalMappingGetBlockIndices()`, `ISLocalToGlobalMappingRestoreBlockIndices()`
1617: @*/
1618: PetscErrorCode ISLocalToGlobalMappingGetIndices(ISLocalToGlobalMapping ltog, const PetscInt *array[])
1619: {
1620: PetscFunctionBegin;
1622: PetscAssertPointer(array, 2);
1623: if (ltog->bs == 1) {
1624: *array = ltog->indices;
1625: } else {
1626: PetscInt *jj, k, i, j, n = ltog->n, bs = ltog->bs;
1627: const PetscInt *ii;
1629: PetscCall(PetscMalloc1(bs * n, &jj));
1630: *array = jj;
1631: k = 0;
1632: ii = ltog->indices;
1633: for (i = 0; i < n; i++)
1634: for (j = 0; j < bs; j++) jj[k++] = bs * ii[i] + j;
1635: }
1636: PetscFunctionReturn(PETSC_SUCCESS);
1637: }
1639: /*@
1640: ISLocalToGlobalMappingRestoreIndices - Restore indices obtained with `ISLocalToGlobalMappingGetIndices()`
1642: Not Collective
1644: Input Parameters:
1645: + ltog - local to global mapping
1646: - array - array of indices
1648: Level: advanced
1650: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingGetIndices()`
1651: @*/
1652: PetscErrorCode ISLocalToGlobalMappingRestoreIndices(ISLocalToGlobalMapping ltog, const PetscInt *array[])
1653: {
1654: PetscFunctionBegin;
1656: PetscAssertPointer(array, 2);
1657: PetscCheck(ltog->bs != 1 || *array == ltog->indices, PETSC_COMM_SELF, PETSC_ERR_ARG_BADPTR, "Trying to return mismatched pointer");
1658: if (ltog->bs > 1) PetscCall(PetscFree(*(void **)array));
1659: PetscFunctionReturn(PETSC_SUCCESS);
1660: }
1662: /*@
1663: ISLocalToGlobalMappingGetBlockIndices - Get global indices for every local block in a `ISLocalToGlobalMapping`
1665: Not Collective
1667: Input Parameter:
1668: . ltog - local to global mapping
1670: Output Parameter:
1671: . array - array of indices
1673: Level: advanced
1675: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingApply()`,
1676: `ISLocalToGlobalMappingRestoreBlockIndices()`
1677: @*/
1678: PetscErrorCode ISLocalToGlobalMappingGetBlockIndices(ISLocalToGlobalMapping ltog, const PetscInt *array[])
1679: {
1680: PetscFunctionBegin;
1682: PetscAssertPointer(array, 2);
1683: *array = ltog->indices;
1684: PetscFunctionReturn(PETSC_SUCCESS);
1685: }
1687: /*@
1688: ISLocalToGlobalMappingRestoreBlockIndices - Restore indices obtained with `ISLocalToGlobalMappingGetBlockIndices()`
1690: Not Collective
1692: Input Parameters:
1693: + ltog - local to global mapping
1694: - array - array of indices
1696: Level: advanced
1698: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingGetIndices()`
1699: @*/
1700: PetscErrorCode ISLocalToGlobalMappingRestoreBlockIndices(ISLocalToGlobalMapping ltog, const PetscInt *array[])
1701: {
1702: PetscFunctionBegin;
1704: PetscAssertPointer(array, 2);
1705: PetscCheck(*array == ltog->indices, PETSC_COMM_SELF, PETSC_ERR_ARG_BADPTR, "Trying to return mismatched pointer");
1706: *array = NULL;
1707: PetscFunctionReturn(PETSC_SUCCESS);
1708: }
1710: /*@
1711: ISLocalToGlobalMappingConcatenate - Create a new mapping that concatenates a list of mappings
1713: Not Collective
1715: Input Parameters:
1716: + comm - communicator for the new mapping, must contain the communicator of every mapping to concatenate
1717: . n - number of mappings to concatenate
1718: - ltogs - local to global mappings
1720: Output Parameter:
1721: . ltogcat - new mapping
1723: Level: advanced
1725: Note:
1726: This currently always returns a mapping with block size of 1
1728: Developer Notes:
1729: If all the input mapping have the same block size we could easily handle that as a special case
1731: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingCreate()`
1732: @*/
1733: PetscErrorCode ISLocalToGlobalMappingConcatenate(MPI_Comm comm, PetscInt n, const ISLocalToGlobalMapping ltogs[], ISLocalToGlobalMapping *ltogcat)
1734: {
1735: PetscInt i, cnt, m, *idx;
1737: PetscFunctionBegin;
1738: PetscCheck(n >= 0, comm, PETSC_ERR_ARG_OUTOFRANGE, "Must have a non-negative number of mappings, given %" PetscInt_FMT, n);
1739: if (n > 0) PetscAssertPointer(ltogs, 3);
1741: PetscAssertPointer(ltogcat, 4);
1742: for (cnt = 0, i = 0; i < n; i++) {
1743: PetscCall(ISLocalToGlobalMappingGetSize(ltogs[i], &m));
1744: cnt += m;
1745: }
1746: PetscCall(PetscMalloc1(cnt, &idx));
1747: for (cnt = 0, i = 0; i < n; i++) {
1748: const PetscInt *subidx;
1749: PetscCall(ISLocalToGlobalMappingGetSize(ltogs[i], &m));
1750: PetscCall(ISLocalToGlobalMappingGetIndices(ltogs[i], &subidx));
1751: PetscCall(PetscArraycpy(&idx[cnt], subidx, m));
1752: PetscCall(ISLocalToGlobalMappingRestoreIndices(ltogs[i], &subidx));
1753: cnt += m;
1754: }
1755: PetscCall(ISLocalToGlobalMappingCreate(comm, 1, cnt, idx, PETSC_OWN_POINTER, ltogcat));
1756: PetscFunctionReturn(PETSC_SUCCESS);
1757: }
1759: /*MC
1760: ISLOCALTOGLOBALMAPPINGBASIC - basic implementation of the `ISLocalToGlobalMapping` object. When `ISGlobalToLocalMappingApply()` is
1761: used this is good for only small and moderate size problems.
1763: Options Database Key:
1764: . -islocaltoglobalmapping_type basic - select this method
1766: Level: beginner
1768: Developer Note:
1769: This stores all the mapping information on each MPI rank.
1771: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingSetType()`, `ISLOCALTOGLOBALMAPPINGHASH`
1772: M*/
1773: PETSC_EXTERN PetscErrorCode ISLocalToGlobalMappingCreate_Basic(ISLocalToGlobalMapping ltog)
1774: {
1775: PetscFunctionBegin;
1776: ltog->ops->globaltolocalmappingapply = ISGlobalToLocalMappingApply_Basic;
1777: ltog->ops->globaltolocalmappingsetup = ISGlobalToLocalMappingSetUp_Basic;
1778: ltog->ops->globaltolocalmappingapplyblock = ISGlobalToLocalMappingApplyBlock_Basic;
1779: ltog->ops->destroy = ISLocalToGlobalMappingDestroy_Basic;
1780: PetscFunctionReturn(PETSC_SUCCESS);
1781: }
1783: /*MC
1784: ISLOCALTOGLOBALMAPPINGHASH - hash implementation of the `ISLocalToGlobalMapping` object. When `ISGlobalToLocalMappingApply()` is
1785: used this is good for large memory problems.
1787: Options Database Key:
1788: . -islocaltoglobalmapping_type hash - select this method
1790: Level: beginner
1792: Note:
1793: This is selected automatically for large problems if the user does not set the type.
1795: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingSetType()`, `ISLOCALTOGLOBALMAPPINGBASIC`
1796: M*/
1797: PETSC_EXTERN PetscErrorCode ISLocalToGlobalMappingCreate_Hash(ISLocalToGlobalMapping ltog)
1798: {
1799: PetscFunctionBegin;
1800: ltog->ops->globaltolocalmappingapply = ISGlobalToLocalMappingApply_Hash;
1801: ltog->ops->globaltolocalmappingsetup = ISGlobalToLocalMappingSetUp_Hash;
1802: ltog->ops->globaltolocalmappingapplyblock = ISGlobalToLocalMappingApplyBlock_Hash;
1803: ltog->ops->destroy = ISLocalToGlobalMappingDestroy_Hash;
1804: PetscFunctionReturn(PETSC_SUCCESS);
1805: }
1807: /*@
1808: ISLocalToGlobalMappingRegister - Registers a method for applying a global to local mapping with an `ISLocalToGlobalMapping`
1810: Not Collective, No Fortran Support
1812: Input Parameters:
1813: + sname - name of a new method
1814: - function - routine to create method context
1816: Example Usage:
1817: .vb
1818: ISLocalToGlobalMappingRegister("my_mapper", MyCreate);
1819: .ve
1821: Then, your mapping can be chosen with the procedural interface via
1822: .vb
1823: ISLocalToGlobalMappingSetType(ltog, "my_mapper")
1824: .ve
1825: or at runtime via the option
1826: .vb
1827: -islocaltoglobalmapping_type my_mapper
1828: .ve
1830: Level: advanced
1832: Note:
1833: `ISLocalToGlobalMappingRegister()` may be called multiple times to add several user-defined mappings.
1835: .seealso: [](sec_scatter), `ISLocalToGlobalMappingRegisterAll()`, `ISLocalToGlobalMappingRegisterDestroy()`, `ISLOCALTOGLOBALMAPPINGBASIC`,
1836: `ISLOCALTOGLOBALMAPPINGHASH`, `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingApply()`
1837: @*/
1838: PetscErrorCode ISLocalToGlobalMappingRegister(const char sname[], PetscErrorCode (*function)(ISLocalToGlobalMapping))
1839: {
1840: PetscFunctionBegin;
1841: PetscCall(ISInitializePackage());
1842: PetscCall(PetscFunctionListAdd(&ISLocalToGlobalMappingList, sname, function));
1843: PetscFunctionReturn(PETSC_SUCCESS);
1844: }
1846: /*@
1847: ISLocalToGlobalMappingSetType - Sets the implementation type `ISLocalToGlobalMapping` will use
1849: Logically Collective
1851: Input Parameters:
1852: + ltog - the `ISLocalToGlobalMapping` object
1853: - type - a known method
1855: Options Database Key:
1856: . -islocaltoglobalmapping_type (basic|hash) - Sets the method used for applying the mapping, see `ISLocalToGlobalMappingType`
1858: Level: intermediate
1860: Notes:
1861: See `ISLocalToGlobalMappingType` for available methods
1863: Normally, it is best to use the `ISLocalToGlobalMappingSetFromOptions()` command and
1864: then set the `ISLocalToGlobalMappingType` from the options database rather than by using
1865: this routine.
1867: Developer Notes:
1868: `ISLocalToGlobalMappingRegister()` is used to add new types to `ISLocalToGlobalMappingList` from which they
1869: are accessed by `ISLocalToGlobalMappingSetType()`.
1871: .seealso: [](sec_scatter), `ISLocalToGlobalMappingType`, `ISLocalToGlobalMappingRegister()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingGetType()`
1872: @*/
1873: PetscErrorCode ISLocalToGlobalMappingSetType(ISLocalToGlobalMapping ltog, ISLocalToGlobalMappingType type)
1874: {
1875: PetscBool match;
1876: PetscErrorCode (*r)(ISLocalToGlobalMapping) = NULL;
1878: PetscFunctionBegin;
1880: if (type) PetscAssertPointer(type, 2);
1882: PetscCall(PetscObjectTypeCompare((PetscObject)ltog, type, &match));
1883: if (match) PetscFunctionReturn(PETSC_SUCCESS);
1885: /* L2G maps defer type setup at globaltolocal calls, allow passing NULL here */
1886: if (type) {
1887: PetscCall(PetscFunctionListFind(ISLocalToGlobalMappingList, type, &r));
1888: PetscCheck(r, PetscObjectComm((PetscObject)ltog), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unable to find requested ISLocalToGlobalMapping type %s", type);
1889: }
1890: /* Destroy the previous private LTOG context */
1891: PetscTryTypeMethod(ltog, destroy);
1892: ltog->ops->destroy = NULL;
1894: PetscCall(PetscObjectChangeTypeName((PetscObject)ltog, type));
1895: if (r) PetscCall((*r)(ltog));
1896: PetscFunctionReturn(PETSC_SUCCESS);
1897: }
1899: /*@
1900: ISLocalToGlobalMappingGetType - Get the type of the `ISLocalToGlobalMapping`
1902: Not Collective
1904: Input Parameter:
1905: . ltog - the `ISLocalToGlobalMapping` object
1907: Output Parameter:
1908: . type - the type
1910: Level: intermediate
1912: Note:
1913: `type` should not be retained for later use as it will be an invalid pointer if the `ISLocalToGlobalMappingType` of `ltog` is changed.
1915: .seealso: [](sec_scatter), `ISLocalToGlobalMappingType`, `ISLocalToGlobalMappingRegister()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingSetType()`,
1916: `PetscObjectTypeCompare()`, `PetscObjectTypeCompareAny()`
1917: @*/
1918: PetscErrorCode ISLocalToGlobalMappingGetType(ISLocalToGlobalMapping ltog, ISLocalToGlobalMappingType *type)
1919: {
1920: PetscFunctionBegin;
1922: PetscAssertPointer(type, 2);
1923: *type = ((PetscObject)ltog)->type_name;
1924: PetscFunctionReturn(PETSC_SUCCESS);
1925: }
1927: PetscBool ISLocalToGlobalMappingRegisterAllCalled = PETSC_FALSE;
1929: /*@
1930: ISLocalToGlobalMappingRegisterAll - Registers all of the local to global mapping components in the `IS` package.
1932: Not Collective
1934: Level: advanced
1936: .seealso: [](sec_scatter), `ISRegister()`, `ISLocalToGlobalRegister()`
1937: @*/
1938: PetscErrorCode ISLocalToGlobalMappingRegisterAll(void)
1939: {
1940: PetscFunctionBegin;
1941: if (ISLocalToGlobalMappingRegisterAllCalled) PetscFunctionReturn(PETSC_SUCCESS);
1942: ISLocalToGlobalMappingRegisterAllCalled = PETSC_TRUE;
1943: PetscCall(ISLocalToGlobalMappingRegister(ISLOCALTOGLOBALMAPPINGBASIC, ISLocalToGlobalMappingCreate_Basic));
1944: PetscCall(ISLocalToGlobalMappingRegister(ISLOCALTOGLOBALMAPPINGHASH, ISLocalToGlobalMappingCreate_Hash));
1945: PetscFunctionReturn(PETSC_SUCCESS);
1946: }