Actual source code: index.c

petsc-master 2016-06-26
Report Typos and Errors
  2: /*
  3:    Defines the abstract operations on index sets, i.e. the public interface.
  4: */
  5: #include <petsc/private/isimpl.h>      /*I "petscis.h" I*/
  6: #include <petscviewer.h>

  8: /* Logging support */
  9: PetscClassId IS_CLASSID;

 13: /*@
 14:    ISIdentity - Determines whether index set is the identity mapping.

 16:    Collective on IS

 18:    Input Parmeters:
 19: .  is - the index set

 21:    Output Parameters:
 22: .  ident - PETSC_TRUE if an identity, else PETSC_FALSE

 24:    Level: intermediate

 26:    Concepts: identity mapping
 27:    Concepts: index sets^is identity

 29: .seealso: ISSetIdentity()
 30: @*/
 31: PetscErrorCode  ISIdentity(IS is,PetscBool  *ident)
 32: {

 38:   *ident = is->isidentity;
 39:   if (*ident) return(0);
 40:   if (is->ops->identity) {
 41:     (*is->ops->identity)(is,ident);
 42:   }
 43:   return(0);
 44: }

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

 51:    Logically Collective on IS

 53:    Input Parmeters:
 54: .  is - the index set

 56:    Level: intermediate

 58:    Concepts: identity mapping
 59:    Concepts: index sets^is identity

 61: .seealso: ISIdentity()
 62: @*/
 63: PetscErrorCode  ISSetIdentity(IS is)
 64: {

 69:   is->isidentity = PETSC_TRUE;
 70:   ISSetPermutation(is);
 71:   return(0);
 72: }

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

 79:    Not Collective

 81:    Input Parmeters:
 82: +  is - the index set
 83: .  gstart - global start
 84: .  gend - global end

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

 90:    Level: developer

 92:    Concepts: index sets^is contiguous

 94: .seealso: ISGetLocalSize(), VecGetOwnershipRange()
 95: @*/
 96: PetscErrorCode  ISContiguousLocal(IS is,PetscInt gstart,PetscInt gend,PetscInt *start,PetscBool *contig)
 97: {

104:   if (is->ops->contiguous) {
105:     (*is->ops->contiguous)(is,gstart,gend,start,contig);
106:   } else {
107:     *start  = -1;
108:     *contig = PETSC_FALSE;
109:   }
110:   return(0);
111: }

115: /*@
116:    ISPermutation - PETSC_TRUE or PETSC_FALSE depending on whether the
117:    index set has been declared to be a permutation.

119:    Logically Collective on IS

121:    Input Parmeters:
122: .  is - the index set

124:    Output Parameters:
125: .  perm - PETSC_TRUE if a permutation, else PETSC_FALSE

127:    Level: intermediate

129:   Concepts: permutation
130:   Concepts: index sets^is permutation

132: .seealso: ISSetPermutation()
133: @*/
134: PetscErrorCode  ISPermutation(IS is,PetscBool  *perm)
135: {
139:   *perm = (PetscBool) is->isperm;
140:   return(0);
141: }

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

148:    Logically Collective on IS

150:    Input Parmeters:
151: .  is - the index set

153:    Level: intermediate

155:   Concepts: permutation
156:   Concepts: index sets^permutation

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

161: .seealso: ISPermutation()
162: @*/
163: PetscErrorCode  ISSetPermutation(IS is)
164: {
167: #if defined(PETSC_USE_DEBUG)
168:   {
169:     PetscMPIInt    size;

172:     MPI_Comm_size(PetscObjectComm((PetscObject)is),&size);
173:     if (size == 1) {
174:       PetscInt       i,n,*idx;
175:       const PetscInt *iidx;

177:       ISGetSize(is,&n);
178:       PetscMalloc1(n,&idx);
179:       ISGetIndices(is,&iidx);
180:       PetscMemcpy(idx,iidx,n*sizeof(PetscInt));
181:       PetscSortInt(n,idx);
182:       for (i=0; i<n; i++) {
183:         if (idx[i] != i) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Index set is not a permutation");
184:       }
185:       PetscFree(idx);
186:       ISRestoreIndices(is,&iidx);
187:     }
188:   }
189: #endif
190:   is->isperm = PETSC_TRUE;
191:   return(0);
192: }

