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 an index set (with multiplicities) in a contiguous way.

 17:    Collective on IS

 19:    Input Parmeters:
 20: +  subset - the index set
 21: -  subset_mult - the multiplcity 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:    Level: intermediate

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

 45:   if (subset_mult) {
 47:   }
 48:   if (!N && !subset_n) return(0);
 49:   ISGetLocalSize(subset,&n);
 50:   if (subset_mult) {
 51:     ISGetLocalSize(subset_mult,&i);
 52:     if (i != n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Local subset and multiplicity sizes don't match! %d != %d",n,i);
 53:   }
 54:   /* create workspace layout for computing global indices of subset */
 55:   ISGetIndices(subset,&idxs);
 56:   lbounds[0] = lbounds[1] = 0;
 57:   for (i=0;i<n;i++) {
 58:     if (idxs[i] < lbounds[0]) lbounds[0] = idxs[i];
 59:     else if (idxs[i] > lbounds[1]) lbounds[1] = idxs[i];
 60:   }
 61:   lbounds[0] = -lbounds[0];
 62:   MPIU_Allreduce(lbounds,gbounds,2,MPIU_INT,MPI_MAX,PetscObjectComm((PetscObject)subset));
 63:   gbounds[0] = -gbounds[0];
 64:   N_n  = gbounds[1] - gbounds[0] + 1;

 66:   PetscLayoutCreate(PetscObjectComm((PetscObject)subset),&map);
 67:   PetscLayoutSetBlockSize(map,1);
 68:   PetscLayoutSetSize(map,N_n);
 69:   PetscLayoutSetUp(map);
 70:   PetscLayoutGetLocalSize(map,&Nl);

 72:   /* create sf : leaf_data == multiplicity of indexes, root data == global index in layout */
 73:   PetscMalloc2(n,&leaf_data,Nl,&root_data);
 74:   if (subset_mult) {
 75:     const PetscInt* idxs_mult;

 77:     ISGetIndices(subset_mult,&idxs_mult);
 78:     PetscArraycpy(leaf_data,idxs_mult,n);
 79:     ISRestoreIndices(subset_mult,&idxs_mult);
 80:   } else {
 81:     for (i=0;i<n;i++) leaf_data[i] = 1;
 82:   }
 83:   /* local size of new subset */
 84:   n_n = 0;
 85:   for (i=0;i<n;i++) n_n += leaf_data[i];

 87:   /* global indexes in layout */
 88:   PetscMalloc1(n_n,&gidxs); /* allocating possibly extra space in gidxs which will be used later */
 89:   for (i=0;i<n;i++) gidxs[i] = idxs[i] - gbounds[0];
 90:   ISRestoreIndices(subset,&idxs);
 91:   PetscSFCreate(PetscObjectComm((PetscObject)subset),&sf);
 92:   PetscSFSetGraphLayout(sf,map,n,NULL,PETSC_COPY_VALUES,gidxs);
 93:   PetscLayoutDestroy(&map);

 95:   /* reduce from leaves to roots */
 96:   PetscArrayzero(root_data,Nl);
 97:   PetscSFReduceBegin(sf,MPIU_INT,leaf_data,root_data,MPI_MAX);
 98:   PetscSFReduceEnd(sf,MPIU_INT,leaf_data,root_data,MPI_MAX);

100:   /* count indexes in local part of layout */
101:   nlocals = 0;
102:   first_index = -1;
103:   first_found = PETSC_FALSE;
104:   for (i=0;i<Nl;i++) {
105:     if (!first_found && root_data[i]) {
106:       first_found = PETSC_TRUE;
107:       first_index = i;
108:     }
109:     nlocals += root_data[i];
110:   }

112:   /* cumulative of number of indexes and size of subset without holes */
113: #if defined(PETSC_HAVE_MPI_EXSCAN)
114:   start = 0;
115:   MPI_Exscan(&nlocals,&start,1,MPIU_INT,MPI_SUM,PetscObjectComm((PetscObject)subset));
116: #else
117:   MPI_Scan(&nlocals,&start,1,MPIU_INT,MPI_SUM,PetscObjectComm((PetscObject)subset));
118:   start = start-nlocals;
119: #endif

121:   if (N) { /* compute total size of new subset if requested */
122:     *N   = start + nlocals;
123:     MPI_Comm_size(PetscObjectComm((PetscObject)subset),&commsize);
124:     MPI_Bcast(N,1,MPIU_INT,commsize-1,PetscObjectComm((PetscObject)subset));
125:   }

127:   if (!subset_n) {
128:     PetscFree(gidxs);
129:     PetscSFDestroy(&sf);
130:     PetscFree2(leaf_data,root_data);
131:     return(0);
132:   }

134:   /* adapt root data with cumulative */
135:   if (first_found) {
136:     PetscInt old_index;

138:     root_data[first_index] += start;
139:     old_index = first_index;
140:     for (i=first_index+1;i<Nl;i++) {
141:       if (root_data[i]) {
142:         root_data[i] += root_data[old_index];
143:         old_index = i;
144:       }
145:     }
146:   }

148:   /* from roots to leaves */
149:   PetscSFBcastBegin(sf,MPIU_INT,root_data,leaf_data,MPI_REPLACE);
150:   PetscSFBcastEnd(sf,MPIU_INT,root_data,leaf_data,MPI_REPLACE);
151:   PetscSFDestroy(&sf);

153:   /* create new IS with global indexes without holes */
154:   if (subset_mult) {
155:     const PetscInt* idxs_mult;
156:     PetscInt        cum;

158:     cum = 0;
159:     ISGetIndices(subset_mult,&idxs_mult);
160:     for (i=0;i<n;i++) {
161:       PetscInt j;
162:       for (j=0;j<idxs_mult[i];j++) gidxs[cum++] = leaf_data[i] - idxs_mult[i] + j;
163:     }
164:     ISRestoreIndices(subset_mult,&idxs_mult);
165:   } else {
166:     for (i=0;i<n;i++) {
167:       gidxs[i] = leaf_data[i]-1;
168:     }
169:   }
170:   ISCreateGeneral(PetscObjectComm((PetscObject)subset),n_n,gidxs,PETSC_OWN_POINTER,subset_n);
171:   PetscFree2(leaf_data,root_data);
172:   return(0);
173: }

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

178:    Collective on IS

180:    Input Parmeters:
181: +  is - the index set
182: -  comps - which components we will extract from is

184:    Output Parameters:
185: .  subis - the new sub index set

187:    Level: intermediate

189:    Example usage:
190:    We have an index set (is) living on 3 processes with the following values:
191:    | 4 9 0 | 2 6 7 | 10 11 1|
192:    and another index set (comps) used to indicate which components of is  we want to take,
193:    | 7 5  | 1 2 | 0 4|
194:    The output index set (subis) should look like:
195:    | 11 7 | 9 0 | 4 6|

