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:       PetscMalloc1(n, &idx);
647:       ISGetIndices(is, &idx);
648:       for (i = 1; i < n; i++) if (idx[i] != idx[i-1] + 1) break;
649:       if (i < n) intervalLocal = PETSC_FALSE;
650:       ISRestoreIndices(is, &idx);
651:     }

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

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

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

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

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

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

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

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

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

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

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

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

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

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

743:    Level: advanced

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

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

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

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

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

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

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

821:    Collective on IS

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

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

829:    Level: intermediate

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

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

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

850:    Logically Collective on IS

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

855:    Level: intermediate

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

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

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

873:    Not Collective

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

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

884:    Level: developer

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

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

905:    Logically Collective on IS

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

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

913:    Level: intermediate

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

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

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

934:    Logically Collective on IS

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

939:    Level: intermediate

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

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

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

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

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

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

980:    Collective on IS

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

985:    Level: beginner

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

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

1015:    Collective on IS

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

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

1025:    Level: intermediate

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

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

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

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

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

1063:    Not Collective

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

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

1071:    Level: beginner

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

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

1085:    Not Collective

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

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

1093:    Level: beginner

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

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

1107:    Not Collective

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

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

1115:    Level: developer

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

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

1130:    Collective

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

1136:    Level: developer

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

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

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

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

1159:    Not Collective

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

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

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

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

1190:    Level: intermediate

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

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

1205:    Not Collective

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

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

1214:    Level: intermediate

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

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

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

1233:   Not Collective

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

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

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

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

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

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

1278:    Not Collective

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

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

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

1301:    Level: intermediate

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

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

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


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

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

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

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

1351:    Collective on IS

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

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

1360:    Level: intermediate

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

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

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

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

1393:    Not Collective.

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

1399:    Level: intermediate

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

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

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

1422:    Collective on IS

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

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

1433:    Level: intermediate

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

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

1449:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
1450:   if (size == 1) *indices = NULL;
1451:   else {
1452:     if (!is->total) {
1453:       ISGatherTotal_Private(is);
1454:     }
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 on IS

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

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

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:   (*is->ops->view)(is,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 on PetscViewer

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:   (*is->ops->load)(is, 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 on IS

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

1645:    Level: intermediate

1647: .seealso: ISSortRemoveDups(), ISSorted()
1648: @*/
1649: PetscErrorCode  ISSort(IS is)
1650: {
1652:   (*is->ops->sort)(is);
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 on IS

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:   (*is->ops->sortremovedups)(is);
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 on IS

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

1687:    Level: intermediate

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

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

1703:    Collective on IS

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

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

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

1717:    Level: intermediate

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

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

1732:    Collective on IS

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

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

1740:    Level: beginner

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

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

1756:    Collective on IS

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

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

1764:    Level: beginner

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

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

1784:    Collective on IS

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

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

1794:    Level: advanced

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

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

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

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

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

1823:    Logicall Collective on IS

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

1829:    Level: intermediate

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

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

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

1866:    Not Collective

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

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

1874:    Level: intermediate

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

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

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

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

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

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

1910:     Not collective

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

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

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

1928:     Level: intermediate

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

1932: M*/

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

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

1941:     Not collective

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

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

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

1959:     Level: intermediate

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

1963: M*/

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

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

1973:     Not collective

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

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

1990:     Level: intermediate

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

1995: M*/

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

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

2004:     Not Collective

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

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

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

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

2025:     Level: intermediate

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

2029: M*/