196: /*@
197:    ISDestroy - Destroys an index set.

199:    Collective on IS

201:    Input Parameters:
202: .  is - the index set

204:    Level: beginner

206: .seealso: ISCreateGeneral(), ISCreateStride(), ISCreateBlocked()
207: @*/
208: PetscErrorCode  ISDestroy(IS *is)
209: {

213:   if (!*is) return(0);
215:   if (--((PetscObject)(*is))->refct > 0) {*is = 0; return(0);}
216:   if ((*is)->complement) {
217:     PetscInt refcnt;
218:     PetscObjectGetReference((PetscObject)((*is)->complement), &refcnt);
219:     if (refcnt > 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Nonlocal IS has not been restored");
220:     ISDestroy(&(*is)->complement);
221:   }
222:   if ((*is)->ops->destroy) {
223:     (*(*is)->ops->destroy)(*is);
224:   }
225:   PetscLayoutDestroy(&(*is)->map);
226:   /* Destroy local representations of offproc data. */
227:   PetscFree((*is)->total);
228:   PetscFree((*is)->nonlocal);
229:   PetscHeaderDestroy(is);
230:   return(0);
231: }

235: /*@
236:    ISInvertPermutation - Creates a new permutation that is the inverse of
237:                          a given permutation.

239:    Collective on IS

241:    Input Parameter:
242: +  is - the index set
243: -  nlocal - number of indices on this processor in result (ignored for 1 proccessor) or
244:             use PETSC_DECIDE

246:    Output Parameter:
247: .  isout - the inverse permutation

249:    Level: intermediate

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

254:    Concepts: inverse permutation
255:    Concepts: permutation^inverse
256:    Concepts: index sets^inverting
257: @*/
258: PetscErrorCode  ISInvertPermutation(IS is,PetscInt nlocal,IS *isout)
259: {

265:   if (!is->isperm) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Not a permutation, must call ISSetPermutation() on the IS first");
266:   if (is->isidentity) {
267:     ISDuplicate(is,isout);
268:   } else {
269:     (*is->ops->invertpermutation)(is,nlocal,isout);
270:     ISSetPermutation(*isout);
271:   }
272:   return(0);
273: }

277: /*@
278:    ISGetSize - Returns the global length of an index set.

280:    Not Collective

282:    Input Parameter:
283: .  is - the index set

285:    Output Parameter:
286: .  size - the global size

288:    Level: beginner

290:    Concepts: size^of index set
291:    Concepts: index sets^size

293: @*/
294: PetscErrorCode  ISGetSize(IS is,PetscInt *size)
295: {

301:   (*is->ops->getsize)(is,size);
302:   return(0);
303: }

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

310:    Not Collective

312:    Input Parameter:
313: .  is - the index set

315:    Output Parameter:
316: .  size - the local size

318:    Level: beginner

320:    Concepts: size^of index set
321:    Concepts: local size^of index set
322:    Concepts: index sets^local size

324: @*/
325: PetscErrorCode  ISGetLocalSize(IS is,PetscInt *size)
326: {

332:   (*is->ops->getlocalsize)(is,size);
333:   return(0);
334: }

338: /*@C
339:    ISGetIndices - Returns a pointer to the indices.  The user should call
340:    ISRestoreIndices() after having looked at the indices.  The user should
341:    NOT change the indices.

343:    Not Collective

345:    Input Parameter:
346: .  is - the index set

348:    Output Parameter:
349: .  ptr - the location to put the pointer to the indices

351:    Fortran Note:
352:    This routine is used differently from Fortran
353: $    IS          is
354: $    integer     is_array(1)
355: $    PetscOffset i_is
356: $    int         ierr
357: $       call ISGetIndices(is,is_array,i_is,ierr)
358: $
359: $   Access first local entry in list
360: $      value = is_array(i_is + 1)
361: $
362: $      ...... other code
363: $       call ISRestoreIndices(is,is_array,i_is,ierr)

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

368:    Level: intermediate

370:    Concepts: index sets^getting indices
371:    Concepts: indices of index set

373: .seealso: ISRestoreIndices(), ISGetIndicesF90()
374: @*/
375: PetscErrorCode  ISGetIndices(IS is,const PetscInt *ptr[])
376: {

382:   (*is->ops->getindices)(is,ptr);
383:   return(0);
384: }

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

391:    Not Collective

393:    Input Parameter:
394: .  is - the index set

396:    Output Parameter:
397: +   min - the minimum value
398: -   max - the maximum value

400:    Level: intermediate

402:    Concepts: index sets^getting indices
403:    Concepts: indices of index set

405: .seealso: ISGetIndices(), ISRestoreIndices(), ISGetIndicesF90()
406: @*/
407: PetscErrorCode  ISGetMinMax(IS is,PetscInt *min,PetscInt *max)
408: {
411:   if (min) *min = is->min;
412:   if (max) *max = is->max;
413:   return(0);
414: }

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

421:   Not Collective

423:   Input Parameter:
424: + is - the index set
425: - key - the search key

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

430:   Level: intermediate
431: @*/
432: PetscErrorCode ISLocate(IS is, PetscInt key, PetscInt *location)
433: {

437:   if (is->ops->locate) {
438:     (*is->ops->locate)(is,key,location);
439:   } else {
440:     PetscInt       numIdx;
441:     PetscBool      sorted;
442:     const PetscInt *idx;

444:     ISGetLocalSize(is,&numIdx);
445:     ISGetIndices(is,&idx);
446:     ISSorted(is,&sorted);
447:     if (sorted) {
448:       PetscFindInt(key,numIdx,idx,location);
449:     } else {
450:       PetscInt i;

452:       *location = -1;
453:       for (i = 0; i < numIdx; i++) {
454:         if (idx[i] == key) {
455:           *location = i;
456:           break;
457:         }
458:       }
459:     }
460:     ISRestoreIndices(is,&idx);
461:   }
462:   return(0);
463: }

467: /*@C
468:    ISRestoreIndices - Restores an index set to a usable state after a call
469:                       to ISGetIndices().

471:    Not Collective

473:    Input Parameters:
474: +  is - the index set
475: -  ptr - the pointer obtained by ISGetIndices()

477:    Fortran Note:
478:    This routine is used differently from Fortran
479: $    IS          is
480: $    integer     is_array(1)
481: $    PetscOffset i_is
482: $    int         ierr
483: $       call ISGetIndices(is,is_array,i_is,ierr)
484: $
485: $   Access first local entry in list
486: $      value = is_array(i_is + 1)
487: $
488: $      ...... other code
489: $       call ISRestoreIndices(is,is_array,i_is,ierr)

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

494:    Level: intermediate

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

499: .seealso: ISGetIndices(), ISRestoreIndicesF90()
500: @*/
501: PetscErrorCode  ISRestoreIndices(IS is,const PetscInt *ptr[])
502: {

508:   if (is->ops->restoreindices) {
509:     (*is->ops->restoreindices)(is,ptr);
510:   }
511:   return(0);
512: }

516: static PetscErrorCode ISGatherTotal_Private(IS is)
517: {
519:   PetscInt       i,n,N;
520:   const PetscInt *lindices;
521:   MPI_Comm       comm;
522:   PetscMPIInt    rank,size,*sizes = NULL,*offsets = NULL,nn;


527:   PetscObjectGetComm((PetscObject)is,&comm);
528:   MPI_Comm_size(comm,&size);
529:   MPI_Comm_rank(comm,&rank);
530:   ISGetLocalSize(is,&n);
531:   PetscMalloc2(size,&sizes,size,&offsets);

533:   PetscMPIIntCast(n,&nn);
534:   MPI_Allgather(&nn,1,MPI_INT,sizes,1,MPI_INT,comm);
535:   offsets[0] = 0;
536:   for (i=1; i<size; ++i) offsets[i] = offsets[i-1] + sizes[i-1];
537:   N = offsets[size-1] + sizes[size-1];

539:   PetscMalloc1(N,&(is->total));
540:   ISGetIndices(is,&lindices);
541:   MPI_Allgatherv((void*)lindices,nn,MPIU_INT,is->total,sizes,offsets,MPIU_INT,comm);
542:   ISRestoreIndices(is,&lindices);
543:   is->local_offset = offsets[rank];
544:   PetscFree2(sizes,offsets);
545:   return(0);
546: }

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

553:    Collective on IS

555:    Input Parameter:
556: .  is - the index set

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

562:    Level: intermediate

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

571:    Concepts: index sets^getting nonlocal indices
572: .seealso: ISRestoreTotalIndices(), ISGetNonlocalIndices(), ISGetSize()
573: @*/
574: PetscErrorCode ISGetTotalIndices(IS is, const PetscInt *indices[])
575: {
577:   PetscMPIInt    size;

582:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
583:   if (size == 1) {
584:     (*is->ops->getindices)(is,indices);
585:   } else {
586:     if (!is->total) {
587:       ISGatherTotal_Private(is);
588:     }
589:     *indices = is->total;
590:   }
591:   return(0);
592: }

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

599:    Not Collective.

601:    Input Parameter:
602: +  is - the index set
603: -  indices - index array; must be the array obtained with ISGetTotalIndices()

605:    Level: intermediate

607:    Concepts: index sets^getting nonlocal indices
608:    Concepts: index sets^restoring nonlocal indices
609: .seealso: ISRestoreTotalIndices(), ISGetNonlocalIndices()
610: @*/
611: PetscErrorCode  ISRestoreTotalIndices(IS is, const PetscInt *indices[])
612: {
614:   PetscMPIInt    size;

619:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
620:   if (size == 1) {
621:     (*is->ops->restoreindices)(is,indices);
622:   } else {
623:     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.");
624:   }
625:   return(0);
626: }
629: /*@C
630:    ISGetNonlocalIndices - Retrieve an array of indices from remote processors
631:                        in this communicator.

633:    Collective on IS

635:    Input Parameter:
636: .  is - the index set

638:    Output Parameter:
639: .  indices - indices with rank 0 indices first, and so on,  omitting
640:              the current rank.  Total number of indices is the difference
641:              total and local, obtained with ISGetSize() and ISGetLocalSize(),
642:              respectively.

644:    Level: intermediate

646:    Notes: restore the indices using ISRestoreNonlocalIndices().
647:           The same scalability considerations as those for ISGetTotalIndices
648:           apply here.

650:    Concepts: index sets^getting nonlocal indices
651: .seealso: ISGetTotalIndices(), ISRestoreNonlocalIndices(), ISGetSize(), ISGetLocalSize().
652: @*/
653: PetscErrorCode  ISGetNonlocalIndices(IS is, const PetscInt *indices[])
654: {
656:   PetscMPIInt    size;
657:   PetscInt       n, N;

662:   MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
663:   if (size == 1) *indices = NULL;
664:   else {
665:     if (!is->total) {
666:       ISGatherTotal_Private(is);
667:     }
668:     ISGetLocalSize(is,&n);
669:     ISGetSize(is,&N);
670:     PetscMalloc1(N-n, &(is->nonlocal));
671:     PetscMemcpy(is->nonlocal, is->total, sizeof(PetscInt)*is->local_offset);
672:     PetscMemcpy(is->nonlocal+is->local_offset, is->total+is->local_offset+n, sizeof(PetscInt)*(N - is->local_offset - n));
673:     *indices = is->nonlocal;
674:   }
675:   return(0);
676: }

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

683:    Not Collective.

685:    Input Parameter:
686: +  is - the index set
687: -  indices - index array; must be the array obtained with ISGetNonlocalIndices()

689:    Level: intermediate

691:    Concepts: index sets^getting nonlocal indices
692:    Concepts: index sets^restoring nonlocal indices
693: .seealso: ISGetTotalIndices(), ISGetNonlocalIndices(), ISRestoreTotalIndices()
694: @*/
695: PetscErrorCode  ISRestoreNonlocalIndices(IS is, const PetscInt *indices[])
696: {
700:   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.");
701:   return(0);
702: }

706: /*@
707:    ISGetNonlocalIS - Gather all nonlocal indices for this IS and present
708:                      them as another sequential index set.


711:    Collective on IS

713:    Input Parameter:
714: .  is - the index set

716:    Output Parameter:
717: .  complement - sequential IS with indices identical to the result of
718:                 ISGetNonlocalIndices()

720:    Level: intermediate

722:    Notes: complement represents the result of ISGetNonlocalIndices as an IS.
723:           Therefore scalability issues similar to ISGetNonlocalIndices apply.
724:           The resulting IS must be restored using ISRestoreNonlocalIS().

726:    Concepts: index sets^getting nonlocal indices
727: .seealso: ISGetNonlocalIndices(), ISRestoreNonlocalIndices(),  ISAllGather(), ISGetSize()
728: @*/
729: PetscErrorCode  ISGetNonlocalIS(IS is, IS *complement)
730: {

736:   /* Check if the complement exists already. */
737:   if (is->complement) {
738:     *complement = is->complement;
739:     PetscObjectReference((PetscObject)(is->complement));
740:   } else {
741:     PetscInt       N, n;
742:     const PetscInt *idx;
743:     ISGetSize(is, &N);
744:     ISGetLocalSize(is,&n);
745:     ISGetNonlocalIndices(is, &idx);
746:     ISCreateGeneral(PETSC_COMM_SELF, N-n,idx, PETSC_USE_POINTER, &(is->complement));
747:     PetscObjectReference((PetscObject)is->complement);
748:     *complement = is->complement;
749:   }
750:   return(0);
751: }


756: /*@
757:    ISRestoreNonlocalIS - Restore the IS obtained with ISGetNonlocalIS().

759:    Not collective.

761:    Input Parameter:
762: +  is         - the index set
763: -  complement - index set of is's nonlocal indices

765:    Level: intermediate


768:    Concepts: index sets^getting nonlocal indices
769:    Concepts: index sets^restoring nonlocal indices
770: .seealso: ISGetNonlocalIS(), ISGetNonlocalIndices(), ISRestoreNonlocalIndices()
771: @*/
772: PetscErrorCode  ISRestoreNonlocalIS(IS is, IS *complement)
773: {
775:   PetscInt       refcnt;

780:   if (*complement != is->complement) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Complement IS being restored was not obtained with ISGetNonlocalIS()");
781:   PetscObjectGetReference((PetscObject)(is->complement), &refcnt);
782:   if (refcnt <= 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Duplicate call to ISRestoreNonlocalIS() detected");
783:   PetscObjectDereference((PetscObject)(is->complement));
784:   return(0);
785: }

789: /*@C
790:    ISView - Displays an index set.

792:    Collective on IS

794:    Input Parameters:
795: +  is - the index set
796: -  viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.

798:    Level: intermediate

800: .seealso: PetscViewerASCIIOpen()
801: @*/
802: PetscErrorCode  ISView(IS is,PetscViewer viewer)
803: {

808:   if (!viewer) {
809:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)is),&viewer);
810:   }

