Actual source code: general.c
2: /*
3: Provides the functions for index sets (IS) defined by a list of integers.
4: */
5: #include <../src/vec/is/is/impls/general/general.h>
6: #include <petsc/private/viewerhdf5impl.h>
8: static PetscErrorCode ISDuplicate_General(IS is, IS *newIS)
9: {
10: IS_General *sub = (IS_General *)is->data;
11: PetscInt n;
13: PetscLayoutGetLocalSize(is->map, &n);
14: ISCreateGeneral(PetscObjectComm((PetscObject)is), n, sub->idx, PETSC_COPY_VALUES, newIS);
15: return 0;
16: }
18: static PetscErrorCode ISDestroy_General(IS is)
19: {
20: IS_General *is_general = (IS_General *)is->data;
22: if (is_general->allocated) PetscFree(is_general->idx);
23: PetscObjectComposeFunction((PetscObject)is, "ISGeneralSetIndices_C", NULL);
24: PetscObjectComposeFunction((PetscObject)is, "ISGeneralFilter_C", NULL);
25: PetscObjectComposeFunction((PetscObject)is, "ISGeneralSetIndicesFromMask_C", NULL);
26: PetscObjectComposeFunction((PetscObject)is, "ISShift_C", NULL);
27: PetscFree(is->data);
28: return 0;
29: }
31: static PetscErrorCode ISCopy_General(IS is, IS isy)
32: {
33: IS_General *is_general = (IS_General *)is->data, *isy_general = (IS_General *)isy->data;
34: PetscInt n;
36: PetscLayoutGetLocalSize(is->map, &n);
37: PetscArraycpy(isy_general->idx, is_general->idx, n);
38: return 0;
39: }
41: PetscErrorCode ISShift_General(IS is, PetscInt shift, IS isy)
42: {
43: IS_General *is_general = (IS_General *)is->data, *isy_general = (IS_General *)isy->data;
44: PetscInt i, n;
46: PetscLayoutGetLocalSize(is->map, &n);
47: for (i = 0; i < n; i++) isy_general->idx[i] = is_general->idx[i] + shift;
48: return 0;
49: }
51: static PetscErrorCode ISOnComm_General(IS is, MPI_Comm comm, PetscCopyMode mode, IS *newis)
52: {
53: IS_General *sub = (IS_General *)is->data;
54: PetscInt n;
57: PetscLayoutGetLocalSize(is->map, &n);
58: ISCreateGeneral(comm, n, sub->idx, mode, newis);
59: return 0;
60: }
62: static PetscErrorCode ISSetBlockSize_General(IS is, PetscInt bs)
63: {
64: PetscLayoutSetBlockSize(is->map, bs);
65: return 0;
66: }
68: static PetscErrorCode ISContiguousLocal_General(IS is, PetscInt gstart, PetscInt gend, PetscInt *start, PetscBool *contig)
69: {
70: IS_General *sub = (IS_General *)is->data;
71: PetscInt n, i, p;
73: *start = 0;
74: *contig = PETSC_TRUE;
75: PetscLayoutGetLocalSize(is->map, &n);
76: if (!n) return 0;
77: p = sub->idx[0];
78: if (p < gstart) goto nomatch;
79: *start = p - gstart;
80: if (n > gend - p) goto nomatch;
81: for (i = 1; i < n; i++, p++) {
82: if (sub->idx[i] != p + 1) goto nomatch;
83: }
84: return 0;
85: nomatch:
86: *start = -1;
87: *contig = PETSC_FALSE;
88: return 0;
89: }
91: static PetscErrorCode ISLocate_General(IS is, PetscInt key, PetscInt *location)
92: {
93: IS_General *sub = (IS_General *)is->data;
94: PetscInt numIdx, i;
95: PetscBool sorted;
97: PetscLayoutGetLocalSize(is->map, &numIdx);
98: ISGetInfo(is, IS_SORTED, IS_LOCAL, PETSC_TRUE, &sorted);
99: if (sorted) PetscFindInt(key, numIdx, sub->idx, location);
100: else {
101: const PetscInt *idx = sub->idx;
103: *location = -1;
104: for (i = 0; i < numIdx; i++) {
105: if (idx[i] == key) {
106: *location = i;
107: return 0;
108: }
109: }
110: }
111: return 0;
112: }
114: static PetscErrorCode ISGetIndices_General(IS in, const PetscInt *idx[])
115: {
116: IS_General *sub = (IS_General *)in->data;
118: *idx = sub->idx;
119: return 0;
120: }
122: static PetscErrorCode ISRestoreIndices_General(IS in, const PetscInt *idx[])
123: {
124: IS_General *sub = (IS_General *)in->data;
126: /* F90Array1dCreate() inside ISRestoreArrayF90() does not keep array when zero length array */
128: return 0;
129: }
131: static PetscErrorCode ISInvertPermutation_General(IS is, PetscInt nlocal, IS *isout)
132: {
133: IS_General *sub = (IS_General *)is->data;
134: PetscInt i, *ii, n, nstart;
135: const PetscInt *idx = sub->idx;
136: PetscMPIInt size;
137: IS istmp, nistmp;
139: PetscLayoutGetLocalSize(is->map, &n);
140: MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
141: if (size == 1) {
142: PetscMalloc1(n, &ii);
143: for (i = 0; i < n; i++) ii[idx[i]] = i;
144: ISCreateGeneral(PETSC_COMM_SELF, n, ii, PETSC_OWN_POINTER, isout);
145: ISSetPermutation(*isout);
146: } else {
147: /* crude, nonscalable get entire IS on each processor */
148: ISAllGather(is, &istmp);
149: ISSetPermutation(istmp);
150: ISInvertPermutation(istmp, PETSC_DECIDE, &nistmp);
151: ISDestroy(&istmp);
152: /* get the part we need */
153: if (nlocal == PETSC_DECIDE) nlocal = n;
154: MPI_Scan(&nlocal, &nstart, 1, MPIU_INT, MPI_SUM, PetscObjectComm((PetscObject)is));
155: if (PetscDefined(USE_DEBUG)) {
156: PetscInt N;
157: PetscMPIInt rank;
158: MPI_Comm_rank(PetscObjectComm((PetscObject)is), &rank);
159: PetscLayoutGetSize(is->map, &N);
161: }
162: nstart -= nlocal;
163: ISGetIndices(nistmp, &idx);
164: ISCreateGeneral(PetscObjectComm((PetscObject)is), nlocal, idx + nstart, PETSC_COPY_VALUES, isout);
165: ISRestoreIndices(nistmp, &idx);
166: ISDestroy(&nistmp);
167: }
168: return 0;
169: }
171: #if defined(PETSC_HAVE_HDF5)
172: static PetscErrorCode ISView_General_HDF5(IS is, PetscViewer viewer)
173: {
174: PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *)viewer->data;
175: hid_t filespace; /* file dataspace identifier */
176: hid_t chunkspace; /* chunk dataset property identifier */
177: hid_t dset_id; /* dataset identifier */
178: hid_t memspace; /* memory dataspace identifier */
179: hid_t inttype; /* int type (H5T_NATIVE_INT or H5T_NATIVE_LLONG) */
180: hid_t file_id, group;
181: hsize_t dim, maxDims[3], dims[3], chunkDims[3], count[3], offset[3];
182: PetscBool timestepping;
183: PetscInt bs, N, n, timestep = PETSC_MIN_INT, low;
184: hsize_t chunksize;
185: const PetscInt *ind;
186: const char *isname;
188: ISGetBlockSize(is, &bs);
189: bs = PetscMax(bs, 1); /* If N = 0, bs = 0 as well */
190: PetscViewerHDF5OpenGroup(viewer, NULL, &file_id, &group);
191: PetscViewerHDF5IsTimestepping(viewer, ×tepping);
192: if (timestepping) PetscViewerHDF5GetTimestep(viewer, ×tep);
194: /* Create the dataspace for the dataset.
195: *
196: * dims - holds the current dimensions of the dataset
197: *
198: * maxDims - holds the maximum dimensions of the dataset (unlimited
199: * for the number of time steps with the current dimensions for the
200: * other dimensions; so only additional time steps can be added).
201: *
202: * chunkDims - holds the size of a single time step (required to
203: * permit extending dataset).
204: */
205: dim = 0;
206: chunksize = 1;
207: if (timestep >= 0) {
208: dims[dim] = timestep + 1;
209: maxDims[dim] = H5S_UNLIMITED;
210: chunkDims[dim] = 1;
211: ++dim;
212: }
213: ISGetSize(is, &N);
214: ISGetLocalSize(is, &n);
215: PetscHDF5IntCast(N / bs, dims + dim);
217: maxDims[dim] = dims[dim];
218: chunkDims[dim] = PetscMax(1, dims[dim]);
219: chunksize *= chunkDims[dim];
220: ++dim;
221: if (bs >= 1) {
222: dims[dim] = bs;
223: maxDims[dim] = dims[dim];
224: chunkDims[dim] = dims[dim];
225: chunksize *= chunkDims[dim];
226: ++dim;
227: }
228: /* hdf5 chunks must be less than 4GB */
229: if (chunksize > PETSC_HDF5_MAX_CHUNKSIZE / 64) {
230: if (bs >= 1) {
231: if (chunkDims[dim - 2] > (hsize_t)PetscSqrtReal((PetscReal)(PETSC_HDF5_MAX_CHUNKSIZE / 64))) chunkDims[dim - 2] = (hsize_t)PetscSqrtReal((PetscReal)(PETSC_HDF5_MAX_CHUNKSIZE / 64));
232: if (chunkDims[dim - 1] > (hsize_t)PetscSqrtReal((PetscReal)(PETSC_HDF5_MAX_CHUNKSIZE / 64))) chunkDims[dim - 1] = (hsize_t)PetscSqrtReal((PetscReal)(PETSC_HDF5_MAX_CHUNKSIZE / 64));
233: } else {
234: chunkDims[dim - 1] = PETSC_HDF5_MAX_CHUNKSIZE / 64;
235: }
236: }
237: PetscCallHDF5Return(filespace, H5Screate_simple, (dim, dims, maxDims));
239: #if defined(PETSC_USE_64BIT_INDICES)
240: inttype = H5T_NATIVE_LLONG;
241: #else
242: inttype = H5T_NATIVE_INT;
243: #endif
245: /* Create the dataset with default properties and close filespace */
246: PetscObjectGetName((PetscObject)is, &isname);
247: if (!H5Lexists(group, isname, H5P_DEFAULT)) {
248: /* Create chunk */
249: PetscCallHDF5Return(chunkspace, H5Pcreate, (H5P_DATASET_CREATE));
250: PetscCallHDF5(H5Pset_chunk, (chunkspace, dim, chunkDims));
252: PetscCallHDF5Return(dset_id, H5Dcreate2, (group, isname, inttype, filespace, H5P_DEFAULT, chunkspace, H5P_DEFAULT));
253: PetscCallHDF5(H5Pclose, (chunkspace));
254: } else {
255: PetscCallHDF5Return(dset_id, H5Dopen2, (group, isname, H5P_DEFAULT));
256: PetscCallHDF5(H5Dset_extent, (dset_id, dims));
257: }
258: PetscCallHDF5(H5Sclose, (filespace));
260: /* Each process defines a dataset and writes it to the hyperslab in the file */
261: dim = 0;
262: if (timestep >= 0) {
263: count[dim] = 1;
264: ++dim;
265: }
266: PetscHDF5IntCast(n / bs, count + dim);
267: ++dim;
268: if (bs >= 1) {
269: count[dim] = bs;
270: ++dim;
271: }
272: if (n > 0 || H5_VERSION_GE(1, 10, 0)) {
273: PetscCallHDF5Return(memspace, H5Screate_simple, (dim, count, NULL));
274: } else {
275: /* Can't create dataspace with zero for any dimension, so create null dataspace. */
276: PetscCallHDF5Return(memspace, H5Screate, (H5S_NULL));
277: }
279: /* Select hyperslab in the file */
280: PetscLayoutGetRange(is->map, &low, NULL);
281: dim = 0;
282: if (timestep >= 0) {
283: offset[dim] = timestep;
284: ++dim;
285: }
286: PetscHDF5IntCast(low / bs, offset + dim);
287: ++dim;
288: if (bs >= 1) {
289: offset[dim] = 0;
290: ++dim;
291: }
292: if (n > 0 || H5_VERSION_GE(1, 10, 0)) {
293: PetscCallHDF5Return(filespace, H5Dget_space, (dset_id));
294: PetscCallHDF5(H5Sselect_hyperslab, (filespace, H5S_SELECT_SET, offset, NULL, count, NULL));
295: } else {
296: /* Create null filespace to match null memspace. */
297: PetscCallHDF5Return(filespace, H5Screate, (H5S_NULL));
298: }
300: ISGetIndices(is, &ind);
301: PetscCallHDF5(H5Dwrite, (dset_id, inttype, memspace, filespace, hdf5->dxpl_id, ind));
302: PetscCallHDF5(H5Fflush, (file_id, H5F_SCOPE_GLOBAL));
303: ISRestoreIndices(is, &ind);
305: /* Close/release resources */
306: PetscCallHDF5(H5Gclose, (group));
307: PetscCallHDF5(H5Sclose, (filespace));
308: PetscCallHDF5(H5Sclose, (memspace));
309: PetscCallHDF5(H5Dclose, (dset_id));
311: if (timestepping) PetscViewerHDF5WriteObjectAttribute(viewer, (PetscObject)is, "timestepping", PETSC_BOOL, ×tepping);
312: PetscInfo(is, "Wrote IS object with name %s\n", isname);
313: return 0;
314: }
315: #endif
317: static PetscErrorCode ISView_General(IS is, PetscViewer viewer)
318: {
319: IS_General *sub = (IS_General *)is->data;
320: PetscInt i, n, *idx = sub->idx;
321: PetscBool iascii, isbinary, ishdf5;
323: PetscLayoutGetLocalSize(is->map, &n);
324: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
325: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary);
326: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERHDF5, &ishdf5);
327: if (iascii) {
328: MPI_Comm comm;
329: PetscMPIInt rank, size;
330: PetscViewerFormat fmt;
331: PetscBool isperm;
333: PetscObjectGetComm((PetscObject)viewer, &comm);
334: MPI_Comm_rank(comm, &rank);
335: MPI_Comm_size(comm, &size);
337: PetscViewerGetFormat(viewer, &fmt);
338: ISGetInfo(is, IS_PERMUTATION, IS_LOCAL, PETSC_FALSE, &isperm);
339: if (isperm && fmt != PETSC_VIEWER_ASCII_MATLAB) PetscViewerASCIIPrintf(viewer, "Index set is permutation\n");
340: PetscViewerASCIIPushSynchronized(viewer);
341: if (size > 1) {
342: if (fmt == PETSC_VIEWER_ASCII_MATLAB) {
343: const char *name;
345: PetscObjectGetName((PetscObject)is, &name);
346: PetscViewerASCIISynchronizedPrintf(viewer, "%s_%d = [...\n", name, rank);
347: for (i = 0; i < n; i++) PetscViewerASCIISynchronizedPrintf(viewer, "%" PetscInt_FMT "\n", idx[i] + 1);
348: PetscViewerASCIISynchronizedPrintf(viewer, "];\n");
349: } else {
350: PetscInt st = 0;
352: if (fmt == PETSC_VIEWER_ASCII_INDEX) st = is->map->rstart;
353: PetscViewerASCIISynchronizedPrintf(viewer, "[%d] Number of indices in set %" PetscInt_FMT "\n", rank, n);
354: for (i = 0; i < n; i++) PetscViewerASCIISynchronizedPrintf(viewer, "[%d] %" PetscInt_FMT " %" PetscInt_FMT "\n", rank, i + st, idx[i]);
355: }
356: } else {
357: if (fmt == PETSC_VIEWER_ASCII_MATLAB) {
358: const char *name;
360: PetscObjectGetName((PetscObject)is, &name);
361: PetscViewerASCIISynchronizedPrintf(viewer, "%s = [...\n", name);
362: for (i = 0; i < n; i++) PetscViewerASCIISynchronizedPrintf(viewer, "%" PetscInt_FMT "\n", idx[i] + 1);
363: PetscViewerASCIISynchronizedPrintf(viewer, "];\n");
364: } else {
365: PetscInt st = 0;
367: if (fmt == PETSC_VIEWER_ASCII_INDEX) st = is->map->rstart;
368: PetscViewerASCIISynchronizedPrintf(viewer, "Number of indices in set %" PetscInt_FMT "\n", n);
369: for (i = 0; i < n; i++) PetscViewerASCIISynchronizedPrintf(viewer, "%" PetscInt_FMT " %" PetscInt_FMT "\n", i + st, idx[i]);
370: }
371: }
372: PetscViewerFlush(viewer);
373: PetscViewerASCIIPopSynchronized(viewer);
374: } else if (isbinary) {
375: ISView_Binary(is, viewer);
376: } else if (ishdf5) {
377: #if defined(PETSC_HAVE_HDF5)
378: ISView_General_HDF5(is, viewer);
379: #endif
380: }
381: return 0;
382: }
384: static PetscErrorCode ISSort_General(IS is)
385: {
386: IS_General *sub = (IS_General *)is->data;
387: PetscInt n;
389: PetscLayoutGetLocalSize(is->map, &n);
390: PetscIntSortSemiOrdered(n, sub->idx);
391: return 0;
392: }
394: static PetscErrorCode ISSortRemoveDups_General(IS is)
395: {
396: IS_General *sub = (IS_General *)is->data;
397: PetscLayout map;
398: PetscInt n;
399: PetscBool sorted;
401: PetscLayoutGetLocalSize(is->map, &n);
402: ISGetInfo(is, IS_SORTED, IS_LOCAL, PETSC_TRUE, &sorted);
403: if (sorted) {
404: PetscSortedRemoveDupsInt(&n, sub->idx);
405: } else {
406: PetscSortRemoveDupsInt(&n, sub->idx);
407: }
408: PetscLayoutCreateFromSizes(PetscObjectComm((PetscObject)is), n, PETSC_DECIDE, is->map->bs, &map);
409: PetscLayoutDestroy(&is->map);
410: is->map = map;
411: return 0;
412: }
414: static PetscErrorCode ISSorted_General(IS is, PetscBool *flg)
415: {
416: ISGetInfo(is, IS_SORTED, IS_LOCAL, PETSC_TRUE, flg);
417: return 0;
418: }
420: PetscErrorCode ISToGeneral_General(IS is)
421: {
422: return 0;
423: }
425: static struct _ISOps myops = {ISGetIndices_General, ISRestoreIndices_General, ISInvertPermutation_General, ISSort_General, ISSortRemoveDups_General, ISSorted_General, ISDuplicate_General, ISDestroy_General, ISView_General, ISLoad_Default, ISCopy_General, ISToGeneral_General, ISOnComm_General, ISSetBlockSize_General, ISContiguousLocal_General, ISLocate_General,
426: /* no specializations of {sorted,unique,perm,interval}{local,global}
427: * because the default checks in ISGetInfo_XXX in index.c are exactly
428: * what we would do for ISGeneral */
429: NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL};
431: PETSC_INTERN PetscErrorCode ISSetUp_General(IS);
433: PetscErrorCode ISSetUp_General(IS is)
434: {
435: IS_General *sub = (IS_General *)is->data;
436: const PetscInt *idx = sub->idx;
437: PetscInt n, i, min, max;
439: PetscLayoutGetLocalSize(is->map, &n);
441: if (n) {
442: min = max = idx[0];
443: for (i = 1; i < n; i++) {
444: if (idx[i] < min) min = idx[i];
445: if (idx[i] > max) max = idx[i];
446: }
447: is->min = min;
448: is->max = max;
449: } else {
450: is->min = PETSC_MAX_INT;
451: is->max = PETSC_MIN_INT;
452: }
453: return 0;
454: }
456: /*@
457: ISCreateGeneral - Creates a data structure for an index set containing a list of integers.
459: Collective
461: Input Parameters:
462: + comm - the MPI communicator
463: . n - the length of the index set
464: . idx - the list of integers
465: - mode - `PETSC_COPY_VALUES`, `PETSC_OWN_POINTER`, or `PETSC_USE_POINTER`; see `PetscCopyMode` for meaning of this flag.
467: Output Parameter:
468: . is - the new index set
470: Level: beginner
472: Notes:
473: When the communicator is not `MPI_COMM_SELF`, the operations on IS are NOT
474: conceptually the same as `MPI_Group` operations. The `IS` are then
475: distributed sets of indices and thus certain operations on them are
476: collective.
478: Use `ISGeneralSetIndices()` to provide indices to an already existing `IS` of `ISType` `ISGENERAL`
480: .seealso: [](sec_scatter), `IS`, `ISGENERAL`, `ISCreateStride()`, `ISCreateBlock()`, `ISAllGather()`, `PETSC_COPY_VALUES`, `PETSC_OWN_POINTER`,
481: `PETSC_USE_POINTER`, `PetscCopyMode`, `ISGeneralSetIndicesFromMask()`
482: @*/
483: PetscErrorCode ISCreateGeneral(MPI_Comm comm, PetscInt n, const PetscInt idx[], PetscCopyMode mode, IS *is)
484: {
485: ISCreate(comm, is);
486: ISSetType(*is, ISGENERAL);
487: ISGeneralSetIndices(*is, n, idx, mode);
488: return 0;
489: }
491: /*@
492: ISGeneralSetIndices - Sets the indices for an `ISGENERAL` index set
494: Logically Collective
496: Input Parameters:
497: + is - the index set
498: . n - the length of the index set
499: . idx - the list of integers
500: - mode - see `PetscCopyMode` for meaning of this flag.
502: Level: beginner
504: Note:
505: Use `ISCreateGeneral()` to create the `IS` and set its indices in a single function call
507: .seealso: [](sec_scatter), `IS`, `ISBLOCK`, `ISCreateGeneral()`, `ISGeneralSetIndicesFromMask()`, `ISBlockSetIndices()`, `ISGENERAL`, `PetscCopyMode`
508: @*/
509: PetscErrorCode ISGeneralSetIndices(IS is, PetscInt n, const PetscInt idx[], PetscCopyMode mode)
510: {
513: ISClearInfoCache(is, PETSC_FALSE);
514: PetscUseMethod(is, "ISGeneralSetIndices_C", (IS, PetscInt, const PetscInt[], PetscCopyMode), (is, n, idx, mode));
515: return 0;
516: }
518: PetscErrorCode ISGeneralSetIndices_General(IS is, PetscInt n, const PetscInt idx[], PetscCopyMode mode)
519: {
520: PetscLayout map;
521: IS_General *sub = (IS_General *)is->data;
526: PetscLayoutCreateFromSizes(PetscObjectComm((PetscObject)is), n, PETSC_DECIDE, is->map->bs, &map);
527: PetscLayoutDestroy(&is->map);
528: is->map = map;
530: if (sub->allocated) PetscFree(sub->idx);
531: if (mode == PETSC_COPY_VALUES) {
532: PetscMalloc1(n, &sub->idx);
533: PetscArraycpy(sub->idx, idx, n);
534: sub->allocated = PETSC_TRUE;
535: } else if (mode == PETSC_OWN_POINTER) {
536: sub->idx = (PetscInt *)idx;
537: sub->allocated = PETSC_TRUE;
538: } else {
539: sub->idx = (PetscInt *)idx;
540: sub->allocated = PETSC_FALSE;
541: }
543: ISSetUp_General(is);
544: ISViewFromOptions(is, NULL, "-is_view");
545: return 0;
546: }
548: /*@
549: ISGeneralSetIndicesFromMask - Sets the indices for an `ISGENERAL` index set using a boolean mask
551: Collective
553: Input Parameters:
554: + is - the index set
555: . rstart - the range start index (inclusive)
556: . rend - the range end index (exclusive)
557: - mask - the boolean mask array of length rend-rstart, indices will be set for each `PETSC_TRUE` value in the array
559: Level: beginner
561: Note:
562: The mask array may be freed by the user after this call.
564: Example:
565: .vb
566: PetscBool mask[] = {PETSC_FALSE, PETSC_TRUE, PETSC_FALSE, PETSC_FALSE, PETSC_TRUE};
567: ISGeneralSetIndicesFromMask(is,10,15,mask);
568: .ve
569: will feed the `IS` with indices
570: .vb
571: {11, 14}
572: .ve
573: locally.
575: .seealso: [](sec_scatter), `IS`, `ISCreateGeneral()`, `ISGeneralSetIndices()`, `ISGENERAL`
576: @*/
577: PetscErrorCode ISGeneralSetIndicesFromMask(IS is, PetscInt rstart, PetscInt rend, const PetscBool mask[])
578: {
581: ISClearInfoCache(is, PETSC_FALSE);
582: PetscUseMethod(is, "ISGeneralSetIndicesFromMask_C", (IS, PetscInt, PetscInt, const PetscBool[]), (is, rstart, rend, mask));
583: return 0;
584: }
586: PetscErrorCode ISGeneralSetIndicesFromMask_General(IS is, PetscInt rstart, PetscInt rend, const PetscBool mask[])
587: {
588: PetscInt i, nidx;
589: PetscInt *idx;
591: for (i = 0, nidx = 0; i < rend - rstart; i++)
592: if (mask[i]) nidx++;
593: PetscMalloc1(nidx, &idx);
594: for (i = 0, nidx = 0; i < rend - rstart; i++) {
595: if (mask[i]) {
596: idx[nidx] = i + rstart;
597: nidx++;
598: }
599: }
600: ISGeneralSetIndices_General(is, nidx, idx, PETSC_OWN_POINTER);
601: return 0;
602: }
604: static PetscErrorCode ISGeneralFilter_General(IS is, PetscInt start, PetscInt end)
605: {
606: IS_General *sub = (IS_General *)is->data;
607: PetscInt *idx = sub->idx, *idxnew;
608: PetscInt i, n = is->map->n, nnew = 0, o;
610: for (i = 0; i < n; ++i)
611: if (idx[i] >= start && idx[i] < end) nnew++;
612: PetscMalloc1(nnew, &idxnew);
613: for (o = 0, i = 0; i < n; i++) {
614: if (idx[i] >= start && idx[i] < end) idxnew[o++] = idx[i];
615: }
616: ISGeneralSetIndices_General(is, nnew, idxnew, PETSC_OWN_POINTER);
617: return 0;
618: }
620: /*@
621: ISGeneralFilter - Remove all indices outside of [start, end) from an `ISGENERAL`
623: Collective
625: Input Parameters:
626: + is - the index set
627: . start - the lowest index kept
628: - end - one more than the highest index kept
630: Level: beginner
632: .seealso: [](sec_scatter), `IS`, `ISGENERAL`, `ISCreateGeneral()`, `ISGeneralSetIndices()`
633: @*/
634: PetscErrorCode ISGeneralFilter(IS is, PetscInt start, PetscInt end)
635: {
637: ISClearInfoCache(is, PETSC_FALSE);
638: PetscUseMethod(is, "ISGeneralFilter_C", (IS, PetscInt, PetscInt), (is, start, end));
639: return 0;
640: }
642: PETSC_EXTERN PetscErrorCode ISCreate_General(IS is)
643: {
644: IS_General *sub;
646: PetscNew(&sub);
647: is->data = (void *)sub;
648: PetscMemcpy(is->ops, &myops, sizeof(myops));
649: PetscObjectComposeFunction((PetscObject)is, "ISGeneralSetIndices_C", ISGeneralSetIndices_General);
650: PetscObjectComposeFunction((PetscObject)is, "ISGeneralSetIndicesFromMask_C", ISGeneralSetIndicesFromMask_General);
651: PetscObjectComposeFunction((PetscObject)is, "ISGeneralFilter_C", ISGeneralFilter_General);
652: PetscObjectComposeFunction((PetscObject)is, "ISShift_C", ISShift_General);
653: return 0;
654: }