Actual source code: index.c

petsc-master 2020-08-09
Report Typos and Errors
  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);
150:   PetscSFBcastEnd(sf,MPIU_INT,root_data,leaf_data);
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: }


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

179:    Collective on IS

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

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

188:    Level: intermediate

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

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


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

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

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

254:    Not collective

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

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

263:    Level: developer

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

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

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

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

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

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

389:    Logically Collective on IS if type is IS_GLOBAL

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

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


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

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

413:    Level: advanced

415: .seealso:  ISInfo, ISInfoType, IS

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

427:   comm = PetscObjectComm((PetscObject)is);
428:   if (type == IS_GLOBAL) {
432:     errcomm = comm;
433:   } else {
434:     errcomm = PETSC_COMM_SELF;
435:   }

437:   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);

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

721: /*@
722:    ISGetInfo - Determine whether an index set satisfies a given property

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

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

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

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

737:    Level: advanced

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

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

754:   comm = PetscObjectComm((PetscObject)is);
755:   if (type == IS_GLOBAL) {
757:     errcomm = comm;
758:   } else {
759:     errcomm = PETSC_COMM_SELF;
760:   }

762:   MPI_Comm_size(comm, &size);
763:   MPI_Comm_rank(comm, &rank);

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

807: static PetscErrorCode ISCopyInfo(IS source, IS dest)
808: {

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

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

820:    Collective on IS

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

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

828:    Level: intermediate

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

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

845:   ISGetInfo(is,IS_IDENTITY,IS_GLOBAL,PETSC_TRUE,ident);
846:   return(0);
847: }

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

852:    Logically Collective on IS

854:    Input Parmeters:
855: .  is - the index set

857:    Level: intermediate

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

863: .seealso: ISIdentity(), ISSetInfo(), ISClearInfoCache()
864: @*/
865: PetscErrorCode  ISSetIdentity(IS is)
866: {

871:   ISSetInfo(is,IS_IDENTITY,IS_GLOBAL,PETSC_TRUE,PETSC_TRUE);
872:   return(0);
873: }

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

878:    Not Collective

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

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

889:    Level: developer

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

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

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

914:    Logically Collective on IS

916:    Input Parmeters:
917: .  is - the index set

919:    Output Parameters:
920: .  perm - PETSC_TRUE if a permutation, else PETSC_FALSE

922:    Level: intermediate

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

930: .seealso: ISSetPermutation(), ISGetInfo()
931: @*/
932: PetscErrorCode  ISPermutation(IS is,PetscBool  *perm)
933: {

939:   ISGetInfo(is,IS_PERMUTATION,IS_GLOBAL,PETSC_FALSE,perm);
940:   return(0);
941: }

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

946:    Logically Collective on IS

948:    Input Parmeters:
949: .  is - the index set

951:    Level: intermediate


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

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

