Actual source code: index.c

  1: /*
  2:    Defines the abstract operations on index sets, i.e. the public interface.
  3: */
  4: #include <petsc/private/isimpl.h>
  5: #include <petscviewer.h>
  6: #include <petscsf.h>

  8: /* Logging support */
  9: PetscClassId IS_CLASSID;
 10: /* TODO: Much more events are missing! */
 11: PetscLogEvent IS_View;
 12: PetscLogEvent IS_Load;

 14: /*@
 15:    ISRenumber - Renumbers the non-negative entries of an index set in a contiguous way, starting from 0.

 17:    Collective

 19:    Input Parameters:
 20: +  subset - the index set
 21: -  subset_mult - the multiplicity of each entry in subset (optional, can be NULL)

 23:    Output Parameters:
 24: +  N - one past the largest entry of the new IS
 25: -  subset_n - the new IS

 27:    Notes: All negative entries are mapped to -1. Indices with non positive multiplicities are skipped.

 29:    Level: intermediate

 31: .seealso:
 32: @*/
 33: PetscErrorCode ISRenumber(IS subset, IS subset_mult, PetscInt *N, IS *subset_n)
 34: {
 35:   PetscSF         sf;
 36:   PetscLayout     map;
 37:   const PetscInt *idxs, *idxs_mult = NULL;
 38:   PetscInt       *leaf_data, *root_data, *gidxs, *ilocal, *ilocalneg;
 39:   PetscInt        N_n, n, i, lbounds[2], gbounds[2], Nl, ibs;
 40:   PetscInt        n_n, nlocals, start, first_index, npos, nneg;
 41:   PetscMPIInt     commsize;
 42:   PetscBool       first_found, isblock;

 47:   else if (!subset_n) return 0;
 48:   ISGetLocalSize(subset, &n);
 49:   if (subset_mult) {
 50:     ISGetLocalSize(subset_mult, &i);
 52:   }
 53:   /* create workspace layout for computing global indices of subset */
 54:   PetscMalloc1(n, &ilocal);
 55:   PetscMalloc1(n, &ilocalneg);
 56:   ISGetIndices(subset, &idxs);
 57:   ISGetBlockSize(subset, &ibs);
 58:   PetscObjectTypeCompare((PetscObject)subset, ISBLOCK, &isblock);
 59:   if (subset_mult) ISGetIndices(subset_mult, &idxs_mult);
 60:   lbounds[0] = PETSC_MAX_INT;
 61:   lbounds[1] = PETSC_MIN_INT;
 62:   for (i = 0, npos = 0, nneg = 0; i < n; i++) {
 63:     if (idxs[i] < 0) {
 64:       ilocalneg[nneg++] = i;
 65:       continue;
 66:     }
 67:     if (idxs[i] < lbounds[0]) lbounds[0] = idxs[i];
 68:     if (idxs[i] > lbounds[1]) lbounds[1] = idxs[i];
 69:     ilocal[npos++] = i;
 70:   }
 71:   if (npos == n) {
 72:     PetscFree(ilocal);
 73:     PetscFree(ilocalneg);
 74:   }

 76:   /* create sf : leaf_data == multiplicity of indexes, root data == global index in layout */
 77:   PetscMalloc1(n, &leaf_data);
 78:   for (i = 0; i < n; i++) leaf_data[i] = idxs_mult ? PetscMax(idxs_mult[i], 0) : 1;

 80:   /* local size of new subset */
 81:   n_n = 0;
 82:   for (i = 0; i < n; i++) n_n += leaf_data[i];
 83:   if (ilocalneg)
 84:     for (i = 0; i < nneg; i++) leaf_data[ilocalneg[i]] = 0;
 85:   PetscFree(ilocalneg);
 86:   PetscMalloc1(PetscMax(n_n, n), &gidxs); /* allocating extra space to reuse gidxs */
 87:   /* check for early termination (all negative) */
 88:   PetscGlobalMinMaxInt(PetscObjectComm((PetscObject)subset), lbounds, gbounds);
 89:   if (gbounds[1] < gbounds[0]) {
 90:     if (N) *N = 0;
 91:     if (subset_n) { /* all negative */
 92:       for (i = 0; i < n_n; i++) gidxs[i] = -1;
 93:       ISCreateGeneral(PetscObjectComm((PetscObject)subset), n_n, gidxs, PETSC_COPY_VALUES, subset_n);
 94:     }
 95:     PetscFree(leaf_data);
 96:     PetscFree(gidxs);
 97:     ISRestoreIndices(subset, &idxs);
 98:     if (subset_mult) ISRestoreIndices(subset_mult, &idxs_mult);
 99:     PetscFree(ilocal);
100:     PetscFree(ilocalneg);
101:     return 0;
102:   }

104:   /* split work */
105:   N_n = gbounds[1] - gbounds[0] + 1;
106:   PetscLayoutCreate(PetscObjectComm((PetscObject)subset), &map);
107:   PetscLayoutSetBlockSize(map, 1);
108:   PetscLayoutSetSize(map, N_n);
109:   PetscLayoutSetUp(map);
110:   PetscLayoutGetLocalSize(map, &Nl);

112:   /* global indexes in layout */
113:   for (i = 0; i < npos; i++) gidxs[i] = (ilocal ? idxs[ilocal[i]] : idxs[i]) - gbounds[0];
114:   ISRestoreIndices(subset, &idxs);
115:   PetscSFCreate(PetscObjectComm((PetscObject)subset), &sf);
116:   PetscSFSetGraphLayout(sf, map, npos, ilocal, PETSC_USE_POINTER, gidxs);
117:   PetscLayoutDestroy(&map);

119:   /* reduce from leaves to roots */
120:   PetscCalloc1(Nl, &root_data);
121:   PetscSFReduceBegin(sf, MPIU_INT, leaf_data, root_data, MPI_MAX);
122:   PetscSFReduceEnd(sf, MPIU_INT, leaf_data, root_data, MPI_MAX);

124:   /* count indexes in local part of layout */
125:   nlocals     = 0;
126:   first_index = -1;
127:   first_found = PETSC_FALSE;
128:   for (i = 0; i < Nl; i++) {
129:     if (!first_found && root_data[i]) {
130:       first_found = PETSC_TRUE;
131:       first_index = i;
132:     }
133:     nlocals += root_data[i];
134:   }

136:   /* cumulative of number of indexes and size of subset without holes */
137: #if defined(PETSC_HAVE_MPI_EXSCAN)
138:   start = 0;
139:   MPI_Exscan(&nlocals, &start, 1, MPIU_INT, MPI_SUM, PetscObjectComm((PetscObject)subset));
140: #else
141:   MPI_Scan(&nlocals, &start, 1, MPIU_INT, MPI_SUM, PetscObjectComm((PetscObject)subset));
142:   start = start - nlocals;
143: #endif

145:   if (N) { /* compute total size of new subset if requested */
146:     *N = start + nlocals;
147:     MPI_Comm_size(PetscObjectComm((PetscObject)subset), &commsize);
148:     MPI_Bcast(N, 1, MPIU_INT, commsize - 1, PetscObjectComm((PetscObject)subset));
149:   }

151:   if (!subset_n) {
152:     PetscFree(gidxs);
153:     PetscSFDestroy(&sf);
154:     PetscFree(leaf_data);
155:     PetscFree(root_data);
156:     PetscFree(ilocal);
157:     if (subset_mult) ISRestoreIndices(subset_mult, &idxs_mult);
158:     return 0;
159:   }

161:   /* adapt root data with cumulative */
162:   if (first_found) {
163:     PetscInt old_index;

165:     root_data[first_index] += start;
166:     old_index = first_index;
167:     for (i = first_index + 1; i < Nl; i++) {
168:       if (root_data[i]) {
169:         root_data[i] += root_data[old_index];
170:         old_index = i;
171:       }
172:     }
173:   }

175:   /* from roots to leaves */
176:   PetscSFBcastBegin(sf, MPIU_INT, root_data, leaf_data, MPI_REPLACE);
177:   PetscSFBcastEnd(sf, MPIU_INT, root_data, leaf_data, MPI_REPLACE);
178:   PetscSFDestroy(&sf);

180:   /* create new IS with global indexes without holes */
181:   for (i = 0; i < n_n; i++) gidxs[i] = -1;
182:   if (subset_mult) {
183:     PetscInt cum;

185:     isblock = PETSC_FALSE;
186:     for (i = 0, cum = 0; i < n; i++)
187:       for (PetscInt j = 0; j < idxs_mult[i]; j++) gidxs[cum++] = leaf_data[i] - idxs_mult[i] + j;
188:   } else
189:     for (i = 0; i < n; i++) gidxs[i] = leaf_data[i] - 1;

191:   if (isblock) {
192:     if (ibs > 1)
193:       for (i = 0; i < n_n / ibs; i++) gidxs[i] = gidxs[i * ibs] / ibs;
194:     ISCreateBlock(PetscObjectComm((PetscObject)subset), ibs, n_n / ibs, gidxs, PETSC_COPY_VALUES, subset_n);
195:   } else {
196:     ISCreateGeneral(PetscObjectComm((PetscObject)subset), n_n, gidxs, PETSC_COPY_VALUES, subset_n);
197:   }
198:   if (subset_mult) ISRestoreIndices(subset_mult, &idxs_mult);
199:   PetscFree(gidxs);
200:   PetscFree(leaf_data);
201:   PetscFree(root_data);
202:   PetscFree(ilocal);
203:   return 0;
204: }

