Actual source code: index.c
petsc-3.5.4 2015-05-23
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: /*@C
419: ISRestoreIndices - Restores an index set to a usable state after a call
420: to ISGetIndices().
422: Not Collective
424: Input Parameters:
425: + is - the index set
426: - ptr - the pointer obtained by ISGetIndices()
428: Fortran Note:
429: This routine is used differently from Fortran
430: $ IS is
431: $ integer is_array(1)
432: $ PetscOffset i_is
433: $ int ierr
434: $ call ISGetIndices(is,is_array,i_is,ierr)
435: $
436: $ Access first local entry in list
437: $ value = is_array(i_is + 1)
438: $
439: $ ...... other code
440: $ call ISRestoreIndices(is,is_array,i_is,ierr)
442: See the Fortran chapter of the users manual and
443: petsc/src/is/examples/[tutorials,tests] for details.
445: Level: intermediate
447: Note:
448: This routine zeros out ptr. This is to prevent accidental us of the array after it has been restored.
450: .seealso: ISGetIndices(), ISRestoreIndicesF90()
451: @*/
452: PetscErrorCode ISRestoreIndices(IS is,const PetscInt *ptr[])
453: {
459: if (is->ops->restoreindices) {
460: (*is->ops->restoreindices)(is,ptr);
461: }
462: return(0);
463: }
467: static PetscErrorCode ISGatherTotal_Private(IS is)
468: {
470: PetscInt i,n,N;
471: const PetscInt *lindices;
472: MPI_Comm comm;
473: PetscMPIInt rank,size,*sizes = NULL,*offsets = NULL,nn;
478: PetscObjectGetComm((PetscObject)is,&comm);
479: MPI_Comm_size(comm,&size);
480: MPI_Comm_rank(comm,&rank);
481: ISGetLocalSize(is,&n);
482: PetscMalloc2(size,&sizes,size,&offsets);
484: PetscMPIIntCast(n,&nn);
485: MPI_Allgather(&nn,1,MPI_INT,sizes,1,MPI_INT,comm);
486: offsets[0] = 0;
487: for (i=1; i<size; ++i) offsets[i] = offsets[i-1] + sizes[i-1];
488: N = offsets[size-1] + sizes[size-1];
490: PetscMalloc1(N,&(is->total));
491: ISGetIndices(is,&lindices);
492: MPI_Allgatherv((void*)lindices,nn,MPIU_INT,is->total,sizes,offsets,MPIU_INT,comm);
493: ISRestoreIndices(is,&lindices);
494: is->local_offset = offsets[rank];
495: PetscFree2(sizes,offsets);
496: return(0);
497: }
501: /*@C
502: ISGetTotalIndices - Retrieve an array containing all indices across the communicator.
504: Collective on IS
506: Input Parameter:
507: . is - the index set
509: Output Parameter:
510: . indices - total indices with rank 0 indices first, and so on; total array size is
511: the same as returned with ISGetSize().
513: Level: intermediate
515: Notes: this is potentially nonscalable, but depends on the size of the total index set
516: and the size of the communicator. This may be feasible for index sets defined on
517: subcommunicators, such that the set size does not grow with PETSC_WORLD_COMM.
518: Note also that there is no way to tell where the local part of the indices starts
519: (use ISGetIndices() and ISGetNonlocalIndices() to retrieve just the local and just
520: the nonlocal part (complement), respectively).
522: Concepts: index sets^getting nonlocal indices
523: .seealso: ISRestoreTotalIndices(), ISGetNonlocalIndices(), ISGetSize()
524: @*/
525: PetscErrorCode ISGetTotalIndices(IS is, const PetscInt *indices[])
526: {
528: PetscMPIInt size;
533: MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
534: if (size == 1) {
535: (*is->ops->getindices)(is,indices);
536: } else {
537: if (!is->total) {
538: ISGatherTotal_Private(is);
539: }
540: *indices = is->total;
541: }
542: return(0);
543: }
547: /*@C
548: ISRestoreTotalIndices - Restore the index array obtained with ISGetTotalIndices().
550: Not Collective.
552: Input Parameter:
553: + is - the index set
554: - indices - index array; must be the array obtained with ISGetTotalIndices()
556: Level: intermediate
558: Concepts: index sets^getting nonlocal indices
559: Concepts: index sets^restoring nonlocal indices
560: .seealso: ISRestoreTotalIndices(), ISGetNonlocalIndices()
561: @*/
562: PetscErrorCode ISRestoreTotalIndices(IS is, const PetscInt *indices[])
563: {
565: PetscMPIInt size;
570: MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
571: if (size == 1) {
572: (*is->ops->restoreindices)(is,indices);
573: } else {
574: 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.");
575: }
576: return(0);
577: }
580: /*@C
581: ISGetNonlocalIndices - Retrieve an array of indices from remote processors
582: in this communicator.
584: Collective on IS
586: Input Parameter:
587: . is - the index set
589: Output Parameter:
590: . indices - indices with rank 0 indices first, and so on, omitting
591: the current rank. Total number of indices is the difference
592: total and local, obtained with ISGetSize() and ISGetLocalSize(),
593: respectively.
595: Level: intermediate
597: Notes: restore the indices using ISRestoreNonlocalIndices().
598: The same scalability considerations as those for ISGetTotalIndices
599: apply here.
601: Concepts: index sets^getting nonlocal indices
602: .seealso: ISGetTotalIndices(), ISRestoreNonlocalIndices(), ISGetSize(), ISGetLocalSize().
603: @*/
604: PetscErrorCode ISGetNonlocalIndices(IS is, const PetscInt *indices[])
605: {
607: PetscMPIInt size;
608: PetscInt n, N;
613: MPI_Comm_size(PetscObjectComm((PetscObject)is), &size);
614: if (size == 1) *indices = NULL;
615: else {
616: if (!is->total) {
617: ISGatherTotal_Private(is);
618: }
619: ISGetLocalSize(is,&n);
620: ISGetSize(is,&N);
621: PetscMalloc(sizeof(PetscInt)*(N-n), &(is->nonlocal));
622: PetscMemcpy(is->nonlocal, is->total, sizeof(PetscInt)*is->local_offset);
623: PetscMemcpy(is->nonlocal+is->local_offset, is->total+is->local_offset+n, sizeof(PetscInt)*(N - is->local_offset - n));
624: *indices = is->nonlocal;
625: }
626: return(0);
627: }
631: /*@C
632: ISRestoreTotalIndices - Restore the index array obtained with ISGetNonlocalIndices().
634: Not Collective.
636: Input Parameter:
637: + is - the index set
638: - indices - index array; must be the array obtained with ISGetNonlocalIndices()
640: Level: intermediate
642: Concepts: index sets^getting nonlocal indices
643: Concepts: index sets^restoring nonlocal indices
644: .seealso: ISGetTotalIndices(), ISGetNonlocalIndices(), ISRestoreTotalIndices()
645: @*/
646: PetscErrorCode ISRestoreNonlocalIndices(IS is, const PetscInt *indices[])
647: {
651: 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.");
652: return(0);
653: }
657: /*@
658: ISGetNonlocalIS - Gather all nonlocal indices for this IS and present
659: them as another sequential index set.
662: Collective on IS
664: Input Parameter:
665: . is - the index set
667: Output Parameter:
668: . complement - sequential IS with indices identical to the result of
669: ISGetNonlocalIndices()
671: Level: intermediate
673: Notes: complement represents the result of ISGetNonlocalIndices as an IS.
674: Therefore scalability issues similar to ISGetNonlocalIndices apply.
675: The resulting IS must be restored using ISRestoreNonlocalIS().
677: Concepts: index sets^getting nonlocal indices
678: .seealso: ISGetNonlocalIndices(), ISRestoreNonlocalIndices(), ISAllGather(), ISGetSize()
679: @*/
680: PetscErrorCode ISGetNonlocalIS(IS is, IS *complement)
681: {
687: /* Check if the complement exists already. */
688: if (is->complement) {
689: *complement = is->complement;
690: PetscObjectReference((PetscObject)(is->complement));
691: } else {
692: PetscInt N, n;
693: const PetscInt *idx;
694: ISGetSize(is, &N);
695: ISGetLocalSize(is,&n);
696: ISGetNonlocalIndices(is, &idx);
697: ISCreateGeneral(PETSC_COMM_SELF, N-n,idx, PETSC_USE_POINTER, &(is->complement));
698: PetscObjectReference((PetscObject)is->complement);
699: *complement = is->complement;
700: }
701: return(0);
702: }
707: /*@
708: ISRestoreNonlocalIS - Restore the IS obtained with ISGetNonlocalIS().
710: Not collective.
712: Input Parameter:
713: + is - the index set
714: - complement - index set of is's nonlocal indices
716: Level: intermediate
719: Concepts: index sets^getting nonlocal indices
720: Concepts: index sets^restoring nonlocal indices
721: .seealso: ISGetNonlocalIS(), ISGetNonlocalIndices(), ISRestoreNonlocalIndices()
722: @*/
723: PetscErrorCode ISRestoreNonlocalIS(IS is, IS *complement)
724: {
726: PetscInt refcnt;
731: if (*complement != is->complement) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Complement IS being restored was not obtained with ISGetNonlocalIS()");
732: PetscObjectGetReference((PetscObject)(is->complement), &refcnt);
733: if (refcnt <= 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Duplicate call to ISRestoreNonlocalIS() detected");
734: PetscObjectDereference((PetscObject)(is->complement));
735: return(0);
736: }
740: /*@C
741: ISView - Displays an index set.
743: Collective on IS
745: Input Parameters:
746: + is - the index set
747: - viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.
749: Level: intermediate
751: .seealso: PetscViewerASCIIOpen()
752: @*/
753: PetscErrorCode ISView(IS is,PetscViewer viewer)
754: {
759: if (!viewer) {
760: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)is),&viewer);
761: }
765: PetscObjectPrintClassNamePrefixType((PetscObject)is,viewer);
766: (*is->ops->view)(is,viewer);
767: return(0);
768: }
772: /*@
773: ISLoad - Loads a vector that has been stored in binary or HDF5 format with ISView().
775: Collective on PetscViewer
777: Input Parameters:
778: + is - the newly loaded vector, this needs to have been created with ISCreate() or some related function before a call to ISLoad().
779: - viewer - binary file viewer, obtained from PetscViewerBinaryOpen() or HDF5 file viewer, obtained from PetscViewerHDF5Open()
781: Level: intermediate
783: Notes:
784: IF using HDF5, you must assign the IS the same name as was used in the IS
785: that was stored in the file using PetscObjectSetName(). Otherwise you will
786: get the error message: "Cannot H5DOpen2() with Vec name NAMEOFOBJECT"
788: Concepts: index set^loading from file
790: .seealso: PetscViewerBinaryOpen(), ISView(), MatLoad(), VecLoad()
791: @*/
792: PetscErrorCode ISLoad(IS is, PetscViewer viewer)
793: {
794: PetscBool isbinary, ishdf5;
800: PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERBINARY, &isbinary);
801: PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERHDF5, &ishdf5);
802: if (!isbinary && !ishdf5) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid viewer; open viewer with PetscViewerBinaryOpen()");
803: if (!((PetscObject)is)->type_name) {ISSetType(is, ISGENERAL);}
804: (*is->ops->load)(is, viewer);
805: return(0);
806: }
810: /*@
811: ISSort - Sorts the indices of an index set.
813: Collective on IS
815: Input Parameters:
816: . is - the index set
818: Level: intermediate
820: Concepts: index sets^sorting
821: Concepts: sorting^index set
823: .seealso: ISSortRemoveDups(), ISSorted()
824: @*/
825: PetscErrorCode ISSort(IS is)
826: {
831: (*is->ops->sort)(is);
832: return(0);
833: }
837: /*@
838: ISSortRemoveDups - Sorts the indices of an index set, removing duplicates.
840: Collective on IS
842: Input Parameters:
843: . is - the index set
845: Level: intermediate
847: Concepts: index sets^sorting
848: Concepts: sorting^index set
850: .seealso: ISSort(), ISSorted()
851: @*/
852: PetscErrorCode ISSortRemoveDups(IS is)
853: {
858: (*is->ops->sortremovedups)(is);
859: return(0);
860: }
864: /*@
865: ISToGeneral - Converts an IS object of any type to ISGENERAL type
867: Collective on IS
869: Input Parameters:
870: . is - the index set
872: Level: intermediate
874: Concepts: index sets^sorting
875: Concepts: sorting^index set
877: .seealso: ISSorted()
878: @*/
879: PetscErrorCode ISToGeneral(IS is)
880: {
885: if (is->ops->togeneral) {
886: (*is->ops->togeneral)(is);
887: } else SETERRQ1(PetscObjectComm((PetscObject)is),PETSC_ERR_SUP,"Not written for this type %s",((PetscObject)is)->type_name);
888: return(0);
889: }
893: /*@
894: ISSorted - Checks the indices to determine whether they have been sorted.
896: Collective on IS
898: Input Parameter:
899: . is - the index set
901: Output Parameter:
902: . flg - output flag, either PETSC_TRUE if the index set is sorted,
903: or PETSC_FALSE otherwise.
905: Notes: For parallel IS objects this only indicates if the local part of the IS
906: is sorted. So some processors may return PETSC_TRUE while others may
907: return PETSC_FALSE.
909: Level: intermediate
911: .seealso: ISSort(), ISSortRemoveDups()
912: @*/
913: PetscErrorCode ISSorted(IS is,PetscBool *flg)
914: {
920: (*is->ops->sorted)(is,flg);
921: return(0);
922: }
926: /*@
927: ISDuplicate - Creates a duplicate copy of an index set.
929: Collective on IS
931: Input Parmeters:
932: . is - the index set
934: Output Parameters:
935: . isnew - the copy of the index set
937: Level: beginner
939: Concepts: index sets^duplicating
941: .seealso: ISCreateGeneral(), ISCopy()
942: @*/
943: PetscErrorCode ISDuplicate(IS is,IS *newIS)
944: {
950: (*is->ops->duplicate)(is,newIS);
951: (*newIS)->isidentity = is->isidentity;
952: (*newIS)->isperm = is->isperm;
953: return(0);
954: }
958: /*@
959: ISCopy - Copies an index set.
961: Collective on IS
963: Input Parmeters:
964: . is - the index set
966: Output Parameters:
967: . isy - the copy of the index set
969: Level: beginner
971: Concepts: index sets^copying
973: .seealso: ISDuplicate()
974: @*/
975: PetscErrorCode ISCopy(IS is,IS isy)
976: {
983: if (is == isy) return(0);
984: (*is->ops->copy)(is,isy);
985: isy->isperm = is->isperm;
986: isy->max = is->max;
987: isy->min = is->min;
988: isy->isidentity = is->isidentity;
989: return(0);
990: }
994: /*@
995: ISOnComm - Split a parallel IS on subcomms (usually self) or concatenate index sets on subcomms into a parallel index set
997: Collective on IS and comm
999: Input Arguments:
1000: + is - index set
1001: . comm - communicator for new index set
1002: - mode - copy semantics, PETSC_USE_POINTER for no-copy if possible, otherwise PETSC_COPY_VALUES
1004: Output Arguments:
1005: . newis - new IS on comm
1007: Level: advanced
1009: Notes:
1010: It is usually desirable to create a parallel IS and look at the local part when necessary.
1012: This function is useful if serial ISs must be created independently, or to view many
1013: logically independent serial ISs.
1015: The input IS must have the same type on every process.
1017: .seealso: ISSplit()
1018: @*/
1019: PetscErrorCode ISOnComm(IS is,MPI_Comm comm,PetscCopyMode mode,IS *newis)
1020: {
1022: PetscMPIInt match;
1027: MPI_Comm_compare(PetscObjectComm((PetscObject)is),comm,&match);
1028: if (mode != PETSC_COPY_VALUES && (match == MPI_IDENT || match == MPI_CONGRUENT)) {
1029: PetscObjectReference((PetscObject)is);
1030: *newis = is;
1031: } else {
1032: (*is->ops->oncomm)(is,comm,mode,newis);
1033: }
1034: return(0);
1035: }
1039: /*@
1040: ISSetBlockSize - informs an index set that it has a given block size
1042: Logicall Collective on IS
1044: Input Arguments:
1045: + is - index set
1046: - bs - block size
1048: Level: intermediate
1050: .seealso: ISGetBlockSize(), ISCreateBlock()
1051: @*/
1052: PetscErrorCode ISSetBlockSize(IS is,PetscInt bs)
1053: {
1059: if (bs < 1) SETERRQ1(PetscObjectComm((PetscObject)is),PETSC_ERR_ARG_OUTOFRANGE,"Block size %D, must be positive",bs);
1060: (*is->ops->setblocksize)(is,bs);
1061: return(0);
1062: }
1066: /*@
1067: ISGetBlockSize - Returns the number of elements in a block.
1069: Not Collective
1071: Input Parameter:
1072: . is - the index set
1074: Output Parameter:
1075: . size - the number of elements in a block
1077: Level: intermediate
1079: Concepts: IS^block size
1080: Concepts: index sets^block size
1082: .seealso: ISBlockGetSize(), ISGetSize(), ISCreateBlock(), ISSetBlockSize()
1083: @*/
1084: PetscErrorCode ISGetBlockSize(IS is,PetscInt *size)
1085: {
1089: PetscLayoutGetBlockSize(is->map, size);
1090: return(0);
1091: }
1095: PetscErrorCode ISGetIndicesCopy(IS is, PetscInt idx[])
1096: {
1098: PetscInt len,i;
1099: const PetscInt *ptr;
1102: ISGetSize(is,&len);
1103: ISGetIndices(is,&ptr);
1104: for (i=0; i<len; i++) idx[i] = ptr[i];
1105: ISRestoreIndices(is,&ptr);
1106: return(0);
1107: }
1109: /*MC
1110: ISGetIndicesF90 - Accesses the elements of an index set from Fortran90.
1111: The users should call ISRestoreIndicesF90() after having looked at the
1112: indices. The user should NOT change the indices.
1114: Synopsis:
1115: ISGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
1117: Not collective
1119: Input Parameter:
1120: . x - index set
1122: Output Parameters:
1123: + xx_v - the Fortran90 pointer to the array
1124: - ierr - error code
1126: Example of Usage:
1127: .vb
1128: PetscScalar, pointer xx_v(:)
1129: ....
1130: call ISGetIndicesF90(x,xx_v,ierr)
1131: a = xx_v(3)
1132: call ISRestoreIndicesF90(x,xx_v,ierr)
1133: .ve
1135: Notes:
1136: Not yet supported for all F90 compilers.
1138: Level: intermediate
1140: .seealso: ISRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices()
1142: Concepts: index sets^getting indices in f90
1143: Concepts: indices of index set in f90
1145: M*/
1147: /*MC
1148: ISRestoreIndicesF90 - Restores an index set to a usable state after
1149: a call to ISGetIndicesF90().
1151: Synopsis:
1152: ISRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
1154: Not collective
1156: Input Parameters:
1157: . x - index set
1158: . xx_v - the Fortran90 pointer to the array
1160: Output Parameter:
1161: . ierr - error code
1164: Example of Usage:
1165: .vb
1166: PetscScalar, pointer xx_v(:)
1167: ....
1168: call ISGetIndicesF90(x,xx_v,ierr)
1169: a = xx_v(3)
1170: call ISRestoreIndicesF90(x,xx_v,ierr)
1171: .ve
1173: Notes:
1174: Not yet supported for all F90 compilers.
1176: Level: intermediate
1178: .seealso: ISGetIndicesF90(), ISGetIndices(), ISRestoreIndices()
1180: M*/
1182: /*MC
1183: ISBlockGetIndicesF90 - Accesses the elements of an index set from Fortran90.
1184: The users should call ISBlockRestoreIndicesF90() after having looked at the
1185: indices. The user should NOT change the indices.
1187: Synopsis:
1188: ISBlockGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
1190: Not collective
1192: Input Parameter:
1193: . x - index set
1195: Output Parameters:
1196: + xx_v - the Fortran90 pointer to the array
1197: - ierr - error code
1198: Example of Usage:
1199: .vb
1200: PetscScalar, pointer xx_v(:)
1201: ....
1202: call ISBlockGetIndicesF90(x,xx_v,ierr)
1203: a = xx_v(3)
1204: call ISBlockRestoreIndicesF90(x,xx_v,ierr)
1205: .ve
1207: Notes:
1208: Not yet supported for all F90 compilers
1210: Level: intermediate
1212: .seealso: ISBlockRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices(),
1213: ISRestoreIndices()
1215: Concepts: index sets^getting block indices in f90
1216: Concepts: indices of index set in f90
1217: Concepts: block^ indices of index set in f90
1219: M*/
1221: /*MC
1222: ISBlockRestoreIndicesF90 - Restores an index set to a usable state after
1223: a call to ISBlockGetIndicesF90().
1225: Synopsis:
1226: ISBlockRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
1228: Not Collective
1230: Input Parameters:
1231: + x - index set
1232: - xx_v - the Fortran90 pointer to the array
1234: Output Parameter:
1235: . ierr - error code
1237: Example of Usage:
1238: .vb
1239: PetscScalar, pointer xx_v(:)
1240: ....
1241: call ISBlockGetIndicesF90(x,xx_v,ierr)
1242: a = xx_v(3)
1243: call ISBlockRestoreIndicesF90(x,xx_v,ierr)
1244: .ve
1246: Notes:
1247: Not yet supported for all F90 compilers
1249: Level: intermediate
1251: .seealso: ISBlockGetIndicesF90(), ISGetIndices(), ISRestoreIndices(), ISRestoreIndicesF90()
1253: M*/