197: .seealso: VecGetSubVector(), MatCreateSubMatrix()
198: @*/
199: PetscErrorCode ISCreateSubIS(IS is,IS comps,IS *subis)
200: {
201:   PetscSF         sf;
202:   const PetscInt  *is_indices,*comps_indices;
203:   PetscInt        *subis_indices,nroots,nleaves,*mine,i,lidx;
204:   PetscMPIInt     owner;
205:   PetscSFNode     *remote;
206:   PetscErrorCode  ierr;
207:   MPI_Comm        comm;


214:   PetscObjectGetComm((PetscObject)is, &comm);
215:   ISGetLocalSize(comps,&nleaves);
216:   ISGetLocalSize(is,&nroots);
217:   PetscMalloc1(nleaves,&remote);
218:   PetscMalloc1(nleaves,&mine);
219:   ISGetIndices(comps,&comps_indices);
220:   /*
221:    * Construct a PetscSF in which "is" data serves as roots and "subis" is leaves.
222:    * Root data are sent to leaves using PetscSFBcast().
223:    * */
224:   for (i=0; i<nleaves; i++) {
225:     mine[i] = i;
226:     /* Connect a remote root with the current leaf. The value on the remote root
227:      * will be received by the current local leaf.
228:      * */
229:     owner = -1;
230:     lidx =  -1;
231:     PetscLayoutFindOwnerIndex(is->map,comps_indices[i],&owner,&lidx);
232:     remote[i].rank = owner;
233:     remote[i].index = lidx;
234:   }
235:   ISRestoreIndices(comps,&comps_indices);
236:   PetscSFCreate(comm,&sf);
237:   PetscSFSetFromOptions(sf);\
238:   PetscSFSetGraph(sf,nroots,nleaves,mine,PETSC_OWN_POINTER,remote,PETSC_OWN_POINTER);

240:   PetscMalloc1(nleaves,&subis_indices);
241:   ISGetIndices(is, &is_indices);
242:   PetscSFBcastBegin(sf,MPIU_INT,is_indices,subis_indices,MPI_REPLACE);
243:   PetscSFBcastEnd(sf,MPIU_INT,is_indices,subis_indices,MPI_REPLACE);
244:   ISRestoreIndices(is,&is_indices);
245:   PetscSFDestroy(&sf);
246:   ISCreateGeneral(comm,nleaves,subis_indices,PETSC_OWN_POINTER,subis);
247:   return(0);
248: }

250: /*@
251:    ISClearInfoCache - clear the cache of computed index set properties

253:    Not collective

255:    Input Parameters:
256: +  is - the index set
257: -  clear_permanent_local - whether to remove the permanent status of local properties

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

262:    Level: developer

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

266: @*/
267: PetscErrorCode ISClearInfoCache(IS is, PetscBool clear_permanent_local)
268: {
269:   PetscInt i, j;

274:   for (i = 0; i < IS_INFO_MAX; i++) {
275:     if (clear_permanent_local) is->info_permanent[IS_LOCAL][i] = PETSC_FALSE;
276:     for (j = 0; j < 2; j++) {
277:       if (!is->info_permanent[j][i]) is->info[j][i] = IS_INFO_UNKNOWN;
278:     }
279:   }
280:   return(0);
281: }

283: static PetscErrorCode ISSetInfo_Internal(IS is, ISInfo info, ISInfoType type, ISInfoBool ipermanent, PetscBool flg)
284: {
285:   ISInfoBool     iflg = flg ? IS_INFO_TRUE : IS_INFO_FALSE;
286:   PetscInt       itype = (type == IS_LOCAL) ? 0 : 1;
287:   PetscBool      permanent_set = (ipermanent == IS_INFO_UNKNOWN) ? PETSC_FALSE : PETSC_TRUE;
288:   PetscBool      permanent = (ipermanent == IS_INFO_TRUE) ? PETSC_TRUE : PETSC_FALSE;

291:   /* set this property */
292:   is->info[itype][(int)info] = iflg;
293:   if (permanent_set) is->info_permanent[itype][(int)info] = permanent;
294:   /* set implications */
295:   switch (info) {
296:   case IS_SORTED:
297:     if (flg && type == IS_GLOBAL) { /* an array that is globally sorted is also locally sorted */
298:       is->info[IS_LOCAL][(int)info] = IS_INFO_TRUE;
299:       /* global permanence implies local permanence */
300:       if (permanent_set && permanent) is->info_permanent[IS_LOCAL][(int)info] = PETSC_TRUE;
301:     }
302:     if (!flg) { /* if an array is not sorted, it cannot be an interval or the identity */
303:       is->info[itype][IS_INTERVAL] = IS_INFO_FALSE;
304:       is->info[itype][IS_IDENTITY] = IS_INFO_FALSE;
305:       if (permanent_set) {
306:         is->info_permanent[itype][IS_INTERVAL] = permanent;
307:         is->info_permanent[itype][IS_IDENTITY] = permanent;
308:       }
309:     }
310:     break;
311:   case IS_UNIQUE:
312:     if (flg && type == IS_GLOBAL) { /* an array that is globally unique is also locally unique */
313:       is->info[IS_LOCAL][(int)info] = IS_INFO_TRUE;
314:       /* global permanence implies local permanence */
315:       if (permanent_set && permanent) is->info_permanent[IS_LOCAL][(int)info] = PETSC_TRUE;
316:     }
317:     if (!flg) { /* if an array is not unique, it cannot be a permutation, and interval, or the identity */
318:       is->info[itype][IS_PERMUTATION] = IS_INFO_FALSE;
319:       is->info[itype][IS_INTERVAL]    = IS_INFO_FALSE;
320:       is->info[itype][IS_IDENTITY]    = IS_INFO_FALSE;
321:       if (permanent_set) {
322:         is->info_permanent[itype][IS_PERMUTATION] = permanent;
323:         is->info_permanent[itype][IS_INTERVAL]    = permanent;
324:         is->info_permanent[itype][IS_IDENTITY]    = permanent;
325:       }
326:     }
327:     break;
328:   case IS_PERMUTATION:
329:     if (flg) { /* an array that is a permutation is unique and is unique locally */
330:       is->info[itype][IS_UNIQUE] = IS_INFO_TRUE;
331:       is->info[IS_LOCAL][IS_UNIQUE] = IS_INFO_TRUE;
332:       if (permanent_set && permanent) {
333:         is->info_permanent[itype][IS_UNIQUE] = PETSC_TRUE;
334:         is->info_permanent[IS_LOCAL][IS_UNIQUE] = PETSC_TRUE;
335:       }
336:     } else { /* an array that is not a permutation cannot be the identity */
337:       is->info[itype][IS_IDENTITY] = IS_INFO_FALSE;
338:       if (permanent_set) is->info_permanent[itype][IS_IDENTITY] = permanent;
339:     }
340:     break;
341:   case IS_INTERVAL:
342:     if (flg) { /* an array that is an interval is sorted and unique */
343:       is->info[itype][IS_SORTED]         = IS_INFO_TRUE;
344:       is->info[IS_LOCAL][IS_SORTED]      = IS_INFO_TRUE;
345:       is->info[itype][IS_UNIQUE]         = IS_INFO_TRUE;
346:       is->info[IS_LOCAL][IS_UNIQUE]      = IS_INFO_TRUE;
347:       if (permanent_set && permanent) {
348:         is->info_permanent[itype][IS_SORTED]    = PETSC_TRUE;
349:         is->info_permanent[IS_LOCAL][IS_SORTED] = PETSC_TRUE;
350:         is->info_permanent[itype][IS_UNIQUE]    = PETSC_TRUE;
351:         is->info_permanent[IS_LOCAL][IS_UNIQUE] = PETSC_TRUE;
352:       }
353:     } else { /* an array that is not an interval cannot be the identity */
354:       is->info[itype][IS_IDENTITY] = IS_INFO_FALSE;
355:       if (permanent_set) is->info_permanent[itype][IS_IDENTITY] = permanent;
356:     }
357:     break;
358:   case IS_IDENTITY:
359:     if (flg) { /* an array that is the identity is sorted, unique, an interval, and a permutation */
360:       is->info[itype][IS_SORTED]         = IS_INFO_TRUE;
361:       is->info[IS_LOCAL][IS_SORTED]      = IS_INFO_TRUE;
362:       is->info[itype][IS_UNIQUE]         = IS_INFO_TRUE;
363:       is->info[IS_LOCAL][IS_UNIQUE]      = IS_INFO_TRUE;
364:       is->info[itype][IS_PERMUTATION]    = IS_INFO_TRUE;
365:       is->info[itype][IS_INTERVAL]       = IS_INFO_TRUE;
366:       is->info[IS_LOCAL][IS_INTERVAL]    = 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:         is->info_permanent[itype][IS_PERMUTATION]    = PETSC_TRUE;
373:         is->info_permanent[itype][IS_INTERVAL]       = PETSC_TRUE;
374:         is->info_permanent[IS_LOCAL][IS_INTERVAL]    = PETSC_TRUE;
375:       }
376:     }
377:     break;
378:   default:
379:     if (type == IS_LOCAL) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Unknown IS property");
380:     else SETERRQ(PetscObjectComm((PetscObject)is), PETSC_ERR_ARG_OUTOFRANGE, "Unknown IS property");
381:   }
382:   return(0);
383: }

