Actual source code: index.c
1: /*$Id: index.c,v 1.85 2001/08/06 21:14:32 bsmith Exp $*/
2: /*
3: Defines the abstract operations on index sets, i.e. the public interface.
4: */
5: #include src/vec/is/isimpl.h
7: /* Logging support */
8: int IS_COOKIE;
10: #undef __FUNCT__
12: /*@C
13: ISIdentity - Determines whether index set is the identity mapping.
15: Collective on IS
17: Input Parmeters:
18: . is - the index set
20: Output Parameters:
21: . ident - PETSC_TRUE if an identity, else PETSC_FALSE
23: Level: intermediate
25: Concepts: identity mapping
26: Concepts: index sets^is identity
28: .seealso: ISSetIdentity()
29: @*/
30: int ISIdentity(IS is,PetscTruth *ident)
31: {
37: *ident = is->isidentity;
38: if (*ident) return(0);
39: if (is->ops->identity) {
40: (*is->ops->identity)(is,ident);
41: }
42: return(0);
43: }
45: #undef __FUNCT__
47: /*@
48: ISSetIdentity - Informs the index set that it is an identity.
50: Collective on IS
52: Input Parmeters:
53: . is - the index set
55: Level: intermediate
57: Concepts: identity mapping
58: Concepts: index sets^is identity
60: .seealso: ISIdentity()
61: @*/
62: int ISSetIdentity(IS is)
63: {
66: is->isidentity = PETSC_TRUE;
67: return(0);
68: }
70: #undef __FUNCT__
72: /*@C
73: ISPermutation - PETSC_TRUE or PETSC_FALSE depending on whether the
74: index set has been declared to be a permutation.
76: Collective on IS
78: Input Parmeters:
79: . is - the index set
81: Output Parameters:
82: . perm - PETSC_TRUE if a permutation, else PETSC_FALSE
84: Level: intermediate
86: Concepts: permutation
87: Concepts: index sets^is permutation
89: .seealso: ISSetPermutation()
90: @*/
91: int ISPermutation(IS is,PetscTruth *perm)
92: {
96: *perm = (PetscTruth) is->isperm;
97: return(0);
98: }
100: #undef __FUNCT__
102: /*@
103: ISSetPermutation - Informs the index set that it is a permutation.
105: Collective on IS
107: Input Parmeters:
108: . is - the index set
110: Level: intermediate
112: Concepts: permutation
113: Concepts: index sets^permutation
115: .seealso: ISPermutation()
116: @*/
117: int ISSetPermutation(IS is)
118: {
121: is->isperm = PETSC_TRUE;
122: return(0);
123: }
125: #undef __FUNCT__
127: /*@C
128: ISDestroy - Destroys an index set.
130: Collective on IS
132: Input Parameters:
133: . is - the index set
135: Level: beginner
137: .seealso: ISCreateGeneral(), ISCreateStride(), ISCreateBlocked()
138: @*/
139: int ISDestroy(IS is)
140: {
145: if (--is->refct > 0) return(0);
147: /* if memory was published with AMS then destroy it */
148: PetscObjectDepublish(is);
150: (*is->ops->destroy)(is);
151: return(0);
152: }
154: #undef __FUNCT__
156: /*@C
157: ISInvertPermutation - Creates a new permutation that is the inverse of
158: a given permutation.
160: Collective on IS
162: Input Parameter:
163: + is - the index set
164: - nlocal - number of indices on this processor in result (ignored for 1 proccessor) or
165: use PETSC_DECIDE
167: Output Parameter:
168: . isout - the inverse permutation
170: Level: intermediate
172: Notes: For parallel index sets this does the complete parallel permutation, but the
173: code is not efficient for huge index sets (10,000,000 indices).
175: Concepts: inverse permutation
176: Concepts: permutation^inverse
177: Concepts: index sets^inverting
178: @*/
179: int ISInvertPermutation(IS is,int nlocal,IS *isout)
180: {
185: if (!is->isperm) SETERRQ(PETSC_ERR_ARG_WRONG,"not a permutation");
186: (*is->ops->invertpermutation)(is,nlocal,isout);
187: ISSetPermutation(*isout);
188: return(0);
189: }
191: #if defined(__cplusplus)
192: #undef __FUNCT__
194: int ISGetSizeNew(IS is,int *size)
195: {
201: is->cops->getsize(size);
202: return(0);
203: }
204: #endif
206: #undef __FUNCT__
208: /*@
209: ISGetSize - Returns the global length of an index set.
211: Not Collective
213: Input Parameter:
214: . is - the index set
216: Output Parameter:
217: . size - the global size
219: Level: beginner
221: Concepts: size^of index set
222: Concepts: index sets^size
224: @*/
225: int ISGetSize(IS is,int *size)
226: {
232: (*is->ops->getsize)(is,size);
233: return(0);
234: }
236: #undef __FUNCT__
238: /*@
239: ISGetLocalSize - Returns the local (processor) length of an index set.
241: Not Collective
243: Input Parameter:
244: . is - the index set
246: Output Parameter:
247: . size - the local size
249: Level: beginner
251: Concepts: size^of index set
252: Concepts: local size^of index set
253: Concepts: index sets^local size
254:
255: @*/
256: int ISGetLocalSize(IS is,int *size)
257: {
263: (*is->ops->getlocalsize)(is,size);
264: return(0);
265: }
267: #undef __FUNCT__
269: /*@C
270: ISGetIndices - Returns a pointer to the indices. The user should call
271: ISRestoreIndices() after having looked at the indices. The user should
272: NOT change the indices.
274: Not Collective
276: Input Parameter:
277: . is - the index set
279: Output Parameter:
280: . ptr - the location to put the pointer to the indices
282: Fortran Note:
283: This routine is used differently from Fortran
284: $ IS is
285: $ integer is_array(1)
286: $ PetscOffset i_is
287: $ int ierr
288: $ call ISGetIndices(is,is_array,i_is,ierr)
289: $
290: $ Access first local entry in list
291: $ value = is_array(i_is + 1)
292: $
293: $ ...... other code
294: $ call ISRestoreIndices(is,is_array,i_is,ierr)
296: See the Fortran chapter of the users manual and
297: petsc/src/is/examples/[tutorials,tests] for details.
299: Level: intermediate
301: Concepts: index sets^getting indices
302: Concepts: indices of index set
304: .seealso: ISRestoreIndices(), ISGetIndicesF90()
305: @*/
306: int ISGetIndices(IS is,int *ptr[])
307: {
313: (*is->ops->getindices)(is,ptr);
314: return(0);
315: }
317: #undef __FUNCT__
319: /*@C
320: ISRestoreIndices - Restores an index set to a usable state after a call
321: to ISGetIndices().
323: Not Collective
325: Input Parameters:
326: + is - the index set
327: - ptr - the pointer obtained by ISGetIndices()
329: Fortran Note:
330: This routine is used differently from Fortran
331: $ IS is
332: $ integer is_array(1)
333: $ PetscOffset i_is
334: $ int ierr
335: $ call ISGetIndices(is,is_array,i_is,ierr)
336: $
337: $ Access first local entry in list
338: $ value = is_array(i_is + 1)
339: $
340: $ ...... other code
341: $ call ISRestoreIndices(is,is_array,i_is,ierr)
343: See the Fortran chapter of the users manual and
344: petsc/src/is/examples/[tutorials,tests] for details.
346: Level: intermediate
348: .seealso: ISGetIndices(), ISRestoreIndicesF90()
349: @*/
350: int ISRestoreIndices(IS is,int *ptr[])
351: {
357: if (is->ops->restoreindices) {
358: (*is->ops->restoreindices)(is,ptr);
359: }
360: return(0);
361: }
363: #undef __FUNCT__
365: /*@C
366: ISView - Displays an index set.
368: Collective on IS
370: Input Parameters:
371: + is - the index set
372: - viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.
374: Level: intermediate
376: .seealso: PetscViewerASCIIOpen()
377: @*/
378: int ISView(IS is,PetscViewer viewer)
379: {
384: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(is->comm);
387:
388: (*is->ops->view)(is,viewer);
389: return(0);
390: }
392: #undef __FUNCT__
394: /*@
395: ISSort - Sorts the indices of an index set.
397: Collective on IS
399: Input Parameters:
400: . is - the index set
402: Level: intermediate
404: Concepts: index sets^sorting
405: Concepts: sorting^index set
407: .seealso: ISSorted()
408: @*/
409: int ISSort(IS is)
410: {
415: (*is->ops->sortindices)(is);
416: return(0);
417: }
419: #undef __FUNCT__
421: /*@C
422: ISSorted - Checks the indices to determine whether they have been sorted.
424: Collective on IS
426: Input Parameter:
427: . is - the index set
429: Output Parameter:
430: . flg - output flag, either PETSC_TRUE if the index set is sorted,
431: or PETSC_FALSE otherwise.
433: Level: intermediate
435: .seealso: ISSort()
436: @*/
437: int ISSorted(IS is,PetscTruth *flg)
438: {
444: (*is->ops->sorted)(is,flg);
445: return(0);
446: }
448: #undef __FUNCT__
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: int 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*/