814:   PetscObjectPrintClassNamePrefixType((PetscObject)is,viewer);
815:   (*is->ops->view)(is,viewer);
816:   return(0);
817: }

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

824:   Collective on PetscViewer

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

830:   Level: intermediate

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

837:   Concepts: index set^loading from file

839: .seealso: PetscViewerBinaryOpen(), ISView(), MatLoad(), VecLoad()
840: @*/
841: PetscErrorCode ISLoad(IS is, PetscViewer viewer)
842: {
843:   PetscBool      isbinary, ishdf5;

849:   PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERBINARY, &isbinary);
850:   PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERHDF5, &ishdf5);
851:   if (!isbinary && !ishdf5) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid viewer; open viewer with PetscViewerBinaryOpen()");
852:   if (!((PetscObject)is)->type_name) {ISSetType(is, ISGENERAL);}
853:   (*is->ops->load)(is, viewer);
854:   return(0);
855: }

859: /*@
860:    ISSort - Sorts the indices of an index set.

862:    Collective on IS

864:    Input Parameters:
865: .  is - the index set

867:    Level: intermediate

869:    Concepts: index sets^sorting
870:    Concepts: sorting^index set

872: .seealso: ISSortRemoveDups(), ISSorted()
873: @*/
874: PetscErrorCode  ISSort(IS is)
875: {

880:   (*is->ops->sort)(is);
881:   return(0);
882: }

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

