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: /*@C
11: ISIdentity - Determines whether index set is the identity mapping.
13: Collective on IS
15: Input Parmeters:
16: . is - the index set
18: Output Parameters:
19: . ident - PETSC_TRUE if an identity, else PETSC_FALSE
21: Level: intermediate
23: Concepts: identity mapping
24: Concepts: index sets^is identity
26: .seealso: ISSetIdentity()
27: @*/
28: int ISIdentity(IS is,PetscTruth *ident)
29: {
35: *ident = is->isidentity;
36: if (*ident) return(0);
37: if (is->ops->identity) {
38: (*is->ops->identity)(is,ident);
39: }
40: return(0);
41: }
43: /*@
44: ISSetIdentity - Informs the index set that it is an identity.
46: Collective on IS
48: Input Parmeters:
49: . is - the index set
51: Level: intermediate
53: Concepts: identity mapping
54: Concepts: index sets^is identity
56: .seealso: ISIdentity()
57: @*/
58: int ISSetIdentity(IS is)
59: {
62: is->isidentity = PETSC_TRUE;
63: return(0);
64: }
66: /*@C
67: ISPermutation - PETSC_TRUE or PETSC_FALSE depending on whether the
68: index set has been declared to be a permutation.
70: Collective on IS
72: Input Parmeters:
73: . is - the index set
75: Output Parameters:
76: . perm - PETSC_TRUE if a permutation, else PETSC_FALSE
78: Level: intermediate
80: Concepts: permutation
81: Concepts: index sets^is permutation
83: .seealso: ISSetPermutation()
84: @*/
85: int ISPermutation(IS is,PetscTruth *perm)
86: {
90: *perm = (PetscTruth) is->isperm;
91: return(0);
92: }
94: /*@
95: ISSetPermutation - Informs the index set that it is a permutation.
97: Collective on IS
99: Input Parmeters:
100: . is - the index set
102: Level: intermediate
104: Concepts: permutation
105: Concepts: index sets^permutation
107: .seealso: ISPermutation()
108: @*/
109: int ISSetPermutation(IS is)
110: {
113: is->isperm = PETSC_TRUE;
114: return(0);
115: }
117: /*@C
118: ISDestroy - Destroys an index set.
120: Collective on IS
122: Input Parameters:
123: . is - the index set
125: Level: beginner
127: .seealso: ISCreateGeneral(), ISCreateStride(), ISCreateBlocked()
128: @*/
129: int ISDestroy(IS is)
130: {
135: if (--is->refct > 0) return(0);
137: /* if memory was published with AMS then destroy it */
138: PetscObjectDepublish(is);
140: (*is->ops->destroy)(is);
141: return(0);
142: }
144: /*@C
145: ISInvertPermutation - Creates a new permutation that is the inverse of
146: a given permutation.
148: Collective on IS
150: Input Parameter:
151: + is - the index set
152: - nlocal - number of indices on this processor in result (ignored for 1 proccessor) or
153: use PETSC_DECIDE
155: Output Parameter:
156: . isout - the inverse permutation
158: Level: intermediate
160: Notes: For parallel index sets this does the complete parallel permutation, but the
161: code is not efficient for huge index sets (10,000,000 indices).
163: Concepts: inverse permutation
164: Concepts: permutation^inverse
165: Concepts: index sets^inverting
166: @*/
167: int ISInvertPermutation(IS is,int nlocal,IS *isout)
168: {
173: if (!is->isperm) SETERRQ(PETSC_ERR_ARG_WRONG,"not a permutation");
174: (*is->ops->invertpermutation)(is,nlocal,isout);
175: ISSetPermutation(*isout);
176: return(0);
177: }
179: #if defined(__cplusplus)
180: int ISGetSizeNew(IS is,int *size)
181: {
187: is->cops->getsize(size);
188: return(0);
189: }
190: #endif
192: /*@
193: ISGetSize - Returns the global length of an index set.
195: Not Collective
197: Input Parameter:
198: . is - the index set
200: Output Parameter:
201: . size - the global size
203: Level: beginner
205: Concepts: size^of index set
206: Concepts: index sets^size
208: @*/
209: int ISGetSize(IS is,int *size)
210: {
216: (*is->ops->getsize)(is,size);
217: return(0);
218: }
220: /*@
221: ISGetLocalSize - Returns the local (processor) length of an index set.
223: Not Collective
225: Input Parameter:
226: . is - the index set
228: Output Parameter:
229: . size - the local size
231: Level: beginner
233: Concepts: size^of index set
234: Concepts: local size^of index set
235: Concepts: index sets^local size
236:
237: @*/
238: int ISGetLocalSize(IS is,int *size)
239: {
245: (*is->ops->getlocalsize)(is,size);
246: return(0);
247: }
249: /*@C
250: ISGetIndices - Returns a pointer to the indices. The user should call
251: ISRestoreIndices() after having looked at the indices. The user should
252: NOT change the indices.
254: Not Collective
256: Input Parameter:
257: . is - the index set
259: Output Parameter:
260: . ptr - the location to put the pointer to the indices
262: Fortran Note:
263: This routine is used differently from Fortran
264: $ IS is
265: $ integer is_array(1)
266: $ PetscOffset i_is
267: $ int ierr
268: $ call ISGetIndices(is,is_array,i_is,ierr)
269: $
270: $ Access first local entry in list
271: $ value = is_array(i_is + 1)
272: $
273: $ ...... other code
274: $ call ISRestoreIndices(is,is_array,i_is,ierr)
276: See the Fortran chapter of the users manual and
277: petsc/src/is/examples/[tutorials,tests] for details.
279: Level: intermediate
281: Concepts: index sets^getting indices
282: Concepts: indices of index set
284: .seealso: ISRestoreIndices(), ISGetIndicesF90()
285: @*/
286: int ISGetIndices(IS is,int *ptr[])
287: {
293: (*is->ops->getindices)(is,ptr);
294: return(0);
295: }
297: /*@C
298: ISRestoreIndices - Restores an index set to a usable state after a call
299: to ISGetIndices().
301: Not Collective
303: Input Parameters:
304: + is - the index set
305: - ptr - the pointer obtained by ISGetIndices()
307: Fortran Note:
308: This routine is used differently from Fortran
309: $ IS is
310: $ integer is_array(1)
311: $ PetscOffset i_is
312: $ int ierr
313: $ call ISGetIndices(is,is_array,i_is,ierr)
314: $
315: $ Access first local entry in list
316: $ value = is_array(i_is + 1)
317: $
318: $ ...... other code
319: $ call ISRestoreIndices(is,is_array,i_is,ierr)
321: See the Fortran chapter of the users manual and
322: petsc/src/is/examples/[tutorials,tests] for details.
324: Level: intermediate
326: .seealso: ISGetIndices(), ISRestoreIndicesF90()
327: @*/
328: int ISRestoreIndices(IS is,int *ptr[])
329: {
335: if (is->ops->restoreindices) {
336: (*is->ops->restoreindices)(is,ptr);
337: }
338: return(0);
339: }
341: /*@C
342: ISView - Displays an index set.
344: Collective on IS
346: Input Parameters:
347: + is - the index set
348: - viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.
350: Level: intermediate
352: .seealso: PetscViewerASCIIOpen()
353: @*/
354: int ISView(IS is,PetscViewer viewer)
355: {
360: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(is->comm);
363:
364: (*is->ops->view)(is,viewer);
365: return(0);
366: }
368: /*@
369: ISSort - Sorts the indices of an index set.
371: Collective on IS
373: Input Parameters:
374: . is - the index set
376: Level: intermediate
378: Concepts: index sets^sorting
379: Concepts: sorting^index set
381: .seealso: ISSorted()
382: @*/
383: int ISSort(IS is)
384: {
389: (*is->ops->sortindices)(is);
390: return(0);
391: }
393: /*@C
394: ISSorted - Checks the indices to determine whether they have been sorted.
396: Collective on IS
398: Input Parameter:
399: . is - the index set
401: Output Parameter:
402: . flg - output flag, either PETSC_TRUE if the index set is sorted,
403: or PETSC_FALSE otherwise.
405: Level: intermediate
407: .seealso: ISSort()
408: @*/
409: int ISSorted(IS is,PetscTruth *flg)
410: {
416: (*is->ops->sorted)(is,flg);
417: return(0);
418: }
420: /*@C
421: ISDuplicate - Creates a duplicate copy of an index set.
423: Collective on IS
425: Input Parmeters:
426: . is - the index set
428: Output Parameters:
429: . isnew - the copy of the index set
431: Level: beginner
433: Concepts: index sets^duplicating
435: .seealso: ISCreateGeneral()
436: @*/
437: int ISDuplicate(IS is,IS *newIS)
438: {
444: (*is->ops->duplicate)(is,newIS);
445: return(0);
446: }
449: /*MC
450: ISGetIndicesF90 - Accesses the elements of an index set from Fortran90.
451: The users should call ISRestoreIndicesF90() after having looked at the
452: indices. The user should NOT change the indices.
454: Synopsis:
455: ISGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
457: Not collective
459: Input Parameter:
460: . x - index set
462: Output Parameters:
463: + xx_v - the Fortran90 pointer to the array
464: - ierr - error code
466: Example of Usage:
467: .vb
468: PetscScalar, pointer xx_v(:)
469: ....
470: call ISGetIndicesF90(x,xx_v,ierr)
471: a = xx_v(3)
472: call ISRestoreIndicesF90(x,xx_v,ierr)
473: .ve
475: Notes:
476: Not yet supported for all F90 compilers.
478: Level: intermediate
480: .seealso: ISRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices()
482: Concepts: index sets^getting indices in f90
483: Concepts: indices of index set in f90
485: M*/
487: /*MC
488: ISRestoreIndicesF90 - Restores an index set to a usable state after
489: a call to ISGetIndicesF90().
491: Synopsis:
492: ISRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
494: Not collective
496: Input Parameters:
497: . x - index set
498: . xx_v - the Fortran90 pointer to the array
500: Output Parameter:
501: . ierr - error code
504: Example of Usage:
505: .vb
506: PetscScalar, pointer xx_v(:)
507: ....
508: call ISGetIndicesF90(x,xx_v,ierr)
509: a = xx_v(3)
510: call ISRestoreIndicesF90(x,xx_v,ierr)
511: .ve
512:
513: Notes:
514: Not yet supported for all F90 compilers.
516: Level: intermediate
518: .seealso: ISGetIndicesF90(), ISGetIndices(), ISRestoreIndices()
520: M*/
522: /*MC
523: ISBlockGetIndicesF90 - Accesses the elements of an index set from Fortran90.
524: The users should call ISBlockRestoreIndicesF90() after having looked at the
525: indices. The user should NOT change the indices.
527: Synopsis:
528: ISBlockGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
530: Not collective
532: Input Parameter:
533: . x - index set
535: Output Parameters:
536: + xx_v - the Fortran90 pointer to the array
537: - ierr - error code
538: Example of Usage:
539: .vb
540: PetscScalar, pointer xx_v(:)
541: ....
542: call ISBlockGetIndicesF90(x,xx_v,ierr)
543: a = xx_v(3)
544: call ISBlockRestoreIndicesF90(x,xx_v,ierr)
545: .ve
547: Notes:
548: Not yet supported for all F90 compilers
550: Level: intermediate
552: .seealso: ISBlockRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices(),
553: ISRestoreIndices()
555: Concepts: index sets^getting block indices in f90
556: Concepts: indices of index set in f90
557: Concepts: block^ indices of index set in f90
559: M*/
561: /*MC
562: ISBlockRestoreIndicesF90 - Restores an index set to a usable state after
563: a call to ISBlockGetIndicesF90().
565: Synopsis:
566: ISBlockRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)
568: Input Parameters:
569: + x - index set
570: - xx_v - the Fortran90 pointer to the array
572: Output Parameter:
573: . ierr - error code
575: Example of Usage:
576: .vb
577: PetscScalar, pointer xx_v(:)
578: ....
579: call ISBlockGetIndicesF90(x,xx_v,ierr)
580: a = xx_v(3)
581: call ISBlockRestoreIndicesF90(x,xx_v,ierr)
582: .ve
583:
584: Notes:
585: Not yet supported for all F90 compilers
587: Level: intermediate
589: .seealso: ISBlockGetIndicesF90(), ISGetIndices(), ISRestoreIndices(), ISRestoreIndicesF90()
591: M*/