385: /*@
386:    ISSetInfo - Set known information about an index set.

388:    Logically Collective on IS if type is IS_GLOBAL

390:    Input Parameters:
391: +  is - the index set
392: .  info - describing a property of the index set, one of those listed below,
393: .  type - IS_LOCAL if the information describes the local portion of the index set,
394:           IS_GLOBAL if it describes the whole index set
395: .  permanent - PETSC_TRUE if it is known that the property will persist through changes to the index set, PETSC_FALSE otherwise
396:                If the user sets a property as permanently known, it will bypass computation of that property
397: -  flg - set the described property as true (PETSC_TRUE) or false (PETSC_FALSE)

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

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

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

411:    Level: advanced

413: .seealso:  ISInfo, ISInfoType, IS

415: @*/
416: PetscErrorCode ISSetInfo(IS is, ISInfo info, ISInfoType type, PetscBool permanent, PetscBool flg)
417: {
418:   MPI_Comm       comm, errcomm;
419:   PetscMPIInt    size;

425:   comm = PetscObjectComm((PetscObject)is);
426:   if (type == IS_GLOBAL) {
430:     errcomm = comm;
431:   } else {
432:     errcomm = PETSC_COMM_SELF;
433:   }

435:   if (((int) info) <= IS_INFO_MIN || ((int) info) >= IS_INFO_MAX) SETERRQ1(errcomm,PETSC_ERR_ARG_OUTOFRANGE,"Options %d is out of range",(int)info);

437:   MPI_Comm_size(comm, &size);
438:   /* do not use global values if size == 1: it makes it easier to keep the implications straight */
439:   if (size == 1) type = IS_LOCAL;
440:   ISSetInfo_Internal(is, info, type, permanent ? IS_INFO_TRUE : IS_INFO_FALSE, flg);
441:   return(0);
442: }