889:   Collective on IS

891:   Input Parameters:
892: . is - the index set

894:   Level: intermediate

896:   Concepts: index sets^sorting
897:   Concepts: sorting^index set

899: .seealso: ISSort(), ISSorted()
900: @*/
901: PetscErrorCode ISSortRemoveDups(IS is)
902: {

907:   (*is->ops->sortremovedups)(is);
908:   return(0);
909: }

913: /*@
914:    ISToGeneral - Converts an IS object of any type to ISGENERAL type

916:    Collective on IS

918:    Input Parameters:
919: .  is - the index set

921:    Level: intermediate

923:    Concepts: index sets^sorting
924:    Concepts: sorting^index set

926: .seealso: ISSorted()
927: @*/
928: PetscErrorCode  ISToGeneral(IS is)
929: {

934:   if (is->ops->togeneral) {
935:     (*is->ops->togeneral)(is);
936:   } else SETERRQ1(PetscObjectComm((PetscObject)is),PETSC_ERR_SUP,"Not written for this type %s",((PetscObject)is)->type_name);
937:   return(0);
938: }

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

945:    Collective on IS

947:    Input Parameter:
948: .  is - the index set

950:    Output Parameter:
951: .  flg - output flag, either PETSC_TRUE if the index set is sorted,
952:          or PETSC_FALSE otherwise.

