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, &timestepping);
192:   if (timestepping) PetscViewerHDF5GetTimestep(viewer, &timestep);

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, &timestepping);
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: }