961: .seealso: ISPermutation(), ISSetInfo(), ISClearInfoCache().
962: @*/
963: PetscErrorCode  ISSetPermutation(IS is)
964: {

969:   if (PetscDefined(USE_DEBUG)) {
970:     PetscMPIInt    size;

972:     MPI_Comm_size(PetscObjectComm((PetscObject)is),&size);
973:     if (size == 1) {
974:       PetscInt       i,n,*idx;
975:       const PetscInt *iidx;

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

993: /*@
994:    ISDestroy - Destroys an index set.

996:    Collective on IS

998:    Input Parameters:
999: .  is - the index set

1001:    Level: beginner

1003: .seealso: ISCreateGeneral(), ISCreateStride(), ISCreateBlocked()
1004: @*/
1005: PetscErrorCode  ISDestroy(IS *is)
1006: {

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

1030: /*@
1031:    ISInvertPermutation - Creates a new permutation that is the inverse of
1032:                          a given permutation.

1034:    Collective on IS

1036:    Input Parameter:
1037: +  is - the index set
1038: -  nlocal - number of indices on this processor in result (ignored for 1 proccessor) or
1039:             use PETSC_DECIDE

1041:    Output Parameter:
1042: .  isout - the inverse permutation

1044:    Level: intermediate

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

1050: @*/
1051: PetscErrorCode  ISInvertPermutation(IS is,PetscInt nlocal,IS *isout)
1052: {
1053:   PetscBool      isperm, isidentity, issame;

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

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

1081: /*@
1082:    ISGetSize - Returns the global length of an index set.

1084:    Not Collective

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

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

1092:    Level: beginner


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

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

1108:    Not Collective

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

1113:    Output Parameter:
1114: .  size - the local size

1116:    Level: beginner

1118: @*/
1119: PetscErrorCode  ISGetLocalSize(IS is,PetscInt *size)
1120: {
1124:   *size = is->map->n;
1125:   return(0);
1126: }

1128: /*@
1129:    ISGetLayout - get PetscLayout describing index set layout

1131:    Not Collective

1133:    Input Arguments:
1134: .  is - the index set

1136:    Output Arguments:
1137: .  map - the layout

1139:    Level: developer

1141: .seealso: ISGetSize(), ISGetLocalSize()
1142: @*/
1143: PetscErrorCode ISGetLayout(IS is,PetscLayout *map)
1144: {

1149:   *map = is->map;
1150:   return(0);
1151: }

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

1158:    Not Collective

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

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

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



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

1191:    Level: intermediate


1194: .seealso: ISRestoreIndices(), ISGetIndicesF90()
1195: @*/
1196: PetscErrorCode  ISGetIndices(IS is,const PetscInt *ptr[])
1197: {

1203:   (*is->ops->getindices)(is,ptr);
1204:   return(0);
1205: }

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

1210:    Not Collective

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

1215:    Output Parameter:
1216: +   min - the minimum value
1217: -   max - the maximum value

1219:    Level: intermediate

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


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

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

1240:   Not Collective

1242:   Input Parameter:
1243: + is - the index set
1244: - key - the search key

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

1249:   Level: intermediate
1250: @*/
1251: PetscErrorCode ISLocate(IS is, PetscInt key, PetscInt *location)
1252: {

1256:   if (is->ops->locate) {
1257:     (*is->ops->locate)(is,key,location);
1258:   } else {
1259:     PetscInt       numIdx;
1260:     PetscBool      sorted;
1261:     const PetscInt *idx;

1263:     ISGetLocalSize(is,&numIdx);
1264:     ISGetIndices(is,&idx);
1265:     ISSorted(is,&sorted);
1266:     if (sorted) {
1267:       PetscFindInt(key,numIdx,idx,location);
1268:     } else {
1269:       PetscInt i;

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

1284: /*@C
1285:    ISRestoreIndices - Restores an index set to a usable state after a call
1286:                       to ISGetIndices().

1288:    Not Collective

1290:    Input Parameters:
1291: +  is - the index set
1292: -  ptr - the pointer obtained by ISGetIndices()

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

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

1311:    Level: intermediate

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

1316: .seealso: ISGetIndices(), ISRestoreIndicesF90()
1317: @*/
1318: PetscErrorCode  ISRestoreIndices(IS is,const PetscInt *ptr[])
1319: {

1325:   if (is->ops->restoreindices) {
1326:     (*is->ops->restoreindices)(is,ptr);
1327:   }
1328:   return(0);
1329: }

1331: static PetscErrorCode ISGatherTotal_Private(IS is)
1332: {
1334:   PetscInt       i,n,N;
1335:   const PetscInt *lindices;
1336:   MPI_Comm       comm;
1337:   PetscMPIInt    rank,size,*sizes = NULL,*offsets = NULL,nn;


1342:   PetscObjectGetComm((PetscObject)is,&comm);
1343:   MPI_Comm_size(comm,&size);
1344:   MPI_Comm_rank(comm,&rank);
1345:   ISGetLocalSize(is,&n);
1346:   PetscMalloc2(size,&sizes,size,&offsets);

1348:   PetscMPIIntCast(n,&nn);
1349:   MPI_Allgather(&nn,1,MPI_INT,sizes,1,MPI_INT,comm);
1350:   offsets[0] = 0;
1351:   for (i=1; i<size; ++i) offsets[i] = offsets[i-1] + sizes[i-1];
1352:   N = offsets[size-1] + sizes[size-1];

1354:   PetscMalloc1(N,&(is->total));
1355:   ISGetIndices(is,&lindices);
1356:   MPI_Allgatherv((void*)lindices,nn,MPIU_INT,is->total,sizes,offsets,MPIU_INT,comm);
1357:   ISRestoreIndices(is,&lindices);
1358:   is->local_offset = offsets[rank];
1359:   PetscFree2(sizes,offsets);
1360:   return(0);
1361: }

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

1366:    Collective on IS

1368:    Input Parameter:
1369: .  is - the index set

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

1375:    Level: intermediate

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

1385: .seealso: ISRestoreTotalIndices(), ISGetNonlocalIndices(), ISGetSize()
1386: @*/
1387: PetscErrorCode ISGetTotalIndices(IS is, const PetscInt *indices[])
1388: {
1390:   PetscMPIInt    size;

1395:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
1396:   if (size == 1) {
1397:     (*is->ops->getindices)(is,indices);
1398:   } else {
1399:     if (!is->total) {
1400:       ISGatherTotal_Private(is);
1401:     }
1402:     *indices = is->total;
1403:   }
1404:   return(0);
1405: }

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

1410:    Not Collective.

1412:    Input Parameter:
1413: +  is - the index set
1414: -  indices - index array; must be the array obtained with ISGetTotalIndices()

1416:    Level: intermediate

1418: .seealso: ISRestoreTotalIndices(), ISGetNonlocalIndices()
1419: @*/
1420: PetscErrorCode  ISRestoreTotalIndices(IS is, const PetscInt *indices[])
1421: {
1423:   PetscMPIInt    size;

1428:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
1429:   if (size == 1) {
1430:     (*is->ops->restoreindices)(is,indices);
1431:   } else {
1432:     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.");
1433:   }
1434:   return(0);
1435: }
1436: /*@C
1437:    ISGetNonlocalIndices - Retrieve an array of indices from remote processors
1438:                        in this communicator.

1440:    Collective on IS

1442:    Input Parameter:
1443: .  is - the index set

1445:    Output Parameter:
1446: .  indices - indices with rank 0 indices first, and so on,  omitting
1447:              the current rank.  Total number of indices is the difference
1448:              total and local, obtained with ISGetSize() and ISGetLocalSize(),
1449:              respectively.

1451:    Level: intermediate

1453:    Notes:
1454:     restore the indices using ISRestoreNonlocalIndices().
1455:           The same scalability considerations as those for ISGetTotalIndices
1456:           apply here.

1458: .seealso: ISGetTotalIndices(), ISRestoreNonlocalIndices(), ISGetSize(), ISGetLocalSize().
1459: @*/
1460: PetscErrorCode  ISGetNonlocalIndices(IS is, const PetscInt *indices[])
1461: {
1463:   PetscMPIInt    size;
1464:   PetscInt       n, N;

1469:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
1470:   if (size == 1) *indices = NULL;
1471:   else {
1472:     if (!is->total) {
1473:       ISGatherTotal_Private(is);
1474:     }
1475:     ISGetLocalSize(is,&n);
1476:     ISGetSize(is,&N);
1477:     PetscMalloc1(N-n, &(is->nonlocal));
1478:     PetscArraycpy(is->nonlocal, is->total, is->local_offset);
1479:     PetscArraycpy(is->nonlocal+is->local_offset, is->total+is->local_offset+n,N - is->local_offset - n);
1480:     *indices = is->nonlocal;
1481:   }
1482:   return(0);
1483: }

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

1488:    Not Collective.

1490:    Input Parameter:
1491: +  is - the index set
1492: -  indices - index array; must be the array obtained with ISGetNonlocalIndices()

1494:    Level: intermediate

1496: .seealso: ISGetTotalIndices(), ISGetNonlocalIndices(), ISRestoreTotalIndices()
1497: @*/
1498: PetscErrorCode  ISRestoreNonlocalIndices(IS is, const PetscInt *indices[])
1499: {
1503:   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.");
1504:   return(0);
1505: }

1507: /*@
1508:    ISGetNonlocalIS - Gather all nonlocal indices for this IS and present
1509:                      them as another sequential index set.


1512:    Collective on IS

1514:    Input Parameter:
1515: .  is - the index set

1517:    Output Parameter:
1518: .  complement - sequential IS with indices identical to the result of
1519:                 ISGetNonlocalIndices()

1521:    Level: intermediate

1523:    Notes:
1524:     complement represents the result of ISGetNonlocalIndices as an IS.
1525:           Therefore scalability issues similar to ISGetNonlocalIndices apply.
1526:           The resulting IS must be restored using ISRestoreNonlocalIS().

1528: .seealso: ISGetNonlocalIndices(), ISRestoreNonlocalIndices(),  ISAllGather(), ISGetSize()
1529: @*/
1530: PetscErrorCode  ISGetNonlocalIS(IS is, IS *complement)
1531: {

1537:   /* Check if the complement exists already. */
1538:   if (is->complement) {
1539:     *complement = is->complement;
1540:     PetscObjectReference((PetscObject)(is->complement));
1541:   } else {
1542:     PetscInt       N, n;
1543:     const PetscInt *idx;
1544:     ISGetSize(is, &N);
1545:     ISGetLocalSize(is,&n);
1546:     ISGetNonlocalIndices(is, &idx);
1547:     ISCreateGeneral(PETSC_COMM_SELF, N-n,idx, PETSC_USE_POINTER, &(is->complement));
1548:     PetscObjectReference((PetscObject)is->complement);
1549:     *complement = is->complement;
1550:   }
1551:   return(0);
1552: }


1555: /*@
1556:    ISRestoreNonlocalIS - Restore the IS obtained with ISGetNonlocalIS().

1558:    Not collective.

1560:    Input Parameter:
1561: +  is         - the index set
1562: -  complement - index set of is's nonlocal indices

1564:    Level: intermediate


1567: .seealso: ISGetNonlocalIS(), ISGetNonlocalIndices(), ISRestoreNonlocalIndices()
1568: @*/
1569: PetscErrorCode  ISRestoreNonlocalIS(IS is, IS *complement)
1570: {
1572:   PetscInt       refcnt;

1577:   if (*complement != is->complement) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Complement IS being restored was not obtained with ISGetNonlocalIS()");
1578:   PetscObjectGetReference((PetscObject)(is->complement), &refcnt);
1579:   if (refcnt <= 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Duplicate call to ISRestoreNonlocalIS() detected");
1580:   PetscObjectDereference((PetscObject)(is->complement));
1581:   return(0);
1582: }

1584: /*@C
1585:    ISViewFromOptions - View from Options

1587:    Collective on IS

1589:    Input Parameters:
1590: +  A - the index set
1591: .  obj - Optional object
1592: -  name - command line option

1594:    Level: intermediate
1595: .seealso:  IS, ISView, PetscObjectViewFromOptions(), ISCreate()
1596: @*/
1597: PetscErrorCode  ISViewFromOptions(IS A,PetscObject obj,const char name[])
1598: {

1603:   PetscObjectViewFromOptions((PetscObject)A,obj,name);
1604:   return(0);
1605: }

1607: /*@C
1608:    ISView - Displays an index set.

1610:    Collective on IS

1612:    Input Parameters:
1613: +  is - the index set
1614: -  viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.

1616:    Level: intermediate

1618: .seealso: PetscViewerASCIIOpen()
1619: @*/
1620: PetscErrorCode  ISView(IS is,PetscViewer viewer)
1621: {

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

1630:   PetscObjectPrintClassNamePrefixType((PetscObject)is,viewer);
1631:   PetscLogEventBegin(IS_View,is,viewer,0,0);
1632:   (*is->ops->view)(is,viewer);
1633:   PetscLogEventEnd(IS_View,is,viewer,0,0);
1634:   return(0);
1635: }

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

1640:   Collective on PetscViewer

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

1646:   Level: intermediate

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

1653: .seealso: PetscViewerBinaryOpen(), ISView(), MatLoad(), VecLoad()
1654: @*/
1655: PetscErrorCode ISLoad(IS is, PetscViewer viewer)
1656: {
1657:   PetscBool      isbinary, ishdf5;

1664:   PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERBINARY, &isbinary);
1665:   PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERHDF5, &ishdf5);
1666:   if (!isbinary && !ishdf5) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid viewer; open viewer with PetscViewerBinaryOpen()");
1667:   if (!((PetscObject)is)->type_name) {ISSetType(is, ISGENERAL);}
1668:   PetscLogEventBegin(IS_Load,is,viewer,0,0);
1669:   (*is->ops->load)(is, viewer);
1670:   PetscLogEventEnd(IS_Load,is,viewer,0,0);
1671:   return(0);
1672: }

1674: /*@
1675:    ISSort - Sorts the indices of an index set.

1677:    Collective on IS

1679:    Input Parameters:
1680: .  is - the index set

1682:    Level: intermediate


1685: .seealso: ISSortRemoveDups(), ISSorted()
1686: @*/
1687: PetscErrorCode  ISSort(IS is)
1688: {

1693:   (*is->ops->sort)(is);
1694:   ISSetInfo(is,IS_SORTED,IS_LOCAL,is->info_permanent[IS_LOCAL][IS_SORTED],PETSC_TRUE);
1695:   return(0);
1696: }

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

1701:   Collective on IS

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

1706:   Level: intermediate


1709: .seealso: ISSort(), ISSorted()
1710: @*/
1711: PetscErrorCode ISSortRemoveDups(IS is)
1712: {

1717:   ISClearInfoCache(is,PETSC_FALSE);
1718:   (*is->ops->sortremovedups)(is);
1719:   ISSetInfo(is,IS_SORTED,IS_LOCAL,is->info_permanent[IS_LOCAL][IS_SORTED],PETSC_TRUE);
1720:   ISSetInfo(is,IS_UNIQUE,IS_LOCAL,is->info_permanent[IS_LOCAL][IS_UNIQUE],PETSC_TRUE);
1721:   return(0);
1722: }

1724: /*@
1725:    ISToGeneral - Converts an IS object of any type to ISGENERAL type

1727:    Collective on IS

1729:    Input Parameters:
1730: .  is - the index set

1732:    Level: intermediate


1735: .seealso: ISSorted()
1736: @*/
1737: PetscErrorCode  ISToGeneral(IS is)
1738: {

1743:   if (is->ops->togeneral) {
1744:     (*is->ops->togeneral)(is);
1745:   } else SETERRQ1(PetscObjectComm((PetscObject)is),PETSC_ERR_SUP,"Not written for this type %s",((PetscObject)is)->type_name);
1746:   return(0);
1747: }

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

1752:    Collective on IS

1754:    Input Parameter:
1755: .  is - the index set

1757:    Output Parameter:
1758: .  flg - output flag, either PETSC_TRUE if the index set is sorted,
1759:          or PETSC_FALSE otherwise.

1761:    Notes:
1762:     For parallel IS objects this only indicates if the local part of the IS
1763:           is sorted. So some processors may return PETSC_TRUE while others may
1764:           return PETSC_FALSE.

1766:    Level: intermediate

1768: .seealso: ISSort(), ISSortRemoveDups()
1769: @*/
1770: PetscErrorCode  ISSorted(IS is,PetscBool  *flg)
1771: {

1777:   ISGetInfo(is,IS_SORTED,IS_LOCAL,PETSC_TRUE,flg);
1778:   return(0);
1779: }

1781: /*@
1782:    ISDuplicate - Creates a duplicate copy of an index set.

1784:    Collective on IS

1786:    Input Parmeters:
1787: .  is - the index set

1789:    Output Parameters:
1790: .  isnew - the copy of the index set

1792:    Level: beginner

1794: .seealso: ISCreateGeneral(), ISCopy()
1795: @*/
1796: PetscErrorCode  ISDuplicate(IS is,IS *newIS)
1797: {

1803:   (*is->ops->duplicate)(is,newIS);
1804:   ISCopyInfo(is,*newIS);
1805:   return(0);
1806: }

1808: /*@
1809:    ISCopy - Copies an index set.

1811:    Collective on IS

1813:    Input Parmeters:
1814: .  is - the index set

1816:    Output Parameters:
1817: .  isy - the copy of the index set

1819:    Level: beginner

1821: .seealso: ISDuplicate()
1822: @*/
1823: PetscErrorCode  ISCopy(IS is,IS isy)
1824: {

1831:   if (is == isy) return(0);
1832:   ISCopyInfo(is,isy);
1833:   isy->max        = is->max;
1834:   isy->min        = is->min;
1835:   (*is->ops->copy)(is,isy);
1836:   return(0);
1837: }

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

1842:    Collective on IS

1844:    Input Arguments:
1845: + is - index set
1846: . comm - communicator for new index set
1847: - mode - copy semantics, PETSC_USE_POINTER for no-copy if possible, otherwise PETSC_COPY_VALUES

1849:    Output Arguments:
1850: . newis - new IS on comm

1852:    Level: advanced

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

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

1860:    The input IS must have the same type on every process.

1862: .seealso: ISSplit()
1863: @*/
1864: PetscErrorCode  ISOnComm(IS is,MPI_Comm comm,PetscCopyMode mode,IS *newis)
1865: {
1867:   PetscMPIInt    match;

1872:   MPI_Comm_compare(PetscObjectComm((PetscObject)is),comm,&match);
1873:   if (mode != PETSC_COPY_VALUES && (match == MPI_IDENT || match == MPI_CONGRUENT)) {
1874:     PetscObjectReference((PetscObject)is);
1875:     *newis = is;
1876:   } else {
1877:     (*is->ops->oncomm)(is,comm,mode,newis);
1878:   }
1879:   return(0);
1880: }

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

1885:    Logicall Collective on IS

1887:    Input Arguments:
1888: + is - index set
1889: - bs - block size

1891:    Level: intermediate
1892:    
1893:    Notes: 
1894:    This is much like the block size for Vecs. It indicates that one can think of the indices as 
1895:    being in a collection of equal size blocks. For ISBlock() these collections of blocks are all contiquous
1896:    within a block but this is not the case for other IS.
1897:    ISBlockGetIndices() only works for ISBlock IS, not others.

1899: .seealso: ISGetBlockSize(), ISCreateBlock(), ISBlockGetIndices(), 
1900: @*/
1901: PetscErrorCode  ISSetBlockSize(IS is,PetscInt bs)
1902: {

1908:   if (bs < 1) SETERRQ1(PetscObjectComm((PetscObject)is),PETSC_ERR_ARG_OUTOFRANGE,"Block size %D, must be positive",bs);
1909:   (*is->ops->setblocksize)(is,bs);
1910:   return(0);
1911: }

1913: /*@
1914:    ISGetBlockSize - Returns the number of elements in a block.

1916:    Not Collective

1918:    Input Parameter:
1919: .  is - the index set

1921:    Output Parameter:
1922: .  size - the number of elements in a block

1924:    Level: intermediate

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

1932: .seealso: ISBlockGetSize(), ISGetSize(), ISCreateBlock(), ISSetBlockSize()
1933: @*/
1934: PetscErrorCode  ISGetBlockSize(IS is,PetscInt *size)
1935: {

1939:   PetscLayoutGetBlockSize(is->map, size);
1940:   return(0);
1941: }

1943: PetscErrorCode ISGetIndicesCopy(IS is, PetscInt idx[])
1944: {
1946:   PetscInt       len,i;
1947:   const PetscInt *ptr;

1950:   ISGetLocalSize(is,&len);
1951:   ISGetIndices(is,&ptr);
1952:   for (i=0; i<len; i++) idx[i] = ptr[i];
1953:   ISRestoreIndices(is,&ptr);
1954:   return(0);
1955: }

1957: /*MC
1958:     ISGetIndicesF90 - Accesses the elements of an index set from Fortran90.
1959:     The users should call ISRestoreIndicesF90() after having looked at the
1960:     indices.  The user should NOT change the indices.

1962:     Synopsis:
1963:     ISGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

1965:     Not collective

1967:     Input Parameter:
1968: .   x - index set

1970:     Output Parameters:
1971: +   xx_v - the Fortran90 pointer to the array
1972: -   ierr - error code

1974:     Example of Usage:
1975: .vb
1976:     PetscInt, pointer xx_v(:)
1977:     ....
1978:     call ISGetIndicesF90(x,xx_v,ierr)
1979:     a = xx_v(3)
1980:     call ISRestoreIndicesF90(x,xx_v,ierr)
1981: .ve

1983:     Level: intermediate

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


1988: M*/

1990: /*MC
1991:     ISRestoreIndicesF90 - Restores an index set to a usable state after
1992:     a call to ISGetIndicesF90().

1994:     Synopsis:
1995:     ISRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

1997:     Not collective

1999:     Input Parameters:
2000: +   x - index set
2001: -   xx_v - the Fortran90 pointer to the array

2003:     Output Parameter:
2004: .   ierr - error code


2007:     Example of Usage:
2008: .vb
2009:     PetscInt, pointer xx_v(:)
2010:     ....
2011:     call ISGetIndicesF90(x,xx_v,ierr)
2012:     a = xx_v(3)
2013:     call ISRestoreIndicesF90(x,xx_v,ierr)
2014: .ve

2016:     Level: intermediate

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

2020: M*/

2022: /*MC
2023:     ISBlockGetIndicesF90 - Accesses the elements of an index set from Fortran90.
2024:     The users should call ISBlockRestoreIndicesF90() after having looked at the
2025:     indices.  The user should NOT change the indices.

2027:     Synopsis:
2028:     ISBlockGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

2030:     Not collective

2032:     Input Parameter:
2033: .   x - index set

2035:     Output Parameters:
2036: +   xx_v - the Fortran90 pointer to the array
2037: -   ierr - error code
2038:     Example of Usage:
2039: .vb
2040:     PetscInt, pointer xx_v(:)
2041:     ....
2042:     call ISBlockGetIndicesF90(x,xx_v,ierr)
2043:     a = xx_v(3)
2044:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
2045: .ve

2047:     Level: intermediate

2049: .seealso:  ISBlockRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices(),
2050:            ISRestoreIndices()

2052: M*/

2054: /*MC
2055:     ISBlockRestoreIndicesF90 - Restores an index set to a usable state after
2056:     a call to ISBlockGetIndicesF90().

2058:     Synopsis:
2059:     ISBlockRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

2061:     Not Collective

2063:     Input Parameters:
2064: +   x - index set
2065: -   xx_v - the Fortran90 pointer to the array

2067:     Output Parameter:
2068: .   ierr - error code

2070:     Example of Usage:
2071: .vb
2072:     PetscInt, pointer xx_v(:)
2073:     ....
2074:     call ISBlockGetIndicesF90(x,xx_v,ierr)
2075:     a = xx_v(3)
2076:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
2077: .ve

2079:     Notes:
2080:     Not yet supported for all F90 compilers

2082:     Level: intermediate

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

2086: M*/