Actual source code: index.c
1: #define PETSCVEC_DLL
3: /*
4: Defines the abstract operations on index sets, i.e. the public interface.
5: */
6: #include private/isimpl.h
8: /* Logging support */
9: PetscCookie IS_COOKIE;
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,PetscTruth *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: 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: {
67: is->isidentity = PETSC_TRUE;
68: return(0);
69: }
73: /*@
74: ISPermutation - PETSC_TRUE or PETSC_FALSE depending on whether the
75: index set has been declared to be a permutation.
77: Collective on IS
79: Input Parmeters:
80: . is - the index set
82: Output Parameters:
83: . perm - PETSC_TRUE if a permutation, else PETSC_FALSE
85: Level: intermediate
87: Concepts: permutation
88: Concepts: index sets^is permutation
90: .seealso: ISSetPermutation()
91: @*/
92: PetscErrorCode ISPermutation(IS is,PetscTruth *perm)
93: {
97: *perm = (PetscTruth) is->isperm;
98: return(0);
99: }
103: /*@
104: ISSetPermutation - Informs the index set that it is a permutation.
106: Collective on IS
108: Input Parmeters:
109: . is - the index set
111: Level: intermediate
113: Concepts: permutation
114: Concepts: index sets^permutation
116: The debug version of the libraries (config/configure.py --with-debugging=1) checks if the
117: index set is actually a permutation. The optimized version just believes you.
119: .seealso: ISPermutation()
120: @*/
121: PetscErrorCode ISSetPermutation(IS is)
122: {
125: #if defined(PETSC_USE_DEBUG)
126: {
127: PetscMPIInt size;
130: MPI_Comm_size(((PetscObject)is)->comm,&size);
131: if (size == 1) {
132: PetscInt i,n,*idx;
133: const PetscInt *iidx;
134:
135: ISGetSize(is,&n);
136: PetscMalloc(n*sizeof(PetscInt),&idx);
137: ISGetIndices(is,&iidx);
138: PetscMemcpy(idx,iidx,n*sizeof(PetscInt));
139: PetscSortInt(n,idx);
140: for (i=0; i<n; i++) {
141: if (idx[i] != i) SETERRQ(PETSC_ERR_ARG_WRONG,"Index set is not a permutation");
142: }
143: PetscFree(idx);
144: }
145: }
146: #endif
147: is->isperm = PETSC_TRUE;
148: return(0);
149: }
153: /*@
154: ISDestroy - Destroys an index set.
156: Collective on IS
158: Input Parameters:
159: . is - the index set
161: Level: beginner
163: .seealso: ISCreateGeneral(), ISCreateStride(), ISCreateBlocked()
164: @*/
165: PetscErrorCode ISDestroy(IS is)
166: {
171: if (--((PetscObject)is)->refct > 0) return(0);
172: if (is->ops->destroy) {
173: (*is->ops->destroy)(is);
174: }
175: PetscHeaderDestroy(is);
176: return(0);
177: }
181: /*@
182: ISInvertPermutation - Creates a new permutation that is the inverse of
183: a given permutation.
185: Collective on IS
187: Input Parameter:
188: + is - the index set
189: - nlocal - number of indices on this processor in result (ignored for 1 proccessor) or
190: use PETSC_DECIDE
192: Output Parameter:
193: . isout - the inverse permutation
195: Level: intermediate
197: Notes: For parallel index sets this does the complete parallel permutation, but the
198: code is not efficient for huge index sets (10,000,000 indices).
200: Concepts: inverse permutation
201: Concepts: permutation^inverse
202: Concepts: index sets^inverting
203: @*/
204: PetscErrorCode ISInvertPermutation(IS is,PetscInt nlocal,IS *isout)
205: {
211: if (!is->isperm) SETERRQ(PETSC_ERR_ARG_WRONG,"Not a permutation, must call ISSetPermutation() on the IS first");
212: (*is->ops->invertpermutation)(is,nlocal,isout);
213: ISSetPermutation(*isout);
214: return(0);
215: }
219: /*@
220: ISGetSize - Returns the global length of an index set.
222: Not Collective
224: Input Parameter:
225: . is - the index set
227: Output Parameter:
228: . size - the global size
230: Level: beginner
232: Concepts: size^of index set
233: Concepts: index sets^size
235: @*/
236: PetscErrorCode ISGetSize(IS is,PetscInt *size)
237: {
243: (*is->ops->getsize)(is,size);
244: return(0);
245: }
249: /*@
250: ISGetLocalSize - Returns the local (processor) length of an index set.
252: Not Collective
254: Input Parameter:
255: . is - the index set
257: Output Parameter:
258: . size - the local size
260: Level: beginner
262: Concepts: size^of index set
263: Concepts: local size^of index set
264: Concepts: index sets^local size
265:
266: @*/
267: PetscErrorCode ISGetLocalSize(IS is,PetscInt *size)
268: {
274: (*is->ops->getlocalsize)(is,size);
275: return(0);
276: }
280: /*@C
281: ISGetIndices - Returns a pointer to the indices. The user should call
282: ISRestoreIndices() after having looked at the indices. The user should
283: NOT change the indices.
285: Not Collective
287: Input Parameter:
288: . is - the index set
290: Output Parameter:
291: . ptr - the location to put the pointer to the indices
293: Fortran Note:
294: This routine is used differently from Fortran
295: $ IS is
296: $ integer is_array(1)
297: $ PetscOffset i_is
298: $ int ierr
299: $ call ISGetIndices(is,is_array,i_is,ierr)
300: $
301: $ Access first local entry in list
302: $ value = is_array(i_is + 1)
303: $
304: $ ...... other code
305: $ call ISRestoreIndices(is,is_array,i_is,ierr)
307: See the Fortran chapter of the users manual and
308: petsc/src/is/examples/[tutorials,tests] for details.
310: Level: intermediate
312: Concepts: index sets^getting indices
313: Concepts: indices of index set
315: .seealso: ISRestoreIndices(), ISGetIndicesF90()
316: @*/
317: PetscErrorCode ISGetIndices(IS is,const PetscInt *ptr[])
318: {
324: (*is->ops->getindices)(is,ptr);
325: return(0);
326: }
330: /*@C
331: ISRestoreIndices - Restores an index set to a usable state after a call
332: to ISGetIndices().
334: Not Collective
336: Input Parameters:
337: + is - the index set
338: - ptr - the pointer obtained by ISGetIndices()
340: Fortran Note:
341: This routine is used differently from Fortran
342: $ IS is
343: $ integer is_array(1)
344: $ PetscOffset i_is
345: $ int ierr
346: $ call ISGetIndices(is,is_array,i_is,ierr)
347: $
348: $ Access first local entry in list
349: $ value = is_array(i_is + 1)
350: $
351: $ ...... other code
352: $ call ISRestoreIndices(is,is_array,i_is,ierr)
354: See the Fortran chapter of the users manual and
355: petsc/src/is/examples/[tutorials,tests] for details.
357: Level: intermediate
359: .seealso: ISGetIndices(), ISRestoreIndicesF90()
360: @*/
361: PetscErrorCode ISRestoreIndices(IS is,const PetscInt *ptr[])
362: {
368: if (is->ops->restoreindices) {
369: (*is->ops->restoreindices)(is,ptr);
370: }
371: return(0);
372: }
376: /*@C
377: ISView - Displays an index set.
379: Collective on IS
381: Input Parameters:
382: + is - the index set
383: - viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.
385: Level: intermediate
387: .seealso: PetscViewerASCIIOpen()
388: @*/
389: PetscErrorCode ISView(IS is,PetscViewer viewer)
390: {
395: if (!viewer) {
396: PetscViewerASCIIGetStdout(((PetscObject)is)->comm,&viewer);
397: }
400:
401: (*is->ops->view)(is,viewer);
402: return(0);
403: }
407: /*@
408: ISSort - Sorts the indices of an index set.
410: Collective on IS
412: Input Parameters:
413: . is - the index set
415: Level: intermediate
417: Concepts: index sets^sorting
418: Concepts: sorting^index set
420: .seealso: ISSorted()
421: @*/
422: PetscErrorCode ISSort(IS is)
423: {
428: (*is->ops->sortindices)(is);
429: return(0);
430: }
434: /*@
435: ISSorted - Checks the indices to determine whether they have been sorted.
437: Collective on IS
439: Input Parameter:
440: . is - the index set
442: Output Parameter:
443: . flg - output flag, either PETSC_TRUE if the index set is sorted,
444: or PETSC_FALSE otherwise.
446: Notes: For parallel IS objects this only indicates if the local part of the IS
447: is sorted. So some processors may return PETSC_TRUE while others may
448: return PETSC_FALSE.
450: Level: intermediate
452: .seealso: ISSort()
453: @*/
454: PetscErrorCode ISSorted(IS is,PetscTruth *flg)
455: {
461: (*is->ops->sorted)(is,flg);
462: return(0);
463: }
467: /*@
468: ISDuplicate - Creates a duplicate copy of an index set.
470: Collective on IS
472: Input Parmeters:
473: . is - the index set
475: Output Parameters:
476: . isnew - the copy of the index set
478: Notes:
479: ISDuplicate() does not copy the index set, but rather allocates storage
480: for the new one. Use ISCopy() to copy an index set.
482: Level: beginner
484: Concepts: index sets^duplicating
486: .seealso: ISCreateGeneral(), ISCopy()
487: @*/
488: PetscErrorCode ISDuplicate(IS is,IS *newIS)
489: {
495: (*is->ops->duplicate)(is,newIS);
496: return(0);
497: }
501: /*@
502: ISCopy - Copies an index set.
504: Collective on IS
506: Input Parmeters:
507: . is - the index set
509: Output Parameters:
510: . isy - the copy of the index set
512: Level: beginner
514: Concepts: index sets^copying
516: .seealso: ISDuplicate()
517: @*/
518: PetscErrorCode ISCopy(IS is,IS isy)
519: {
526: if (is == isy) return(0);
527: (*is->ops->copy)(is,isy);
528: isy->isperm = is->isperm;
529: isy->max = is->max;
530: isy->min = is->min;
531: isy->isidentity = is->isidentity;
532: return(0);
533: }
535: /*MC
536: ISGetIndicesF90 - Accesses the elements of an index set from Fortran90.
537: The users should call ISRestoreIndicesF90() after having looked at the
538: indices. The user should NOT change the indices.
540: Synopsis:
541: ISGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
543: Not collective
545: Input Parameter:
546: . x - index set
548: Output Parameters:
549: + xx_v - the Fortran90 pointer to the array
550: - ierr - error code
552: Example of Usage:
553: .vb
554: PetscScalar, pointer xx_v(:)
555: ....
556: call ISGetIndicesF90(x,xx_v,ierr)
557: a = xx_v(3)
558: call ISRestoreIndicesF90(x,xx_v,ierr)
559: .ve
561: Notes:
562: Not yet supported for all F90 compilers.
564: Level: intermediate
566: .seealso: ISRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices()
568: Concepts: index sets^getting indices in f90
569: Concepts: indices of index set in f90
571: M*/
573: /*MC
574: ISRestoreIndicesF90 - Restores an index set to a usable state after
575: a call to ISGetIndicesF90().
577: Synopsis:
578: ISRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
580: Not collective
582: Input Parameters:
583: . x - index set
584: . xx_v - the Fortran90 pointer to the array
586: Output Parameter:
587: . ierr - error code
590: Example of Usage:
591: .vb
592: PetscScalar, pointer xx_v(:)
593: ....
594: call ISGetIndicesF90(x,xx_v,ierr)
595: a = xx_v(3)
596: call ISRestoreIndicesF90(x,xx_v,ierr)
597: .ve
598:
599: Notes:
600: Not yet supported for all F90 compilers.
602: Level: intermediate
604: .seealso: ISGetIndicesF90(), ISGetIndices(), ISRestoreIndices()
606: M*/
608: /*MC
609: ISBlockGetIndicesF90 - Accesses the elements of an index set from Fortran90.
610: The users should call ISBlockRestoreIndicesF90() after having looked at the
611: indices. The user should NOT change the indices.
613: Synopsis:
614: ISBlockGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
616: Not collective
618: Input Parameter:
619: . x - index set
621: Output Parameters:
622: + xx_v - the Fortran90 pointer to the array
623: - ierr - error code
624: Example of Usage:
625: .vb
626: PetscScalar, pointer xx_v(:)
627: ....
628: call ISBlockGetIndicesF90(x,xx_v,ierr)
629: a = xx_v(3)
630: call ISBlockRestoreIndicesF90(x,xx_v,ierr)
631: .ve
633: Notes:
634: Not yet supported for all F90 compilers
636: Level: intermediate
638: .seealso: ISBlockRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices(),
639: ISRestoreIndices()
641: Concepts: index sets^getting block indices in f90
642: Concepts: indices of index set in f90
643: Concepts: block^ indices of index set in f90
645: M*/
647: /*MC
648: ISBlockRestoreIndicesF90 - Restores an index set to a usable state after
649: a call to ISBlockGetIndicesF90().
651: Synopsis:
652: ISBlockRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
654: Input Parameters:
655: + x - index set
656: - xx_v - the Fortran90 pointer to the array
658: Output Parameter:
659: . ierr - error code
661: Example of Usage:
662: .vb
663: PetscScalar, pointer xx_v(:)
664: ....
665: call ISBlockGetIndicesF90(x,xx_v,ierr)
666: a = xx_v(3)
667: call ISBlockRestoreIndicesF90(x,xx_v,ierr)
668: .ve
669:
670: Notes:
671: Not yet supported for all F90 compilers
673: Level: intermediate
675: .seealso: ISBlockGetIndicesF90(), ISGetIndices(), ISRestoreIndices(), ISRestoreIndicesF90()
677: M*/