954:    Notes: For parallel IS objects this only indicates if the local part of the IS
955:           is sorted. So some processors may return PETSC_TRUE while others may
956:           return PETSC_FALSE.

958:    Level: intermediate

960: .seealso: ISSort(), ISSortRemoveDups()
961: @*/
962: PetscErrorCode  ISSorted(IS is,PetscBool  *flg)
963: {

969:   (*is->ops->sorted)(is,flg);
970:   return(0);
971: }

975: /*@
976:    ISDuplicate - Creates a duplicate copy of an index set.

978:    Collective on IS

980:    Input Parmeters:
981: .  is - the index set

983:    Output Parameters:
984: .  isnew - the copy of the index set

986:    Level: beginner

988:    Concepts: index sets^duplicating

990: .seealso: ISCreateGeneral(), ISCopy()
991: @*/
992: PetscErrorCode  ISDuplicate(IS is,IS *newIS)
993: {

999:   (*is->ops->duplicate)(is,newIS);
1000:   (*newIS)->isidentity = is->isidentity;
1001:   (*newIS)->isperm     = is->isperm;
1002:   return(0);
1003: }

1007: /*@
1008:    ISCopy - Copies an index set.

1010:    Collective on IS

1012:    Input Parmeters:
1013: .  is - the index set

1015:    Output Parameters:
1016: .  isy - the copy of the index set

1018:    Level: beginner

1020:    Concepts: index sets^copying

1022: .seealso: ISDuplicate()
1023: @*/
1024: PetscErrorCode  ISCopy(IS is,IS isy)
1025: {

1032:   if (is == isy) return(0);
1033:   (*is->ops->copy)(is,isy);
1034:   isy->isperm     = is->isperm;
1035:   isy->max        = is->max;
1036:   isy->min        = is->min;
1037:   isy->isidentity = is->isidentity;
1038:   return(0);
1039: }

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

