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 on IS

 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 - the maximum 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) { ilocalneg[nneg++] = i; continue; }
 64:     if (idxs[i] < lbounds[0]) lbounds[0] = idxs[i];
 65:     if (idxs[i] > lbounds[1]) lbounds[1] = idxs[i];
 66:     ilocal[npos++] = i;
 67:   }
 68:   if (npos == n) {
 69:     PetscFree(ilocal);
 70:     PetscFree(ilocalneg);
 71:   }

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

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

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

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

115:   /* reduce from leaves to roots */
116:   PetscCalloc1(Nl,&root_data);
117:   PetscSFReduceBegin(sf,MPIU_INT,leaf_data,root_data,MPI_MAX);
118:   PetscSFReduceEnd(sf,MPIU_INT,leaf_data,root_data,MPI_MAX);

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

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

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

147:   if (!subset_n) {
148:     PetscFree(gidxs);
149:     PetscSFDestroy(&sf);
150:     PetscFree(leaf_data);
151:     PetscFree(root_data);
152:     PetscFree(ilocal);
153:     if (subset_mult) ISRestoreIndices(subset_mult,&idxs_mult);
154:     return 0;
155:   }

157:   /* adapt root data with cumulative */
158:   if (first_found) {
159:     PetscInt old_index;

161:     root_data[first_index] += start;
162:     old_index = first_index;
163:     for (i=first_index+1;i<Nl;i++) {
164:       if (root_data[i]) {
165:         root_data[i] += root_data[old_index];
166:         old_index = i;
167:       }
168:     }
169:   }

171:   /* from roots to leaves */
172:   PetscSFBcastBegin(sf,MPIU_INT,root_data,leaf_data,MPI_REPLACE);
173:   PetscSFBcastEnd(sf,MPIU_INT,root_data,leaf_data,MPI_REPLACE);
174:   PetscSFDestroy(&sf);

176:   /* create new IS with global indexes without holes */
177:   for (i=0;i<n_n;i++) gidxs[i] = -1;
178:   if (subset_mult) {
179:     PetscInt cum;

181:     isblock = PETSC_FALSE;
182:     for (i=0,cum=0;i<n;i++) for (PetscInt j=0;j<idxs_mult[i];j++) gidxs[cum++] = leaf_data[i] - idxs_mult[i] + j;
183:   } else for (i=0;i<n;i++) gidxs[i] = leaf_data[i]-1;

185:   if (isblock) {
186:     if (ibs > 1) for (i=0;i<n_n/ibs;i++) gidxs[i] = gidxs[i*ibs]/ibs;
187:     ISCreateBlock(PetscObjectComm((PetscObject)subset),ibs,n_n/ibs,gidxs,PETSC_COPY_VALUES,subset_n);
188:   } else {
189:     ISCreateGeneral(PetscObjectComm((PetscObject)subset),n_n,gidxs,PETSC_COPY_VALUES,subset_n);
190:   }
191:   if (subset_mult) ISRestoreIndices(subset_mult,&idxs_mult);
192:   PetscFree(gidxs);
193:   PetscFree(leaf_data);
194:   PetscFree(root_data);
195:   PetscFree(ilocal);
196:   return 0;
197: }

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

202:    Collective on IS

204:    Input Parameters:
205: +  is - the index set
206: -  comps - which components we will extract from is

208:    Output Parameters:
209: .  subis - the new sub index set

211:    Level: intermediate

213:    Example usage:
214:    We have an index set (is) living on 3 processes with the following values:
215:    | 4 9 0 | 2 6 7 | 10 11 1|
216:    and another index set (comps) used to indicate which components of is  we want to take,
217:    | 7 5  | 1 2 | 0 4|
218:    The output index set (subis) should look like:
219:    | 11 7 | 9 0 | 4 6|

221: .seealso: VecGetSubVector(), MatCreateSubMatrix()
222: @*/
223: PetscErrorCode ISCreateSubIS(IS is,IS comps,IS *subis)
224: {
225:   PetscSF         sf;
226:   const PetscInt  *is_indices,*comps_indices;
227:   PetscInt        *subis_indices,nroots,nleaves,*mine,i,lidx;
228:   PetscMPIInt     owner;
229:   PetscSFNode     *remote;
230:   MPI_Comm        comm;


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

262:   PetscMalloc1(nleaves,&subis_indices);
263:   ISGetIndices(is, &is_indices);
264:   PetscSFBcastBegin(sf,MPIU_INT,is_indices,subis_indices,MPI_REPLACE);
265:   PetscSFBcastEnd(sf,MPIU_INT,is_indices,subis_indices,MPI_REPLACE);
266:   ISRestoreIndices(is,&is_indices);
267:   PetscSFDestroy(&sf);
268:   ISCreateGeneral(comm,nleaves,subis_indices,PETSC_OWN_POINTER,subis);
269:   return 0;
270: }

272: /*@
273:    ISClearInfoCache - clear the cache of computed index set properties

275:    Not collective

277:    Input Parameters:
278: +  is - the index set
279: -  clear_permanent_local - whether to remove the permanent status of local properties

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

284:    Level: developer

286: .seealso:  ISInfo, ISInfoType, ISSetInfo(), ISClearInfoCache()

288: @*/
289: PetscErrorCode ISClearInfoCache(IS is, PetscBool clear_permanent_local)
290: {
291:   PetscInt i, j;

295:   for (i = 0; i < IS_INFO_MAX; i++) {
296:     if (clear_permanent_local) is->info_permanent[IS_LOCAL][i] = PETSC_FALSE;
297:     for (j = 0; j < 2; j++) {
298:       if (!is->info_permanent[j][i]) is->info[j][i] = IS_INFO_UNKNOWN;
299:     }
300:   }
301:   return 0;
302: }