206: /*@
207:    ISCreateSubIS - Create a sub index set from a global index set selecting some components.

209:    Collective

211:    Input Parameters:
212: +  is - the index set
213: -  comps - which components we will extract from is

215:    Output Parameters:
216: .  subis - the new sub index set

218:    Level: intermediate

220:    Example usage:
221:    We have an index set (is) living on 3 processes with the following values:
222:    | 4 9 0 | 2 6 7 | 10 11 1|
223:    and another index set (comps) used to indicate which components of is  we want to take,
224:    | 7 5  | 1 2 | 0 4|
225:    The output index set (subis) should look like:
226:    | 11 7 | 9 0 | 4 6|

228: .seealso: `VecGetSubVector()`, `MatCreateSubMatrix()`
229: @*/
230: PetscErrorCode ISCreateSubIS(IS is, IS comps, IS *subis)
231: {
232:   PetscSF         sf;
233:   const PetscInt *is_indices, *comps_indices;
234:   PetscInt       *subis_indices, nroots, nleaves, *mine, i, lidx;
235:   PetscMPIInt     owner;
236:   PetscSFNode    *remote;
237:   MPI_Comm        comm;


243:   PetscObjectGetComm((PetscObject)is, &comm);
244:   ISGetLocalSize(comps, &nleaves);
245:   ISGetLocalSize(is, &nroots);
246:   PetscMalloc1(nleaves, &remote);
247:   PetscMalloc1(nleaves, &mine);
248:   ISGetIndices(comps, &comps_indices);
249:   /*
250:    * Construct a PetscSF in which "is" data serves as roots and "subis" is leaves.
251:    * Root data are sent to leaves using PetscSFBcast().
252:    * */
253:   for (i = 0; i < nleaves; i++) {
254:     mine[i] = i;
255:     /* Connect a remote root with the current leaf. The value on the remote root
256:      * will be received by the current local leaf.
257:      * */
258:     owner = -1;
259:     lidx  = -1;
260:     PetscLayoutFindOwnerIndex(is->map, comps_indices[i], &owner, &lidx);
261:     remote[i].rank  = owner;
262:     remote[i].index = lidx;
263:   }
264:   ISRestoreIndices(comps, &comps_indices);
265:   PetscSFCreate(comm, &sf);
266:   PetscSFSetFromOptions(sf);
267:   PetscSFSetGraph(sf, nroots, nleaves, mine, PETSC_OWN_POINTER, remote, PETSC_OWN_POINTER);

269:   PetscMalloc1(nleaves, &subis_indices);
270:   ISGetIndices(is, &is_indices);
271:   PetscSFBcastBegin(sf, MPIU_INT, is_indices, subis_indices, MPI_REPLACE);
272:   PetscSFBcastEnd(sf, MPIU_INT, is_indices, subis_indices, MPI_REPLACE);
273:   ISRestoreIndices(is, &is_indices);
274:   PetscSFDestroy(&sf);
275:   ISCreateGeneral(comm, nleaves, subis_indices, PETSC_OWN_POINTER, subis);
276:   return 0;
277: }

279: /*@
280:    ISClearInfoCache - clear the cache of computed index set properties

282:    Not collective

284:    Input Parameters:
285: +  is - the index set
286: -  clear_permanent_local - whether to remove the permanent status of local properties

288:    NOTE: because all processes must agree on the global permanent status of a property,
289:    the permanent status can only be changed with ISSetInfo(), because this routine is not collective

291:    Level: developer

293: .seealso: `ISInfo`, `ISInfoType`, `ISSetInfo()`, `ISClearInfoCache()`

295: @*/
296: PetscErrorCode ISClearInfoCache(IS is, PetscBool clear_permanent_local)
297: {
298:   PetscInt i, j;

302:   for (i = 0; i < IS_INFO_MAX; i++) {
303:     if (clear_permanent_local) is->info_permanent[IS_LOCAL][i] = PETSC_FALSE;
304:     for (j = 0; j < 2; j++) {
305:       if (!is->info_permanent[j][i]) is->info[j][i] = IS_INFO_UNKNOWN;
306:     }
307:   }
308:   return 0;
309: }

311: static PetscErrorCode ISSetInfo_Internal(IS is, ISInfo info, ISInfoType type, ISInfoBool ipermanent, PetscBool flg)
312: {
313:   ISInfoBool iflg          = flg ? IS_INFO_TRUE : IS_INFO_FALSE;
314:   PetscInt   itype         = (type == IS_LOCAL) ? 0 : 1;
315:   PetscBool  permanent_set = (ipermanent == IS_INFO_UNKNOWN) ? PETSC_FALSE : PETSC_TRUE;
316:   PetscBool  permanent     = (ipermanent == IS_INFO_TRUE) ? PETSC_TRUE : PETSC_FALSE;

318:   /* set this property */
319:   is->info[itype][(int)info] = iflg;
320:   if (permanent_set) is->info_permanent[itype][(int)info] = permanent;
321:   /* set implications */
322:   switch (info) {
323:   case IS_SORTED:
324:     if (flg && type == IS_GLOBAL) { /* an array that is globally sorted is also locally sorted */
325:       is->info[IS_LOCAL][(int)info] = IS_INFO_TRUE;
326:       /* global permanence implies local permanence */
327:       if (permanent_set && permanent) is->info_permanent[IS_LOCAL][(int)info] = PETSC_TRUE;
328:     }
329:     if (!flg) { /* if an array is not sorted, it cannot be an interval or the identity */
330:       is->info[itype][IS_INTERVAL] = IS_INFO_FALSE;
331:       is->info[itype][IS_IDENTITY] = IS_INFO_FALSE;
332:       if (permanent_set) {
333:         is->info_permanent[itype][IS_INTERVAL] = permanent;
334:         is->info_permanent[itype][IS_IDENTITY] = permanent;
335:       }
336:     }
337:     break;
338:   case IS_UNIQUE:
339:     if (flg && type == IS_GLOBAL) { /* an array that is globally unique is also locally unique */
340:       is->info[IS_LOCAL][(int)info] = IS_INFO_TRUE;
341:       /* global permanence implies local permanence */
342:       if (permanent_set && permanent) is->info_permanent[IS_LOCAL][(int)info] = PETSC_TRUE;
343:     }
344:     if (!flg) { /* if an array is not unique, it cannot be a permutation, and interval, or the identity */
345:       is->info[itype][IS_PERMUTATION] = IS_INFO_FALSE;
346:       is->info[itype][IS_INTERVAL]    = IS_INFO_FALSE;
347:       is->info[itype][IS_IDENTITY]    = IS_INFO_FALSE;
348:       if (permanent_set) {
349:         is->info_permanent[itype][IS_PERMUTATION] = permanent;
350:         is->info_permanent[itype][IS_INTERVAL]    = permanent;
351:         is->info_permanent[itype][IS_IDENTITY]    = permanent;
352:       }
353:     }
354:     break;
355:   case IS_PERMUTATION:
356:     if (flg) { /* an array that is a permutation is unique and is unique locally */
357:       is->info[itype][IS_UNIQUE]    = IS_INFO_TRUE;
358:       is->info[IS_LOCAL][IS_UNIQUE] = IS_INFO_TRUE;
359:       if (permanent_set && permanent) {
360:         is->info_permanent[itype][IS_UNIQUE]    = PETSC_TRUE;
361:         is->info_permanent[IS_LOCAL][IS_UNIQUE] = PETSC_TRUE;
362:       }
363:     } else { /* an array that is not a permutation cannot be the identity */
364:       is->info[itype][IS_IDENTITY] = IS_INFO_FALSE;
365:       if (permanent_set) is->info_permanent[itype][IS_IDENTITY] = permanent;
366:     }
367:     break;
368:   case IS_INTERVAL:
369:     if (flg) { /* an array that is an interval is sorted and unique */
370:       is->info[itype][IS_SORTED]    = IS_INFO_TRUE;
371:       is->info[IS_LOCAL][IS_SORTED] = IS_INFO_TRUE;
372:       is->info[itype][IS_UNIQUE]    = IS_INFO_TRUE;
373:       is->info[IS_LOCAL][IS_UNIQUE] = IS_INFO_TRUE;
374:       if (permanent_set && permanent) {
375:         is->info_permanent[itype][IS_SORTED]    = PETSC_TRUE;
376:         is->info_permanent[IS_LOCAL][IS_SORTED] = PETSC_TRUE;
377:         is->info_permanent[itype][IS_UNIQUE]    = PETSC_TRUE;
378:         is->info_permanent[IS_LOCAL][IS_UNIQUE] = PETSC_TRUE;
379:       }
380:     } else { /* an array that is not an interval cannot be the identity */
381:       is->info[itype][IS_IDENTITY] = IS_INFO_FALSE;
382:       if (permanent_set) is->info_permanent[itype][IS_IDENTITY] = permanent;
383:     }
384:     break;
385:   case IS_IDENTITY:
386:     if (flg) { /* an array that is the identity is sorted, unique, an interval, and a permutation */
387:       is->info[itype][IS_SORTED]      = IS_INFO_TRUE;
388:       is->info[IS_LOCAL][IS_SORTED]   = IS_INFO_TRUE;
389:       is->info[itype][IS_UNIQUE]      = IS_INFO_TRUE;
390:       is->info[IS_LOCAL][IS_UNIQUE]   = IS_INFO_TRUE;
391:       is->info[itype][IS_PERMUTATION] = IS_INFO_TRUE;
392:       is->info[itype][IS_INTERVAL]    = IS_INFO_TRUE;
393:       is->info[IS_LOCAL][IS_INTERVAL] = IS_INFO_TRUE;
394:       if (permanent_set && permanent) {
395:         is->info_permanent[itype][IS_SORTED]      = PETSC_TRUE;
396:         is->info_permanent[IS_LOCAL][IS_SORTED]   = PETSC_TRUE;
397:         is->info_permanent[itype][IS_UNIQUE]      = PETSC_TRUE;
398:         is->info_permanent[IS_LOCAL][IS_UNIQUE]   = PETSC_TRUE;
399:         is->info_permanent[itype][IS_PERMUTATION] = PETSC_TRUE;
400:         is->info_permanent[itype][IS_INTERVAL]    = PETSC_TRUE;
401:         is->info_permanent[IS_LOCAL][IS_INTERVAL] = PETSC_TRUE;
402:       }
403:     }
404:     break;
405:   default:
407:     SETERRQ(PetscObjectComm((PetscObject)is), PETSC_ERR_ARG_OUTOFRANGE, "Unknown IS property");
408:   }
409:   return 0;
410: }