1046:    Collective on IS and comm

1048:    Input Arguments:
1049: + is - index set
1050: . comm - communicator for new index set
1051: - mode - copy semantics, PETSC_USE_POINTER for no-copy if possible, otherwise PETSC_COPY_VALUES

1053:    Output Arguments:
1054: . newis - new IS on comm

1056:    Level: advanced

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

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

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

1066: .seealso: ISSplit()
1067: @*/
1068: PetscErrorCode  ISOnComm(IS is,MPI_Comm comm,PetscCopyMode mode,IS *newis)
1069: {
1071:   PetscMPIInt    match;

1076:   MPI_Comm_compare(PetscObjectComm((PetscObject)is),comm,&match);
1077:   if (mode != PETSC_COPY_VALUES && (match == MPI_IDENT || match == MPI_CONGRUENT)) {
1078:     PetscObjectReference((PetscObject)is);
1079:     *newis = is;
1080:   } else {
1081:     (*is->ops->oncomm)(is,comm,mode,newis);
1082:   }
1083:   return(0);
1084: }

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

1091:    Logicall Collective on IS

1093:    Input Arguments:
1094: + is - index set
1095: - bs - block size

1097:    Level: intermediate

1099: .seealso: ISGetBlockSize(), ISCreateBlock()
1100: @*/
1101: PetscErrorCode  ISSetBlockSize(IS is,PetscInt bs)
1102: {

1108:   if (bs < 1) SETERRQ1(PetscObjectComm((PetscObject)is),PETSC_ERR_ARG_OUTOFRANGE,"Block size %D, must be positive",bs);
1109:   (*is->ops->setblocksize)(is,bs);
1110:   return(0);
1111: }

1115: /*@
1116:    ISGetBlockSize - Returns the number of elements in a block.

1118:    Not Collective

1120:    Input Parameter:
1121: .  is - the index set

1123:    Output Parameter:
1124: .  size - the number of elements in a block

1126:    Level: intermediate

1128:    Concepts: IS^block size
1129:    Concepts: index sets^block size

1131: .seealso: ISBlockGetSize(), ISGetSize(), ISCreateBlock(), ISSetBlockSize()
1132: @*/
1133: PetscErrorCode  ISGetBlockSize(IS is,PetscInt *size)
1134: {

1138:   PetscLayoutGetBlockSize(is->map, size);
1139:   return(0);
1140: }