304: static PetscErrorCode ISSetInfo_Internal(IS is, ISInfo info, ISInfoType type, ISInfoBool ipermanent, PetscBool flg)
305: {
306:   ISInfoBool     iflg = flg ? IS_INFO_TRUE : IS_INFO_FALSE;
307:   PetscInt       itype = (type == IS_LOCAL) ? 0 : 1;
308:   PetscBool      permanent_set = (ipermanent == IS_INFO_UNKNOWN) ? PETSC_FALSE : PETSC_TRUE;
309:   PetscBool      permanent = (ipermanent == IS_INFO_TRUE) ? PETSC_TRUE : PETSC_FALSE;

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

405: /*@
406:    ISSetInfo - Set known information about an index set.

408:    Logically Collective on IS if type is IS_GLOBAL

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

419:   Info Describing IS Structure:
420: +    IS_SORTED - the [local part of the] index set is sorted in ascending order
421: .    IS_UNIQUE - each entry in the [local part of the] index set is unique
422: .    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
423: .    IS_INTERVAL - the [local part of the] index set is equal to a contiguous range of integers {f, f + 1, ..., f + N-1}
424: -    IS_IDENTITY - the [local part of the] index set is equal to the integers {0, 1, ..., N-1}

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

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

431:    Level: advanced

433: .seealso:  ISInfo, ISInfoType, IS

435: @*/
436: PetscErrorCode ISSetInfo(IS is, ISInfo info, ISInfoType type, PetscBool permanent, PetscBool flg)
437: {
438:   MPI_Comm       comm, errcomm;
439:   PetscMPIInt    size;

443:   comm = PetscObjectComm((PetscObject)is);
444:   if (type == IS_GLOBAL) {
448:     errcomm = comm;
449:   } else {
450:     errcomm = PETSC_COMM_SELF;
451:   }


455:   MPI_Comm_size(comm, &size);
456:   /* do not use global values if size == 1: it makes it easier to keep the implications straight */
457:   if (size == 1) type = IS_LOCAL;
458:   ISSetInfo_Internal(is, info, type, permanent ? IS_INFO_TRUE : IS_INFO_FALSE, flg);
459:   return 0;
460: }

462: static PetscErrorCode ISGetInfo_Sorted(IS is, ISInfoType type, PetscBool *flg)
463: {
464:   MPI_Comm       comm;
465:   PetscMPIInt    size, rank;

467:   comm = PetscObjectComm((PetscObject)is);
468:   MPI_Comm_size(comm, &size);
469:   MPI_Comm_size(comm, &rank);
470:   if (type == IS_GLOBAL && is->ops->sortedglobal) {
471:     (*is->ops->sortedglobal)(is,flg);
472:   } else {
473:     PetscBool sortedLocal = PETSC_FALSE;

475:     /* determine if the array is locally sorted */
476:     if (type == IS_GLOBAL && size > 1) {
477:       /* call ISGetInfo so that a cached value will be used if possible */
478:       ISGetInfo(is, IS_SORTED, IS_LOCAL, PETSC_TRUE, &sortedLocal);
479:     } else if (is->ops->sortedlocal) {
480:       (*is->ops->sortedlocal)(is,&sortedLocal);
481:     } else {
482:       /* default: get the local indices and directly check */
483:       const PetscInt *idx;
484:       PetscInt n;

486:       ISGetIndices(is, &idx);
487:       ISGetLocalSize(is, &n);
488:       PetscSortedInt(n, idx, &sortedLocal);
489:       ISRestoreIndices(is, &idx);
490:     }

492:     if (type == IS_LOCAL || size == 1) {
493:       *flg = sortedLocal;
494:     } else {
495:       MPI_Allreduce(&sortedLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
496:       if (*flg) {
497:         PetscInt  n, min = PETSC_MAX_INT, max = PETSC_MIN_INT;
498:         PetscInt  maxprev;

500:         ISGetLocalSize(is, &n);
501:         if (n) ISGetMinMax(is, &min, &max);
502:         maxprev = PETSC_MIN_INT;
503:         MPI_Exscan(&max, &maxprev, 1, MPIU_INT, MPI_MAX, comm);
504:         if (rank && (maxprev > min)) sortedLocal = PETSC_FALSE;
505:         MPI_Allreduce(&sortedLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
506:       }
507:     }
508:   }
509:   return 0;
510: }

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

514: static PetscErrorCode ISGetInfo_Unique(IS is, ISInfoType type, PetscBool *flg)
515: {
516:   MPI_Comm       comm;
517:   PetscMPIInt    size, rank;
518:   PetscInt       i;

520:   comm = PetscObjectComm((PetscObject)is);
521:   MPI_Comm_size(comm, &size);
522:   MPI_Comm_size(comm, &rank);
523:   if (type == IS_GLOBAL && is->ops->uniqueglobal) {
524:     (*is->ops->uniqueglobal)(is,flg);
525:   } else {
526:     PetscBool uniqueLocal;
527:     PetscInt  n = -1;
528:     PetscInt  *idx = NULL;

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

547:     PetscFree(idx);
548:     if (type == IS_LOCAL || size == 1) {
549:       *flg = uniqueLocal;
550:     } else {
551:       MPI_Allreduce(&uniqueLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
552:       if (*flg) {
553:         PetscInt  min = PETSC_MAX_INT, max = PETSC_MIN_INT, maxprev;

555:         if (!idx) {
556:           ISGetLocalSize(is, &n);
557:           PetscMalloc1(n, &idx);
558:           ISGetIndicesCopy(is, idx);
559:         }
560:         PetscParallelSortInt(is->map, is->map, idx, idx);
561:         if (n) {
562:           min = idx[0];
563:           max = idx[n - 1];
564:         }
565:         for (i = 1; i < n; i++) if (idx[i] == idx[i-1]) break;
566:         if (i < n) uniqueLocal = PETSC_FALSE;
567:         maxprev = PETSC_MIN_INT;
568:         MPI_Exscan(&max, &maxprev, 1, MPIU_INT, MPI_MAX, comm);
569:         if (rank && (maxprev == min)) uniqueLocal = PETSC_FALSE;
570:         MPI_Allreduce(&uniqueLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
571:       }
572:     }
573:     PetscFree(idx);
574:   }
575:   return 0;
576: }

578: static PetscErrorCode ISGetInfo_Permutation(IS is, ISInfoType type, PetscBool *flg)
579: {
580:   MPI_Comm       comm;
581:   PetscMPIInt    size, rank;

583:   comm = PetscObjectComm((PetscObject)is);
584:   MPI_Comm_size(comm, &size);
585:   MPI_Comm_size(comm, &rank);
586:   if (type == IS_GLOBAL && is->ops->permglobal) {
587:     (*is->ops->permglobal)(is,flg);
588:   } else if (type == IS_LOCAL && is->ops->permlocal) {
589:     (*is->ops->permlocal)(is,flg);
590:   } else {
591:     PetscBool permLocal;
592:     PetscInt  n, i, rStart;
593:     PetscInt  *idx;

595:     ISGetLocalSize(is, &n);
596:     PetscMalloc1(n, &idx);
597:     ISGetIndicesCopy(is, idx);
598:     if (type == IS_GLOBAL) {
599:       PetscParallelSortInt(is->map, is->map, idx, idx);
600:       PetscLayoutGetRange(is->map, &rStart, NULL);
601:     } else {
602:       PetscIntSortSemiOrdered(n, idx);
603:       rStart = 0;
604:     }
605:     permLocal = PETSC_TRUE;
606:     for (i = 0; i < n; i++) {
607:       if (idx[i] != rStart + i) break;
608:     }
609:     if (i < n) permLocal = PETSC_FALSE;
610:     if (type == IS_LOCAL || size == 1) {
611:       *flg = permLocal;
612:     } else {
613:       MPI_Allreduce(&permLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
614:     }
615:     PetscFree(idx);
616:   }
617:   return 0;
618: }

620: static PetscErrorCode ISGetInfo_Interval(IS is, ISInfoType type, PetscBool *flg)
621: {
622:   MPI_Comm       comm;
623:   PetscMPIInt    size, rank;
624:   PetscInt       i;

626:   comm = PetscObjectComm((PetscObject)is);
627:   MPI_Comm_size(comm, &size);
628:   MPI_Comm_size(comm, &rank);
629:   if (type == IS_GLOBAL && is->ops->intervalglobal) {
630:     (*is->ops->intervalglobal)(is,flg);
631:   } else {
632:     PetscBool intervalLocal;

634:     /* determine if the array is locally an interval */
635:     if (type == IS_GLOBAL && size > 1) {
636:       /* call ISGetInfo so that a cached value will be used if possible */
637:       ISGetInfo(is, IS_INTERVAL, IS_LOCAL, PETSC_TRUE, &intervalLocal);
638:     } else if (is->ops->intervallocal) {
639:       (*is->ops->intervallocal)(is,&intervalLocal);
640:     } else {
641:       PetscInt        n;
642:       const PetscInt  *idx;
643:       /* default: get the local indices and directly check */
644:       intervalLocal = PETSC_TRUE;
645:       ISGetLocalSize(is, &n);
646:       ISGetIndices(is, &idx);
647:       for (i = 1; i < n; i++) if (idx[i] != idx[i-1] + 1) break;
648:       if (i < n) intervalLocal = PETSC_FALSE;
649:       ISRestoreIndices(is, &idx);
650:     }

652:     if (type == IS_LOCAL || size == 1) {
653:       *flg = intervalLocal;
654:     } else {
655:       MPI_Allreduce(&intervalLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
656:       if (*flg) {
657:         PetscInt  n, min = PETSC_MAX_INT, max = PETSC_MIN_INT;
658:         PetscInt  maxprev;

660:         ISGetLocalSize(is, &n);
661:         if (n) ISGetMinMax(is, &min, &max);
662:         maxprev = PETSC_MIN_INT;
663:         MPI_Exscan(&max, &maxprev, 1, MPIU_INT, MPI_MAX, comm);
664:         if (rank && n && (maxprev != min - 1)) intervalLocal = PETSC_FALSE;
665:         MPI_Allreduce(&intervalLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
666:       }
667:     }
668:   }
669:   return 0;
670: }

672: static PetscErrorCode ISGetInfo_Identity(IS is, ISInfoType type, PetscBool *flg)
673: {
674:   MPI_Comm       comm;
675:   PetscMPIInt    size, rank;

677:   comm = PetscObjectComm((PetscObject)is);
678:   MPI_Comm_size(comm, &size);
679:   MPI_Comm_size(comm, &rank);
680:   if (type == IS_GLOBAL && is->ops->intervalglobal) {
681:     PetscBool isinterval;

683:     (*is->ops->intervalglobal)(is,&isinterval);
684:     *flg = PETSC_FALSE;
685:     if (isinterval) {
686:       PetscInt  min;

688:       ISGetMinMax(is, &min, NULL);
689:       MPI_Bcast(&min, 1, MPIU_INT, 0, comm);
690:       if (min == 0) *flg = PETSC_TRUE;
691:     }
692:   } else if (type == IS_LOCAL && is->ops->intervallocal) {
693:     PetscBool isinterval;

695:     (*is->ops->intervallocal)(is,&isinterval);
696:     *flg = PETSC_FALSE;
697:     if (isinterval) {
698:       PetscInt  min;

700:       ISGetMinMax(is, &min, NULL);
701:       if (min == 0) *flg = PETSC_TRUE;
702:     }
703:   } else {
704:     PetscBool identLocal;
705:     PetscInt  n, i, rStart;
706:     const PetscInt *idx;

708:     ISGetLocalSize(is, &n);
709:     ISGetIndices(is, &idx);
710:     PetscLayoutGetRange(is->map, &rStart, NULL);
711:     identLocal = PETSC_TRUE;
712:     for (i = 0; i < n; i++) {
713:       if (idx[i] != rStart + i) break;
714:     }
715:     if (i < n) identLocal = PETSC_FALSE;
716:     if (type == IS_LOCAL || size == 1) {
717:       *flg = identLocal;
718:     } else {
719:       MPI_Allreduce(&identLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
720:     }
721:     ISRestoreIndices(is, &idx);
722:   }
723:   return 0;
724: }

726: /*@
727:    ISGetInfo - Determine whether an index set satisfies a given property

729:    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())

731:    Input Parameters:
732: +  is - the index set
733: .  info - describing a property of the index set, one of those listed in the documentation of ISSetInfo()
734: .  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
735: -  type - whether the property is local (IS_LOCAL) or global (IS_GLOBAL)

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

740:    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().

742:    Level: advanced

744: .seealso:  ISInfo, ISInfoType, ISSetInfo(), ISClearInfoCache()

746: @*/
747: PetscErrorCode ISGetInfo(IS is, ISInfo info, ISInfoType type, PetscBool compute, PetscBool *flg)
748: {
749:   MPI_Comm       comm, errcomm;
750:   PetscMPIInt    rank, size;
751:   PetscInt       itype;
752:   PetscBool      hasprop;
753:   PetscBool      infer;

757:   comm = PetscObjectComm((PetscObject)is);
758:   if (type == IS_GLOBAL) {
760:     errcomm = comm;
761:   } else {
762:     errcomm = PETSC_COMM_SELF;
763:   }

765:   MPI_Comm_size(comm, &size);
766:   MPI_Comm_rank(comm, &rank);

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

810: static PetscErrorCode ISCopyInfo(IS source, IS dest)
811: {
812:   PetscArraycpy(&dest->info[0], &source->info[0], 2);
813:   PetscArraycpy(&dest->info_permanent[0], &source->info_permanent[0], 2);
814:   return 0;
815: }

817: /*@
818:    ISIdentity - Determines whether index set is the identity mapping.

820:    Collective on IS

822:    Input Parameters:
823: .  is - the index set

825:    Output Parameters:
826: .  ident - PETSC_TRUE if an identity, else PETSC_FALSE

828:    Level: intermediate

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

836: .seealso: ISSetIdentity(), ISGetInfo()
837: @*/
838: PetscErrorCode  ISIdentity(IS is,PetscBool  *ident)
839: {
842:   ISGetInfo(is,IS_IDENTITY,IS_GLOBAL,PETSC_TRUE,ident);
843:   return 0;
844: }

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

849:    Logically Collective on IS

851:    Input Parameter:
852: .  is - the index set

854:    Level: intermediate

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

860: .seealso: ISIdentity(), ISSetInfo(), ISClearInfoCache()
861: @*/
862: PetscErrorCode  ISSetIdentity(IS is)
863: {
865:   ISSetInfo(is,IS_IDENTITY,IS_GLOBAL,PETSC_TRUE,PETSC_TRUE);
866:   return 0;
867: }

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

872:    Not Collective

874:    Input Parameters:
875: +  is - the index set
876: .  gstart - global start
877: -  gend - global end

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

883:    Level: developer

885: .seealso: ISGetLocalSize(), VecGetOwnershipRange()
886: @*/
887: PetscErrorCode  ISContiguousLocal(IS is,PetscInt gstart,PetscInt gend,PetscInt *start,PetscBool *contig)
888: {
892:   *start  = -1;
893:   *contig = PETSC_FALSE;
894:   if (is->ops->contiguous) {
895:     (*is->ops->contiguous)(is,gstart,gend,start,contig);
896:   }
897:   return 0;
898: }

900: /*@
901:    ISPermutation - PETSC_TRUE or PETSC_FALSE depending on whether the
902:    index set has been declared to be a permutation.

904:    Logically Collective on IS

906:    Input Parameter:
907: .  is - the index set

909:    Output Parameter:
910: .  perm - PETSC_TRUE if a permutation, else PETSC_FALSE

912:    Level: intermediate

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

920: .seealso: ISSetPermutation(), ISGetInfo()
921: @*/
922: PetscErrorCode  ISPermutation(IS is,PetscBool  *perm)
923: {
926:   ISGetInfo(is,IS_PERMUTATION,IS_GLOBAL,PETSC_FALSE,perm);
927:   return 0;
928: }

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

933:    Logically Collective on IS

935:    Input Parameter:
936: .  is - the index set

938:    Level: intermediate

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

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

947: .seealso: ISPermutation(), ISSetInfo(), ISClearInfoCache().
948: @*/
949: PetscErrorCode  ISSetPermutation(IS is)
950: {
952:   if (PetscDefined(USE_DEBUG)) {
953:     PetscMPIInt    size;

955:     MPI_Comm_size(PetscObjectComm((PetscObject)is),&size);
956:     if (size == 1) {
957:       PetscInt       i,n,*idx;
958:       const PetscInt *iidx;

960:       ISGetSize(is,&n);
961:       PetscMalloc1(n,&idx);
962:       ISGetIndices(is,&iidx);
963:       PetscArraycpy(idx,iidx,n);
964:       PetscIntSortSemiOrdered(n,idx);
965:       for (i=0; i<n; i++) {
967:       }
968:       PetscFree(idx);
969:       ISRestoreIndices(is,&iidx);
970:     }
971:   }
972:   ISSetInfo(is,IS_PERMUTATION,IS_GLOBAL,PETSC_TRUE,PETSC_TRUE);
973:   return 0;
974: }

976: /*@C
977:    ISDestroy - Destroys an index set.

979:    Collective on IS

981:    Input Parameters:
982: .  is - the index set

984:    Level: beginner

986: .seealso: ISCreateGeneral(), ISCreateStride(), ISCreateBlocked()
987: @*/
988: PetscErrorCode  ISDestroy(IS *is)
989: {
990:   if (!*is) return 0;
992:   if (--((PetscObject)(*is))->refct > 0) {*is = NULL; return 0;}
993:   if ((*is)->complement) {
994:     PetscInt refcnt;
995:     PetscObjectGetReference((PetscObject)((*is)->complement), &refcnt);
997:     ISDestroy(&(*is)->complement);
998:   }
999:   if ((*is)->ops->destroy) {
1000:     (*(*is)->ops->destroy)(*is);
1001:   }
1002:   PetscLayoutDestroy(&(*is)->map);
1003:   /* Destroy local representations of offproc data. */
1004:   PetscFree((*is)->total);
1005:   PetscFree((*is)->nonlocal);
1006:   PetscHeaderDestroy(is);
1007:   return 0;
1008: }

1010: /*@
1011:    ISInvertPermutation - Creates a new permutation that is the inverse of
1012:                          a given permutation.

1014:    Collective on IS

1016:    Input Parameters:
1017: +  is - the index set
1018: -  nlocal - number of indices on this processor in result (ignored for 1 proccessor) or
1019:             use PETSC_DECIDE

1021:    Output Parameter:
1022: .  isout - the inverse permutation

1024:    Level: intermediate

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

1030: @*/
1031: PetscErrorCode  ISInvertPermutation(IS is,PetscInt nlocal,IS *isout)
1032: {
1033:   PetscBool      isperm, isidentity, issame;

1037:   ISGetInfo(is,IS_PERMUTATION,IS_GLOBAL,PETSC_TRUE,&isperm);
1039:   ISGetInfo(is,IS_IDENTITY,IS_GLOBAL,PETSC_TRUE,&isidentity);
1040:   issame = PETSC_FALSE;
1041:   if (isidentity) {
1042:     PetscInt n;
1043:     PetscBool isallsame;

1045:     ISGetLocalSize(is, &n);
1046:     issame = (PetscBool) (n == nlocal);
1047:     MPI_Allreduce(&issame, &isallsame, 1, MPIU_BOOL, MPI_LAND, PetscObjectComm((PetscObject)is));
1048:     issame = isallsame;
1049:   }
1050:   if (issame) {
1051:     ISDuplicate(is,isout);
1052:   } else {
1053:     (*is->ops->invertpermutation)(is,nlocal,isout);
1054:     ISSetPermutation(*isout);
1055:   }
1056:   return 0;
1057: }

1059: /*@
1060:    ISGetSize - Returns the global length of an index set.

1062:    Not Collective

1064:    Input Parameter:
1065: .  is - the index set

1067:    Output Parameter:
1068: .  size - the global size

1070:    Level: beginner

1072: @*/
1073: PetscErrorCode  ISGetSize(IS is,PetscInt *size)
1074: {
1077:   *size = is->map->N;
1078:   return 0;
1079: }

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

1084:    Not Collective

1086:    Input Parameter:
1087: .  is - the index set

1089:    Output Parameter:
1090: .  size - the local size

1092:    Level: beginner

1094: @*/
1095: PetscErrorCode  ISGetLocalSize(IS is,PetscInt *size)
1096: {
1099:   *size = is->map->n;
1100:   return 0;
1101: }

1103: /*@
1104:    ISGetLayout - get PetscLayout describing index set layout

1106:    Not Collective

1108:    Input Parameter:
1109: .  is - the index set

1111:    Output Parameter:
1112: .  map - the layout

1114:    Level: developer

1116: .seealso: ISSetLayout(), ISGetSize(), ISGetLocalSize()
1117: @*/
1118: PetscErrorCode ISGetLayout(IS is,PetscLayout *map)
1119: {
1122:   *map = is->map;
1123:   return 0;
1124: }

1126: /*@
1127:    ISSetLayout - set PetscLayout describing index set layout

1129:    Collective

1131:    Input Arguments:
1132: +  is - the index set
1133: -  map - the layout

1135:    Level: developer

1137:    Notes:
1138:    Users should typically use higher level functions such as ISCreateGeneral().

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

1143: .seealso: ISCreate(), ISGetLayout(), ISGetSize(), ISGetLocalSize()
1144: @*/
1145: PetscErrorCode ISSetLayout(IS is,PetscLayout map)
1146: {
1149:   PetscLayoutReference(map,&is->map);
1150:   return 0;
1151: }

1153: /*@C
1154:    ISGetIndices - Returns a pointer to the indices.  The user should call
1155:    ISRestoreIndices() after having looked at the indices.  The user should
1156:    NOT change the indices.

1158:    Not Collective

1160:    Input Parameter:
1161: .  is - the index set

1163:    Output Parameter:
1164: .  ptr - the location to put the pointer to the indices

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

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

1189:    Level: intermediate

1191: .seealso: ISRestoreIndices(), ISGetIndicesF90()
1192: @*/
1193: PetscErrorCode  ISGetIndices(IS is,const PetscInt *ptr[])
1194: {
1197:   (*is->ops->getindices)(is,ptr);
1198:   return 0;
1199: }

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

1204:    Not Collective

1206:    Input Parameter:
1207: .  is - the index set

1209:    Output Parameters:
1210: +   min - the minimum value
1211: -   max - the maximum value

1213:    Level: intermediate

1215:    Notes:
1216:     Empty index sets return min=PETSC_MAX_INT and max=PETSC_MIN_INT.
1217:     In parallel, it returns the min and max of the local portion of the IS

1219: .seealso: ISGetIndices(), ISRestoreIndices(), ISGetIndicesF90()
1220: @*/
1221: PetscErrorCode  ISGetMinMax(IS is,PetscInt *min,PetscInt *max)
1222: {
1224:   if (min) *min = is->min;
1225:   if (max) *max = is->max;
1226:   return 0;
1227: }

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

1232:   Not Collective

1234:   Input Parameters:
1235: + is - the index set
1236: - key - the search key

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

1241:   Level: intermediate
1242: @*/
1243: PetscErrorCode ISLocate(IS is, PetscInt key, PetscInt *location)
1244: {
1245:   if (is->ops->locate) {
1246:     (*is->ops->locate)(is,key,location);
1247:   } else {
1248:     PetscInt       numIdx;
1249:     PetscBool      sorted;
1250:     const PetscInt *idx;

1252:     ISGetLocalSize(is,&numIdx);
1253:     ISGetIndices(is,&idx);
1254:     ISSorted(is,&sorted);
1255:     if (sorted) {
1256:       PetscFindInt(key,numIdx,idx,location);
1257:     } else {
1258:       PetscInt i;

1260:       *location = -1;
1261:       for (i = 0; i < numIdx; i++) {
1262:         if (idx[i] == key) {
1263:           *location = i;
1264:           break;
1265:         }
1266:       }
1267:     }
1268:     ISRestoreIndices(is,&idx);
1269:   }
1270:   return 0;
1271: }

1273: /*@C
1274:    ISRestoreIndices - Restores an index set to a usable state after a call
1275:                       to ISGetIndices().

1277:    Not Collective

1279:    Input Parameters:
1280: +  is - the index set
1281: -  ptr - the pointer obtained by ISGetIndices()

1283:    Fortran Note:
1284:    This routine is used differently from Fortran
1285: $    IS          is
1286: $    integer     is_array(1)
1287: $    PetscOffset i_is
1288: $    int         ierr
1289: $       call ISGetIndices(is,is_array,i_is,ierr)
1290: $
1291: $   Access first local entry in list
1292: $      value = is_array(i_is + 1)
1293: $
1294: $      ...... other code
1295: $       call ISRestoreIndices(is,is_array,i_is,ierr)

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

1300:    Level: intermediate

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

1305: .seealso: ISGetIndices(), ISRestoreIndicesF90()
1306: @*/
1307: PetscErrorCode  ISRestoreIndices(IS is,const PetscInt *ptr[])
1308: {
1311:   if (is->ops->restoreindices) {
1312:     (*is->ops->restoreindices)(is,ptr);
1313:   }
1314:   return 0;
1315: }

1317: static PetscErrorCode ISGatherTotal_Private(IS is)
1318: {
1319:   PetscInt       i,n,N;
1320:   const PetscInt *lindices;
1321:   MPI_Comm       comm;
1322:   PetscMPIInt    rank,size,*sizes = NULL,*offsets = NULL,nn;


1326:   PetscObjectGetComm((PetscObject)is,&comm);
1327:   MPI_Comm_size(comm,&size);
1328:   MPI_Comm_rank(comm,&rank);
1329:   ISGetLocalSize(is,&n);
1330:   PetscMalloc2(size,&sizes,size,&offsets);

1332:   PetscMPIIntCast(n,&nn);
1333:   MPI_Allgather(&nn,1,MPI_INT,sizes,1,MPI_INT,comm);
1334:   offsets[0] = 0;
1335:   for (i=1; i<size; ++i) offsets[i] = offsets[i-1] + sizes[i-1];
1336:   N = offsets[size-1] + sizes[size-1];

1338:   PetscMalloc1(N,&(is->total));
1339:   ISGetIndices(is,&lindices);
1340:   MPI_Allgatherv((void*)lindices,nn,MPIU_INT,is->total,sizes,offsets,MPIU_INT,comm);
1341:   ISRestoreIndices(is,&lindices);
1342:   is->local_offset = offsets[rank];
1343:   PetscFree2(sizes,offsets);
1344:   return 0;
1345: }

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

1350:    Collective on IS

1352:    Input Parameter:
1353: .  is - the index set

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

1359:    Level: intermediate

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

1369: .seealso: ISRestoreTotalIndices(), ISGetNonlocalIndices(), ISGetSize()
1370: @*/
1371: PetscErrorCode ISGetTotalIndices(IS is, const PetscInt *indices[])
1372: {
1373:   PetscMPIInt    size;

1377:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
1378:   if (size == 1) {
1379:     (*is->ops->getindices)(is,indices);
1380:   } else {
1381:     if (!is->total) {
1382:       ISGatherTotal_Private(is);
1383:     }
1384:     *indices = is->total;
1385:   }
1386:   return 0;
1387: }

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

1392:    Not Collective.

1394:    Input Parameters:
1395: +  is - the index set
1396: -  indices - index array; must be the array obtained with ISGetTotalIndices()

1398:    Level: intermediate

1400: .seealso: ISRestoreTotalIndices(), ISGetNonlocalIndices()
1401: @*/
1402: PetscErrorCode  ISRestoreTotalIndices(IS is, const PetscInt *indices[])
1403: {
1404:   PetscMPIInt    size;

1408:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
1409:   if (size == 1) {
1410:     (*is->ops->restoreindices)(is,indices);
1411:   } else {
1413:   }
1414:   return 0;
1415: }

1417: /*@C
1418:    ISGetNonlocalIndices - Retrieve an array of indices from remote processors
1419:                        in this communicator.

1421:    Collective on IS

1423:    Input Parameter:
1424: .  is - the index set

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

1432:    Level: intermediate

1434:    Notes:
1435:     restore the indices using ISRestoreNonlocalIndices().
1436:           The same scalability considerations as those for ISGetTotalIndices
1437:           apply here.

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

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

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

1467:    Not Collective.

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

1473:    Level: intermediate

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

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

1489:    Collective on IS

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

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

1498:    Level: intermediate

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

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

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

1531:    Not collective.

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

1537:    Level: intermediate

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

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

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

1557:    Collective on IS

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

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

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

1577:    Collective on IS

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

1583:    Level: intermediate

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

1594:   PetscObjectPrintClassNamePrefixType((PetscObject)is,viewer);
1595:   PetscLogEventBegin(IS_View,is,viewer,0,0);
1596:   (*is->ops->view)(is,viewer);
1597:   PetscLogEventEnd(IS_View,is,viewer,0,0);
1598:   return 0;
1599: }

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

1604:   Collective on PetscViewer

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

1610:   Level: intermediate

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

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

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

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

1639:    Collective on IS

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

1644:    Level: intermediate

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

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

1659:   Collective on IS

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

1664:   Level: intermediate

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

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

1681:    Collective on IS

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

1686:    Level: intermediate

1688: .seealso: ISSorted()
1689: @*/
1690: PetscErrorCode  ISToGeneral(IS is)
1691: {
1693:   if (is->ops->togeneral) {
1694:     (*is->ops->togeneral)(is);
1695:   } else SETERRQ(PetscObjectComm((PetscObject)is),PETSC_ERR_SUP,"Not written for this type %s",((PetscObject)is)->type_name);
1696:   return 0;
1697: }

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

1702:    Collective on IS

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

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

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

1716:    Level: intermediate

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

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

1731:    Collective on IS

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

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

1739:    Level: beginner

1741: .seealso: ISCreateGeneral(), ISCopy()
1742: @*/
1743: PetscErrorCode  ISDuplicate(IS is,IS *newIS)
1744: {
1747:   (*is->ops->duplicate)(is,newIS);
1748:   ISCopyInfo(is,*newIS);
1749:   return 0;
1750: }

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

1755:    Collective on IS

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

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

1763:    Level: beginner

1765: .seealso: ISDuplicate()
1766: @*/
1767: PetscErrorCode  ISCopy(IS is,IS isy)
1768: {
1772:   if (is == isy) return 0;
1773:   ISCopyInfo(is,isy);
1774:   isy->max        = is->max;
1775:   isy->min        = is->min;
1776:   (*is->ops->copy)(is,isy);
1777:   return 0;
1778: }

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

1783:    Collective on IS

1785:    Input Parameters:
1786: + is - index set
1787: . comm - communicator for new index set
1788: - mode - copy semantics, PETSC_USE_POINTER for no-copy if possible, otherwise PETSC_COPY_VALUES

1790:    Output Parameter:
1791: . newis - new IS on comm

1793:    Level: advanced

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

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

1801:    The input IS must have the same type on every process.
1802: @*/
1803: PetscErrorCode  ISOnComm(IS is,MPI_Comm comm,PetscCopyMode mode,IS *newis)
1804: {
1805:   PetscMPIInt    match;

1809:   MPI_Comm_compare(PetscObjectComm((PetscObject)is),comm,&match);
1810:   if (mode != PETSC_COPY_VALUES && (match == MPI_IDENT || match == MPI_CONGRUENT)) {
1811:     PetscObjectReference((PetscObject)is);
1812:     *newis = is;
1813:   } else {
1814:     (*is->ops->oncomm)(is,comm,mode,newis);
1815:   }
1816:   return 0;
1817: }

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

1822:    Logicall Collective on IS

1824:    Input Parameters:
1825: + is - index set
1826: - bs - block size

1828:    Level: intermediate

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

1836: .seealso: ISGetBlockSize(), ISCreateBlock(), ISBlockGetIndices(),
1837: @*/
1838: PetscErrorCode  ISSetBlockSize(IS is,PetscInt bs)
1839: {
1843:   if (PetscDefined(USE_DEBUG)) {
1844:     const PetscInt *indices;
1845:     PetscInt       length,i,j;
1846:     ISGetIndices(is,&indices);
1847:     if (indices) {
1848:       ISGetLocalSize(is,&length);
1850:       for (i=0;i<length/bs;i+=bs) {
1851:         for (j=0;j<bs-1;j++) {
1853:         }
1854:       }
1855:     }
1856:     ISRestoreIndices(is,&indices);
1857:   }
1858:   (*is->ops->setblocksize)(is,bs);
1859:   return 0;
1860: }

1862: /*@
1863:    ISGetBlockSize - Returns the number of elements in a block.

1865:    Not Collective

1867:    Input Parameter:
1868: .  is - the index set

1870:    Output Parameter:
1871: .  size - the number of elements in a block

1873:    Level: intermediate

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

1881: .seealso: ISBlockGetSize(), ISGetSize(), ISCreateBlock(), ISSetBlockSize()
1882: @*/
1883: PetscErrorCode  ISGetBlockSize(IS is,PetscInt *size)
1884: {
1885:   PetscLayoutGetBlockSize(is->map, size);
1886:   return 0;
1887: }

1889: PetscErrorCode ISGetIndicesCopy(IS is, PetscInt idx[])
1890: {
1891:   PetscInt       len,i;
1892:   const PetscInt *ptr;

1894:   ISGetLocalSize(is,&len);
1895:   ISGetIndices(is,&ptr);
1896:   for (i=0; i<len; i++) idx[i] = ptr[i];
1897:   ISRestoreIndices(is,&ptr);
1898:   return 0;
1899: }

1901: /*MC
1902:     ISGetIndicesF90 - Accesses the elements of an index set from Fortran90.
1903:     The users should call ISRestoreIndicesF90() after having looked at the
1904:     indices.  The user should NOT change the indices.

1906:     Synopsis:
1907:     ISGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

1909:     Not collective

1911:     Input Parameter:
1912: .   x - index set

1914:     Output Parameters:
1915: +   xx_v - the Fortran90 pointer to the array
1916: -   ierr - error code

1918:     Example of Usage:
1919: .vb
1920:     PetscInt, pointer xx_v(:)
1921:     ....
1922:     call ISGetIndicesF90(x,xx_v,ierr)
1923:     a = xx_v(3)
1924:     call ISRestoreIndicesF90(x,xx_v,ierr)
1925: .ve

1927:     Level: intermediate

1929: .seealso:  ISRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices()

1931: M*/

1933: /*MC
1934:     ISRestoreIndicesF90 - Restores an index set to a usable state after
1935:     a call to ISGetIndicesF90().

1937:     Synopsis:
1938:     ISRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

1940:     Not collective

1942:     Input Parameters:
1943: +   x - index set
1944: -   xx_v - the Fortran90 pointer to the array

1946:     Output Parameter:
1947: .   ierr - error code

1949:     Example of Usage:
1950: .vb
1951:     PetscInt, pointer xx_v(:)
1952:     ....
1953:     call ISGetIndicesF90(x,xx_v,ierr)
1954:     a = xx_v(3)
1955:     call ISRestoreIndicesF90(x,xx_v,ierr)
1956: .ve

1958:     Level: intermediate

1960: .seealso:  ISGetIndicesF90(), ISGetIndices(), ISRestoreIndices()

1962: M*/

1964: /*MC
1965:     ISBlockGetIndicesF90 - Accesses the elements of an index set from Fortran90.
1966:     The users should call ISBlockRestoreIndicesF90() after having looked at the
1967:     indices.  The user should NOT change the indices.

1969:     Synopsis:
1970:     ISBlockGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

1972:     Not collective

1974:     Input Parameter:
1975: .   x - index set

1977:     Output Parameters:
1978: +   xx_v - the Fortran90 pointer to the array
1979: -   ierr - error code
1980:     Example of Usage:
1981: .vb
1982:     PetscInt, pointer xx_v(:)
1983:     ....
1984:     call ISBlockGetIndicesF90(x,xx_v,ierr)
1985:     a = xx_v(3)
1986:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
1987: .ve

1989:     Level: intermediate

1991: .seealso:  ISBlockRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices(),
1992:            ISRestoreIndices()

1994: M*/

1996: /*MC
1997:     ISBlockRestoreIndicesF90 - Restores an index set to a usable state after
1998:     a call to ISBlockGetIndicesF90().

2000:     Synopsis:
2001:     ISBlockRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

2003:     Not Collective

2005:     Input Parameters:
2006: +   x - index set
2007: -   xx_v - the Fortran90 pointer to the array

2009:     Output Parameter:
2010: .   ierr - error code

2012:     Example of Usage:
2013: .vb
2014:     PetscInt, pointer xx_v(:)
2015:     ....
2016:     call ISBlockGetIndicesF90(x,xx_v,ierr)
2017:     a = xx_v(3)
2018:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
2019: .ve

2021:     Notes:
2022:     Not yet supported for all F90 compilers

2024:     Level: intermediate

2026: .seealso:  ISBlockGetIndicesF90(), ISGetIndices(), ISRestoreIndices(), ISRestoreIndicesF90()

2028: M*/