412: /*@
413:    ISSetInfo - Set known information about an index set.

415:    Logically Collective on IS if type is IS_GLOBAL

417:    Input Parameters:
418: +  is - the index set
419: .  info - describing a property of the index set, one of those listed below,
420: .  type - IS_LOCAL if the information describes the local portion of the index set,
421:           IS_GLOBAL if it describes the whole index set
422: .  permanent - PETSC_TRUE if it is known that the property will persist through changes to the index set, PETSC_FALSE otherwise
423:                If the user sets a property as permanently known, it will bypass computation of that property
424: -  flg - set the described property as true (PETSC_TRUE) or false (PETSC_FALSE)

426:   Info Describing IS Structure:
427: +    IS_SORTED - the [local part of the] index set is sorted in ascending order
428: .    IS_UNIQUE - each entry in the [local part of the] index set is unique
429: .    IS_PERMUTATION - the [local part of the] index set is a permutation of the integers {0, 1, ..., N-1}, where N is the size of the [local part of the] index set
430: .    IS_INTERVAL - the [local part of the] index set is equal to a contiguous range of integers {f, f + 1, ..., f + N-1}
431: -    IS_IDENTITY - the [local part of the] index set is equal to the integers {0, 1, ..., N-1}

433:    Notes:
434:    If type is IS_GLOBAL, all processes that share the index set must pass the same value in flg

436:    It is possible to set a property with ISSetInfo() that contradicts what would be previously computed with ISGetInfo()

438:    Level: advanced

440: .seealso: `ISInfo`, `ISInfoType`, `IS`

442: @*/
443: PetscErrorCode ISSetInfo(IS is, ISInfo info, ISInfoType type, PetscBool permanent, PetscBool flg)
444: {
445:   MPI_Comm    comm, errcomm;
446:   PetscMPIInt size;

450:   comm = PetscObjectComm((PetscObject)is);
451:   if (type == IS_GLOBAL) {
455:     errcomm = comm;
456:   } else {
457:     errcomm = PETSC_COMM_SELF;
458:   }


462:   MPI_Comm_size(comm, &size);
463:   /* do not use global values if size == 1: it makes it easier to keep the implications straight */
464:   if (size == 1) type = IS_LOCAL;
465:   ISSetInfo_Internal(is, info, type, permanent ? IS_INFO_TRUE : IS_INFO_FALSE, flg);
466:   return 0;
467: }