1144: PetscErrorCode ISGetIndicesCopy(IS is, PetscInt idx[])
1145: {
1147:   PetscInt       len,i;
1148:   const PetscInt *ptr;

1151:   ISGetSize(is,&len);
1152:   ISGetIndices(is,&ptr);
1153:   for (i=0; i<len; i++) idx[i] = ptr[i];
1154:   ISRestoreIndices(is,&ptr);
1155:   return(0);
1156: }

1158: /*MC
1159:     ISGetIndicesF90 - Accesses the elements of an index set from Fortran90.
1160:     The users should call ISRestoreIndicesF90() after having looked at the
1161:     indices.  The user should NOT change the indices.

1163:     Synopsis:
1164:     ISGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

1166:     Not collective

1168:     Input Parameter:
1169: .   x - index set

1171:     Output Parameters:
1172: +   xx_v - the Fortran90 pointer to the array
1173: -   ierr - error code

1175:     Example of Usage:
1176: .vb
1177:     PetscScalar, pointer xx_v(:)
1178:     ....
1179:     call ISGetIndicesF90(x,xx_v,ierr)
1180:     a = xx_v(3)
1181:     call ISRestoreIndicesF90(x,xx_v,ierr)
1182: .ve

1184:     Notes:
1185:     Not yet supported for all F90 compilers.

1187:     Level: intermediate

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

1191:   Concepts: index sets^getting indices in f90
1192:   Concepts: indices of index set in f90

1194: M*/

1196: /*MC
1197:     ISRestoreIndicesF90 - Restores an index set to a usable state after
1198:     a call to ISGetIndicesF90().

1200:     Synopsis:
1201:     ISRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

1203:     Not collective

1205:     Input Parameters:
1206: .   x - index set
1207: .   xx_v - the Fortran90 pointer to the array

1209:     Output Parameter:
1210: .   ierr - error code


1213:     Example of Usage:
1214: .vb
1215:     PetscScalar, pointer xx_v(:)
1216:     ....
1217:     call ISGetIndicesF90(x,xx_v,ierr)
1218:     a = xx_v(3)
1219:     call ISRestoreIndicesF90(x,xx_v,ierr)
1220: .ve

1222:     Notes:
1223:     Not yet supported for all F90 compilers.

1225:     Level: intermediate

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

1229: M*/

1231: /*MC
1232:     ISBlockGetIndicesF90 - Accesses the elements of an index set from Fortran90.
1233:     The users should call ISBlockRestoreIndicesF90() after having looked at the
1234:     indices.  The user should NOT change the indices.

1236:     Synopsis:
1237:     ISBlockGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

1239:     Not collective

1241:     Input Parameter:
1242: .   x - index set

1244:     Output Parameters:
1245: +   xx_v - the Fortran90 pointer to the array
1246: -   ierr - error code
1247:     Example of Usage:
1248: .vb
1249:     PetscScalar, pointer xx_v(:)
1250:     ....
1251:     call ISBlockGetIndicesF90(x,xx_v,ierr)
1252:     a = xx_v(3)
1253:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
1254: .ve

1256:     Notes:
1257:     Not yet supported for all F90 compilers

1259:     Level: intermediate

1261: .seealso:  ISBlockRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices(),
1262:            ISRestoreIndices()

1264:   Concepts: index sets^getting block indices in f90
1265:   Concepts: indices of index set in f90
1266:   Concepts: block^ indices of index set in f90

1268: M*/

1270: /*MC
1271:     ISBlockRestoreIndicesF90 - Restores an index set to a usable state after
1272:     a call to ISBlockGetIndicesF90().

1274:     Synopsis:
1275:     ISBlockRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

1277:     Not Collective

1279:     Input Parameters:
1280: +   x - index set
1281: -   xx_v - the Fortran90 pointer to the array

1283:     Output Parameter:
1284: .   ierr - error code

1286:     Example of Usage:
1287: .vb
1288:     PetscScalar, pointer xx_v(:)
1289:     ....
1290:     call ISBlockGetIndicesF90(x,xx_v,ierr)
1291:     a = xx_v(3)
1292:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
1293: .ve

1295:     Notes:
1296:     Not yet supported for all F90 compilers

1298:     Level: intermediate

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

1302: M*/