Actual source code: isltog.c
2: #include <petsc/private/isimpl.h>
3: #include <petsc/private/hashmapi.h>
4: #include <petscsf.h>
5: #include <petscviewer.h>
7: PetscClassId IS_LTOGM_CLASSID;
8: static PetscErrorCode ISLocalToGlobalMappingGetBlockInfo_Private(ISLocalToGlobalMapping, PetscInt *, PetscInt **, PetscInt **, PetscInt ***);
10: typedef struct {
11: PetscInt *globals;
12: } ISLocalToGlobalMapping_Basic;
14: typedef struct {
15: PetscHMapI globalht;
16: } ISLocalToGlobalMapping_Hash;
18: /*@C
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: /*@C
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: /*@C
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: /* -----------------------------------------------------------------------------------------*/
127: /*
128: Creates the global mapping information in the ISLocalToGlobalMapping structure
130: If the user has not selected how to handle the global to local mapping then use HASH for "large" problems
131: */
132: static PetscErrorCode ISGlobalToLocalMappingSetUp(ISLocalToGlobalMapping mapping)
133: {
134: PetscInt i, *idx = mapping->indices, n = mapping->n, end, start;
136: PetscFunctionBegin;
137: if (mapping->data) PetscFunctionReturn(PETSC_SUCCESS);
138: end = 0;
139: start = PETSC_MAX_INT;
141: for (i = 0; i < n; i++) {
142: if (idx[i] < 0) continue;
143: if (idx[i] < start) start = idx[i];
144: if (idx[i] > end) end = idx[i];
145: }
146: if (start > end) {
147: start = 0;
148: end = -1;
149: }
150: mapping->globalstart = start;
151: mapping->globalend = end;
152: if (!((PetscObject)mapping)->type_name) {
153: if ((end - start) > PetscMax(4 * n, 1000000)) {
154: PetscCall(ISLocalToGlobalMappingSetType(mapping, ISLOCALTOGLOBALMAPPINGHASH));
155: } else {
156: PetscCall(ISLocalToGlobalMappingSetType(mapping, ISLOCALTOGLOBALMAPPINGBASIC));
157: }
158: }
159: PetscTryTypeMethod(mapping, globaltolocalmappingsetup);
160: PetscFunctionReturn(PETSC_SUCCESS);
161: }
163: static PetscErrorCode ISGlobalToLocalMappingSetUp_Basic(ISLocalToGlobalMapping mapping)
164: {
165: PetscInt i, *idx = mapping->indices, n = mapping->n, end, start, *globals;
166: ISLocalToGlobalMapping_Basic *map;
168: PetscFunctionBegin;
169: start = mapping->globalstart;
170: end = mapping->globalend;
171: PetscCall(PetscNew(&map));
172: PetscCall(PetscMalloc1(end - start + 2, &globals));
173: map->globals = globals;
174: for (i = 0; i < end - start + 1; i++) globals[i] = -1;
175: for (i = 0; i < n; i++) {
176: if (idx[i] < 0) continue;
177: globals[idx[i] - start] = i;
178: }
179: mapping->data = (void *)map;
180: PetscFunctionReturn(PETSC_SUCCESS);
181: }
183: static PetscErrorCode ISGlobalToLocalMappingSetUp_Hash(ISLocalToGlobalMapping mapping)
184: {
185: PetscInt i, *idx = mapping->indices, n = mapping->n;
186: ISLocalToGlobalMapping_Hash *map;
188: PetscFunctionBegin;
189: PetscCall(PetscNew(&map));
190: PetscCall(PetscHMapICreate(&map->globalht));
191: for (i = 0; i < n; i++) {
192: if (idx[i] < 0) continue;
193: PetscCall(PetscHMapISet(map->globalht, idx[i], i));
194: }
195: mapping->data = (void *)map;
196: PetscFunctionReturn(PETSC_SUCCESS);
197: }
199: static PetscErrorCode ISLocalToGlobalMappingDestroy_Basic(ISLocalToGlobalMapping mapping)
200: {
201: ISLocalToGlobalMapping_Basic *map = (ISLocalToGlobalMapping_Basic *)mapping->data;
203: PetscFunctionBegin;
204: if (!map) PetscFunctionReturn(PETSC_SUCCESS);
205: PetscCall(PetscFree(map->globals));
206: PetscCall(PetscFree(mapping->data));
207: PetscFunctionReturn(PETSC_SUCCESS);
208: }
210: static PetscErrorCode ISLocalToGlobalMappingDestroy_Hash(ISLocalToGlobalMapping mapping)
211: {
212: ISLocalToGlobalMapping_Hash *map = (ISLocalToGlobalMapping_Hash *)mapping->data;
214: PetscFunctionBegin;
215: if (!map) PetscFunctionReturn(PETSC_SUCCESS);
216: PetscCall(PetscHMapIDestroy(&map->globalht));
217: PetscCall(PetscFree(mapping->data));
218: PetscFunctionReturn(PETSC_SUCCESS);
219: }
221: #define GTOLTYPE _Basic
222: #define GTOLNAME _Basic
223: #define GTOLBS mapping->bs
224: #define GTOL(g, local) \
225: do { \
226: local = map->globals[g / bs - start]; \
227: if (local >= 0) local = bs * local + (g % bs); \
228: } while (0)
230: #include <../src/vec/is/utils/isltog.h>
232: #define GTOLTYPE _Basic
233: #define GTOLNAME Block_Basic
234: #define GTOLBS 1
235: #define GTOL(g, local) \
236: do { \
237: local = map->globals[g - start]; \
238: } while (0)
239: #include <../src/vec/is/utils/isltog.h>
241: #define GTOLTYPE _Hash
242: #define GTOLNAME _Hash
243: #define GTOLBS mapping->bs
244: #define GTOL(g, local) \
245: do { \
246: (void)PetscHMapIGet(map->globalht, g / bs, &local); \
247: if (local >= 0) local = bs * local + (g % bs); \
248: } while (0)
249: #include <../src/vec/is/utils/isltog.h>
251: #define GTOLTYPE _Hash
252: #define GTOLNAME Block_Hash
253: #define GTOLBS 1
254: #define GTOL(g, local) \
255: do { \
256: (void)PetscHMapIGet(map->globalht, g, &local); \
257: } while (0)
258: #include <../src/vec/is/utils/isltog.h>
260: /*@
261: ISLocalToGlobalMappingDuplicate - Duplicates the local to global mapping object
263: Not Collective
265: Input Parameter:
266: . ltog - local to global mapping
268: Output Parameter:
269: . nltog - the duplicated local to global mapping
271: Level: advanced
273: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`
274: @*/
275: PetscErrorCode ISLocalToGlobalMappingDuplicate(ISLocalToGlobalMapping ltog, ISLocalToGlobalMapping *nltog)
276: {
277: ISLocalToGlobalMappingType l2gtype;
279: PetscFunctionBegin;
281: PetscCall(ISLocalToGlobalMappingCreate(PetscObjectComm((PetscObject)ltog), ltog->bs, ltog->n, ltog->indices, PETSC_COPY_VALUES, nltog));
282: PetscCall(ISLocalToGlobalMappingGetType(ltog, &l2gtype));
283: PetscCall(ISLocalToGlobalMappingSetType(*nltog, l2gtype));
284: PetscFunctionReturn(PETSC_SUCCESS);
285: }
287: /*@
288: ISLocalToGlobalMappingGetSize - Gets the local size of a local to global mapping
290: Not Collective
292: Input Parameter:
293: . ltog - local to global mapping
295: Output Parameter:
296: . n - the number of entries in the local mapping, `ISLocalToGlobalMappingGetIndices()` returns an array of this length
298: Level: advanced
300: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`
301: @*/
302: PetscErrorCode ISLocalToGlobalMappingGetSize(ISLocalToGlobalMapping mapping, PetscInt *n)
303: {
304: PetscFunctionBegin;
307: *n = mapping->bs * mapping->n;
308: PetscFunctionReturn(PETSC_SUCCESS);
309: }
311: /*@C
312: ISLocalToGlobalMappingViewFromOptions - View an `ISLocalToGlobalMapping` based on values in the options database
314: Collective
316: Input Parameters:
317: + A - the local to global mapping object
318: . obj - Optional object that provides the options prefix used for the options database query
319: - name - command line option
321: Level: intermediate
323: Note:
324: See `PetscObjectViewFromOptions()` for the available `PetscViewer` and `PetscViewerFormat`
326: .seealso: [](sec_scatter), `PetscViewer`, ``ISLocalToGlobalMapping`, `ISLocalToGlobalMappingView`, `PetscObjectViewFromOptions()`, `ISLocalToGlobalMappingCreate()`
327: @*/
328: PetscErrorCode ISLocalToGlobalMappingViewFromOptions(ISLocalToGlobalMapping A, PetscObject obj, const char name[])
329: {
330: PetscFunctionBegin;
332: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
333: PetscFunctionReturn(PETSC_SUCCESS);
334: }
336: /*@C
337: ISLocalToGlobalMappingView - View a local to global mapping
339: Not Collective
341: Input Parameters:
342: + ltog - local to global mapping
343: - viewer - viewer
345: Level: advanced
347: .seealso: [](sec_scatter), `PetscViewer`, `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`
348: @*/
349: PetscErrorCode ISLocalToGlobalMappingView(ISLocalToGlobalMapping mapping, PetscViewer viewer)
350: {
351: PetscInt i;
352: PetscMPIInt rank;
353: PetscBool iascii;
355: PetscFunctionBegin;
357: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)mapping), &viewer));
360: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)mapping), &rank));
361: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
362: if (iascii) {
363: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)mapping, viewer));
364: PetscCall(PetscViewerASCIIPushSynchronized(viewer));
365: for (i = 0; i < mapping->n; i++) {
366: PetscInt bs = mapping->bs, g = mapping->indices[i];
367: if (bs == 1) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d] %" PetscInt_FMT " %" PetscInt_FMT "\n", rank, i, g));
368: 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));
369: }
370: PetscCall(PetscViewerFlush(viewer));
371: PetscCall(PetscViewerASCIIPopSynchronized(viewer));
372: }
373: PetscFunctionReturn(PETSC_SUCCESS);
374: }
376: /*@
377: ISLocalToGlobalMappingCreateIS - Creates a mapping between a local (0 to n)
378: ordering and a global parallel ordering.
380: Not Collective
382: Input Parameter:
383: . is - index set containing the global numbers for each local number
385: Output Parameter:
386: . mapping - new mapping data structure
388: Level: advanced
390: Note:
391: the block size of the `IS` determines the block size of the mapping
393: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingSetFromOptions()`
394: @*/
395: PetscErrorCode ISLocalToGlobalMappingCreateIS(IS is, ISLocalToGlobalMapping *mapping)
396: {
397: PetscInt n, bs;
398: const PetscInt *indices;
399: MPI_Comm comm;
400: PetscBool isblock;
402: PetscFunctionBegin;
406: PetscCall(PetscObjectGetComm((PetscObject)is, &comm));
407: PetscCall(ISGetLocalSize(is, &n));
408: PetscCall(PetscObjectTypeCompare((PetscObject)is, ISBLOCK, &isblock));
409: if (!isblock) {
410: PetscCall(ISGetIndices(is, &indices));
411: PetscCall(ISLocalToGlobalMappingCreate(comm, 1, n, indices, PETSC_COPY_VALUES, mapping));
412: PetscCall(ISRestoreIndices(is, &indices));
413: } else {
414: PetscCall(ISGetBlockSize(is, &bs));
415: PetscCall(ISBlockGetIndices(is, &indices));
416: PetscCall(ISLocalToGlobalMappingCreate(comm, bs, n / bs, indices, PETSC_COPY_VALUES, mapping));
417: PetscCall(ISBlockRestoreIndices(is, &indices));
418: }
419: PetscFunctionReturn(PETSC_SUCCESS);
420: }
422: /*@C
423: ISLocalToGlobalMappingCreateSF - Creates a mapping between a local (0 to n)
424: ordering and a global parallel ordering.
426: Collective
428: Input Parameters:
429: + sf - star forest mapping contiguous local indices to (rank, offset)
430: - start - first global index on this process, or `PETSC_DECIDE` to compute contiguous global numbering automatically
432: Output Parameter:
433: . mapping - new mapping data structure
435: Level: advanced
437: Note:
438: If any processor calls this with `start` = `PETSC_DECIDE` then all processors must, otherwise the program will hang.
440: .seealso: [](sec_scatter), `PetscSF`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingSetFromOptions()`
441: @*/
442: PetscErrorCode ISLocalToGlobalMappingCreateSF(PetscSF sf, PetscInt start, ISLocalToGlobalMapping *mapping)
443: {
444: PetscInt i, maxlocal, nroots, nleaves, *globals, *ltog;
445: MPI_Comm comm;
447: PetscFunctionBegin;
450: PetscCall(PetscObjectGetComm((PetscObject)sf, &comm));
451: PetscCall(PetscSFGetGraph(sf, &nroots, &nleaves, NULL, NULL));
452: if (start == PETSC_DECIDE) {
453: start = 0;
454: PetscCallMPI(MPI_Exscan(&nroots, &start, 1, MPIU_INT, MPI_SUM, comm));
455: } else PetscCheck(start >= 0, comm, PETSC_ERR_ARG_OUTOFRANGE, "start must be nonnegative or PETSC_DECIDE");
456: PetscCall(PetscSFGetLeafRange(sf, NULL, &maxlocal));
457: ++maxlocal;
458: PetscCall(PetscMalloc1(nroots, &globals));
459: PetscCall(PetscMalloc1(maxlocal, <og));
460: for (i = 0; i < nroots; i++) globals[i] = start + i;
461: for (i = 0; i < maxlocal; i++) ltog[i] = -1;
462: PetscCall(PetscSFBcastBegin(sf, MPIU_INT, globals, ltog, MPI_REPLACE));
463: PetscCall(PetscSFBcastEnd(sf, MPIU_INT, globals, ltog, MPI_REPLACE));
464: PetscCall(ISLocalToGlobalMappingCreate(comm, 1, maxlocal, ltog, PETSC_OWN_POINTER, mapping));
465: PetscCall(PetscFree(globals));
466: PetscFunctionReturn(PETSC_SUCCESS);
467: }
469: /*@
470: ISLocalToGlobalMappingSetBlockSize - Sets the blocksize of the mapping
472: Not Collective
474: Input Parameters:
475: + mapping - mapping data structure
476: - bs - the blocksize
478: Level: advanced
480: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`
481: @*/
482: PetscErrorCode ISLocalToGlobalMappingSetBlockSize(ISLocalToGlobalMapping mapping, PetscInt bs)
483: {
484: PetscInt *nid;
485: const PetscInt *oid;
486: PetscInt i, cn, on, obs, nn;
488: PetscFunctionBegin;
490: PetscCheck(bs >= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid block size %" PetscInt_FMT, bs);
491: if (bs == mapping->bs) PetscFunctionReturn(PETSC_SUCCESS);
492: on = mapping->n;
493: obs = mapping->bs;
494: oid = mapping->indices;
495: nn = (on * obs) / bs;
496: 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);
498: PetscCall(PetscMalloc1(nn, &nid));
499: PetscCall(ISLocalToGlobalMappingGetIndices(mapping, &oid));
500: for (i = 0; i < nn; i++) {
501: PetscInt j;
502: for (j = 0, cn = 0; j < bs - 1; j++) {
503: if (oid[i * bs + j] < 0) {
504: cn++;
505: continue;
506: }
507: 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]);
508: }
509: if (oid[i * bs + j] < 0) cn++;
510: if (cn) {
511: 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);
512: nid[i] = -1;
513: } else {
514: nid[i] = oid[i * bs] / bs;
515: }
516: }
517: PetscCall(ISLocalToGlobalMappingRestoreIndices(mapping, &oid));
519: mapping->n = nn;
520: mapping->bs = bs;
521: PetscCall(PetscFree(mapping->indices));
522: mapping->indices = nid;
523: mapping->globalstart = 0;
524: mapping->globalend = 0;
526: /* reset the cached information */
527: PetscCall(PetscFree(mapping->info_procs));
528: PetscCall(PetscFree(mapping->info_numprocs));
529: if (mapping->info_indices) {
530: PetscInt i;
532: PetscCall(PetscFree((mapping->info_indices)[0]));
533: for (i = 1; i < mapping->info_nproc; i++) PetscCall(PetscFree(mapping->info_indices[i]));
534: PetscCall(PetscFree(mapping->info_indices));
535: }
536: mapping->info_cached = PETSC_FALSE;
538: PetscTryTypeMethod(mapping, destroy);
539: PetscFunctionReturn(PETSC_SUCCESS);
540: }
542: /*@
543: ISLocalToGlobalMappingGetBlockSize - Gets the blocksize of the mapping
544: ordering and a global parallel ordering.
546: Not Collective
548: Input Parameter:
549: . mapping - mapping data structure
551: Output Parameter:
552: . bs - the blocksize
554: Level: advanced
556: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`
557: @*/
558: PetscErrorCode ISLocalToGlobalMappingGetBlockSize(ISLocalToGlobalMapping mapping, PetscInt *bs)
559: {
560: PetscFunctionBegin;
562: *bs = mapping->bs;
563: PetscFunctionReturn(PETSC_SUCCESS);
564: }
566: /*@
567: ISLocalToGlobalMappingCreate - Creates a mapping between a local (0 to n)
568: ordering and a global parallel ordering.
570: Not Collective, but communicator may have more than one process
572: Input Parameters:
573: + comm - MPI communicator
574: . bs - the block size
575: . n - the number of local elements divided by the block size, or equivalently the number of block indices
576: . 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
577: - mode - see PetscCopyMode
579: Output Parameter:
580: . mapping - new mapping data structure
582: Level: advanced
584: Notes:
585: There is one integer value in indices per block and it represents the actual indices bs*idx + j, where j=0,..,bs-1
587: For "small" problems when using `ISGlobalToLocalMappingApply()` and `ISGlobalToLocalMappingApplyBlock()`, the `ISLocalToGlobalMappingType`
588: of `ISLOCALTOGLOBALMAPPINGBASIC` will be used; this uses more memory but is faster; this approach is not scalable for extremely large mappings.
590: For large problems `ISLOCALTOGLOBALMAPPINGHASH` is used, this is scalable.
591: Use `ISLocalToGlobalMappingSetType()` or call `ISLocalToGlobalMappingSetFromOptions()` with the option
592: `-islocaltoglobalmapping_type` <`basic`,`hash`> to control which is used.
594: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingSetFromOptions()`,
595: `ISLOCALTOGLOBALMAPPINGBASIC`, `ISLOCALTOGLOBALMAPPINGHASH`
596: `ISLocalToGlobalMappingSetType()`, `ISLocalToGlobalMappingType`
597: @*/
598: PetscErrorCode ISLocalToGlobalMappingCreate(MPI_Comm comm, PetscInt bs, PetscInt n, const PetscInt indices[], PetscCopyMode mode, ISLocalToGlobalMapping *mapping)
599: {
600: PetscInt *in;
602: PetscFunctionBegin;
606: *mapping = NULL;
607: PetscCall(ISInitializePackage());
609: PetscCall(PetscHeaderCreate(*mapping, IS_LTOGM_CLASSID, "ISLocalToGlobalMapping", "Local to global mapping", "IS", comm, ISLocalToGlobalMappingDestroy, ISLocalToGlobalMappingView));
610: (*mapping)->n = n;
611: (*mapping)->bs = bs;
612: if (mode == PETSC_COPY_VALUES) {
613: PetscCall(PetscMalloc1(n, &in));
614: PetscCall(PetscArraycpy(in, indices, n));
615: (*mapping)->indices = in;
616: (*mapping)->dealloc_indices = PETSC_TRUE;
617: } else if (mode == PETSC_OWN_POINTER) {
618: (*mapping)->indices = (PetscInt *)indices;
619: (*mapping)->dealloc_indices = PETSC_TRUE;
620: } else if (mode == PETSC_USE_POINTER) {
621: (*mapping)->indices = (PetscInt *)indices;
622: } else SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "Invalid mode %d", mode);
623: PetscFunctionReturn(PETSC_SUCCESS);
624: }
626: PetscFunctionList ISLocalToGlobalMappingList = NULL;
628: /*@
629: ISLocalToGlobalMappingSetFromOptions - Set mapping options from the options database.
631: Not Collective
633: Input Parameter:
634: . mapping - mapping data structure
636: Options Database Key:
637: . -islocaltoglobalmapping_type - <basic,hash> nonscalable and scalable versions
639: Level: advanced
641: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingSetFromOptions()`,
642: `ISLOCALTOGLOBALMAPPINGBASIC`, `ISLOCALTOGLOBALMAPPINGHASH`
643: `ISLocalToGlobalMappingSetType()`, `ISLocalToGlobalMappingType`
644: @*/
645: PetscErrorCode ISLocalToGlobalMappingSetFromOptions(ISLocalToGlobalMapping mapping)
646: {
647: char type[256];
648: ISLocalToGlobalMappingType defaulttype = "Not set";
649: PetscBool flg;
651: PetscFunctionBegin;
653: PetscCall(ISLocalToGlobalMappingRegisterAll());
654: PetscObjectOptionsBegin((PetscObject)mapping);
655: PetscCall(PetscOptionsFList("-islocaltoglobalmapping_type", "ISLocalToGlobalMapping method", "ISLocalToGlobalMappingSetType", ISLocalToGlobalMappingList, (char *)(((PetscObject)mapping)->type_name) ? ((PetscObject)mapping)->type_name : defaulttype, type, 256, &flg));
656: if (flg) PetscCall(ISLocalToGlobalMappingSetType(mapping, type));
657: PetscOptionsEnd();
658: PetscFunctionReturn(PETSC_SUCCESS);
659: }
661: /*@
662: ISLocalToGlobalMappingDestroy - Destroys a mapping between a local (0 to n)
663: ordering and a global parallel ordering.
665: Not Collective
667: Input Parameter:
668: . mapping - mapping data structure
670: Level: advanced
672: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingCreate()`
673: @*/
674: PetscErrorCode ISLocalToGlobalMappingDestroy(ISLocalToGlobalMapping *mapping)
675: {
676: PetscFunctionBegin;
677: if (!*mapping) PetscFunctionReturn(PETSC_SUCCESS);
679: if (--((PetscObject)(*mapping))->refct > 0) {
680: *mapping = NULL;
681: PetscFunctionReturn(PETSC_SUCCESS);
682: }
683: if ((*mapping)->dealloc_indices) PetscCall(PetscFree((*mapping)->indices));
684: PetscCall(PetscFree((*mapping)->info_procs));
685: PetscCall(PetscFree((*mapping)->info_numprocs));
686: if ((*mapping)->info_indices) {
687: PetscInt i;
689: PetscCall(PetscFree(((*mapping)->info_indices)[0]));
690: for (i = 1; i < (*mapping)->info_nproc; i++) PetscCall(PetscFree(((*mapping)->info_indices)[i]));
691: PetscCall(PetscFree((*mapping)->info_indices));
692: }
693: if ((*mapping)->info_nodei) PetscCall(PetscFree(((*mapping)->info_nodei)[0]));
694: PetscCall(PetscFree2((*mapping)->info_nodec, (*mapping)->info_nodei));
695: if ((*mapping)->ops->destroy) PetscCall((*(*mapping)->ops->destroy)(*mapping));
696: PetscCall(PetscHeaderDestroy(mapping));
697: *mapping = NULL;
698: PetscFunctionReturn(PETSC_SUCCESS);
699: }
701: /*@
702: ISLocalToGlobalMappingApplyIS - Creates from an `IS` in the local numbering
703: a new index set using the global numbering defined in an `ISLocalToGlobalMapping`
704: context.
706: Collective
708: Input Parameters:
709: + mapping - mapping between local and global numbering
710: - is - index set in local numbering
712: Output Parameter:
713: . newis - index set in global numbering
715: Level: advanced
717: Note:
718: The output `IS` will have the same communicator as the input `IS`.
720: .seealso: [](sec_scatter), `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingCreate()`,
721: `ISLocalToGlobalMappingDestroy()`, `ISGlobalToLocalMappingApply()`
722: @*/
723: PetscErrorCode ISLocalToGlobalMappingApplyIS(ISLocalToGlobalMapping mapping, IS is, IS *newis)
724: {
725: PetscInt n, *idxout;
726: const PetscInt *idxin;
728: PetscFunctionBegin;
733: PetscCall(ISGetLocalSize(is, &n));
734: PetscCall(ISGetIndices(is, &idxin));
735: PetscCall(PetscMalloc1(n, &idxout));
736: PetscCall(ISLocalToGlobalMappingApply(mapping, n, idxin, idxout));
737: PetscCall(ISRestoreIndices(is, &idxin));
738: PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)is), n, idxout, PETSC_OWN_POINTER, newis));
739: PetscFunctionReturn(PETSC_SUCCESS);
740: }
742: /*@
743: ISLocalToGlobalMappingApply - Takes a list of integers in a local numbering
744: and converts them to the global numbering.
746: Not Collective
748: Input Parameters:
749: + mapping - the local to global mapping context
750: . N - number of integers
751: - in - input indices in local numbering
753: Output Parameter:
754: . out - indices in global numbering
756: Level: advanced
758: Note:
759: The `in` and `out` array parameters may be identical.
761: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingApplyBlock()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingDestroy()`,
762: `ISLocalToGlobalMappingApplyIS()`, `AOCreateBasic()`, `AOApplicationToPetsc()`,
763: `AOPetscToApplication()`, `ISGlobalToLocalMappingApply()`
764: @*/
765: PetscErrorCode ISLocalToGlobalMappingApply(ISLocalToGlobalMapping mapping, PetscInt N, const PetscInt in[], PetscInt out[])
766: {
767: PetscInt i, bs, Nmax;
769: PetscFunctionBegin;
771: bs = mapping->bs;
772: Nmax = bs * mapping->n;
773: if (bs == 1) {
774: const PetscInt *idx = mapping->indices;
775: for (i = 0; i < N; i++) {
776: if (in[i] < 0) {
777: out[i] = in[i];
778: continue;
779: }
780: 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);
781: out[i] = idx[in[i]];
782: }
783: } else {
784: const PetscInt *idx = mapping->indices;
785: for (i = 0; i < N; i++) {
786: if (in[i] < 0) {
787: out[i] = in[i];
788: continue;
789: }
790: 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);
791: out[i] = idx[in[i] / bs] * bs + (in[i] % bs);
792: }
793: }
794: PetscFunctionReturn(PETSC_SUCCESS);
795: }
797: /*@
798: ISLocalToGlobalMappingApplyBlock - Takes a list of integers in a local block numbering and converts them to the global block numbering
800: Not Collective
802: Input Parameters:
803: + mapping - the local to global mapping context
804: . N - number of integers
805: - in - input indices in local block numbering
807: Output Parameter:
808: . out - indices in global block numbering
810: Example:
811: 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
812: (the first block) would produce 0 and the mapping applied to 1 (the second block) would produce 3.
814: Level: advanced
816: Note:
817: The `in` and `out` array parameters may be identical.
819: .seealso: [](sec_scatter), `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingDestroy()`,
820: `ISLocalToGlobalMappingApplyIS()`, `AOCreateBasic()`, `AOApplicationToPetsc()`,
821: `AOPetscToApplication()`, `ISGlobalToLocalMappingApply()`
822: @*/
823: PetscErrorCode ISLocalToGlobalMappingApplyBlock(ISLocalToGlobalMapping mapping, PetscInt N, const PetscInt in[], PetscInt out[])
824: {
825: PetscInt i, Nmax;
826: const PetscInt *idx;
828: PetscFunctionBegin;
830: Nmax = mapping->n;
831: idx = mapping->indices;
832: for (i = 0; i < N; i++) {
833: if (in[i] < 0) {
834: out[i] = in[i];
835: continue;
836: }
837: 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);
838: out[i] = idx[in[i]];
839: }
840: PetscFunctionReturn(PETSC_SUCCESS);
841: }
843: /*@
844: ISGlobalToLocalMappingApply - Provides the local numbering for a list of integers
845: specified with a global numbering.
847: Not Collective
849: Input Parameters:
850: + mapping - mapping between local and global numbering
851: . type - `IS_GTOLM_MASK` - maps global indices with no local value to -1 in the output list (i.e., mask them)
852: `IS_GTOLM_DROP` - drops the indices with no local value from the output list
853: . n - number of global indices to map
854: - idx - global indices to map
856: Output Parameters:
857: + nout - number of indices in output array (if type == `IS_GTOLM_MASK` then nout = n)
858: - idxout - local index of each global index, one must pass in an array long enough
859: to hold all the indices. You can call `ISGlobalToLocalMappingApply()` with
860: idxout == NULL to determine the required length (returned in nout)
861: and then allocate the required space and call `ISGlobalToLocalMappingApply()`
862: a second time to set the values.
864: Level: advanced
866: Notes:
867: Either `nout` or `idxout` may be `NULL`. `idx` and `idxout` may be identical.
869: For "small" problems when using `ISGlobalToLocalMappingApply()` and `ISGlobalToLocalMappingApplyBlock()`, the `ISLocalToGlobalMappingType` of
870: `ISLOCALTOGLOBALMAPPINGBASIC` will be used;
871: 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.
872: Use `ISLocalToGlobalMappingSetType()` or call `ISLocalToGlobalMappingSetFromOptions()` with the option -islocaltoglobalmapping_type <basic,hash> to control which is used.
874: Developer Note:
875: The manual page states that `idx` and `idxout` may be identical but the calling
876: sequence declares `idx` as const so it cannot be the same as `idxout`.
878: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingApply()`, `ISGlobalToLocalMappingApplyBlock()`, `ISLocalToGlobalMappingCreate()`,
879: `ISLocalToGlobalMappingDestroy()`
880: @*/
881: PetscErrorCode ISGlobalToLocalMappingApply(ISLocalToGlobalMapping mapping, ISGlobalToLocalMappingMode type, PetscInt n, const PetscInt idx[], PetscInt *nout, PetscInt idxout[])
882: {
883: PetscFunctionBegin;
885: if (!mapping->data) PetscCall(ISGlobalToLocalMappingSetUp(mapping));
886: PetscUseTypeMethod(mapping, globaltolocalmappingapply, type, n, idx, nout, idxout);
887: PetscFunctionReturn(PETSC_SUCCESS);
888: }
890: /*@
891: ISGlobalToLocalMappingApplyIS - Creates from an `IS` in the global numbering
892: a new index set using the local numbering defined in an `ISLocalToGlobalMapping`
893: context.
895: Not Collective
897: Input Parameters:
898: + mapping - mapping between local and global numbering
899: . type - `IS_GTOLM_MASK` - maps global indices with no local value to -1 in the output list (i.e., mask them)
900: `IS_GTOLM_DROP` - drops the indices with no local value from the output list
901: - is - index set in global numbering
903: Output Parameter:
904: . newis - index set in local numbering
906: Level: advanced
908: Note:
909: The output `IS` will be sequential, as it encodes a purely local operation
911: .seealso: [](sec_scatter), `ISGlobalToLocalMapping`, `ISGlobalToLocalMappingApply()`, `ISLocalToGlobalMappingCreate()`,
912: `ISLocalToGlobalMappingDestroy()`
913: @*/
914: PetscErrorCode ISGlobalToLocalMappingApplyIS(ISLocalToGlobalMapping mapping, ISGlobalToLocalMappingMode type, IS is, IS *newis)
915: {
916: PetscInt n, nout, *idxout;
917: const PetscInt *idxin;
919: PetscFunctionBegin;
924: PetscCall(ISGetLocalSize(is, &n));
925: PetscCall(ISGetIndices(is, &idxin));
926: if (type == IS_GTOLM_MASK) {
927: PetscCall(PetscMalloc1(n, &idxout));
928: } else {
929: PetscCall(ISGlobalToLocalMappingApply(mapping, type, n, idxin, &nout, NULL));
930: PetscCall(PetscMalloc1(nout, &idxout));
931: }
932: PetscCall(ISGlobalToLocalMappingApply(mapping, type, n, idxin, &nout, idxout));
933: PetscCall(ISRestoreIndices(is, &idxin));
934: PetscCall(ISCreateGeneral(PETSC_COMM_SELF, nout, idxout, PETSC_OWN_POINTER, newis));
935: PetscFunctionReturn(PETSC_SUCCESS);
936: }
938: /*@
939: ISGlobalToLocalMappingApplyBlock - Provides the local block numbering for a list of integers
940: specified with a block global numbering.
942: Not Collective
944: Input Parameters:
945: + mapping - mapping between local and global numbering
946: . type - `IS_GTOLM_MASK` - maps global indices with no local value to -1 in the output list (i.e., mask them)
947: `IS_GTOLM_DROP` - drops the indices with no local value from the output list
948: . n - number of global indices to map
949: - idx - global indices to map
951: Output Parameters:
952: + nout - number of indices in output array (if type == `IS_GTOLM_MASK` then nout = n)
953: - idxout - local index of each global index, one must pass in an array long enough
954: to hold all the indices. You can call `ISGlobalToLocalMappingApplyBlock()` with
955: idxout == NULL to determine the required length (returned in nout)
956: and then allocate the required space and call `ISGlobalToLocalMappingApplyBlock()`
957: a second time to set the values.
959: Level: advanced
961: Notes:
962: Either `nout` or `idxout` may be `NULL`. `idx` and `idxout` may be identical.
964: For "small" problems when using `ISGlobalToLocalMappingApply()` and `ISGlobalToLocalMappingApplyBlock()`, the `ISLocalToGlobalMappingType` of
965: `ISLOCALTOGLOBALMAPPINGBASIC` will be used;
966: 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.
967: Use `ISLocalToGlobalMappingSetType()` or call `ISLocalToGlobalMappingSetFromOptions()` with the option -islocaltoglobalmapping_type <basic,hash> to control which is used.
969: Developer Note:
970: The manual page states that `idx` and `idxout` may be identical but the calling
971: sequence declares `idx` as const so it cannot be the same as `idxout`.
973: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingApply()`, `ISGlobalToLocalMappingApply()`, `ISLocalToGlobalMappingCreate()`,
974: `ISLocalToGlobalMappingDestroy()`
975: @*/
976: PetscErrorCode ISGlobalToLocalMappingApplyBlock(ISLocalToGlobalMapping mapping, ISGlobalToLocalMappingMode type, PetscInt n, const PetscInt idx[], PetscInt *nout, PetscInt idxout[])
977: {
978: PetscFunctionBegin;
980: if (!mapping->data) PetscCall(ISGlobalToLocalMappingSetUp(mapping));
981: PetscUseTypeMethod(mapping, globaltolocalmappingapplyblock, type, n, idx, nout, idxout);
982: PetscFunctionReturn(PETSC_SUCCESS);
983: }
985: /*@C
986: ISLocalToGlobalMappingGetBlockInfo - Gets the neighbor information for each processor and
987: each index shared by more than one processor
989: Collective
991: Input Parameter:
992: . mapping - the mapping from local to global indexing
994: Output Parameters:
995: + nproc - number of processors that are connected to this one
996: . proc - neighboring processors
997: . numproc - number of indices for each subdomain (processor)
998: - indices - indices of nodes (in local numbering) shared with neighbors (sorted by global numbering)
1000: Level: advanced
1002: Fortran Usage:
1003: .vb
1004: PetscInt indices[nproc][numprocmax],ierr)
1005: ISLocalToGlobalMpngGetInfoSize(ISLocalToGlobalMapping,PetscInt nproc,PetscInt numprocmax,ierr) followed by
1006: ISLocalToGlobalMappingGetInfo(ISLocalToGlobalMapping,PetscInt nproc, PetscInt procs[nproc],PetscInt numprocs[nproc],
1007: .ve
1009: Fortran Note:
1010: There is no `ISLocalToGlobalMappingRestoreInfo()` in Fortran. You must make sure that `procs`[], `numprocs`[] and
1011: `indices`[][] are large enough arrays, either by allocating them dynamically or defining static ones large enough.
1013: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1014: `ISLocalToGlobalMappingRestoreInfo()`
1015: @*/
1016: PetscErrorCode ISLocalToGlobalMappingGetBlockInfo(ISLocalToGlobalMapping mapping, PetscInt *nproc, PetscInt *procs[], PetscInt *numprocs[], PetscInt **indices[])
1017: {
1018: PetscFunctionBegin;
1020: if (mapping->info_cached) {
1021: *nproc = mapping->info_nproc;
1022: *procs = mapping->info_procs;
1023: *numprocs = mapping->info_numprocs;
1024: *indices = mapping->info_indices;
1025: } else {
1026: PetscCall(ISLocalToGlobalMappingGetBlockInfo_Private(mapping, nproc, procs, numprocs, indices));
1027: }
1028: PetscFunctionReturn(PETSC_SUCCESS);
1029: }
1031: static PetscErrorCode ISLocalToGlobalMappingGetBlockInfo_Private(ISLocalToGlobalMapping mapping, PetscInt *nproc, PetscInt *procs[], PetscInt *numprocs[], PetscInt **indices[])
1032: {
1033: PetscMPIInt size, rank, tag1, tag2, tag3, *len, *source, imdex;
1034: PetscInt i, n = mapping->n, Ng, ng, max = 0, *lindices = mapping->indices;
1035: PetscInt *nprocs, *owner, nsends, *sends, j, *starts, nmax, nrecvs, *recvs, proc;
1036: PetscInt cnt, scale, *ownedsenders, *nownedsenders, rstart;
1037: PetscInt node, nownedm, nt, *sends2, nsends2, *starts2, *lens2, *dest, nrecvs2, *starts3, *recvs2, k, *bprocs, *tmp;
1038: PetscInt first_procs, first_numprocs, *first_indices;
1039: MPI_Request *recv_waits, *send_waits;
1040: MPI_Status recv_status, *send_status, *recv_statuses;
1041: MPI_Comm comm;
1042: PetscBool debug = PETSC_FALSE;
1044: PetscFunctionBegin;
1045: PetscCall(PetscObjectGetComm((PetscObject)mapping, &comm));
1046: PetscCallMPI(MPI_Comm_size(comm, &size));
1047: PetscCallMPI(MPI_Comm_rank(comm, &rank));
1048: if (size == 1) {
1049: *nproc = 0;
1050: *procs = NULL;
1051: PetscCall(PetscNew(numprocs));
1052: (*numprocs)[0] = 0;
1053: PetscCall(PetscNew(indices));
1054: (*indices)[0] = NULL;
1055: /* save info for reuse */
1056: mapping->info_nproc = *nproc;
1057: mapping->info_procs = *procs;
1058: mapping->info_numprocs = *numprocs;
1059: mapping->info_indices = *indices;
1060: mapping->info_cached = PETSC_TRUE;
1061: PetscFunctionReturn(PETSC_SUCCESS);
1062: }
1064: PetscCall(PetscOptionsGetBool(((PetscObject)mapping)->options, NULL, "-islocaltoglobalmappinggetinfo_debug", &debug, NULL));
1066: /*
1067: Notes on ISLocalToGlobalMappingGetBlockInfo
1069: globally owned node - the nodes that have been assigned to this processor in global
1070: numbering, just for this routine.
1072: nontrivial globally owned node - node assigned to this processor that is on a subdomain
1073: boundary (i.e. is has more than one local owner)
1075: locally owned node - node that exists on this processors subdomain
1077: nontrivial locally owned node - node that is not in the interior (i.e. has more than one
1078: local subdomain
1079: */
1080: PetscCall(PetscObjectGetNewTag((PetscObject)mapping, &tag1));
1081: PetscCall(PetscObjectGetNewTag((PetscObject)mapping, &tag2));
1082: PetscCall(PetscObjectGetNewTag((PetscObject)mapping, &tag3));
1084: for (i = 0; i < n; i++) {
1085: if (lindices[i] > max) max = lindices[i];
1086: }
1087: PetscCall(MPIU_Allreduce(&max, &Ng, 1, MPIU_INT, MPI_MAX, comm));
1088: Ng++;
1089: PetscCallMPI(MPI_Comm_size(comm, &size));
1090: PetscCallMPI(MPI_Comm_rank(comm, &rank));
1091: scale = Ng / size + 1;
1092: ng = scale;
1093: if (rank == size - 1) ng = Ng - scale * (size - 1);
1094: ng = PetscMax(1, ng);
1095: rstart = scale * rank;
1097: /* determine ownership ranges of global indices */
1098: PetscCall(PetscMalloc1(2 * size, &nprocs));
1099: PetscCall(PetscArrayzero(nprocs, 2 * size));
1101: /* determine owners of each local node */
1102: PetscCall(PetscMalloc1(n, &owner));
1103: for (i = 0; i < n; i++) {
1104: proc = lindices[i] / scale; /* processor that globally owns this index */
1105: nprocs[2 * proc + 1] = 1; /* processor globally owns at least one of ours */
1106: owner[i] = proc;
1107: nprocs[2 * proc]++; /* count of how many that processor globally owns of ours */
1108: }
1109: nsends = 0;
1110: for (i = 0; i < size; i++) nsends += nprocs[2 * i + 1];
1111: PetscCall(PetscInfo(mapping, "Number of global owners for my local data %" PetscInt_FMT "\n", nsends));
1113: /* inform other processors of number of messages and max length*/
1114: PetscCall(PetscMaxSum(comm, nprocs, &nmax, &nrecvs));
1115: PetscCall(PetscInfo(mapping, "Number of local owners for my global data %" PetscInt_FMT "\n", nrecvs));
1117: /* post receives for owned rows */
1118: PetscCall(PetscMalloc1((2 * nrecvs + 1) * (nmax + 1), &recvs));
1119: PetscCall(PetscMalloc1(nrecvs + 1, &recv_waits));
1120: for (i = 0; i < nrecvs; i++) PetscCallMPI(MPI_Irecv(recvs + 2 * nmax * i, 2 * nmax, MPIU_INT, MPI_ANY_SOURCE, tag1, comm, recv_waits + i));
1122: /* pack messages containing lists of local nodes to owners */
1123: PetscCall(PetscMalloc1(2 * n + 1, &sends));
1124: PetscCall(PetscMalloc1(size + 1, &starts));
1125: starts[0] = 0;
1126: for (i = 1; i < size; i++) starts[i] = starts[i - 1] + 2 * nprocs[2 * i - 2];
1127: for (i = 0; i < n; i++) {
1128: sends[starts[owner[i]]++] = lindices[i];
1129: sends[starts[owner[i]]++] = i;
1130: }
1131: PetscCall(PetscFree(owner));
1132: starts[0] = 0;
1133: for (i = 1; i < size; i++) starts[i] = starts[i - 1] + 2 * nprocs[2 * i - 2];
1135: /* send the messages */
1136: PetscCall(PetscMalloc1(nsends + 1, &send_waits));
1137: PetscCall(PetscMalloc1(nsends + 1, &dest));
1138: cnt = 0;
1139: for (i = 0; i < size; i++) {
1140: if (nprocs[2 * i]) {
1141: PetscCallMPI(MPI_Isend(sends + starts[i], 2 * nprocs[2 * i], MPIU_INT, i, tag1, comm, send_waits + cnt));
1142: dest[cnt] = i;
1143: cnt++;
1144: }
1145: }
1146: PetscCall(PetscFree(starts));
1148: /* wait on receives */
1149: PetscCall(PetscMalloc1(nrecvs + 1, &source));
1150: PetscCall(PetscMalloc1(nrecvs + 1, &len));
1151: cnt = nrecvs;
1152: PetscCall(PetscCalloc1(ng + 1, &nownedsenders));
1153: while (cnt) {
1154: PetscCallMPI(MPI_Waitany(nrecvs, recv_waits, &imdex, &recv_status));
1155: /* unpack receives into our local space */
1156: PetscCallMPI(MPI_Get_count(&recv_status, MPIU_INT, &len[imdex]));
1157: source[imdex] = recv_status.MPI_SOURCE;
1158: len[imdex] = len[imdex] / 2;
1159: /* count how many local owners for each of my global owned indices */
1160: for (i = 0; i < len[imdex]; i++) nownedsenders[recvs[2 * imdex * nmax + 2 * i] - rstart]++;
1161: cnt--;
1162: }
1163: PetscCall(PetscFree(recv_waits));
1165: /* count how many globally owned indices are on an edge multiplied by how many processors own them. */
1166: nownedm = 0;
1167: for (i = 0; i < ng; i++) {
1168: if (nownedsenders[i] > 1) nownedm += nownedsenders[i];
1169: }
1171: /* create single array to contain rank of all local owners of each globally owned index */
1172: PetscCall(PetscMalloc1(nownedm + 1, &ownedsenders));
1173: PetscCall(PetscMalloc1(ng + 1, &starts));
1174: starts[0] = 0;
1175: for (i = 1; i < ng; i++) {
1176: if (nownedsenders[i - 1] > 1) starts[i] = starts[i - 1] + nownedsenders[i - 1];
1177: else starts[i] = starts[i - 1];
1178: }
1180: /* for each nontrivial globally owned node list all arriving processors */
1181: for (i = 0; i < nrecvs; i++) {
1182: for (j = 0; j < len[i]; j++) {
1183: node = recvs[2 * i * nmax + 2 * j] - rstart;
1184: if (nownedsenders[node] > 1) ownedsenders[starts[node]++] = source[i];
1185: }
1186: }
1188: if (debug) { /* ----------------------------------- */
1189: starts[0] = 0;
1190: for (i = 1; i < ng; i++) {
1191: if (nownedsenders[i - 1] > 1) starts[i] = starts[i - 1] + nownedsenders[i - 1];
1192: else starts[i] = starts[i - 1];
1193: }
1194: for (i = 0; i < ng; i++) {
1195: if (nownedsenders[i] > 1) {
1196: PetscCall(PetscSynchronizedPrintf(comm, "[%d] global node %" PetscInt_FMT " local owner processors: ", rank, i + rstart));
1197: for (j = 0; j < nownedsenders[i]; j++) PetscCall(PetscSynchronizedPrintf(comm, "%" PetscInt_FMT " ", ownedsenders[starts[i] + j]));
1198: PetscCall(PetscSynchronizedPrintf(comm, "\n"));
1199: }
1200: }
1201: PetscCall(PetscSynchronizedFlush(comm, PETSC_STDOUT));
1202: } /* ----------------------------------- */
1204: /* wait on original sends */
1205: if (nsends) {
1206: PetscCall(PetscMalloc1(nsends, &send_status));
1207: PetscCallMPI(MPI_Waitall(nsends, send_waits, send_status));
1208: PetscCall(PetscFree(send_status));
1209: }
1210: PetscCall(PetscFree(send_waits));
1211: PetscCall(PetscFree(sends));
1212: PetscCall(PetscFree(nprocs));
1214: /* pack messages to send back to local owners */
1215: starts[0] = 0;
1216: for (i = 1; i < ng; i++) {
1217: if (nownedsenders[i - 1] > 1) starts[i] = starts[i - 1] + nownedsenders[i - 1];
1218: else starts[i] = starts[i - 1];
1219: }
1220: nsends2 = nrecvs;
1221: PetscCall(PetscMalloc1(nsends2 + 1, &nprocs)); /* length of each message */
1222: for (i = 0; i < nrecvs; i++) {
1223: nprocs[i] = 1;
1224: for (j = 0; j < len[i]; j++) {
1225: node = recvs[2 * i * nmax + 2 * j] - rstart;
1226: if (nownedsenders[node] > 1) nprocs[i] += 2 + nownedsenders[node];
1227: }
1228: }
1229: nt = 0;
1230: for (i = 0; i < nsends2; i++) nt += nprocs[i];
1232: PetscCall(PetscMalloc1(nt + 1, &sends2));
1233: PetscCall(PetscMalloc1(nsends2 + 1, &starts2));
1235: starts2[0] = 0;
1236: for (i = 1; i < nsends2; i++) starts2[i] = starts2[i - 1] + nprocs[i - 1];
1237: /*
1238: Each message is 1 + nprocs[i] long, and consists of
1239: (0) the number of nodes being sent back
1240: (1) the local node number,
1241: (2) the number of processors sharing it,
1242: (3) the processors sharing it
1243: */
1244: for (i = 0; i < nsends2; i++) {
1245: cnt = 1;
1246: sends2[starts2[i]] = 0;
1247: for (j = 0; j < len[i]; j++) {
1248: node = recvs[2 * i * nmax + 2 * j] - rstart;
1249: if (nownedsenders[node] > 1) {
1250: sends2[starts2[i]]++;
1251: sends2[starts2[i] + cnt++] = recvs[2 * i * nmax + 2 * j + 1];
1252: sends2[starts2[i] + cnt++] = nownedsenders[node];
1253: PetscCall(PetscArraycpy(&sends2[starts2[i] + cnt], &ownedsenders[starts[node]], nownedsenders[node]));
1254: cnt += nownedsenders[node];
1255: }
1256: }
1257: }
1259: /* receive the message lengths */
1260: nrecvs2 = nsends;
1261: PetscCall(PetscMalloc1(nrecvs2 + 1, &lens2));
1262: PetscCall(PetscMalloc1(nrecvs2 + 1, &starts3));
1263: PetscCall(PetscMalloc1(nrecvs2 + 1, &recv_waits));
1264: for (i = 0; i < nrecvs2; i++) PetscCallMPI(MPI_Irecv(&lens2[i], 1, MPIU_INT, dest[i], tag2, comm, recv_waits + i));
1266: /* send the message lengths */
1267: for (i = 0; i < nsends2; i++) PetscCallMPI(MPI_Send(&nprocs[i], 1, MPIU_INT, source[i], tag2, comm));
1269: /* wait on receives of lens */
1270: if (nrecvs2) {
1271: PetscCall(PetscMalloc1(nrecvs2, &recv_statuses));
1272: PetscCallMPI(MPI_Waitall(nrecvs2, recv_waits, recv_statuses));
1273: PetscCall(PetscFree(recv_statuses));
1274: }
1275: PetscCall(PetscFree(recv_waits));
1277: starts3[0] = 0;
1278: nt = 0;
1279: for (i = 0; i < nrecvs2 - 1; i++) {
1280: starts3[i + 1] = starts3[i] + lens2[i];
1281: nt += lens2[i];
1282: }
1283: if (nrecvs2) nt += lens2[nrecvs2 - 1];
1285: PetscCall(PetscMalloc1(nt + 1, &recvs2));
1286: PetscCall(PetscMalloc1(nrecvs2 + 1, &recv_waits));
1287: for (i = 0; i < nrecvs2; i++) PetscCallMPI(MPI_Irecv(recvs2 + starts3[i], lens2[i], MPIU_INT, dest[i], tag3, comm, recv_waits + i));
1289: /* send the messages */
1290: PetscCall(PetscMalloc1(nsends2 + 1, &send_waits));
1291: for (i = 0; i < nsends2; i++) PetscCallMPI(MPI_Isend(sends2 + starts2[i], nprocs[i], MPIU_INT, source[i], tag3, comm, send_waits + i));
1293: /* wait on receives */
1294: if (nrecvs2) {
1295: PetscCall(PetscMalloc1(nrecvs2, &recv_statuses));
1296: PetscCallMPI(MPI_Waitall(nrecvs2, recv_waits, recv_statuses));
1297: PetscCall(PetscFree(recv_statuses));
1298: }
1299: PetscCall(PetscFree(recv_waits));
1300: PetscCall(PetscFree(nprocs));
1302: if (debug) { /* ----------------------------------- */
1303: cnt = 0;
1304: for (i = 0; i < nrecvs2; i++) {
1305: nt = recvs2[cnt++];
1306: for (j = 0; j < nt; j++) {
1307: PetscCall(PetscSynchronizedPrintf(comm, "[%d] local node %" PetscInt_FMT " number of subdomains %" PetscInt_FMT ": ", rank, recvs2[cnt], recvs2[cnt + 1]));
1308: for (k = 0; k < recvs2[cnt + 1]; k++) PetscCall(PetscSynchronizedPrintf(comm, "%" PetscInt_FMT " ", recvs2[cnt + 2 + k]));
1309: cnt += 2 + recvs2[cnt + 1];
1310: PetscCall(PetscSynchronizedPrintf(comm, "\n"));
1311: }
1312: }
1313: PetscCall(PetscSynchronizedFlush(comm, PETSC_STDOUT));
1314: } /* ----------------------------------- */
1316: /* count number subdomains for each local node */
1317: PetscCall(PetscCalloc1(size, &nprocs));
1318: cnt = 0;
1319: for (i = 0; i < nrecvs2; i++) {
1320: nt = recvs2[cnt++];
1321: for (j = 0; j < nt; j++) {
1322: for (k = 0; k < recvs2[cnt + 1]; k++) nprocs[recvs2[cnt + 2 + k]]++;
1323: cnt += 2 + recvs2[cnt + 1];
1324: }
1325: }
1326: nt = 0;
1327: for (i = 0; i < size; i++) nt += (nprocs[i] > 0);
1328: *nproc = nt;
1329: PetscCall(PetscMalloc1(nt + 1, procs));
1330: PetscCall(PetscMalloc1(nt + 1, numprocs));
1331: PetscCall(PetscMalloc1(nt + 1, indices));
1332: for (i = 0; i < nt + 1; i++) (*indices)[i] = NULL;
1333: PetscCall(PetscMalloc1(size, &bprocs));
1334: cnt = 0;
1335: for (i = 0; i < size; i++) {
1336: if (nprocs[i] > 0) {
1337: bprocs[i] = cnt;
1338: (*procs)[cnt] = i;
1339: (*numprocs)[cnt] = nprocs[i];
1340: PetscCall(PetscMalloc1(nprocs[i], &(*indices)[cnt]));
1341: cnt++;
1342: }
1343: }
1345: /* make the list of subdomains for each nontrivial local node */
1346: PetscCall(PetscArrayzero(*numprocs, nt));
1347: cnt = 0;
1348: for (i = 0; i < nrecvs2; i++) {
1349: nt = recvs2[cnt++];
1350: for (j = 0; j < nt; j++) {
1351: for (k = 0; k < recvs2[cnt + 1]; k++) (*indices)[bprocs[recvs2[cnt + 2 + k]]][(*numprocs)[bprocs[recvs2[cnt + 2 + k]]]++] = recvs2[cnt];
1352: cnt += 2 + recvs2[cnt + 1];
1353: }
1354: }
1355: PetscCall(PetscFree(bprocs));
1356: PetscCall(PetscFree(recvs2));
1358: /* sort the node indexing by their global numbers */
1359: nt = *nproc;
1360: for (i = 0; i < nt; i++) {
1361: PetscCall(PetscMalloc1((*numprocs)[i], &tmp));
1362: for (j = 0; j < (*numprocs)[i]; j++) tmp[j] = lindices[(*indices)[i][j]];
1363: PetscCall(PetscSortIntWithArray((*numprocs)[i], tmp, (*indices)[i]));
1364: PetscCall(PetscFree(tmp));
1365: }
1367: if (debug) { /* ----------------------------------- */
1368: nt = *nproc;
1369: for (i = 0; i < nt; i++) {
1370: PetscCall(PetscSynchronizedPrintf(comm, "[%d] subdomain %" PetscInt_FMT " number of indices %" PetscInt_FMT ": ", rank, (*procs)[i], (*numprocs)[i]));
1371: for (j = 0; j < (*numprocs)[i]; j++) PetscCall(PetscSynchronizedPrintf(comm, "%" PetscInt_FMT " ", (*indices)[i][j]));
1372: PetscCall(PetscSynchronizedPrintf(comm, "\n"));
1373: }
1374: PetscCall(PetscSynchronizedFlush(comm, PETSC_STDOUT));
1375: } /* ----------------------------------- */
1377: /* wait on sends */
1378: if (nsends2) {
1379: PetscCall(PetscMalloc1(nsends2, &send_status));
1380: PetscCallMPI(MPI_Waitall(nsends2, send_waits, send_status));
1381: PetscCall(PetscFree(send_status));
1382: }
1384: PetscCall(PetscFree(starts3));
1385: PetscCall(PetscFree(dest));
1386: PetscCall(PetscFree(send_waits));
1388: PetscCall(PetscFree(nownedsenders));
1389: PetscCall(PetscFree(ownedsenders));
1390: PetscCall(PetscFree(starts));
1391: PetscCall(PetscFree(starts2));
1392: PetscCall(PetscFree(lens2));
1394: PetscCall(PetscFree(source));
1395: PetscCall(PetscFree(len));
1396: PetscCall(PetscFree(recvs));
1397: PetscCall(PetscFree(nprocs));
1398: PetscCall(PetscFree(sends2));
1400: /* put the information about myself as the first entry in the list */
1401: first_procs = (*procs)[0];
1402: first_numprocs = (*numprocs)[0];
1403: first_indices = (*indices)[0];
1404: for (i = 0; i < *nproc; i++) {
1405: if ((*procs)[i] == rank) {
1406: (*procs)[0] = (*procs)[i];
1407: (*numprocs)[0] = (*numprocs)[i];
1408: (*indices)[0] = (*indices)[i];
1409: (*procs)[i] = first_procs;
1410: (*numprocs)[i] = first_numprocs;
1411: (*indices)[i] = first_indices;
1412: break;
1413: }
1414: }
1416: /* save info for reuse */
1417: mapping->info_nproc = *nproc;
1418: mapping->info_procs = *procs;
1419: mapping->info_numprocs = *numprocs;
1420: mapping->info_indices = *indices;
1421: mapping->info_cached = PETSC_TRUE;
1422: PetscFunctionReturn(PETSC_SUCCESS);
1423: }
1425: /*@C
1426: ISLocalToGlobalMappingRestoreBlockInfo - Frees the memory allocated by `ISLocalToGlobalMappingGetBlockInfo()`
1428: Collective
1430: Input Parameter:
1431: . mapping - the mapping from local to global indexing
1433: Output Parameters:
1434: + nproc - number of processors that are connected to this one
1435: . proc - neighboring processors
1436: . numproc - number of indices for each processor
1437: - indices - indices of local nodes shared with neighbor (sorted by global numbering)
1439: Level: advanced
1441: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1442: `ISLocalToGlobalMappingGetInfo()`
1443: @*/
1444: PetscErrorCode ISLocalToGlobalMappingRestoreBlockInfo(ISLocalToGlobalMapping mapping, PetscInt *nproc, PetscInt *procs[], PetscInt *numprocs[], PetscInt **indices[])
1445: {
1446: PetscFunctionBegin;
1448: if (mapping->info_free) {
1449: PetscCall(PetscFree(*numprocs));
1450: if (*indices) {
1451: PetscInt i;
1453: PetscCall(PetscFree((*indices)[0]));
1454: for (i = 1; i < *nproc; i++) PetscCall(PetscFree((*indices)[i]));
1455: PetscCall(PetscFree(*indices));
1456: }
1457: }
1458: *nproc = 0;
1459: *procs = NULL;
1460: *numprocs = NULL;
1461: *indices = NULL;
1462: PetscFunctionReturn(PETSC_SUCCESS);
1463: }
1465: /*@C
1466: ISLocalToGlobalMappingGetInfo - Gets the neighbor information for each processor and
1467: each index shared by more than one processor
1469: Collective
1471: Input Parameter:
1472: . mapping - the mapping from local to global indexing
1474: Output Parameters:
1475: + nproc - number of processors that are connected to this one
1476: . proc - neighboring processors
1477: . numproc - number of indices for each subdomain (processor)
1478: - indices - indices of nodes (in local numbering) shared with neighbors (sorted by global numbering)
1480: Level: advanced
1482: Note:
1483: The user needs to call `ISLocalToGlobalMappingRestoreInfo()` when the data is no longer needed.
1485: Fortran Usage:
1486: .vb
1487: PetscInt indices[nproc][numprocmax],ierr)
1488: ISLocalToGlobalMpngGetInfoSize(ISLocalToGlobalMapping,PetscInt nproc,PetscInt numprocmax,ierr) followed by
1489: ISLocalToGlobalMappingGetInfo(ISLocalToGlobalMapping,PetscInt nproc, PetscInt procs[nproc],PetscInt numprocs[nproc],
1490: .ve
1492: Fortran Note:
1493: There is no `ISLocalToGlobalMappingRestoreInfo()` in Fortran. You must make sure that `procs`[], `numprocs`[] and
1494: `indices`[][] are large enough arrays, either by allocating them dynamically or defining static ones large enough.
1496: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1497: `ISLocalToGlobalMappingRestoreInfo()`
1498: @*/
1499: PetscErrorCode ISLocalToGlobalMappingGetInfo(ISLocalToGlobalMapping mapping, PetscInt *nproc, PetscInt *procs[], PetscInt *numprocs[], PetscInt **indices[])
1500: {
1501: PetscInt **bindices = NULL, *bnumprocs = NULL, bs, i, j, k;
1503: PetscFunctionBegin;
1505: bs = mapping->bs;
1506: PetscCall(ISLocalToGlobalMappingGetBlockInfo(mapping, nproc, procs, &bnumprocs, &bindices));
1507: if (bs > 1) { /* we need to expand the cached info */
1508: PetscCall(PetscCalloc1(*nproc, &*indices));
1509: PetscCall(PetscCalloc1(*nproc, &*numprocs));
1510: for (i = 0; i < *nproc; i++) {
1511: PetscCall(PetscMalloc1(bs * bnumprocs[i], &(*indices)[i]));
1512: for (j = 0; j < bnumprocs[i]; j++) {
1513: for (k = 0; k < bs; k++) (*indices)[i][j * bs + k] = bs * bindices[i][j] + k;
1514: }
1515: (*numprocs)[i] = bnumprocs[i] * bs;
1516: }
1517: mapping->info_free = PETSC_TRUE;
1518: } else {
1519: *numprocs = bnumprocs;
1520: *indices = bindices;
1521: }
1522: PetscFunctionReturn(PETSC_SUCCESS);
1523: }
1525: /*@C
1526: ISLocalToGlobalMappingRestoreInfo - Frees the memory allocated by `ISLocalToGlobalMappingGetInfo()`
1528: Collective
1530: Input Parameter:
1531: . mapping - the mapping from local to global indexing
1533: Output Parameters:
1534: + nproc - number of processors that are connected to this one
1535: . proc - neighboring processors
1536: . numproc - number of indices for each processor
1537: - indices - indices of local nodes shared with neighbor (sorted by global numbering)
1539: Level: advanced
1541: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1542: `ISLocalToGlobalMappingGetInfo()`
1543: @*/
1544: PetscErrorCode ISLocalToGlobalMappingRestoreInfo(ISLocalToGlobalMapping mapping, PetscInt *nproc, PetscInt *procs[], PetscInt *numprocs[], PetscInt **indices[])
1545: {
1546: PetscFunctionBegin;
1547: PetscCall(ISLocalToGlobalMappingRestoreBlockInfo(mapping, nproc, procs, numprocs, indices));
1548: PetscFunctionReturn(PETSC_SUCCESS);
1549: }
1551: /*@C
1552: ISLocalToGlobalMappingGetNodeInfo - Gets the neighbor information for each MPI rank
1554: Collective
1556: Input Parameter:
1557: . mapping - the mapping from local to global indexing
1559: Output Parameters:
1560: + nnodes - number of local nodes (same `ISLocalToGlobalMappingGetSize()`)
1561: . count - number of neighboring processors per node
1562: - indices - indices of processes sharing the node (sorted)
1564: Level: advanced
1566: Note:
1567: The user needs to call `ISLocalToGlobalMappingRestoreInfo()` when the data is no longer needed.
1569: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1570: `ISLocalToGlobalMappingGetInfo()`, `ISLocalToGlobalMappingRestoreNodeInfo()`
1571: @*/
1572: PetscErrorCode ISLocalToGlobalMappingGetNodeInfo(ISLocalToGlobalMapping mapping, PetscInt *nnodes, PetscInt *count[], PetscInt **indices[])
1573: {
1574: PetscInt n;
1576: PetscFunctionBegin;
1578: PetscCall(ISLocalToGlobalMappingGetSize(mapping, &n));
1579: if (!mapping->info_nodec) {
1580: PetscInt i, m, n_neigh, *neigh, *n_shared, **shared;
1582: PetscCall(PetscMalloc2(n + 1, &mapping->info_nodec, n, &mapping->info_nodei));
1583: PetscCall(ISLocalToGlobalMappingGetInfo(mapping, &n_neigh, &neigh, &n_shared, &shared));
1584: for (i = 0; i < n; i++) mapping->info_nodec[i] = 1;
1585: m = n;
1586: mapping->info_nodec[n] = 0;
1587: for (i = 1; i < n_neigh; i++) {
1588: PetscInt j;
1590: m += n_shared[i];
1591: for (j = 0; j < n_shared[i]; j++) mapping->info_nodec[shared[i][j]] += 1;
1592: }
1593: if (n) PetscCall(PetscMalloc1(m, &mapping->info_nodei[0]));
1594: for (i = 1; i < n; i++) mapping->info_nodei[i] = mapping->info_nodei[i - 1] + mapping->info_nodec[i - 1];
1595: PetscCall(PetscArrayzero(mapping->info_nodec, n));
1596: for (i = 0; i < n; i++) {
1597: mapping->info_nodec[i] = 1;
1598: mapping->info_nodei[i][0] = neigh[0];
1599: }
1600: for (i = 1; i < n_neigh; i++) {
1601: PetscInt j;
1603: for (j = 0; j < n_shared[i]; j++) {
1604: PetscInt k = shared[i][j];
1606: mapping->info_nodei[k][mapping->info_nodec[k]] = neigh[i];
1607: mapping->info_nodec[k] += 1;
1608: }
1609: }
1610: for (i = 0; i < n; i++) PetscCall(PetscSortRemoveDupsInt(&mapping->info_nodec[i], mapping->info_nodei[i]));
1611: PetscCall(ISLocalToGlobalMappingRestoreInfo(mapping, &n_neigh, &neigh, &n_shared, &shared));
1612: }
1613: if (nnodes) *nnodes = n;
1614: if (count) *count = mapping->info_nodec;
1615: if (indices) *indices = mapping->info_nodei;
1616: PetscFunctionReturn(PETSC_SUCCESS);
1617: }
1619: /*@C
1620: ISLocalToGlobalMappingRestoreNodeInfo - Frees the memory allocated by `ISLocalToGlobalMappingGetNodeInfo()`
1622: Collective
1624: Input Parameter:
1625: . mapping - the mapping from local to global indexing
1627: Output Parameters:
1628: + nnodes - number of local nodes
1629: . count - number of neighboring processors per node
1630: - indices - indices of processes sharing the node (sorted)
1632: Level: advanced
1634: .seealso: [](sec_scatter), `ISLocalToGlobalMappingDestroy()`, `ISLocalToGlobalMappingCreateIS()`, `ISLocalToGlobalMappingCreate()`,
1635: `ISLocalToGlobalMappingGetInfo()`
1636: @*/
1637: PetscErrorCode ISLocalToGlobalMappingRestoreNodeInfo(ISLocalToGlobalMapping mapping, PetscInt *nnodes, PetscInt *count[], PetscInt **indices[])
1638: {
1639: PetscFunctionBegin;
1641: if (nnodes) *nnodes = 0;
1642: if (count) *count = NULL;
1643: if (indices) *indices = NULL;
1644: PetscFunctionReturn(PETSC_SUCCESS);
1645: }
1647: /*MC
1648: ISLocalToGlobalMappingGetIndicesF90 - Get global indices for every local point that is mapped
1650: Synopsis:
1651: ISLocalToGlobalMappingGetIndicesF90(ISLocalToGlobalMapping ltog, PetscInt, pointer :: array(:)}, integer ierr)
1653: Not Collective
1655: Input Parameter:
1656: . A - the matrix
1658: Output Parameter:
1659: . array - array of indices, the length of this array may be obtained with `ISLocalToGlobalMappingGetSize()`
1661: Level: advanced
1663: Note:
1664: Use `ISLocalToGlobalMappingGetIndicesF90()` when you no longer need access to the data
1666: .seealso: [](sec_fortranarrays), [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingGetIndices()`, `ISLocalToGlobalMappingRestoreIndices()`,
1667: `ISLocalToGlobalMappingRestoreIndicesF90()`
1668: M*/
1670: /*MC
1671: ISLocalToGlobalMappingRestoreIndicesF90 - restores the global indices for every local point that is mapped obtained with `ISLocalToGlobalMappingGetIndicesF90()`
1673: Synopsis:
1674: ISLocalToGlobalMappingRestoreIndicesF90(ISLocalToGlobalMapping ltog, PetscInt, pointer :: array(:)}, integer ierr)
1676: Not Collective
1678: Input Parameters:
1679: + A - the matrix
1680: - array - array of indices, the length of this array may be obtained with `ISLocalToGlobalMappingGetSize()`
1682: Level: advanced
1684: .seealso: [](sec_fortranarrays), [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingGetIndices()`, `ISLocalToGlobalMappingRestoreIndices()`,
1685: `ISLocalToGlobalMappingGetIndicesF90()`
1686: M*/
1688: /*@C
1689: ISLocalToGlobalMappingGetIndices - Get global indices for every local point that is mapped
1691: Not Collective
1693: Input Parameter:
1694: . ltog - local to global mapping
1696: Output Parameter:
1697: . array - array of indices, the length of this array may be obtained with `ISLocalToGlobalMappingGetSize()`
1699: Level: advanced
1701: Note:
1702: `ISLocalToGlobalMappingGetSize()` returns the length the this array
1704: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingRestoreIndices()`,
1705: `ISLocalToGlobalMappingGetBlockIndices()`, `ISLocalToGlobalMappingRestoreBlockIndices()`
1706: @*/
1707: PetscErrorCode ISLocalToGlobalMappingGetIndices(ISLocalToGlobalMapping ltog, const PetscInt **array)
1708: {
1709: PetscFunctionBegin;
1712: if (ltog->bs == 1) {
1713: *array = ltog->indices;
1714: } else {
1715: PetscInt *jj, k, i, j, n = ltog->n, bs = ltog->bs;
1716: const PetscInt *ii;
1718: PetscCall(PetscMalloc1(bs * n, &jj));
1719: *array = jj;
1720: k = 0;
1721: ii = ltog->indices;
1722: for (i = 0; i < n; i++)
1723: for (j = 0; j < bs; j++) jj[k++] = bs * ii[i] + j;
1724: }
1725: PetscFunctionReturn(PETSC_SUCCESS);
1726: }
1728: /*@C
1729: ISLocalToGlobalMappingRestoreIndices - Restore indices obtained with `ISLocalToGlobalMappingGetIndices()`
1731: Not Collective
1733: Input Parameters:
1734: + ltog - local to global mapping
1735: - array - array of indices
1737: Level: advanced
1739: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingGetIndices()`
1740: @*/
1741: PetscErrorCode ISLocalToGlobalMappingRestoreIndices(ISLocalToGlobalMapping ltog, const PetscInt **array)
1742: {
1743: PetscFunctionBegin;
1746: PetscCheck(ltog->bs != 1 || *array == ltog->indices, PETSC_COMM_SELF, PETSC_ERR_ARG_BADPTR, "Trying to return mismatched pointer");
1748: if (ltog->bs > 1) PetscCall(PetscFree(*(void **)array));
1749: PetscFunctionReturn(PETSC_SUCCESS);
1750: }
1752: /*@C
1753: ISLocalToGlobalMappingGetBlockIndices - Get global indices for every local block
1755: Not Collective
1757: Input Parameter:
1758: . ltog - local to global mapping
1760: Output Parameter:
1761: . array - array of indices
1763: Level: advanced
1765: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingRestoreBlockIndices()`
1766: @*/
1767: PetscErrorCode ISLocalToGlobalMappingGetBlockIndices(ISLocalToGlobalMapping ltog, const PetscInt **array)
1768: {
1769: PetscFunctionBegin;
1772: *array = ltog->indices;
1773: PetscFunctionReturn(PETSC_SUCCESS);
1774: }
1776: /*@C
1777: ISLocalToGlobalMappingRestoreBlockIndices - Restore indices obtained with `ISLocalToGlobalMappingGetBlockIndices()`
1779: Not Collective
1781: Input Parameters:
1782: + ltog - local to global mapping
1783: - array - array of indices
1785: Level: advanced
1787: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingApply()`, `ISLocalToGlobalMappingGetIndices()`
1788: @*/
1789: PetscErrorCode ISLocalToGlobalMappingRestoreBlockIndices(ISLocalToGlobalMapping ltog, const PetscInt **array)
1790: {
1791: PetscFunctionBegin;
1794: PetscCheck(*array == ltog->indices, PETSC_COMM_SELF, PETSC_ERR_ARG_BADPTR, "Trying to return mismatched pointer");
1795: *array = NULL;
1796: PetscFunctionReturn(PETSC_SUCCESS);
1797: }
1799: /*@C
1800: ISLocalToGlobalMappingConcatenate - Create a new mapping that concatenates a list of mappings
1802: Not Collective
1804: Input Parameters:
1805: + comm - communicator for the new mapping, must contain the communicator of every mapping to concatenate
1806: . n - number of mappings to concatenate
1807: - ltogs - local to global mappings
1809: Output Parameter:
1810: . ltogcat - new mapping
1812: Level: advanced
1814: Note:
1815: This currently always returns a mapping with block size of 1
1817: Developer Note:
1818: If all the input mapping have the same block size we could easily handle that as a special case
1820: .seealso: [](sec_scatter), `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingCreate()`
1821: @*/
1822: PetscErrorCode ISLocalToGlobalMappingConcatenate(MPI_Comm comm, PetscInt n, const ISLocalToGlobalMapping ltogs[], ISLocalToGlobalMapping *ltogcat)
1823: {
1824: PetscInt i, cnt, m, *idx;
1826: PetscFunctionBegin;
1827: PetscCheck(n >= 0, comm, PETSC_ERR_ARG_OUTOFRANGE, "Must have a non-negative number of mappings, given %" PetscInt_FMT, n);
1831: for (cnt = 0, i = 0; i < n; i++) {
1832: PetscCall(ISLocalToGlobalMappingGetSize(ltogs[i], &m));
1833: cnt += m;
1834: }
1835: PetscCall(PetscMalloc1(cnt, &idx));
1836: for (cnt = 0, i = 0; i < n; i++) {
1837: const PetscInt *subidx;
1838: PetscCall(ISLocalToGlobalMappingGetSize(ltogs[i], &m));
1839: PetscCall(ISLocalToGlobalMappingGetIndices(ltogs[i], &subidx));
1840: PetscCall(PetscArraycpy(&idx[cnt], subidx, m));
1841: PetscCall(ISLocalToGlobalMappingRestoreIndices(ltogs[i], &subidx));
1842: cnt += m;
1843: }
1844: PetscCall(ISLocalToGlobalMappingCreate(comm, 1, cnt, idx, PETSC_OWN_POINTER, ltogcat));
1845: PetscFunctionReturn(PETSC_SUCCESS);
1846: }
1848: /*MC
1849: ISLOCALTOGLOBALMAPPINGBASIC - basic implementation of the `ISLocalToGlobalMapping` object. When `ISGlobalToLocalMappingApply()` is
1850: used this is good for only small and moderate size problems.
1852: Options Database Key:
1853: . -islocaltoglobalmapping_type basic - select this method
1855: Level: beginner
1857: Developer Note:
1858: This stores all the mapping information on each MPI rank.
1860: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingSetType()`, `ISLOCALTOGLOBALMAPPINGHASH`
1861: M*/
1862: PETSC_EXTERN PetscErrorCode ISLocalToGlobalMappingCreate_Basic(ISLocalToGlobalMapping ltog)
1863: {
1864: PetscFunctionBegin;
1865: ltog->ops->globaltolocalmappingapply = ISGlobalToLocalMappingApply_Basic;
1866: ltog->ops->globaltolocalmappingsetup = ISGlobalToLocalMappingSetUp_Basic;
1867: ltog->ops->globaltolocalmappingapplyblock = ISGlobalToLocalMappingApplyBlock_Basic;
1868: ltog->ops->destroy = ISLocalToGlobalMappingDestroy_Basic;
1869: PetscFunctionReturn(PETSC_SUCCESS);
1870: }
1872: /*MC
1873: ISLOCALTOGLOBALMAPPINGHASH - hash implementation of the `ISLocalToGlobalMapping` object. When `ISGlobalToLocalMappingApply()` is
1874: used this is good for large memory problems.
1876: Options Database Key:
1877: . -islocaltoglobalmapping_type hash - select this method
1879: Level: beginner
1881: Note:
1882: This is selected automatically for large problems if the user does not set the type.
1884: .seealso: [](sec_scatter), `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingSetType()`, `ISLOCALTOGLOBALMAPPINGBASIC`
1885: M*/
1886: PETSC_EXTERN PetscErrorCode ISLocalToGlobalMappingCreate_Hash(ISLocalToGlobalMapping ltog)
1887: {
1888: PetscFunctionBegin;
1889: ltog->ops->globaltolocalmappingapply = ISGlobalToLocalMappingApply_Hash;
1890: ltog->ops->globaltolocalmappingsetup = ISGlobalToLocalMappingSetUp_Hash;
1891: ltog->ops->globaltolocalmappingapplyblock = ISGlobalToLocalMappingApplyBlock_Hash;
1892: ltog->ops->destroy = ISLocalToGlobalMappingDestroy_Hash;
1893: PetscFunctionReturn(PETSC_SUCCESS);
1894: }
1896: /*@C
1897: ISLocalToGlobalMappingRegister - Registers a method for applying a global to local mapping with an `ISLocalToGlobalMapping`
1899: Not Collective
1901: Input Parameters:
1902: + sname - name of a new method
1903: - function - routine to create method context
1905: Sample usage:
1906: .vb
1907: ISLocalToGlobalMappingRegister("my_mapper", MyCreate);
1908: .ve
1910: Then, your mapping can be chosen with the procedural interface via
1911: $ ISLocalToGlobalMappingSetType(ltog, "my_mapper")
1912: or at runtime via the option
1913: $ -islocaltoglobalmapping_type my_mapper
1915: Level: advanced
1917: Note:
1918: `ISLocalToGlobalMappingRegister()` may be called multiple times to add several user-defined mappings.
1920: .seealso: [](sec_scatter), `ISLocalToGlobalMappingRegisterAll()`, `ISLocalToGlobalMappingRegisterDestroy()`, `ISLOCALTOGLOBALMAPPINGBASIC`,
1921: `ISLOCALTOGLOBALMAPPINGHASH`, `ISLocalToGlobalMapping`, `ISLocalToGlobalMappingApply()`
1922: @*/
1923: PetscErrorCode ISLocalToGlobalMappingRegister(const char sname[], PetscErrorCode (*function)(ISLocalToGlobalMapping))
1924: {
1925: PetscFunctionBegin;
1926: PetscCall(ISInitializePackage());
1927: PetscCall(PetscFunctionListAdd(&ISLocalToGlobalMappingList, sname, function));
1928: PetscFunctionReturn(PETSC_SUCCESS);
1929: }
1931: /*@C
1932: ISLocalToGlobalMappingSetType - Sets the implementation type `ISLocalToGlobalMapping` will use
1934: Logically Collective
1936: Input Parameters:
1937: + ltog - the `ISLocalToGlobalMapping` object
1938: - type - a known method
1940: Options Database Key:
1941: . -islocaltoglobalmapping_type <method> - Sets the method; use -help for a list of available methods (for instance, basic or hash)
1943: Level: intermediate
1945: Notes:
1946: See `ISLocalToGlobalMappingType` for available methods
1948: Normally, it is best to use the `ISLocalToGlobalMappingSetFromOptions()` command and
1949: then set the `ISLocalToGlobalMappingType` from the options database rather than by using
1950: this routine.
1952: Developer Note:
1953: `ISLocalToGlobalMappingRegister()` is used to add new types to `ISLocalToGlobalMappingList` from which they
1954: are accessed by `ISLocalToGlobalMappingSetType()`.
1956: .seealso: [](sec_scatter), `ISLocalToGlobalMappingType`, `ISLocalToGlobalMappingType`, `ISLocalToGlobalMappingRegister()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingGetType()`
1957: @*/
1958: PetscErrorCode ISLocalToGlobalMappingSetType(ISLocalToGlobalMapping ltog, ISLocalToGlobalMappingType type)
1959: {
1960: PetscBool match;
1961: PetscErrorCode (*r)(ISLocalToGlobalMapping) = NULL;
1963: PetscFunctionBegin;
1967: PetscCall(PetscObjectTypeCompare((PetscObject)ltog, type, &match));
1968: if (match) PetscFunctionReturn(PETSC_SUCCESS);
1970: /* L2G maps defer type setup at globaltolocal calls, allow passing NULL here */
1971: if (type) {
1972: PetscCall(PetscFunctionListFind(ISLocalToGlobalMappingList, type, &r));
1973: PetscCheck(r, PetscObjectComm((PetscObject)ltog), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unable to find requested ISLocalToGlobalMapping type %s", type);
1974: }
1975: /* Destroy the previous private LTOG context */
1976: PetscTryTypeMethod(ltog, destroy);
1977: ltog->ops->destroy = NULL;
1979: PetscCall(PetscObjectChangeTypeName((PetscObject)ltog, type));
1980: if (r) PetscCall((*r)(ltog));
1981: PetscFunctionReturn(PETSC_SUCCESS);
1982: }
1984: /*@C
1985: ISLocalToGlobalMappingGetType - Get the type of the `ISLocalToGlobalMapping`
1987: Not Collective
1989: Input Parameter:
1990: . ltog - the `ISLocalToGlobalMapping` object
1992: Output Parameter:
1993: . type - the type
1995: Level: intermediate
1997: .seealso: [](sec_scatter), `ISLocalToGlobalMappingType`, `ISLocalToGlobalMappingType`, `ISLocalToGlobalMappingRegister()`, `ISLocalToGlobalMappingCreate()`, `ISLocalToGlobalMappingSetType()`
1998: @*/
1999: PetscErrorCode ISLocalToGlobalMappingGetType(ISLocalToGlobalMapping ltog, ISLocalToGlobalMappingType *type)
2000: {
2001: PetscFunctionBegin;
2004: *type = ((PetscObject)ltog)->type_name;
2005: PetscFunctionReturn(PETSC_SUCCESS);
2006: }
2008: PetscBool ISLocalToGlobalMappingRegisterAllCalled = PETSC_FALSE;
2010: /*@C
2011: ISLocalToGlobalMappingRegisterAll - Registers all of the local to global mapping components in the `IS` package.
2013: Not Collective
2015: Level: advanced
2017: .seealso: [](sec_scatter), `ISRegister()`, `ISLocalToGlobalRegister()`
2018: @*/
2019: PetscErrorCode ISLocalToGlobalMappingRegisterAll(void)
2020: {
2021: PetscFunctionBegin;
2022: if (ISLocalToGlobalMappingRegisterAllCalled) PetscFunctionReturn(PETSC_SUCCESS);
2023: ISLocalToGlobalMappingRegisterAllCalled = PETSC_TRUE;
2024: PetscCall(ISLocalToGlobalMappingRegister(ISLOCALTOGLOBALMAPPINGBASIC, ISLocalToGlobalMappingCreate_Basic));
2025: PetscCall(ISLocalToGlobalMappingRegister(ISLOCALTOGLOBALMAPPINGHASH, ISLocalToGlobalMappingCreate_Hash));
2026: PetscFunctionReturn(PETSC_SUCCESS);
2027: }