469: static PetscErrorCode ISGetInfo_Sorted(IS is, ISInfoType type, PetscBool *flg)
470: {
471:   MPI_Comm    comm;
472:   PetscMPIInt size, rank;

474:   comm = PetscObjectComm((PetscObject)is);
475:   MPI_Comm_size(comm, &size);
476:   MPI_Comm_size(comm, &rank);
477:   if (type == IS_GLOBAL && is->ops->sortedglobal) {
478:     PetscUseTypeMethod(is, sortedglobal, flg);
479:   } else {
480:     PetscBool sortedLocal = PETSC_FALSE;

482:     /* determine if the array is locally sorted */
483:     if (type == IS_GLOBAL && size > 1) {
484:       /* call ISGetInfo so that a cached value will be used if possible */
485:       ISGetInfo(is, IS_SORTED, IS_LOCAL, PETSC_TRUE, &sortedLocal);
486:     } else if (is->ops->sortedlocal) {
487:       PetscUseTypeMethod(is, sortedlocal, &sortedLocal);
488:     } else {
489:       /* default: get the local indices and directly check */
490:       const PetscInt *idx;
491:       PetscInt        n;

493:       ISGetIndices(is, &idx);
494:       ISGetLocalSize(is, &n);
495:       PetscSortedInt(n, idx, &sortedLocal);
496:       ISRestoreIndices(is, &idx);
497:     }

499:     if (type == IS_LOCAL || size == 1) {
500:       *flg = sortedLocal;
501:     } else {
502:       MPI_Allreduce(&sortedLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
503:       if (*flg) {
504:         PetscInt n, min = PETSC_MAX_INT, max = PETSC_MIN_INT;
505:         PetscInt maxprev;

507:         ISGetLocalSize(is, &n);
508:         if (n) ISGetMinMax(is, &min, &max);
509:         maxprev = PETSC_MIN_INT;
510:         MPI_Exscan(&max, &maxprev, 1, MPIU_INT, MPI_MAX, comm);
511:         if (rank && (maxprev > min)) sortedLocal = PETSC_FALSE;
512:         MPI_Allreduce(&sortedLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
513:       }
514:     }
515:   }
516:   return 0;
517: }

519: PetscErrorCode ISGetIndicesCopy(IS is, PetscInt idx[]);

521: static PetscErrorCode ISGetInfo_Unique(IS is, ISInfoType type, PetscBool *flg)
522: {
523:   MPI_Comm    comm;
524:   PetscMPIInt size, rank;
525:   PetscInt    i;

527:   comm = PetscObjectComm((PetscObject)is);
528:   MPI_Comm_size(comm, &size);
529:   MPI_Comm_size(comm, &rank);
530:   if (type == IS_GLOBAL && is->ops->uniqueglobal) {
531:     PetscUseTypeMethod(is, uniqueglobal, flg);
532:   } else {
533:     PetscBool uniqueLocal;
534:     PetscInt  n   = -1;
535:     PetscInt *idx = NULL;

537:     /* determine if the array is locally unique */
538:     if (type == IS_GLOBAL && size > 1) {
539:       /* call ISGetInfo so that a cached value will be used if possible */
540:       ISGetInfo(is, IS_UNIQUE, IS_LOCAL, PETSC_TRUE, &uniqueLocal);
541:     } else if (is->ops->uniquelocal) {
542:       PetscUseTypeMethod(is, uniquelocal, &uniqueLocal);
543:     } else {
544:       /* default: get the local indices and directly check */
545:       uniqueLocal = PETSC_TRUE;
546:       ISGetLocalSize(is, &n);
547:       PetscMalloc1(n, &idx);
548:       ISGetIndicesCopy(is, idx);
549:       PetscIntSortSemiOrdered(n, idx);
550:       for (i = 1; i < n; i++)
551:         if (idx[i] == idx[i - 1]) break;
552:       if (i < n) uniqueLocal = PETSC_FALSE;
553:     }

555:     PetscFree(idx);
556:     if (type == IS_LOCAL || size == 1) {
557:       *flg = uniqueLocal;
558:     } else {
559:       MPI_Allreduce(&uniqueLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
560:       if (*flg) {
561:         PetscInt min = PETSC_MAX_INT, max = PETSC_MIN_INT, maxprev;

563:         if (!idx) {
564:           ISGetLocalSize(is, &n);
565:           PetscMalloc1(n, &idx);
566:           ISGetIndicesCopy(is, idx);
567:         }
568:         PetscParallelSortInt(is->map, is->map, idx, idx);
569:         if (n) {
570:           min = idx[0];
571:           max = idx[n - 1];
572:         }
573:         for (i = 1; i < n; i++)
574:           if (idx[i] == idx[i - 1]) break;
575:         if (i < n) uniqueLocal = PETSC_FALSE;
576:         maxprev = PETSC_MIN_INT;
577:         MPI_Exscan(&max, &maxprev, 1, MPIU_INT, MPI_MAX, comm);
578:         if (rank && (maxprev == min)) uniqueLocal = PETSC_FALSE;
579:         MPI_Allreduce(&uniqueLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
580:       }
581:     }
582:     PetscFree(idx);
583:   }
584:   return 0;
585: }

587: static PetscErrorCode ISGetInfo_Permutation(IS is, ISInfoType type, PetscBool *flg)
588: {
589:   MPI_Comm    comm;
590:   PetscMPIInt size, rank;

592:   comm = PetscObjectComm((PetscObject)is);
593:   MPI_Comm_size(comm, &size);
594:   MPI_Comm_size(comm, &rank);
595:   if (type == IS_GLOBAL && is->ops->permglobal) {
596:     PetscUseTypeMethod(is, permglobal, flg);
597:   } else if (type == IS_LOCAL && is->ops->permlocal) {
598:     PetscUseTypeMethod(is, permlocal, flg);
599:   } else {
600:     PetscBool permLocal;
601:     PetscInt  n, i, rStart;
602:     PetscInt *idx;

604:     ISGetLocalSize(is, &n);
605:     PetscMalloc1(n, &idx);
606:     ISGetIndicesCopy(is, idx);
607:     if (type == IS_GLOBAL) {
608:       PetscParallelSortInt(is->map, is->map, idx, idx);
609:       PetscLayoutGetRange(is->map, &rStart, NULL);
610:     } else {
611:       PetscIntSortSemiOrdered(n, idx);
612:       rStart = 0;
613:     }
614:     permLocal = PETSC_TRUE;
615:     for (i = 0; i < n; i++) {
616:       if (idx[i] != rStart + i) break;
617:     }
618:     if (i < n) permLocal = PETSC_FALSE;
619:     if (type == IS_LOCAL || size == 1) {
620:       *flg = permLocal;
621:     } else {
622:       MPI_Allreduce(&permLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
623:     }
624:     PetscFree(idx);
625:   }
626:   return 0;
627: }

629: static PetscErrorCode ISGetInfo_Interval(IS is, ISInfoType type, PetscBool *flg)
630: {
631:   MPI_Comm    comm;
632:   PetscMPIInt size, rank;
633:   PetscInt    i;

635:   comm = PetscObjectComm((PetscObject)is);
636:   MPI_Comm_size(comm, &size);
637:   MPI_Comm_size(comm, &rank);
638:   if (type == IS_GLOBAL && is->ops->intervalglobal) {
639:     PetscUseTypeMethod(is, intervalglobal, flg);
640:   } else {
641:     PetscBool intervalLocal;

643:     /* determine if the array is locally an interval */
644:     if (type == IS_GLOBAL && size > 1) {
645:       /* call ISGetInfo so that a cached value will be used if possible */
646:       ISGetInfo(is, IS_INTERVAL, IS_LOCAL, PETSC_TRUE, &intervalLocal);
647:     } else if (is->ops->intervallocal) {
648:       PetscUseTypeMethod(is, intervallocal, &intervalLocal);
649:     } else {
650:       PetscInt        n;
651:       const PetscInt *idx;
652:       /* default: get the local indices and directly check */
653:       intervalLocal = PETSC_TRUE;
654:       ISGetLocalSize(is, &n);
655:       ISGetIndices(is, &idx);
656:       for (i = 1; i < n; i++)
657:         if (idx[i] != idx[i - 1] + 1) break;
658:       if (i < n) intervalLocal = PETSC_FALSE;
659:       ISRestoreIndices(is, &idx);
660:     }

662:     if (type == IS_LOCAL || size == 1) {
663:       *flg = intervalLocal;
664:     } else {
665:       MPI_Allreduce(&intervalLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
666:       if (*flg) {
667:         PetscInt n, min = PETSC_MAX_INT, max = PETSC_MIN_INT;
668:         PetscInt maxprev;

670:         ISGetLocalSize(is, &n);
671:         if (n) ISGetMinMax(is, &min, &max);
672:         maxprev = PETSC_MIN_INT;
673:         MPI_Exscan(&max, &maxprev, 1, MPIU_INT, MPI_MAX, comm);
674:         if (rank && n && (maxprev != min - 1)) intervalLocal = PETSC_FALSE;
675:         MPI_Allreduce(&intervalLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
676:       }
677:     }
678:   }
679:   return 0;
680: }

682: static PetscErrorCode ISGetInfo_Identity(IS is, ISInfoType type, PetscBool *flg)
683: {
684:   MPI_Comm    comm;
685:   PetscMPIInt size, rank;

687:   comm = PetscObjectComm((PetscObject)is);
688:   MPI_Comm_size(comm, &size);
689:   MPI_Comm_size(comm, &rank);
690:   if (type == IS_GLOBAL && is->ops->intervalglobal) {
691:     PetscBool isinterval;

693:     PetscUseTypeMethod(is, intervalglobal, &isinterval);
694:     *flg = PETSC_FALSE;
695:     if (isinterval) {
696:       PetscInt min;

698:       ISGetMinMax(is, &min, NULL);
699:       MPI_Bcast(&min, 1, MPIU_INT, 0, comm);
700:       if (min == 0) *flg = PETSC_TRUE;
701:     }
702:   } else if (type == IS_LOCAL && is->ops->intervallocal) {
703:     PetscBool isinterval;

705:     PetscUseTypeMethod(is, intervallocal, &isinterval);
706:     *flg = PETSC_FALSE;
707:     if (isinterval) {
708:       PetscInt min;

710:       ISGetMinMax(is, &min, NULL);
711:       if (min == 0) *flg = PETSC_TRUE;
712:     }
713:   } else {
714:     PetscBool       identLocal;
715:     PetscInt        n, i, rStart;
716:     const PetscInt *idx;

718:     ISGetLocalSize(is, &n);
719:     ISGetIndices(is, &idx);
720:     PetscLayoutGetRange(is->map, &rStart, NULL);
721:     identLocal = PETSC_TRUE;
722:     for (i = 0; i < n; i++) {
723:       if (idx[i] != rStart + i) break;
724:     }
725:     if (i < n) identLocal = PETSC_FALSE;
726:     if (type == IS_LOCAL || size == 1) {
727:       *flg = identLocal;
728:     } else {
729:       MPI_Allreduce(&identLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
730:     }
731:     ISRestoreIndices(is, &idx);
732:   }
733:   return 0;
734: }

736: /*@
737:    ISGetInfo - Determine whether an index set satisfies a given property

739:    Collective or logically collective on IS if the type is IS_GLOBAL (logically collective if the value of the property has been permanently set with ISSetInfo())

741:    Input Parameters:
742: +  is - the index set
743: .  info - describing a property of the index set, one of those listed in the documentation of ISSetInfo()
744: .  compute - if PETSC_FALSE, the property will not be computed if it is not already known and the property will be assumed to be false
745: -  type - whether the property is local (IS_LOCAL) or global (IS_GLOBAL)

747:    Output Parameter:
748: .  flg - wheter the property is true (PETSC_TRUE) or false (PETSC_FALSE)

750:    Note: ISGetInfo uses cached values when possible, which will be incorrect if ISSetInfo() has been called with incorrect information.  To clear cached values, use ISClearInfoCache().

752:    Level: advanced

754: .seealso: `ISInfo`, `ISInfoType`, `ISSetInfo()`, `ISClearInfoCache()`

756: @*/
757: PetscErrorCode ISGetInfo(IS is, ISInfo info, ISInfoType type, PetscBool compute, PetscBool *flg)
758: {
759:   MPI_Comm    comm, errcomm;
760:   PetscMPIInt rank, size;
761:   PetscInt    itype;
762:   PetscBool   hasprop;
763:   PetscBool   infer;

767:   comm = PetscObjectComm((PetscObject)is);
768:   if (type == IS_GLOBAL) {
770:     errcomm = comm;
771:   } else {
772:     errcomm = PETSC_COMM_SELF;
773:   }

775:   MPI_Comm_size(comm, &size);
776:   MPI_Comm_rank(comm, &rank);

779:   if (size == 1) type = IS_LOCAL;
780:   itype   = (type == IS_LOCAL) ? 0 : 1;
781:   hasprop = PETSC_FALSE;
782:   infer   = PETSC_FALSE;
783:   if (is->info_permanent[itype][(int)info]) {
784:     hasprop = (is->info[itype][(int)info] == IS_INFO_TRUE) ? PETSC_TRUE : PETSC_FALSE;
785:     infer   = PETSC_TRUE;
786:   } else if ((itype == IS_LOCAL) && (is->info[IS_LOCAL][info] != IS_INFO_UNKNOWN)) {
787:     /* we can cache local properties as long as we clear them when the IS changes */
788:     /* NOTE: we only cache local values because there is no ISAssemblyBegin()/ISAssemblyEnd(),
789:      so we have no way of knowing when a cached value has been invalidated by changes on a different process */
790:     hasprop = (is->info[itype][(int)info] == IS_INFO_TRUE) ? PETSC_TRUE : PETSC_FALSE;
791:     infer   = PETSC_TRUE;
792:   } else if (compute) {
793:     switch (info) {
794:     case IS_SORTED:
795:       ISGetInfo_Sorted(is, type, &hasprop);
796:       break;
797:     case IS_UNIQUE:
798:       ISGetInfo_Unique(is, type, &hasprop);
799:       break;
800:     case IS_PERMUTATION:
801:       ISGetInfo_Permutation(is, type, &hasprop);
802:       break;
803:     case IS_INTERVAL:
804:       ISGetInfo_Interval(is, type, &hasprop);
805:       break;
806:     case IS_IDENTITY:
807:       ISGetInfo_Identity(is, type, &hasprop);
808:       break;
809:     default:
810:       SETERRQ(errcomm, PETSC_ERR_ARG_OUTOFRANGE, "Unknown IS property");
811:     }
812:     infer = PETSC_TRUE;
813:   }
814:   /* call ISSetInfo_Internal to keep all of the implications straight */
815:   if (infer) ISSetInfo_Internal(is, info, type, IS_INFO_UNKNOWN, hasprop);
816:   *flg = hasprop;
817:   return 0;
818: }

820: static PetscErrorCode ISCopyInfo(IS source, IS dest)
821: {
822:   PetscArraycpy(&dest->info[0], &source->info[0], 2);
823:   PetscArraycpy(&dest->info_permanent[0], &source->info_permanent[0], 2);
824:   return 0;
825: }

827: /*@
828:    ISIdentity - Determines whether index set is the identity mapping.

830:    Collective

832:    Input Parameters:
833: .  is - the index set

835:    Output Parameters:
836: .  ident - PETSC_TRUE if an identity, else PETSC_FALSE

838:    Level: intermediate

840:    Note: If ISSetIdentity() (or ISSetInfo() for a permanent property) has been called,
841:    ISIdentity() will return its answer without communication between processes, but
842:    otherwise the output ident will be computed from ISGetInfo(),
843:    which may require synchronization on the communicator of IS.  To avoid this computation,
844:    call ISGetInfo() directly with the compute flag set to PETSC_FALSE, and ident will be assumed false.

846: .seealso: `ISSetIdentity()`, `ISGetInfo()`
847: @*/
848: PetscErrorCode ISIdentity(IS is, PetscBool *ident)
849: {
852:   ISGetInfo(is, IS_IDENTITY, IS_GLOBAL, PETSC_TRUE, ident);
853:   return 0;
854: }

856: /*@
857:    ISSetIdentity - Informs the index set that it is an identity.

859:    Logically Collective

861:    Input Parameter:
862: .  is - the index set

864:    Level: intermediate

866:    Note: The IS will be considered the identity permanently, even if indices have been changes (for example, with
867:    ISGeneralSetIndices()).  It's a good idea to only set this property if the IS will not change in the future.
868:    To clear this property, use ISClearInfoCache().

870: .seealso: `ISIdentity()`, `ISSetInfo()`, `ISClearInfoCache()`
871: @*/
872: PetscErrorCode ISSetIdentity(IS is)
873: {
875:   ISSetInfo(is, IS_IDENTITY, IS_GLOBAL, PETSC_TRUE, PETSC_TRUE);
876:   return 0;
877: }

879: /*@
880:    ISContiguousLocal - Locates an index set with contiguous range within a global range, if possible

882:    Not Collective

884:    Input Parameters:
885: +  is - the index set
886: .  gstart - global start
887: -  gend - global end

889:    Output Parameters:
890: +  start - start of contiguous block, as an offset from gstart
891: -  contig - PETSC_TRUE if the index set refers to contiguous entries on this process, else PETSC_FALSE

893:    Level: developer

895: .seealso: `ISGetLocalSize()`, `VecGetOwnershipRange()`
896: @*/
897: PetscErrorCode ISContiguousLocal(IS is, PetscInt gstart, PetscInt gend, PetscInt *start, PetscBool *contig)
898: {
902:   *start  = -1;
903:   *contig = PETSC_FALSE;
904:   PetscTryTypeMethod(is, contiguous, gstart, gend, start, contig);
905:   return 0;
906: }

908: /*@
909:    ISPermutation - PETSC_TRUE or PETSC_FALSE depending on whether the
910:    index set has been declared to be a permutation.

912:    Logically Collective

914:    Input Parameter:
915: .  is - the index set

917:    Output Parameter:
918: .  perm - PETSC_TRUE if a permutation, else PETSC_FALSE

920:    Level: intermediate

922:    Note: If it is not already known that the IS is a permutation (if ISSetPermutation()
923:    or ISSetInfo() has not been called), this routine will not attempt to compute
924:    whether the index set is a permutation and will assume perm is PETSC_FALSE.
925:    To compute the value when it is not already known, use ISGetInfo() with
926:    the compute flag set to PETSC_TRUE.

928: .seealso: `ISSetPermutation()`, `ISGetInfo()`
929: @*/
930: PetscErrorCode ISPermutation(IS is, PetscBool *perm)
931: {
934:   ISGetInfo(is, IS_PERMUTATION, IS_GLOBAL, PETSC_FALSE, perm);
935:   return 0;
936: }

938: /*@
939:    ISSetPermutation - Informs the index set that it is a permutation.

941:    Logically Collective

943:    Input Parameter:
944: .  is - the index set

946:    Level: intermediate

948:    The debug version of the libraries (./configure --with-debugging=1) checks if the
949:   index set is actually a permutation. The optimized version just believes you.

951:    Note: The IS will be considered a permutation permanently, even if indices have been changes (for example, with
952:    ISGeneralSetIndices()).  It's a good idea to only set this property if the IS will not change in the future.
953:    To clear this property, use ISClearInfoCache().

955: .seealso: `ISPermutation()`, `ISSetInfo()`, `ISClearInfoCache().`
956: @*/
957: PetscErrorCode ISSetPermutation(IS is)
958: {
960:   if (PetscDefined(USE_DEBUG)) {
961:     PetscMPIInt size;

963:     MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
964:     if (size == 1) {
965:       PetscInt        i, n, *idx;
966:       const PetscInt *iidx;

968:       ISGetSize(is, &n);
969:       PetscMalloc1(n, &idx);
970:       ISGetIndices(is, &iidx);
971:       PetscArraycpy(idx, iidx, n);
972:       PetscIntSortSemiOrdered(n, idx);
974:       PetscFree(idx);
975:       ISRestoreIndices(is, &iidx);
976:     }
977:   }
978:   ISSetInfo(is, IS_PERMUTATION, IS_GLOBAL, PETSC_TRUE, PETSC_TRUE);
979:   return 0;
980: }

982: /*@C
983:    ISDestroy - Destroys an index set.

985:    Collective

987:    Input Parameters:
988: .  is - the index set

990:    Level: beginner

992: .seealso: `ISCreateGeneral()`, `ISCreateStride()`, `ISCreateBlocked()`
993: @*/
994: PetscErrorCode ISDestroy(IS *is)
995: {
996:   if (!*is) return 0;
998:   if (--((PetscObject)(*is))->refct > 0) {
999:     *is = NULL;
1000:     return 0;
1001:   }
1002:   if ((*is)->complement) {
1003:     PetscInt refcnt;
1004:     PetscObjectGetReference((PetscObject)((*is)->complement), &refcnt);
1006:     ISDestroy(&(*is)->complement);
1007:   }
1008:   if ((*is)->ops->destroy) (*(*is)->ops->destroy)(*is);
1009:   PetscLayoutDestroy(&(*is)->map);
1010:   /* Destroy local representations of offproc data. */
1011:   PetscFree((*is)->total);
1012:   PetscFree((*is)->nonlocal);
1013:   PetscHeaderDestroy(is);
1014:   return 0;
1015: }

1017: /*@
1018:    ISInvertPermutation - Creates a new permutation that is the inverse of
1019:                          a given permutation.

1021:    Collective

1023:    Input Parameters:
1024: +  is - the index set
1025: -  nlocal - number of indices on this processor in result (ignored for 1 proccessor) or
1026:             use PETSC_DECIDE

1028:    Output Parameter:
1029: .  isout - the inverse permutation

1031:    Level: intermediate

1033:    Notes:
1034:     For parallel index sets this does the complete parallel permutation, but the
1035:     code is not efficient for huge index sets (10,000,000 indices).

1037: @*/
1038: PetscErrorCode ISInvertPermutation(IS is, PetscInt nlocal, IS *isout)
1039: {
1040:   PetscBool isperm, isidentity, issame;

1044:   ISGetInfo(is, IS_PERMUTATION, IS_GLOBAL, PETSC_TRUE, &isperm);
1046:   ISGetInfo(is, IS_IDENTITY, IS_GLOBAL, PETSC_TRUE, &isidentity);
1047:   issame = PETSC_FALSE;
1048:   if (isidentity) {
1049:     PetscInt  n;
1050:     PetscBool isallsame;

1052:     ISGetLocalSize(is, &n);
1053:     issame = (PetscBool)(n == nlocal);
1054:     MPI_Allreduce(&issame, &isallsame, 1, MPIU_BOOL, MPI_LAND, PetscObjectComm((PetscObject)is));
1055:     issame = isallsame;
1056:   }
1057:   if (issame) {
1058:     ISDuplicate(is, isout);
1059:   } else {
1060:     PetscUseTypeMethod(is, invertpermutation, nlocal, isout);
1061:     ISSetPermutation(*isout);
1062:   }
1063:   return 0;
1064: }

1066: /*@
1067:    ISGetSize - Returns the global length of an index set.

1069:    Not Collective

1071:    Input Parameter:
1072: .  is - the index set

1074:    Output Parameter:
1075: .  size - the global size

1077:    Level: beginner

1079: @*/
1080: PetscErrorCode ISGetSize(IS is, PetscInt *size)
1081: {
1084:   *size = is->map->N;
1085:   return 0;
1086: }

1088: /*@
1089:    ISGetLocalSize - Returns the local (processor) length of an index set.

1091:    Not Collective

1093:    Input Parameter:
1094: .  is - the index set

1096:    Output Parameter:
1097: .  size - the local size

1099:    Level: beginner

1101: @*/
1102: PetscErrorCode ISGetLocalSize(IS is, PetscInt *size)
1103: {
1106:   *size = is->map->n;
1107:   return 0;
1108: }

1110: /*@
1111:    ISGetLayout - get PetscLayout describing index set layout

1113:    Not Collective

1115:    Input Parameter:
1116: .  is - the index set

1118:    Output Parameter:
1119: .  map - the layout

1121:    Level: developer

1123: .seealso: `ISSetLayout()`, `ISGetSize()`, `ISGetLocalSize()`
1124: @*/
1125: PetscErrorCode ISGetLayout(IS is, PetscLayout *map)
1126: {
1129:   *map = is->map;
1130:   return 0;
1131: }

1133: /*@
1134:    ISSetLayout - set PetscLayout describing index set layout

1136:    Collective

1138:    Input Arguments:
1139: +  is - the index set
1140: -  map - the layout

1142:    Level: developer

1144:    Notes:
1145:    Users should typically use higher level functions such as ISCreateGeneral().

1147:    This function can be useful in some special cases of constructing a new IS, e.g. after ISCreate() and before ISLoad().
1148:    Otherwise, it is only valid to replace the layout with a layout known to be equivalent.

1150: .seealso: `ISCreate()`, `ISGetLayout()`, `ISGetSize()`, `ISGetLocalSize()`
1151: @*/
1152: PetscErrorCode ISSetLayout(IS is, PetscLayout map)
1153: {
1156:   PetscLayoutReference(map, &is->map);
1157:   return 0;
1158: }

1160: /*@C
1161:    ISGetIndices - Returns a pointer to the indices.  The user should call
1162:    ISRestoreIndices() after having looked at the indices.  The user should
1163:    NOT change the indices.

1165:    Not Collective

1167:    Input Parameter:
1168: .  is - the index set

1170:    Output Parameter:
1171: .  ptr - the location to put the pointer to the indices

1173:    Fortran Note:
1174:    This routine has two different interfaces from Fortran; the first is not recommend, it does not require Fortran 90
1175: $    IS          is
1176: $    integer     is_array(1)
1177: $    PetscOffset i_is
1178: $    int         ierr
1179: $       call ISGetIndices(is,is_array,i_is,ierr)
1180: $
1181: $   Access first local entry in list
1182: $      value = is_array(i_is + 1)
1183: $
1184: $      ...... other code
1185: $       call ISRestoreIndices(is,is_array,i_is,ierr)
1186:    The second Fortran interface is recommended.
1187: $          use petscisdef
1188: $          PetscInt, pointer :: array(:)
1189: $          PetscErrorCode  ierr
1190: $          IS       i
1191: $          call ISGetIndicesF90(i,array,ierr)

1193:    See the Fortran chapter of the users manual and
1194:    petsc/src/is/[tutorials,tests] for details.

1196:    Level: intermediate

1198: .seealso: `ISRestoreIndices()`, `ISGetIndicesF90()`
1199: @*/
1200: PetscErrorCode ISGetIndices(IS is, const PetscInt *ptr[])
1201: {
1204:   PetscUseTypeMethod(is, getindices, ptr);
1205:   return 0;
1206: }

1208: /*@C
1209:    ISGetMinMax - Gets the minimum and maximum values in an IS

1211:    Not Collective

1213:    Input Parameter:
1214: .  is - the index set

1216:    Output Parameters:
1217: +   min - the minimum value
1218: -   max - the maximum value

1220:    Level: intermediate

1222:    Notes:
1223:     Empty index sets return min=PETSC_MAX_INT and max=PETSC_MIN_INT.
1224:     In parallel, it returns the min and max of the local portion of the IS

1226: .seealso: `ISGetIndices()`, `ISRestoreIndices()`, `ISGetIndicesF90()`
1227: @*/
1228: PetscErrorCode ISGetMinMax(IS is, PetscInt *min, PetscInt *max)
1229: {
1231:   if (min) *min = is->min;
1232:   if (max) *max = is->max;
1233:   return 0;
1234: }

1236: /*@
1237:   ISLocate - determine the location of an index within the local component of an index set

1239:   Not Collective

1241:   Input Parameters:
1242: + is - the index set
1243: - key - the search key

1245:   Output Parameter:
1246: . location - if >= 0, a location within the index set that is equal to the key, otherwise the key is not in the index set

1248:   Level: intermediate
1249: @*/
1250: PetscErrorCode ISLocate(IS is, PetscInt key, PetscInt *location)
1251: {
1252:   if (is->ops->locate) {
1253:     PetscUseTypeMethod(is, locate, key, location);
1254:   } else {
1255:     PetscInt        numIdx;
1256:     PetscBool       sorted;
1257:     const PetscInt *idx;

1259:     ISGetLocalSize(is, &numIdx);
1260:     ISGetIndices(is, &idx);
1261:     ISSorted(is, &sorted);
1262:     if (sorted) {
1263:       PetscFindInt(key, numIdx, idx, location);
1264:     } else {
1265:       PetscInt i;

1267:       *location = -1;
1268:       for (i = 0; i < numIdx; i++) {
1269:         if (idx[i] == key) {
1270:           *location = i;
1271:           break;
1272:         }
1273:       }
1274:     }
1275:     ISRestoreIndices(is, &idx);
1276:   }
1277:   return 0;
1278: }

1280: /*@C
1281:    ISRestoreIndices - Restores an index set to a usable state after a call
1282:                       to ISGetIndices().

1284:    Not Collective

1286:    Input Parameters:
1287: +  is - the index set
1288: -  ptr - the pointer obtained by ISGetIndices()

1290:    Fortran Note:
1291:    This routine is used differently from Fortran
1292: $    IS          is
1293: $    integer     is_array(1)
1294: $    PetscOffset i_is
1295: $    int         ierr
1296: $       call ISGetIndices(is,is_array,i_is,ierr)
1297: $
1298: $   Access first local entry in list
1299: $      value = is_array(i_is + 1)
1300: $
1301: $      ...... other code
1302: $       call ISRestoreIndices(is,is_array,i_is,ierr)

1304:    See the Fortran chapter of the users manual and
1305:    petsc/src/vec/is/tests for details.

1307:    Level: intermediate

1309:    Note:
1310:    This routine zeros out ptr. This is to prevent accidental us of the array after it has been restored.

1312: .seealso: `ISGetIndices()`, `ISRestoreIndicesF90()`
1313: @*/
1314: PetscErrorCode ISRestoreIndices(IS is, const PetscInt *ptr[])
1315: {
1318:   PetscTryTypeMethod(is, restoreindices, ptr);
1319:   return 0;
1320: }

1322: static PetscErrorCode ISGatherTotal_Private(IS is)
1323: {
1324:   PetscInt        i, n, N;
1325:   const PetscInt *lindices;
1326:   MPI_Comm        comm;
1327:   PetscMPIInt     rank, size, *sizes = NULL, *offsets = NULL, nn;


1331:   PetscObjectGetComm((PetscObject)is, &comm);
1332:   MPI_Comm_size(comm, &size);
1333:   MPI_Comm_rank(comm, &rank);
1334:   ISGetLocalSize(is, &n);
1335:   PetscMalloc2(size, &sizes, size, &offsets);

1337:   PetscMPIIntCast(n, &nn);
1338:   MPI_Allgather(&nn, 1, MPI_INT, sizes, 1, MPI_INT, comm);
1339:   offsets[0] = 0;
1340:   for (i = 1; i < size; ++i) offsets[i] = offsets[i - 1] + sizes[i - 1];
1341:   N = offsets[size - 1] + sizes[size - 1];

1343:   PetscMalloc1(N, &(is->total));
1344:   ISGetIndices(is, &lindices);
1345:   MPI_Allgatherv((void *)lindices, nn, MPIU_INT, is->total, sizes, offsets, MPIU_INT, comm);
1346:   ISRestoreIndices(is, &lindices);
1347:   is->local_offset = offsets[rank];
1348:   PetscFree2(sizes, offsets);
1349:   return 0;
1350: }

1352: /*@C
1353:    ISGetTotalIndices - Retrieve an array containing all indices across the communicator.

1355:    Collective

1357:    Input Parameter:
1358: .  is - the index set

1360:    Output Parameter:
1361: .  indices - total indices with rank 0 indices first, and so on; total array size is
1362:              the same as returned with ISGetSize().

1364:    Level: intermediate

1366:    Notes:
1367:     this is potentially nonscalable, but depends on the size of the total index set
1368:      and the size of the communicator. This may be feasible for index sets defined on
1369:      subcommunicators, such that the set size does not grow with PETSC_WORLD_COMM.
1370:      Note also that there is no way to tell where the local part of the indices starts
1371:      (use ISGetIndices() and ISGetNonlocalIndices() to retrieve just the local and just
1372:       the nonlocal part (complement), respectively).

1374: .seealso: `ISRestoreTotalIndices()`, `ISGetNonlocalIndices()`, `ISGetSize()`
1375: @*/
1376: PetscErrorCode ISGetTotalIndices(IS is, const PetscInt *indices[])
1377: {
1378:   PetscMPIInt size;

1382:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
1383:   if (size == 1) {
1384:     PetscUseTypeMethod(is, getindices, indices);
1385:   } else {
1386:     if (!is->total) ISGatherTotal_Private(is);
1387:     *indices = is->total;
1388:   }
1389:   return 0;
1390: }

1392: /*@C
1393:    ISRestoreTotalIndices - Restore the index array obtained with ISGetTotalIndices().

1395:    Not Collective.

1397:    Input Parameters:
1398: +  is - the index set
1399: -  indices - index array; must be the array obtained with ISGetTotalIndices()

1401:    Level: intermediate

1403: .seealso: `ISRestoreTotalIndices()`, `ISGetNonlocalIndices()`
1404: @*/
1405: PetscErrorCode ISRestoreTotalIndices(IS is, const PetscInt *indices[])
1406: {
1407:   PetscMPIInt size;

1411:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
1412:   if (size == 1) {
1413:     PetscUseTypeMethod(is, restoreindices, indices);
1414:   } else {
1416:   }
1417:   return 0;
1418: }

1420: /*@C
1421:    ISGetNonlocalIndices - Retrieve an array of indices from remote processors
1422:                        in this communicator.

1424:    Collective

1426:    Input Parameter:
1427: .  is - the index set

1429:    Output Parameter:
1430: .  indices - indices with rank 0 indices first, and so on,  omitting
1431:              the current rank.  Total number of indices is the difference
1432:              total and local, obtained with ISGetSize() and ISGetLocalSize(),
1433:              respectively.

1435:    Level: intermediate

1437:    Notes:
1438:     restore the indices using ISRestoreNonlocalIndices().
1439:           The same scalability considerations as those for ISGetTotalIndices
1440:           apply here.

1442: .seealso: `ISGetTotalIndices()`, `ISRestoreNonlocalIndices()`, `ISGetSize()`, `ISGetLocalSize().`
1443: @*/
1444: PetscErrorCode ISGetNonlocalIndices(IS is, const PetscInt *indices[])
1445: {
1446:   PetscMPIInt size;
1447:   PetscInt    n, N;

1451:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
1452:   if (size == 1) *indices = NULL;
1453:   else {
1454:     if (!is->total) ISGatherTotal_Private(is);
1455:     ISGetLocalSize(is, &n);
1456:     ISGetSize(is, &N);
1457:     PetscMalloc1(N - n, &(is->nonlocal));
1458:     PetscArraycpy(is->nonlocal, is->total, is->local_offset);
1459:     PetscArraycpy(is->nonlocal + is->local_offset, is->total + is->local_offset + n, N - is->local_offset - n);
1460:     *indices = is->nonlocal;
1461:   }
1462:   return 0;
1463: }

1465: /*@C
1466:    ISRestoreNonlocalIndices - Restore the index array obtained with ISGetNonlocalIndices().

1468:    Not Collective.

1470:    Input Parameters:
1471: +  is - the index set
1472: -  indices - index array; must be the array obtained with ISGetNonlocalIndices()

1474:    Level: intermediate

1476: .seealso: `ISGetTotalIndices()`, `ISGetNonlocalIndices()`, `ISRestoreTotalIndices()`
1477: @*/
1478: PetscErrorCode ISRestoreNonlocalIndices(IS is, const PetscInt *indices[])
1479: {
1483:   return 0;
1484: }

1486: /*@
1487:    ISGetNonlocalIS - Gather all nonlocal indices for this IS and present
1488:                      them as another sequential index set.

1490:    Collective

1492:    Input Parameter:
1493: .  is - the index set

1495:    Output Parameter:
1496: .  complement - sequential IS with indices identical to the result of
1497:                 ISGetNonlocalIndices()

1499:    Level: intermediate

1501:    Notes:
1502:     complement represents the result of ISGetNonlocalIndices as an IS.
1503:           Therefore scalability issues similar to ISGetNonlocalIndices apply.
1504:           The resulting IS must be restored using ISRestoreNonlocalIS().

1506: .seealso: `ISGetNonlocalIndices()`, `ISRestoreNonlocalIndices()`, `ISAllGather()`, `ISGetSize()`
1507: @*/
1508: PetscErrorCode ISGetNonlocalIS(IS is, IS *complement)
1509: {
1512:   /* Check if the complement exists already. */
1513:   if (is->complement) {
1514:     *complement = is->complement;
1515:     PetscObjectReference((PetscObject)(is->complement));
1516:   } else {
1517:     PetscInt        N, n;
1518:     const PetscInt *idx;
1519:     ISGetSize(is, &N);
1520:     ISGetLocalSize(is, &n);
1521:     ISGetNonlocalIndices(is, &idx);
1522:     ISCreateGeneral(PETSC_COMM_SELF, N - n, idx, PETSC_USE_POINTER, &(is->complement));
1523:     PetscObjectReference((PetscObject)is->complement);
1524:     *complement = is->complement;
1525:   }
1526:   return 0;
1527: }

1529: /*@
1530:    ISRestoreNonlocalIS - Restore the IS obtained with ISGetNonlocalIS().

1532:    Not collective.

1534:    Input Parameters:
1535: +  is         - the index set
1536: -  complement - index set of is's nonlocal indices

1538:    Level: intermediate

1540: .seealso: `ISGetNonlocalIS()`, `ISGetNonlocalIndices()`, `ISRestoreNonlocalIndices()`
1541: @*/
1542: PetscErrorCode ISRestoreNonlocalIS(IS is, IS *complement)
1543: {
1544:   PetscInt refcnt;

1549:   PetscObjectGetReference((PetscObject)(is->complement), &refcnt);
1551:   PetscObjectDereference((PetscObject)(is->complement));
1552:   return 0;
1553: }

1555: /*@C
1556:    ISViewFromOptions - View from Options

1558:    Collective

1560:    Input Parameters:
1561: +  A - the index set
1562: .  obj - Optional object
1563: -  name - command line option

1565:    Level: intermediate
1566: .seealso: `IS`, `ISView`, `PetscObjectViewFromOptions()`, `ISCreate()`
1567: @*/
1568: PetscErrorCode ISViewFromOptions(IS A, PetscObject obj, const char name[])
1569: {
1571:   PetscObjectViewFromOptions((PetscObject)A, obj, name);
1572:   return 0;
1573: }

1575: /*@C
1576:    ISView - Displays an index set.

1578:    Collective

1580:    Input Parameters:
1581: +  is - the index set
1582: -  viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.

1584:    Level: intermediate

1586: .seealso: `PetscViewerASCIIOpen()`
1587: @*/
1588: PetscErrorCode ISView(IS is, PetscViewer viewer)
1589: {
1591:   if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)is), &viewer);

1595:   PetscObjectPrintClassNamePrefixType((PetscObject)is, viewer);
1596:   PetscLogEventBegin(IS_View, is, viewer, 0, 0);
1597:   PetscUseTypeMethod(is, view, viewer);
1598:   PetscLogEventEnd(IS_View, is, viewer, 0, 0);
1599:   return 0;
1600: }

1602: /*@
1603:   ISLoad - Loads a vector that has been stored in binary or HDF5 format with ISView().

1605:   Collective

1607:   Input Parameters:
1608: + is - the newly loaded vector, this needs to have been created with ISCreate() or some related function before a call to ISLoad().
1609: - viewer - binary file viewer, obtained from PetscViewerBinaryOpen() or HDF5 file viewer, obtained from PetscViewerHDF5Open()

1611:   Level: intermediate

1613:   Notes:
1614:   IF using HDF5, you must assign the IS the same name as was used in the IS
1615:   that was stored in the file using PetscObjectSetName(). Otherwise you will
1616:   get the error message: "Cannot H5DOpen2() with Vec name NAMEOFOBJECT"

1618: .seealso: `PetscViewerBinaryOpen()`, `ISView()`, `MatLoad()`, `VecLoad()`
1619: @*/
1620: PetscErrorCode ISLoad(IS is, PetscViewer viewer)
1621: {
1622:   PetscBool isbinary, ishdf5;

1627:   PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary);
1628:   PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERHDF5, &ishdf5);
1630:   if (!((PetscObject)is)->type_name) ISSetType(is, ISGENERAL);
1631:   PetscLogEventBegin(IS_Load, is, viewer, 0, 0);
1632:   PetscUseTypeMethod(is, load, viewer);
1633:   PetscLogEventEnd(IS_Load, is, viewer, 0, 0);
1634:   return 0;
1635: }

1637: /*@
1638:    ISSort - Sorts the indices of an index set.

1640:    Collective

1642:    Input Parameters:
1643: .  is - the index set

1645:    Level: intermediate

1647: .seealso: `ISSortRemoveDups()`, `ISSorted()`
1648: @*/
1649: PetscErrorCode ISSort(IS is)
1650: {
1652:   PetscUseTypeMethod(is, sort);
1653:   ISSetInfo(is, IS_SORTED, IS_LOCAL, is->info_permanent[IS_LOCAL][IS_SORTED], PETSC_TRUE);
1654:   return 0;
1655: }

1657: /*@
1658:   ISSortRemoveDups - Sorts the indices of an index set, removing duplicates.

1660:   Collective

1662:   Input Parameters:
1663: . is - the index set

1665:   Level: intermediate

1667: .seealso: `ISSort()`, `ISSorted()`
1668: @*/
1669: PetscErrorCode ISSortRemoveDups(IS is)
1670: {
1672:   ISClearInfoCache(is, PETSC_FALSE);
1673:   PetscUseTypeMethod(is, sortremovedups);
1674:   ISSetInfo(is, IS_SORTED, IS_LOCAL, is->info_permanent[IS_LOCAL][IS_SORTED], PETSC_TRUE);
1675:   ISSetInfo(is, IS_UNIQUE, IS_LOCAL, is->info_permanent[IS_LOCAL][IS_UNIQUE], PETSC_TRUE);
1676:   return 0;
1677: }

1679: /*@
1680:    ISToGeneral - Converts an IS object of any type to ISGENERAL type

1682:    Collective

1684:    Input Parameters:
1685: .  is - the index set

1687:    Level: intermediate

1689: .seealso: `ISSorted()`
1690: @*/
1691: PetscErrorCode ISToGeneral(IS is)
1692: {
1694:   PetscUseTypeMethod(is, togeneral);
1695:   return 0;
1696: }

1698: /*@
1699:    ISSorted - Checks the indices to determine whether they have been sorted.

1701:    Collective

1703:    Input Parameter:
1704: .  is - the index set

1706:    Output Parameter:
1707: .  flg - output flag, either PETSC_TRUE if the index set is sorted,
1708:          or PETSC_FALSE otherwise.

1710:    Notes:
1711:     For parallel IS objects this only indicates if the local part of the IS
1712:           is sorted. So some processors may return PETSC_TRUE while others may
1713:           return PETSC_FALSE.

1715:    Level: intermediate

1717: .seealso: `ISSort()`, `ISSortRemoveDups()`
1718: @*/
1719: PetscErrorCode ISSorted(IS is, PetscBool *flg)
1720: {
1723:   ISGetInfo(is, IS_SORTED, IS_LOCAL, PETSC_TRUE, flg);
1724:   return 0;
1725: }

1727: /*@
1728:    ISDuplicate - Creates a duplicate copy of an index set.

1730:    Collective

1732:    Input Parameter:
1733: .  is - the index set

1735:    Output Parameter:
1736: .  isnew - the copy of the index set

1738:    Level: beginner

1740: .seealso: `ISCreateGeneral()`, `ISCopy()`
1741: @*/
1742: PetscErrorCode ISDuplicate(IS is, IS *newIS)
1743: {
1746:   PetscUseTypeMethod(is, duplicate, newIS);
1747:   ISCopyInfo(is, *newIS);
1748:   return 0;
1749: }

1751: /*@
1752:    ISCopy - Copies an index set.

1754:    Collective

1756:    Input Parameter:
1757: .  is - the index set

1759:    Output Parameter:
1760: .  isy - the copy of the index set

1762:    Level: beginner

1764: .seealso: `ISDuplicate()`, `ISShift()`
1765: @*/
1766: PetscErrorCode ISCopy(IS is, IS isy)
1767: {
1768:   PetscInt bs, bsy;

1773:   if (is == isy) return 0;
1774:   PetscLayoutGetBlockSize(is->map, &bs);
1775:   PetscLayoutGetBlockSize(isy->map, &bsy);
1779:   ISCopyInfo(is, isy);
1780:   isy->max = is->max;
1781:   isy->min = is->min;
1782:   PetscUseTypeMethod(is, copy, isy);
1783:   return 0;
1784: }

1786: /*@
1787:    ISShift - Shift all indices by given offset

1789:    Collective

1791:    Input Parameters:
1792: +  is - the index set
1793: -  offset - the offset

1795:    Output Parameter:
1796: .  isy - the shifted copy of the input index set

1798:    Notes:
1799:    The offset can be different across processes.
1800:    IS is and isy can be the same.

1802:    Level: beginner

1804: .seealso: `ISDuplicate()`, `ISCopy()`
1805: @*/
1806: PetscErrorCode ISShift(IS is, PetscInt offset, IS isy)
1807: {
1811:   if (!offset) {
1812:     ISCopy(is, isy);
1813:     return 0;
1814:   }
1818:   ISCopyInfo(is, isy);
1819:   isy->max = is->max + offset;
1820:   isy->min = is->min + offset;
1821:   PetscUseMethod(is, "ISShift_C", (IS, PetscInt, IS), (is, offset, isy));
1822:   return 0;
1823: }

1825: /*@
1826:    ISOnComm - Split a parallel IS on subcomms (usually self) or concatenate index sets on subcomms into a parallel index set

1828:    Collective

1830:    Input Parameters:
1831: + is - index set
1832: . comm - communicator for new index set
1833: - mode - copy semantics, PETSC_USE_POINTER for no-copy if possible, otherwise PETSC_COPY_VALUES

1835:    Output Parameter:
1836: . newis - new IS on comm

1838:    Level: advanced

1840:    Notes:
1841:    It is usually desirable to create a parallel IS and look at the local part when necessary.

1843:    This function is useful if serial ISs must be created independently, or to view many
1844:    logically independent serial ISs.

1846:    The input IS must have the same type on every process.
1847: @*/
1848: PetscErrorCode ISOnComm(IS is, MPI_Comm comm, PetscCopyMode mode, IS *newis)
1849: {
1850:   PetscMPIInt match;

1854:   MPI_Comm_compare(PetscObjectComm((PetscObject)is), comm, &match);
1855:   if (mode != PETSC_COPY_VALUES && (match == MPI_IDENT || match == MPI_CONGRUENT)) {
1856:     PetscObjectReference((PetscObject)is);
1857:     *newis = is;
1858:   } else PetscUseTypeMethod(is, oncomm, comm, mode, newis);
1859:   return 0;
1860: }

1862: /*@
1863:    ISSetBlockSize - informs an index set that it has a given block size

1865:    Logicall Collective

1867:    Input Parameters:
1868: + is - index set
1869: - bs - block size

1871:    Level: intermediate

1873:    Notes:
1874:    This is much like the block size for Vecs. It indicates that one can think of the indices as
1875:    being in a collection of equal size blocks. For ISBlock() these collections of blocks are all contiquous
1876:    within a block but this is not the case for other IS.
1877:    ISBlockGetIndices() only works for ISBlock IS, not others.

1879: .seealso: `ISGetBlockSize()`, `ISCreateBlock()`, `ISBlockGetIndices()`,
1880: @*/
1881: PetscErrorCode ISSetBlockSize(IS is, PetscInt bs)
1882: {
1886:   if (PetscDefined(USE_DEBUG)) {
1887:     const PetscInt *indices;
1888:     PetscInt        length, i, j;
1889:     ISGetIndices(is, &indices);
1890:     if (indices) {
1891:       ISGetLocalSize(is, &length);
1893:       for (i = 0; i < length / bs; i += bs) {
1894:         for (j = 0; j < bs - 1; j++) {
1896:         }
1897:       }
1898:     }
1899:     ISRestoreIndices(is, &indices);
1900:   }
1901:   PetscUseTypeMethod(is, setblocksize, bs);
1902:   return 0;
1903: }

1905: /*@
1906:    ISGetBlockSize - Returns the number of elements in a block.

1908:    Not Collective

1910:    Input Parameter:
1911: .  is - the index set

1913:    Output Parameter:
1914: .  size - the number of elements in a block

1916:    Level: intermediate

1918: Notes:
1919:    This is much like the block size for Vecs. It indicates that one can think of the indices as
1920:    being in a collection of equal size blocks. For ISBlock() these collections of blocks are all contiquous
1921:    within a block but this is not the case for other IS.
1922:    ISBlockGetIndices() only works for ISBlock IS, not others.

1924: .seealso: `ISBlockGetSize()`, `ISGetSize()`, `ISCreateBlock()`, `ISSetBlockSize()`
1925: @*/
1926: PetscErrorCode ISGetBlockSize(IS is, PetscInt *size)
1927: {
1928:   PetscLayoutGetBlockSize(is->map, size);
1929:   return 0;
1930: }

1932: PetscErrorCode ISGetIndicesCopy(IS is, PetscInt idx[])
1933: {
1934:   PetscInt        len, i;
1935:   const PetscInt *ptr;

1937:   ISGetLocalSize(is, &len);
1938:   ISGetIndices(is, &ptr);
1939:   for (i = 0; i < len; i++) idx[i] = ptr[i];
1940:   ISRestoreIndices(is, &ptr);
1941:   return 0;
1942: }

1944: /*MC
1945:     ISGetIndicesF90 - Accesses the elements of an index set from Fortran90.
1946:     The users should call ISRestoreIndicesF90() after having looked at the
1947:     indices.  The user should NOT change the indices.

1949:     Synopsis:
1950:     ISGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

1952:     Not collective

1954:     Input Parameter:
1955: .   x - index set

1957:     Output Parameters:
1958: +   xx_v - the Fortran90 pointer to the array
1959: -   ierr - error code

1961:     Example of Usage:
1962: .vb
1963:     PetscInt, pointer xx_v(:)
1964:     ....
1965:     call ISGetIndicesF90(x,xx_v,ierr)
1966:     a = xx_v(3)
1967:     call ISRestoreIndicesF90(x,xx_v,ierr)
1968: .ve

1970:     Level: intermediate

1972: .seealso: `ISRestoreIndicesF90()`, `ISGetIndices()`, `ISRestoreIndices()`

1974: M*/

1976: /*MC
1977:     ISRestoreIndicesF90 - Restores an index set to a usable state after
1978:     a call to ISGetIndicesF90().

1980:     Synopsis:
1981:     ISRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

1983:     Not collective

1985:     Input Parameters:
1986: +   x - index set
1987: -   xx_v - the Fortran90 pointer to the array

1989:     Output Parameter:
1990: .   ierr - error code

1992:     Example of Usage:
1993: .vb
1994:     PetscInt, pointer xx_v(:)
1995:     ....
1996:     call ISGetIndicesF90(x,xx_v,ierr)
1997:     a = xx_v(3)
1998:     call ISRestoreIndicesF90(x,xx_v,ierr)
1999: .ve

2001:     Level: intermediate

2003: .seealso: `ISGetIndicesF90()`, `ISGetIndices()`, `ISRestoreIndices()`

2005: M*/

2007: /*MC
2008:     ISBlockGetIndicesF90 - Accesses the elements of an index set from Fortran90.
2009:     The users should call ISBlockRestoreIndicesF90() after having looked at the
2010:     indices.  The user should NOT change the indices.

2012:     Synopsis:
2013:     ISBlockGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

2015:     Not collective

2017:     Input Parameter:
2018: .   x - index set

2020:     Output Parameters:
2021: +   xx_v - the Fortran90 pointer to the array
2022: -   ierr - error code
2023:     Example of Usage:
2024: .vb
2025:     PetscInt, pointer xx_v(:)
2026:     ....
2027:     call ISBlockGetIndicesF90(x,xx_v,ierr)
2028:     a = xx_v(3)
2029:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
2030: .ve

2032:     Level: intermediate

2034: .seealso: `ISBlockRestoreIndicesF90()`, `ISGetIndices()`, `ISRestoreIndices()`,
2035:           `ISRestoreIndices()`

2037: M*/

2039: /*MC
2040:     ISBlockRestoreIndicesF90 - Restores an index set to a usable state after
2041:     a call to ISBlockGetIndicesF90().

2043:     Synopsis:
2044:     ISBlockRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

2046:     Not Collective

2048:     Input Parameters:
2049: +   x - index set
2050: -   xx_v - the Fortran90 pointer to the array

2052:     Output Parameter:
2053: .   ierr - error code

2055:     Example of Usage:
2056: .vb
2057:     PetscInt, pointer xx_v(:)
2058:     ....
2059:     call ISBlockGetIndicesF90(x,xx_v,ierr)
2060:     a = xx_v(3)
2061:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
2062: .ve

2064:     Notes:
2065:     Not yet supported for all F90 compilers

2067:     Level: intermediate

2069: .seealso: `ISBlockGetIndicesF90()`, `ISGetIndices()`, `ISRestoreIndices()`, `ISRestoreIndicesF90()`

2071: M*/