444: static PetscErrorCode ISGetInfo_Sorted(IS is, ISInfoType type, PetscBool *flg)
445: {
446:   MPI_Comm       comm;
447:   PetscMPIInt    size, rank;

451:   comm = PetscObjectComm((PetscObject)is);
452:   MPI_Comm_size(comm, &size);
453:   MPI_Comm_size(comm, &rank);
454:   if (type == IS_GLOBAL && is->ops->sortedglobal) {
455:     (*is->ops->sortedglobal)(is,flg);
456:   } else {
457:     PetscBool sortedLocal = PETSC_FALSE;

459:     /* determine if the array is locally sorted */
460:     if (type == IS_GLOBAL && size > 1) {
461:       /* call ISGetInfo so that a cached value will be used if possible */
462:       ISGetInfo(is, IS_SORTED, IS_LOCAL, PETSC_TRUE, &sortedLocal);
463:     } else if (is->ops->sortedlocal) {
464:       (*is->ops->sortedlocal)(is,&sortedLocal);
465:     } else {
466:       /* default: get the local indices and directly check */
467:       const PetscInt *idx;
468:       PetscInt n;

470:       ISGetIndices(is, &idx);
471:       ISGetLocalSize(is, &n);
472:       PetscSortedInt(n, idx, &sortedLocal);
473:       ISRestoreIndices(is, &idx);
474:     }

476:     if (type == IS_LOCAL || size == 1) {
477:       *flg = sortedLocal;
478:     } else {
479:       MPI_Allreduce(&sortedLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
480:       if (*flg) {
481:         PetscInt  n, min = PETSC_MAX_INT, max = PETSC_MIN_INT;
482:         PetscInt  maxprev;

484:         ISGetLocalSize(is, &n);
485:         if (n) {ISGetMinMax(is, &min, &max);}
486:         maxprev = PETSC_MIN_INT;
487:         MPI_Exscan(&max, &maxprev, 1, MPIU_INT, MPI_MAX, comm);
488:         if (rank && (maxprev > min)) sortedLocal = PETSC_FALSE;
489:         MPI_Allreduce(&sortedLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
490:       }
491:     }
492:   }
493:   return(0);
494: }

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

498: static PetscErrorCode ISGetInfo_Unique(IS is, ISInfoType type, PetscBool *flg)
499: {
500:   MPI_Comm       comm;
501:   PetscMPIInt    size, rank;
502:   PetscInt       i;

506:   comm = PetscObjectComm((PetscObject)is);
507:   MPI_Comm_size(comm, &size);
508:   MPI_Comm_size(comm, &rank);
509:   if (type == IS_GLOBAL && is->ops->uniqueglobal) {
510:     (*is->ops->uniqueglobal)(is,flg);
511:   } else {
512:     PetscBool uniqueLocal;
513:     PetscInt  n = -1;
514:     PetscInt  *idx = NULL;

516:     /* determine if the array is locally unique */
517:     if (type == IS_GLOBAL && size > 1) {
518:       /* call ISGetInfo so that a cached value will be used if possible */
519:       ISGetInfo(is, IS_UNIQUE, IS_LOCAL, PETSC_TRUE, &uniqueLocal);
520:     } else if (is->ops->uniquelocal) {
521:       (*is->ops->uniquelocal)(is,&uniqueLocal);
522:     } else {
523:       /* default: get the local indices and directly check */
524:       uniqueLocal = PETSC_TRUE;
525:       ISGetLocalSize(is, &n);
526:       PetscMalloc1(n, &idx);
527:       ISGetIndicesCopy(is, idx);
528:       PetscIntSortSemiOrdered(n, idx);
529:       for (i = 1; i < n; i++) if (idx[i] == idx[i-1]) break;
530:       if (i < n) uniqueLocal = PETSC_FALSE;
531:     }

533:     PetscFree(idx);
534:     if (type == IS_LOCAL || size == 1) {
535:       *flg = uniqueLocal;
536:     } else {
537:       MPI_Allreduce(&uniqueLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
538:       if (*flg) {
539:         PetscInt  min = PETSC_MAX_INT, max = PETSC_MIN_INT, maxprev;

541:         if (!idx) {
542:           ISGetLocalSize(is, &n);
543:           PetscMalloc1(n, &idx);
544:           ISGetIndicesCopy(is, idx);
545:         }
546:         PetscParallelSortInt(is->map, is->map, idx, idx);
547:         if (n) {
548:           min = idx[0];
549:           max = idx[n - 1];
550:         }
551:         for (i = 1; i < n; i++) if (idx[i] == idx[i-1]) break;
552:         if (i < n) uniqueLocal = PETSC_FALSE;
553:         maxprev = PETSC_MIN_INT;
554:         MPI_Exscan(&max, &maxprev, 1, MPIU_INT, MPI_MAX, comm);
555:         if (rank && (maxprev == min)) uniqueLocal = PETSC_FALSE;
556:         MPI_Allreduce(&uniqueLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
557:       }
558:     }
559:     PetscFree(idx);
560:   }
561:   return(0);
562: }

564: static PetscErrorCode ISGetInfo_Permutation(IS is, ISInfoType type, PetscBool *flg)
565: {
566:   MPI_Comm       comm;
567:   PetscMPIInt    size, rank;

571:   comm = PetscObjectComm((PetscObject)is);
572:   MPI_Comm_size(comm, &size);
573:   MPI_Comm_size(comm, &rank);
574:   if (type == IS_GLOBAL && is->ops->permglobal) {
575:     (*is->ops->permglobal)(is,flg);
576:   } else if (type == IS_LOCAL && is->ops->permlocal) {
577:     (*is->ops->permlocal)(is,flg);
578:   } else {
579:     PetscBool permLocal;
580:     PetscInt  n, i, rStart;
581:     PetscInt  *idx;

583:     ISGetLocalSize(is, &n);
584:     PetscMalloc1(n, &idx);
585:     ISGetIndicesCopy(is, idx);
586:     if (type == IS_GLOBAL) {
587:       PetscParallelSortInt(is->map, is->map, idx, idx);
588:       PetscLayoutGetRange(is->map, &rStart, NULL);
589:     } else {
590:       PetscIntSortSemiOrdered(n, idx);
591:       rStart = 0;
592:     }
593:     permLocal = PETSC_TRUE;
594:     for (i = 0; i < n; i++) {
595:       if (idx[i] != rStart + i) break;
596:     }
597:     if (i < n) permLocal = PETSC_FALSE;
598:     if (type == IS_LOCAL || size == 1) {
599:       *flg = permLocal;
600:     } else {
601:       MPI_Allreduce(&permLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
602:     }
603:     PetscFree(idx);
604:   }
605:   return(0);
606: }

608: static PetscErrorCode ISGetInfo_Interval(IS is, ISInfoType type, PetscBool *flg)
609: {
610:   MPI_Comm       comm;
611:   PetscMPIInt    size, rank;
612:   PetscInt       i;

616:   comm = PetscObjectComm((PetscObject)is);
617:   MPI_Comm_size(comm, &size);
618:   MPI_Comm_size(comm, &rank);
619:   if (type == IS_GLOBAL && is->ops->intervalglobal) {
620:     (*is->ops->intervalglobal)(is,flg);
621:   } else {
622:     PetscBool intervalLocal;

624:     /* determine if the array is locally an interval */
625:     if (type == IS_GLOBAL && size > 1) {
626:       /* call ISGetInfo so that a cached value will be used if possible */
627:       ISGetInfo(is, IS_INTERVAL, IS_LOCAL, PETSC_TRUE, &intervalLocal);
628:     } else if (is->ops->intervallocal) {
629:       (*is->ops->intervallocal)(is,&intervalLocal);
630:     } else {
631:       PetscInt        n;
632:       const PetscInt  *idx;
633:       /* default: get the local indices and directly check */
634:       intervalLocal = PETSC_TRUE;
635:       ISGetLocalSize(is, &n);
636:       PetscMalloc1(n, &idx);
637:       ISGetIndices(is, &idx);
638:       for (i = 1; i < n; i++) if (idx[i] != idx[i-1] + 1) break;
639:       if (i < n) intervalLocal = PETSC_FALSE;
640:       ISRestoreIndices(is, &idx);
641:     }

643:     if (type == IS_LOCAL || size == 1) {
644:       *flg = intervalLocal;
645:     } else {
646:       MPI_Allreduce(&intervalLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
647:       if (*flg) {
648:         PetscInt  n, min = PETSC_MAX_INT, max = PETSC_MIN_INT;
649:         PetscInt  maxprev;

651:         ISGetLocalSize(is, &n);
652:         if (n) {ISGetMinMax(is, &min, &max);}
653:         maxprev = PETSC_MIN_INT;
654:         MPI_Exscan(&max, &maxprev, 1, MPIU_INT, MPI_MAX, comm);
655:         if (rank && n && (maxprev != min - 1)) intervalLocal = PETSC_FALSE;
656:         MPI_Allreduce(&intervalLocal, flg, 1, MPIU_BOOL, MPI_LAND, comm);
657:       }
658:     }
659:   }
660:   return(0);
661: }

663: static PetscErrorCode ISGetInfo_Identity(IS is, ISInfoType type, PetscBool *flg)
664: {
665:   MPI_Comm       comm;
666:   PetscMPIInt    size, rank;

670:   comm = PetscObjectComm((PetscObject)is);
671:   MPI_Comm_size(comm, &size);
672:   MPI_Comm_size(comm, &rank);
673:   if (type == IS_GLOBAL && is->ops->intervalglobal) {
674:     PetscBool isinterval;

676:     (*is->ops->intervalglobal)(is,&isinterval);
677:     *flg = PETSC_FALSE;
678:     if (isinterval) {
679:       PetscInt  min;

681:       ISGetMinMax(is, &min, NULL);
682:       MPI_Bcast(&min, 1, MPIU_INT, 0, comm);
683:       if (min == 0) *flg = PETSC_TRUE;
684:     }
685:   } else if (type == IS_LOCAL && is->ops->intervallocal) {
686:     PetscBool isinterval;

688:     (*is->ops->intervallocal)(is,&isinterval);
689:     *flg = PETSC_FALSE;
690:     if (isinterval) {
691:       PetscInt  min;

693:       ISGetMinMax(is, &min, NULL);
694:       if (min == 0) *flg = PETSC_TRUE;
695:     }
696:   } else {
697:     PetscBool identLocal;
698:     PetscInt  n, i, rStart;
699:     const PetscInt *idx;

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

719: /*@
720:    ISGetInfo - Determine whether an index set satisfies a given property

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

724:    Input Parameters:
725: +  is - the index set
726: .  info - describing a property of the index set, one of those listed in the documentation of ISSetInfo()
727: .  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
728: -  type - whether the property is local (IS_LOCAL) or global (IS_GLOBAL)

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

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

735:    Level: advanced

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

739: @*/
740: PetscErrorCode ISGetInfo(IS is, ISInfo info, ISInfoType type, PetscBool compute, PetscBool *flg)
741: {
742:   MPI_Comm       comm, errcomm;
743:   PetscMPIInt    rank, size;
744:   PetscInt       itype;
745:   PetscBool      hasprop;
746:   PetscBool      infer;

752:   comm = PetscObjectComm((PetscObject)is);
753:   if (type == IS_GLOBAL) {
755:     errcomm = comm;
756:   } else {
757:     errcomm = PETSC_COMM_SELF;
758:   }

760:   MPI_Comm_size(comm, &size);
761:   MPI_Comm_rank(comm, &rank);

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

805: static PetscErrorCode ISCopyInfo(IS source, IS dest)
806: {

810:   PetscArraycpy(&dest->info[0], &source->info[0], 2);
811:   PetscArraycpy(&dest->info_permanent[0], &source->info_permanent[0], 2);
812:   return(0);
813: }

815: /*@
816:    ISIdentity - Determines whether index set is the identity mapping.

818:    Collective on IS

820:    Input Parmeters:
821: .  is - the index set

823:    Output Parameters:
824: .  ident - PETSC_TRUE if an identity, else PETSC_FALSE

826:    Level: intermediate

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

834: .seealso: ISSetIdentity(), ISGetInfo()
835: @*/
836: PetscErrorCode  ISIdentity(IS is,PetscBool  *ident)
837: {

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 Parmeters:
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: {

869:   ISSetInfo(is,IS_IDENTITY,IS_GLOBAL,PETSC_TRUE,PETSC_TRUE);
870:   return(0);
871: }

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

876:    Not Collective

878:    Input Parmeters:
879: +  is - the index set
880: .  gstart - global start
881: -  gend - global end

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

887:    Level: developer

889: .seealso: ISGetLocalSize(), VecGetOwnershipRange()
890: @*/
891: PetscErrorCode  ISContiguousLocal(IS is,PetscInt gstart,PetscInt gend,PetscInt *start,PetscBool *contig)
892: {

899:   *start  = -1;
900:   *contig = PETSC_FALSE;
901:   if (is->ops->contiguous) {
902:     (*is->ops->contiguous)(is,gstart,gend,start,contig);
903:   }
904:   return(0);
905: }

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

911:    Logically Collective on IS

913:    Input Parmeters:
914: .  is - the index set

916:    Output Parameters:
917: .  perm - PETSC_TRUE if a permutation, else PETSC_FALSE

919:    Level: intermediate

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

927: .seealso: ISSetPermutation(), ISGetInfo()
928: @*/
929: PetscErrorCode  ISPermutation(IS is,PetscBool  *perm)
930: {

936:   ISGetInfo(is,IS_PERMUTATION,IS_GLOBAL,PETSC_FALSE,perm);
937:   return(0);
938: }

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

943:    Logically Collective on IS

945:    Input Parmeters:
946: .  is - the index set

948:    Level: intermediate

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

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

957: .seealso: ISPermutation(), ISSetInfo(), ISClearInfoCache().
958: @*/
959: PetscErrorCode  ISSetPermutation(IS is)
960: {

965:   if (PetscDefined(USE_DEBUG)) {
966:     PetscMPIInt    size;

968:     MPI_Comm_size(PetscObjectComm((PetscObject)is),&size);
969:     if (size == 1) {
970:       PetscInt       i,n,*idx;
971:       const PetscInt *iidx;

973:       ISGetSize(is,&n);
974:       PetscMalloc1(n,&idx);
975:       ISGetIndices(is,&iidx);
976:       PetscArraycpy(idx,iidx,n);
977:       PetscIntSortSemiOrdered(n,idx);
978:       for (i=0; i<n; i++) {
979:         if (idx[i] != i) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Index set is not a permutation");
980:       }
981:       PetscFree(idx);
982:       ISRestoreIndices(is,&iidx);
983:     }
984:   }
985:   ISSetInfo(is,IS_PERMUTATION,IS_GLOBAL,PETSC_TRUE,PETSC_TRUE);
986:   return(0);
987: }

989: /*@C
990:    ISDestroy - Destroys an index set.

992:    Collective on IS

994:    Input Parameters:
995: .  is - the index set

997:    Level: beginner

999: .seealso: ISCreateGeneral(), ISCreateStride(), ISCreateBlocked()
1000: @*/
1001: PetscErrorCode  ISDestroy(IS *is)
1002: {

1006:   if (!*is) return(0);
1008:   if (--((PetscObject)(*is))->refct > 0) {*is = NULL; return(0);}
1009:   if ((*is)->complement) {
1010:     PetscInt refcnt;
1011:     PetscObjectGetReference((PetscObject)((*is)->complement), &refcnt);
1012:     if (refcnt > 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Nonlocal IS has not been restored");
1013:     ISDestroy(&(*is)->complement);
1014:   }
1015:   if ((*is)->ops->destroy) {
1016:     (*(*is)->ops->destroy)(*is);
1017:   }
1018:   PetscLayoutDestroy(&(*is)->map);
1019:   /* Destroy local representations of offproc data. */
1020:   PetscFree((*is)->total);
1021:   PetscFree((*is)->nonlocal);
1022:   PetscHeaderDestroy(is);
1023:   return(0);
1024: }

1026: /*@
1027:    ISInvertPermutation - Creates a new permutation that is the inverse of
1028:                          a given permutation.

1030:    Collective on IS

1032:    Input Parameter:
1033: +  is - the index set
1034: -  nlocal - number of indices on this processor in result (ignored for 1 proccessor) or
1035:             use PETSC_DECIDE

1037:    Output Parameter:
1038: .  isout - the inverse permutation

1040:    Level: intermediate

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

1046: @*/
1047: PetscErrorCode  ISInvertPermutation(IS is,PetscInt nlocal,IS *isout)
1048: {
1049:   PetscBool      isperm, isidentity, issame;

1055:   ISGetInfo(is,IS_PERMUTATION,IS_GLOBAL,PETSC_TRUE,&isperm);
1056:   if (!isperm) SETERRQ(PetscObjectComm((PetscObject)is),PETSC_ERR_ARG_WRONG,"Not a permutation");
1057:   ISGetInfo(is,IS_IDENTITY,IS_GLOBAL,PETSC_TRUE,&isidentity);
1058:   issame = PETSC_FALSE;
1059:   if (isidentity) {
1060:     PetscInt n;
1061:     PetscBool isallsame;

1063:     ISGetLocalSize(is, &n);
1064:     issame = (PetscBool) (n == nlocal);
1065:     MPI_Allreduce(&issame, &isallsame, 1, MPIU_BOOL, MPI_LAND, PetscObjectComm((PetscObject)is));
1066:     issame = isallsame;
1067:   }
1068:   if (issame) {
1069:     ISDuplicate(is,isout);
1070:   } else {
1071:     (*is->ops->invertpermutation)(is,nlocal,isout);
1072:     ISSetPermutation(*isout);
1073:   }
1074:   return(0);
1075: }

1077: /*@
1078:    ISGetSize - Returns the global length of an index set.

1080:    Not Collective

1082:    Input Parameter:
1083: .  is - the index set

1085:    Output Parameter:
1086: .  size - the global size

1088:    Level: beginner

1090: @*/
1091: PetscErrorCode  ISGetSize(IS is,PetscInt *size)
1092: {
1096:   *size = is->map->N;
1097:   return(0);
1098: }

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

1103:    Not Collective

1105:    Input Parameter:
1106: .  is - the index set

1108:    Output Parameter:
1109: .  size - the local size

1111:    Level: beginner

1113: @*/
1114: PetscErrorCode  ISGetLocalSize(IS is,PetscInt *size)
1115: {
1119:   *size = is->map->n;
1120:   return(0);
1121: }

1123: /*@
1124:    ISGetLayout - get PetscLayout describing index set layout

1126:    Not Collective

1128:    Input Arguments:
1129: .  is - the index set

1131:    Output Arguments:
1132: .  map - the layout

1134:    Level: developer

1136: .seealso: ISGetSize(), ISGetLocalSize()
1137: @*/
1138: PetscErrorCode ISGetLayout(IS is,PetscLayout *map)
1139: {

1144:   *map = is->map;
1145:   return(0);
1146: }

1148: /*@C
1149:    ISGetIndices - Returns a pointer to the indices.  The user should call
1150:    ISRestoreIndices() after having looked at the indices.  The user should
1151:    NOT change the indices.

1153:    Not Collective

1155:    Input Parameter:
1156: .  is - the index set

1158:    Output Parameter:
1159: .  ptr - the location to put the pointer to the indices

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

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

1184:    Level: intermediate

1186: .seealso: ISRestoreIndices(), ISGetIndicesF90()
1187: @*/
1188: PetscErrorCode  ISGetIndices(IS is,const PetscInt *ptr[])
1189: {

1195:   (*is->ops->getindices)(is,ptr);
1196:   return(0);
1197: }

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

1202:    Not Collective

1204:    Input Parameter:
1205: .  is - the index set

1207:    Output Parameter:
1208: +   min - the minimum value
1209: -   max - the maximum value

1211:    Level: intermediate

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

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

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

1231:   Not Collective

1233:   Input Parameter:
1234: + is - the index set
1235: - key - the search key

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

1240:   Level: intermediate
1241: @*/
1242: PetscErrorCode ISLocate(IS is, PetscInt key, PetscInt *location)
1243: {

1247:   if (is->ops->locate) {
1248:     (*is->ops->locate)(is,key,location);
1249:   } else {
1250:     PetscInt       numIdx;
1251:     PetscBool      sorted;
1252:     const PetscInt *idx;

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

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

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

1279:    Not Collective

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

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

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

1302:    Level: intermediate

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

1307: .seealso: ISGetIndices(), ISRestoreIndicesF90()
1308: @*/
1309: PetscErrorCode  ISRestoreIndices(IS is,const PetscInt *ptr[])
1310: {

1316:   if (is->ops->restoreindices) {
1317:     (*is->ops->restoreindices)(is,ptr);
1318:   }
1319:   return(0);
1320: }

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


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

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

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

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

1357:    Collective on IS

1359:    Input Parameter:
1360: .  is - the index set

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

1366:    Level: intermediate

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

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

1386:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
1387:   if (size == 1) {
1388:     (*is->ops->getindices)(is,indices);
1389:   } else {
1390:     if (!is->total) {
1391:       ISGatherTotal_Private(is);
1392:     }
1393:     *indices = is->total;
1394:   }
1395:   return(0);
1396: }

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

1401:    Not Collective.

1403:    Input Parameter:
1404: +  is - the index set
1405: -  indices - index array; must be the array obtained with ISGetTotalIndices()

1407:    Level: intermediate

1409: .seealso: ISRestoreTotalIndices(), ISGetNonlocalIndices()
1410: @*/
1411: PetscErrorCode  ISRestoreTotalIndices(IS is, const PetscInt *indices[])
1412: {
1414:   PetscMPIInt    size;

1419:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
1420:   if (size == 1) {
1421:     (*is->ops->restoreindices)(is,indices);
1422:   } else {
1423:     if (is->total != *indices) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Index array pointer being restored does not point to the array obtained from the IS.");
1424:   }
1425:   return(0);
1426: }
1427: /*@C
1428:    ISGetNonlocalIndices - Retrieve an array of indices from remote processors
1429:                        in this communicator.

1431:    Collective on IS

1433:    Input Parameter:
1434: .  is - the index set

1436:    Output Parameter:
1437: .  indices - indices with rank 0 indices first, and so on,  omitting
1438:              the current rank.  Total number of indices is the difference
1439:              total and local, obtained with ISGetSize() and ISGetLocalSize(),
1440:              respectively.

1442:    Level: intermediate

1444:    Notes:
1445:     restore the indices using ISRestoreNonlocalIndices().
1446:           The same scalability considerations as those for ISGetTotalIndices
1447:           apply here.

1449: .seealso: ISGetTotalIndices(), ISRestoreNonlocalIndices(), ISGetSize(), ISGetLocalSize().
1450: @*/
1451: PetscErrorCode  ISGetNonlocalIndices(IS is, const PetscInt *indices[])
1452: {
1454:   PetscMPIInt    size;
1455:   PetscInt       n, N;

1460:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
1461:   if (size == 1) *indices = NULL;
1462:   else {
1463:     if (!is->total) {
1464:       ISGatherTotal_Private(is);
1465:     }
1466:     ISGetLocalSize(is,&n);
1467:     ISGetSize(is,&N);
1468:     PetscMalloc1(N-n, &(is->nonlocal));
1469:     PetscArraycpy(is->nonlocal, is->total, is->local_offset);
1470:     PetscArraycpy(is->nonlocal+is->local_offset, is->total+is->local_offset+n,N - is->local_offset - n);
1471:     *indices = is->nonlocal;
1472:   }
1473:   return(0);
1474: }

1476: /*@C
1477:    ISRestoreTotalIndices - Restore the index array obtained with ISGetNonlocalIndices().

1479:    Not Collective.

1481:    Input Parameter:
1482: +  is - the index set
1483: -  indices - index array; must be the array obtained with ISGetNonlocalIndices()

1485:    Level: intermediate

1487: .seealso: ISGetTotalIndices(), ISGetNonlocalIndices(), ISRestoreTotalIndices()
1488: @*/
1489: PetscErrorCode  ISRestoreNonlocalIndices(IS is, const PetscInt *indices[])
1490: {
1494:   if (is->nonlocal != *indices) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Index array pointer being restored does not point to the array obtained from the IS.");
1495:   return(0);
1496: }

1498: /*@
1499:    ISGetNonlocalIS - Gather all nonlocal indices for this IS and present
1500:                      them as another sequential index set.

1502:    Collective on IS

1504:    Input Parameter:
1505: .  is - the index set

1507:    Output Parameter:
1508: .  complement - sequential IS with indices identical to the result of
1509:                 ISGetNonlocalIndices()

1511:    Level: intermediate

1513:    Notes:
1514:     complement represents the result of ISGetNonlocalIndices as an IS.
1515:           Therefore scalability issues similar to ISGetNonlocalIndices apply.
1516:           The resulting IS must be restored using ISRestoreNonlocalIS().

1518: .seealso: ISGetNonlocalIndices(), ISRestoreNonlocalIndices(),  ISAllGather(), ISGetSize()
1519: @*/
1520: PetscErrorCode  ISGetNonlocalIS(IS is, IS *complement)
1521: {

1527:   /* Check if the complement exists already. */
1528:   if (is->complement) {
1529:     *complement = is->complement;
1530:     PetscObjectReference((PetscObject)(is->complement));
1531:   } else {
1532:     PetscInt       N, n;
1533:     const PetscInt *idx;
1534:     ISGetSize(is, &N);
1535:     ISGetLocalSize(is,&n);
1536:     ISGetNonlocalIndices(is, &idx);
1537:     ISCreateGeneral(PETSC_COMM_SELF, N-n,idx, PETSC_USE_POINTER, &(is->complement));
1538:     PetscObjectReference((PetscObject)is->complement);
1539:     *complement = is->complement;
1540:   }
1541:   return(0);
1542: }

1544: /*@
1545:    ISRestoreNonlocalIS - Restore the IS obtained with ISGetNonlocalIS().

1547:    Not collective.

1549:    Input Parameter:
1550: +  is         - the index set
1551: -  complement - index set of is's nonlocal indices

1553:    Level: intermediate

1555: .seealso: ISGetNonlocalIS(), ISGetNonlocalIndices(), ISRestoreNonlocalIndices()
1556: @*/
1557: PetscErrorCode  ISRestoreNonlocalIS(IS is, IS *complement)
1558: {
1560:   PetscInt       refcnt;

1565:   if (*complement != is->complement) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Complement IS being restored was not obtained with ISGetNonlocalIS()");
1566:   PetscObjectGetReference((PetscObject)(is->complement), &refcnt);
1567:   if (refcnt <= 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Duplicate call to ISRestoreNonlocalIS() detected");
1568:   PetscObjectDereference((PetscObject)(is->complement));
1569:   return(0);
1570: }

1572: /*@C
1573:    ISViewFromOptions - View from Options

1575:    Collective on IS

1577:    Input Parameters:
1578: +  A - the index set
1579: .  obj - Optional object
1580: -  name - command line option

1582:    Level: intermediate
1583: .seealso:  IS, ISView, PetscObjectViewFromOptions(), ISCreate()
1584: @*/
1585: PetscErrorCode  ISViewFromOptions(IS A,PetscObject obj,const char name[])
1586: {

1591:   PetscObjectViewFromOptions((PetscObject)A,obj,name);
1592:   return(0);
1593: }

1595: /*@C
1596:    ISView - Displays an index set.

1598:    Collective on IS

1600:    Input Parameters:
1601: +  is - the index set
1602: -  viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.

1604:    Level: intermediate

1606: .seealso: PetscViewerASCIIOpen()
1607: @*/
1608: PetscErrorCode  ISView(IS is,PetscViewer viewer)
1609: {

1614:   if (!viewer) {PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)is),&viewer);}

1618:   PetscObjectPrintClassNamePrefixType((PetscObject)is,viewer);
1619:   PetscLogEventBegin(IS_View,is,viewer,0,0);
1620:   (*is->ops->view)(is,viewer);
1621:   PetscLogEventEnd(IS_View,is,viewer,0,0);
1622:   return(0);
1623: }

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

1628:   Collective on PetscViewer

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

1634:   Level: intermediate

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

1641: .seealso: PetscViewerBinaryOpen(), ISView(), MatLoad(), VecLoad()
1642: @*/
1643: PetscErrorCode ISLoad(IS is, PetscViewer viewer)
1644: {
1645:   PetscBool      isbinary, ishdf5;

1652:   PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERBINARY, &isbinary);
1653:   PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERHDF5, &ishdf5);
1654:   if (!isbinary && !ishdf5) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid viewer; open viewer with PetscViewerBinaryOpen()");
1655:   if (!((PetscObject)is)->type_name) {ISSetType(is, ISGENERAL);}
1656:   PetscLogEventBegin(IS_Load,is,viewer,0,0);
1657:   (*is->ops->load)(is, viewer);
1658:   PetscLogEventEnd(IS_Load,is,viewer,0,0);
1659:   return(0);
1660: }

1662: /*@
1663:    ISSort - Sorts the indices of an index set.

1665:    Collective on IS

1667:    Input Parameters:
1668: .  is - the index set

1670:    Level: intermediate

1672: .seealso: ISSortRemoveDups(), ISSorted()
1673: @*/
1674: PetscErrorCode  ISSort(IS is)
1675: {

1680:   (*is->ops->sort)(is);
1681:   ISSetInfo(is,IS_SORTED,IS_LOCAL,is->info_permanent[IS_LOCAL][IS_SORTED],PETSC_TRUE);
1682:   return(0);
1683: }

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

1688:   Collective on IS

1690:   Input Parameters:
1691: . is - the index set

1693:   Level: intermediate

1695: .seealso: ISSort(), ISSorted()
1696: @*/
1697: PetscErrorCode ISSortRemoveDups(IS is)
1698: {

1703:   ISClearInfoCache(is,PETSC_FALSE);
1704:   (*is->ops->sortremovedups)(is);
1705:   ISSetInfo(is,IS_SORTED,IS_LOCAL,is->info_permanent[IS_LOCAL][IS_SORTED],PETSC_TRUE);
1706:   ISSetInfo(is,IS_UNIQUE,IS_LOCAL,is->info_permanent[IS_LOCAL][IS_UNIQUE],PETSC_TRUE);
1707:   return(0);
1708: }

1710: /*@
1711:    ISToGeneral - Converts an IS object of any type to ISGENERAL type

1713:    Collective on IS

1715:    Input Parameters:
1716: .  is - the index set

1718:    Level: intermediate

1720: .seealso: ISSorted()
1721: @*/
1722: PetscErrorCode  ISToGeneral(IS is)
1723: {

1728:   if (is->ops->togeneral) {
1729:     (*is->ops->togeneral)(is);
1730:   } else SETERRQ1(PetscObjectComm((PetscObject)is),PETSC_ERR_SUP,"Not written for this type %s",((PetscObject)is)->type_name);
1731:   return(0);
1732: }

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

1737:    Collective on IS

1739:    Input Parameter:
1740: .  is - the index set

1742:    Output Parameter:
1743: .  flg - output flag, either PETSC_TRUE if the index set is sorted,
1744:          or PETSC_FALSE otherwise.

1746:    Notes:
1747:     For parallel IS objects this only indicates if the local part of the IS
1748:           is sorted. So some processors may return PETSC_TRUE while others may
1749:           return PETSC_FALSE.

1751:    Level: intermediate

1753: .seealso: ISSort(), ISSortRemoveDups()
1754: @*/
1755: PetscErrorCode  ISSorted(IS is,PetscBool  *flg)
1756: {

1762:   ISGetInfo(is,IS_SORTED,IS_LOCAL,PETSC_TRUE,flg);
1763:   return(0);
1764: }

1766: /*@
1767:    ISDuplicate - Creates a duplicate copy of an index set.

1769:    Collective on IS

1771:    Input Parmeters:
1772: .  is - the index set

1774:    Output Parameters:
1775: .  isnew - the copy of the index set

1777:    Level: beginner

1779: .seealso: ISCreateGeneral(), ISCopy()
1780: @*/
1781: PetscErrorCode  ISDuplicate(IS is,IS *newIS)
1782: {

1788:   (*is->ops->duplicate)(is,newIS);
1789:   ISCopyInfo(is,*newIS);
1790:   return(0);
1791: }

1793: /*@
1794:    ISCopy - Copies an index set.

1796:    Collective on IS

1798:    Input Parmeters:
1799: .  is - the index set

1801:    Output Parameters:
1802: .  isy - the copy of the index set

1804:    Level: beginner

1806: .seealso: ISDuplicate()
1807: @*/
1808: PetscErrorCode  ISCopy(IS is,IS isy)
1809: {

1816:   if (is == isy) return(0);
1817:   ISCopyInfo(is,isy);
1818:   isy->max        = is->max;
1819:   isy->min        = is->min;
1820:   (*is->ops->copy)(is,isy);
1821:   return(0);
1822: }

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

1827:    Collective on IS

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

1834:    Output Arguments:
1835: . newis - new IS on comm

1837:    Level: advanced

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

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

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

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

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

1868:    Logicall Collective on IS

1870:    Input Arguments:
1871: + is - index set
1872: - bs - block size

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: ISGetBlockSize(), ISCreateBlock(), ISBlockGetIndices(),
1883: @*/
1884: PetscErrorCode  ISSetBlockSize(IS is,PetscInt bs)
1885: {

1891:   if (bs < 1) SETERRQ1(PetscObjectComm((PetscObject)is),PETSC_ERR_ARG_OUTOFRANGE,"Block size %D, must be positive",bs);
1892:   (*is->ops->setblocksize)(is,bs);
1893:   return(0);
1894: }

1896: /*@
1897:    ISGetBlockSize - Returns the number of elements in a block.

1899:    Not Collective

1901:    Input Parameter:
1902: .  is - the index set

1904:    Output Parameter:
1905: .  size - the number of elements in a block

1907:    Level: intermediate

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

1915: .seealso: ISBlockGetSize(), ISGetSize(), ISCreateBlock(), ISSetBlockSize()
1916: @*/
1917: PetscErrorCode  ISGetBlockSize(IS is,PetscInt *size)
1918: {

1922:   PetscLayoutGetBlockSize(is->map, size);
1923:   return(0);
1924: }

1926: PetscErrorCode ISGetIndicesCopy(IS is, PetscInt idx[])
1927: {
1929:   PetscInt       len,i;
1930:   const PetscInt *ptr;

1933:   ISGetLocalSize(is,&len);
1934:   ISGetIndices(is,&ptr);
1935:   for (i=0; i<len; i++) idx[i] = ptr[i];
1936:   ISRestoreIndices(is,&ptr);
1937:   return(0);
1938: }

1940: /*MC
1941:     ISGetIndicesF90 - Accesses the elements of an index set from Fortran90.
1942:     The users should call ISRestoreIndicesF90() after having looked at the
1943:     indices.  The user should NOT change the indices.

1945:     Synopsis:
1946:     ISGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

1948:     Not collective

1950:     Input Parameter:
1951: .   x - index set

1953:     Output Parameters:
1954: +   xx_v - the Fortran90 pointer to the array
1955: -   ierr - error code

1957:     Example of Usage:
1958: .vb
1959:     PetscInt, pointer xx_v(:)
1960:     ....
1961:     call ISGetIndicesF90(x,xx_v,ierr)
1962:     a = xx_v(3)
1963:     call ISRestoreIndicesF90(x,xx_v,ierr)
1964: .ve

1966:     Level: intermediate

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

1970: M*/

1972: /*MC
1973:     ISRestoreIndicesF90 - Restores an index set to a usable state after
1974:     a call to ISGetIndicesF90().

1976:     Synopsis:
1977:     ISRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

1979:     Not collective

1981:     Input Parameters:
1982: +   x - index set
1983: -   xx_v - the Fortran90 pointer to the array

1985:     Output Parameter:
1986: .   ierr - error code

1988:     Example of Usage:
1989: .vb
1990:     PetscInt, pointer xx_v(:)
1991:     ....
1992:     call ISGetIndicesF90(x,xx_v,ierr)
1993:     a = xx_v(3)
1994:     call ISRestoreIndicesF90(x,xx_v,ierr)
1995: .ve

1997:     Level: intermediate

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

2001: M*/

2003: /*MC
2004:     ISBlockGetIndicesF90 - Accesses the elements of an index set from Fortran90.
2005:     The users should call ISBlockRestoreIndicesF90() after having looked at the
2006:     indices.  The user should NOT change the indices.

2008:     Synopsis:
2009:     ISBlockGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

2011:     Not collective

2013:     Input Parameter:
2014: .   x - index set

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

2028:     Level: intermediate

2030: .seealso:  ISBlockRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices(),
2031:            ISRestoreIndices()

2033: M*/

2035: /*MC
2036:     ISBlockRestoreIndicesF90 - Restores an index set to a usable state after
2037:     a call to ISBlockGetIndicesF90().

2039:     Synopsis:
2040:     ISBlockRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

2042:     Not Collective

2044:     Input Parameters:
2045: +   x - index set
2046: -   xx_v - the Fortran90 pointer to the array

2048:     Output Parameter:
2049: .   ierr - error code

2051:     Example of Usage:
2052: .vb
2053:     PetscInt, pointer xx_v(:)
2054:     ....
2055:     call ISBlockGetIndicesF90(x,xx_v,ierr)
2056:     a = xx_v(3)
2057:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
2058: .ve

2060:     Notes:
2061:     Not yet supported for all F90 compilers

2063:     Level: intermediate

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

2067: M*/