Actual source code: index.c
1: /*
2: Defines the abstract operations on index sets, i.e. the public interface.
3: */
4: #include src/vec/is/isimpl.h
6: /* Logging support */
7: PetscCookie IS_COOKIE = 0;
11: /*@C
12: ISIdentity - Determines whether index set is the identity mapping.
14: Collective on IS
16: Input Parmeters:
17: . is - the index set
19: Output Parameters:
20: . ident - PETSC_TRUE if an identity, else PETSC_FALSE
22: Level: intermediate
24: Concepts: identity mapping
25: Concepts: index sets^is identity
27: .seealso: ISSetIdentity()
28: @*/
29: PetscErrorCode ISIdentity(IS is,PetscTruth *ident)
30: {
36: *ident = is->isidentity;
37: if (*ident) return(0);
38: if (is->ops->identity) {
39: (*is->ops->identity)(is,ident);
40: }
41: return(0);
42: }
46: /*@
47: ISSetIdentity - Informs the index set that it is an identity.
49: Collective on IS
51: Input Parmeters:
52: . is - the index set
54: Level: intermediate
56: Concepts: identity mapping
57: Concepts: index sets^is identity
59: .seealso: ISIdentity()
60: @*/
61: PetscErrorCode ISSetIdentity(IS is)
62: {
65: is->isidentity = PETSC_TRUE;
66: return(0);
67: }
71: /*@C
72: ISPermutation - PETSC_TRUE or PETSC_FALSE depending on whether the
73: index set has been declared to be a permutation.
75: Collective on IS
77: Input Parmeters:
78: . is - the index set
80: Output Parameters:
81: . perm - PETSC_TRUE if a permutation, else PETSC_FALSE
83: Level: intermediate
85: Concepts: permutation
86: Concepts: index sets^is permutation
88: .seealso: ISSetPermutation()
89: @*/
90: PetscErrorCode ISPermutation(IS is,PetscTruth *perm)
91: {
95: *perm = (PetscTruth) is->isperm;
96: return(0);
97: }
101: /*@
102: ISSetPermutation - Informs the index set that it is a permutation.
104: Collective on IS
106: Input Parmeters:
107: . is - the index set
109: Level: intermediate
111: Concepts: permutation
112: Concepts: index sets^permutation
114: .seealso: ISPermutation()
115: @*/
116: PetscErrorCode ISSetPermutation(IS is)
117: {
120: is->isperm = PETSC_TRUE;
121: return(0);
122: }
126: /*@C
127: ISDestroy - Destroys an index set.
129: Collective on IS
131: Input Parameters:
132: . is - the index set
134: Level: beginner
136: .seealso: ISCreateGeneral(), ISCreateStride(), ISCreateBlocked()
137: @*/
138: PetscErrorCode ISDestroy(IS is)
139: {
144: if (--is->refct > 0) return(0);
145: (*is->ops->destroy)(is);
146: return(0);
147: }
151: /*@C
152: ISInvertPermutation - Creates a new permutation that is the inverse of
153: a given permutation.
155: Collective on IS
157: Input Parameter:
158: + is - the index set
159: - nlocal - number of indices on this processor in result (ignored for 1 proccessor) or
160: use PETSC_DECIDE
162: Output Parameter:
163: . isout - the inverse permutation
165: Level: intermediate
167: Notes: For parallel index sets this does the complete parallel permutation, but the
168: code is not efficient for huge index sets (10,000,000 indices).
170: Concepts: inverse permutation
171: Concepts: permutation^inverse
172: Concepts: index sets^inverting
173: @*/
174: PetscErrorCode ISInvertPermutation(IS is,PetscInt nlocal,IS *isout)
175: {
181: if (!is->isperm) SETERRQ(PETSC_ERR_ARG_WRONG,"not a permutation");
182: (*is->ops->invertpermutation)(is,nlocal,isout);
183: ISSetPermutation(*isout);
184: return(0);
185: }
187: #if defined(__cplusplus)
190: PetscErrorCode ISGetSizeNew(IS is,PetscInt *size)
191: {
197: is->cops->getsize(size);
198: return(0);
199: }
200: #endif
204: /*@
205: ISGetSize - Returns the global length of an index set.
207: Not Collective
209: Input Parameter:
210: . is - the index set
212: Output Parameter:
213: . size - the global size
215: Level: beginner
217: Concepts: size^of index set
218: Concepts: index sets^size
220: @*/
221: PetscErrorCode ISGetSize(IS is,PetscInt *size)
222: {
228: (*is->ops->getsize)(is,size);
229: return(0);
230: }
234: /*@
235: ISGetLocalSize - Returns the local (processor) length of an index set.
237: Not Collective
239: Input Parameter:
240: . is - the index set
242: Output Parameter:
243: . size - the local size
245: Level: beginner
247: Concepts: size^of index set
248: Concepts: local size^of index set
249: Concepts: index sets^local size
250:
251: @*/
252: PetscErrorCode ISGetLocalSize(IS is,PetscInt *size)
253: {
259: (*is->ops->getlocalsize)(is,size);
260: return(0);
261: }
265: /*@C
266: ISGetIndices - Returns a pointer to the indices. The user should call
267: ISRestoreIndices() after having looked at the indices. The user should
268: NOT change the indices.
270: Not Collective
272: Input Parameter:
273: . is - the index set
275: Output Parameter:
276: . ptr - the location to put the pointer to the indices
278: Fortran Note:
279: This routine is used differently from Fortran
280: $ IS is
281: $ integer is_array(1)
282: $ PetscOffset i_is
283: $ int ierr
284: $ call ISGetIndices(is,is_array,i_is,ierr)
285: $
286: $ Access first local entry in list
287: $ value = is_array(i_is + 1)
288: $
289: $ ...... other code
290: $ call ISRestoreIndices(is,is_array,i_is,ierr)
292: See the Fortran chapter of the users manual and
293: petsc/src/is/examples/[tutorials,tests] for details.
295: Level: intermediate
297: Concepts: index sets^getting indices
298: Concepts: indices of index set
300: .seealso: ISRestoreIndices(), ISGetIndicesF90()
301: @*/
302: PetscErrorCode ISGetIndices(IS is,PetscInt *ptr[])
303: {
309: (*is->ops->getindices)(is,ptr);
310: return(0);
311: }
315: /*@C
316: ISRestoreIndices - Restores an index set to a usable state after a call
317: to ISGetIndices().
319: Not Collective
321: Input Parameters:
322: + is - the index set
323: - ptr - the pointer obtained by ISGetIndices()
325: Fortran Note:
326: This routine is used differently from Fortran
327: $ IS is
328: $ integer is_array(1)
329: $ PetscOffset i_is
330: $ int ierr
331: $ call ISGetIndices(is,is_array,i_is,ierr)
332: $
333: $ Access first local entry in list
334: $ value = is_array(i_is + 1)
335: $
336: $ ...... other code
337: $ call ISRestoreIndices(is,is_array,i_is,ierr)
339: See the Fortran chapter of the users manual and
340: petsc/src/is/examples/[tutorials,tests] for details.
342: Level: intermediate
344: .seealso: ISGetIndices(), ISRestoreIndicesF90()
345: @*/
346: PetscErrorCode ISRestoreIndices(IS is,PetscInt *ptr[])
347: {
353: if (is->ops->restoreindices) {
354: (*is->ops->restoreindices)(is,ptr);
355: }
356: return(0);
357: }
361: /*@C
362: ISView - Displays an index set.
364: Collective on IS
366: Input Parameters:
367: + is - the index set
368: - viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.
370: Level: intermediate
372: .seealso: PetscViewerASCIIOpen()
373: @*/
374: PetscErrorCode ISView(IS is,PetscViewer viewer)
375: {
380: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(is->comm);
383:
384: (*is->ops->view)(is,viewer);
385: return(0);
386: }
390: /*@
391: ISSort - Sorts the indices of an index set.
393: Collective on IS
395: Input Parameters:
396: . is - the index set
398: Level: intermediate
400: Concepts: index sets^sorting
401: Concepts: sorting^index set
403: .seealso: ISSorted()
404: @*/
405: PetscErrorCode ISSort(IS is)
406: {
411: (*is->ops->sortindices)(is);
412: return(0);
413: }
417: /*@C
418: ISSorted - Checks the indices to determine whether they have been sorted.
420: Collective on IS
422: Input Parameter:
423: . is - the index set
425: Output Parameter:
426: . flg - output flag, either PETSC_TRUE if the index set is sorted,
427: or PETSC_FALSE otherwise.
429: Notes: For parallel IS objects this only indicates if the local part of the IS
430: is sorted. So some processors may return PETSC_TRUE while others may
431: return PETSC_FALSE.
433: Level: intermediate
435: .seealso: ISSort()
436: @*/
437: PetscErrorCode ISSorted(IS is,PetscTruth *flg)
438: {
444: (*is->ops->sorted)(is,flg);
445: return(0);
446: }
450: /*@C
451: ISDuplicate - Creates a duplicate copy of an index set.
453: Collective on IS
455: Input Parmeters:
456: . is - the index set
458: Output Parameters:
459: . isnew - the copy of the index set
461: Level: beginner
463: Concepts: index sets^duplicating
465: .seealso: ISCreateGeneral()
466: @*/
467: PetscErrorCode ISDuplicate(IS is,IS *newIS)
468: {
474: (*is->ops->duplicate)(is,newIS);
475: return(0);
476: }
479: /*MC
480: ISGetIndicesF90 - Accesses the elements of an index set from Fortran90.
481: The users should call ISRestoreIndicesF90() after having looked at the
482: indices. The user should NOT change the indices.
484: Synopsis:
485: ISGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
487: Not collective
489: Input Parameter:
490: . x - index set
492: Output Parameters:
493: + xx_v - the Fortran90 pointer to the array
494: - ierr - error code
496: Example of Usage:
497: .vb
498: PetscScalar, pointer xx_v(:)
499: ....
500: call ISGetIndicesF90(x,xx_v,ierr)
501: a = xx_v(3)
502: call ISRestoreIndicesF90(x,xx_v,ierr)
503: .ve
505: Notes:
506: Not yet supported for all F90 compilers.
508: Level: intermediate
510: .seealso: ISRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices()
512: Concepts: index sets^getting indices in f90
513: Concepts: indices of index set in f90
515: M*/
517: /*MC
518: ISRestoreIndicesF90 - Restores an index set to a usable state after
519: a call to ISGetIndicesF90().
521: Synopsis:
522: ISRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
524: Not collective
526: Input Parameters:
527: . x - index set
528: . xx_v - the Fortran90 pointer to the array
530: Output Parameter:
531: . ierr - error code
534: Example of Usage:
535: .vb
536: PetscScalar, pointer xx_v(:)
537: ....
538: call ISGetIndicesF90(x,xx_v,ierr)
539: a = xx_v(3)
540: call ISRestoreIndicesF90(x,xx_v,ierr)
541: .ve
542:
543: Notes:
544: Not yet supported for all F90 compilers.
546: Level: intermediate
548: .seealso: ISGetIndicesF90(), ISGetIndices(), ISRestoreIndices()
550: M*/
552: /*MC
553: ISBlockGetIndicesF90 - Accesses the elements of an index set from Fortran90.
554: The users should call ISBlockRestoreIndicesF90() after having looked at the
555: indices. The user should NOT change the indices.
557: Synopsis:
558: ISBlockGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
560: Not collective
562: Input Parameter:
563: . x - index set
565: Output Parameters:
566: + xx_v - the Fortran90 pointer to the array
567: - ierr - error code
568: Example of Usage:
569: .vb
570: PetscScalar, pointer xx_v(:)
571: ....
572: call ISBlockGetIndicesF90(x,xx_v,ierr)
573: a = xx_v(3)
574: call ISBlockRestoreIndicesF90(x,xx_v,ierr)
575: .ve
577: Notes:
578: Not yet supported for all F90 compilers
580: Level: intermediate
582: .seealso: ISBlockRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices(),
583: ISRestoreIndices()
585: Concepts: index sets^getting block indices in f90
586: Concepts: indices of index set in f90
587: Concepts: block^ indices of index set in f90
589: M*/
591: /*MC
592: ISBlockRestoreIndicesF90 - Restores an index set to a usable state after
593: a call to ISBlockGetIndicesF90().
595: Synopsis:
596: ISBlockRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
598: Input Parameters:
599: + x - index set
600: - xx_v - the Fortran90 pointer to the array
602: Output Parameter:
603: . ierr - error code
605: Example of Usage:
606: .vb
607: PetscScalar, pointer xx_v(:)
608: ....
609: call ISBlockGetIndicesF90(x,xx_v,ierr)
610: a = xx_v(3)
611: call ISBlockRestoreIndicesF90(x,xx_v,ierr)
612: .ve
613:
614: Notes:
615: Not yet supported for all F90 compilers
617: Level: intermediate
619: .seealso: ISBlockGetIndicesF90(), ISGetIndices(), ISRestoreIndices(), ISRestoreIndicesF90()
621: M*/