Actual source code: index.c

  1: /*$Id: index.c,v 1.82 2001/03/23 23:21:06 balay Exp $*/
  2: /*  
  3:    Defines the abstract operations on index sets, i.e. the public interface. 
  4: */
  5: #include "src/vec/is/isimpl.h" 
  6:      /*I "petscis.h" I*/

  8: /*@C
  9:    ISIdentity - Determines whether index set is the identity mapping.

 11:    Collective on IS

 13:    Input Parmeters:
 14: .  is - the index set

 16:    Output Parameters:
 17: .  ident - PETSC_TRUE if an identity, else PETSC_FALSE

 19:    Level: intermediate

 21:    Concepts: identity mapping
 22:    Concepts: index sets^is identity

 24: .seealso: ISSetIdentity()
 25: @*/
 26: int ISIdentity(IS is,PetscTruth *ident)
 27: {

 33:   *ident = is->isidentity;
 34:   if (*ident) return(0);
 35:   if (is->ops->identity) {
 36:     (*is->ops->identity)(is,ident);
 37:   }
 38:   return(0);
 39: }

 41: /*@
 42:    ISSetIdentity - Informs the index set that it is an identity.

 44:    Collective on IS

 46:    Input Parmeters:
 47: .  is - the index set

 49:    Level: intermediate

 51:    Concepts: identity mapping
 52:    Concepts: index sets^is identity

 54: .seealso: ISIdentity()
 55: @*/
 56: int ISSetIdentity(IS is)
 57: {
 60:   is->isidentity = PETSC_TRUE;
 61:   return(0);
 62: }

 64: /*@C
 65:    ISPermutation - PETSC_TRUE or PETSC_FALSE depending on whether the 
 66:    index set has been declared to be a permutation.

 68:    Collective on IS

 70:    Input Parmeters:
 71: .  is - the index set

 73:    Output Parameters:
 74: .  perm - PETSC_TRUE if a permutation, else PETSC_FALSE

 76:    Level: intermediate

 78:   Concepts: permutation
 79:   Concepts: index sets^is permutation

 81: .seealso: ISSetPermutation()
 82: @*/
 83: int ISPermutation(IS is,PetscTruth *perm)
 84: {
 88:   *perm = (PetscTruth) is->isperm;
 89:   return(0);
 90: }

 92: /*@
 93:    ISSetPermutation - Informs the index set that it is a permutation.

 95:    Collective on IS

 97:    Input Parmeters:
 98: .  is - the index set

100:    Level: intermediate

102:   Concepts: permutation
103:   Concepts: index sets^permutation

105: .seealso: ISPermutation()
106: @*/
107: int ISSetPermutation(IS is)
108: {
111:   is->isperm = PETSC_TRUE;
112:   return(0);
113: }

115: /*@C
116:    ISDestroy - Destroys an index set.

118:    Collective on IS

120:    Input Parameters:
121: .  is - the index set

123:    Level: beginner

125: .seealso: ISCreateGeneral(), ISCreateStride(), ISCreateBlocked()
126: @*/
127: int ISDestroy(IS is)
128: {

133:   if (--is->refct > 0) return(0);

135:   /* if memory was published with AMS then destroy it */
136:   PetscObjectDepublish(is);

138:   (*is->ops->destroy)(is);
139:   return(0);
140: }

142: /*@C
143:    ISInvertPermutation - Creates a new permutation that is the inverse of 
144:                          a given permutation.

146:    Collective on IS

148:    Input Parameter:
149: +  is - the index set
150: -  nlocal - number of indices on this processor in result (ignored for 1 proccessor) or
151:             use PETSC_DECIDE

153:    Output Parameter:
154: .  isout - the inverse permutation

156:    Level: intermediate

158:    Notes: For parallel index sets this does the complete parallel permutation, but the 
159:     code is not efficient for huge index sets (10,000,000 indices).

161:    Concepts: inverse permutation
162:    Concepts: permutation^inverse
163:    Concepts: index sets^inverting
164: @*/
165: int ISInvertPermutation(IS is,int nlocal,IS *isout)
166: {

171:   if (!is->isperm) SETERRQ(PETSC_ERR_ARG_WRONG,"not a permutation");
172:   (*is->ops->invertpermutation)(is,nlocal,isout);
173:   return(0);
174: }

176: /*@
177:    ISGetSize - Returns the global length of an index set. 

179:    Not Collective

181:    Input Parameter:
182: .  is - the index set

184:    Output Parameter:
185: .  size - the global size

187:    Level: beginner

189:    Concepts: size^of index set
190:    Concepts: index sets^size

192: @*/
193: int ISGetSize(IS is,int *size)
194: {

200:   (*is->ops->getsize)(is,size);
201:   return(0);
202: }

204: /*@
205:    ISGetLocalSize - Returns the local (processor) length of an index set. 

207:    Not Collective

209:    Input Parameter:
210: .  is - the index set

212:    Output Parameter:
213: .  size - the local size

215:    Level: beginner

217:    Concepts: size^of index set
218:    Concepts: local size^of index set
219:    Concepts: index sets^local size
220:   
221: @*/
222: int ISGetLocalSize(IS is,int *size)
223: {

229:   (*is->ops->getlocalsize)(is,size);
230:   return(0);
231: }

233: /*@C
234:    ISGetIndices - Returns a pointer to the indices.  The user should call 
235:    ISRestoreIndices() after having looked at the indices.  The user should 
236:    NOT change the indices.

238:    Not Collective

240:    Input Parameter:
241: .  is - the index set

243:    Output Parameter:
244: .  ptr - the location to put the pointer to the indices

246:    Fortran Note:
247:    This routine is used differently from Fortran
248: $    IS          is
249: $    integer     is_array(1)
250: $    PetscOffset i_is
251: $    int         ierr
252: $       call ISGetIndices(is,is_array,i_is,ierr)
253: $
254: $   Access first local entry in list
255: $      value = is_array(i_is + 1)
256: $
257: $      ...... other code
258: $       call ISRestoreIndices(is,is_array,i_is,ierr)

260:    See the Fortran chapter of the users manual and 
261:    petsc/src/is/examples/[tutorials,tests] for details.

263:    Level: intermediate

265:    Concepts: index sets^getting indices
266:    Concepts: indices of index set

268: .seealso: ISRestoreIndices(), ISGetIndicesF90()
269: @*/
270: int ISGetIndices(IS is,int *ptr[])
271: {

277:   (*is->ops->getindices)(is,ptr);
278:   return(0);
279: }

281: /*@C
282:    ISRestoreIndices - Restores an index set to a usable state after a call 
283:                       to ISGetIndices().

285:    Not Collective

287:    Input Parameters:
288: +  is - the index set
289: -  ptr - the pointer obtained by ISGetIndices()

291:    Fortran Note:
292:    This routine is used differently from Fortran
293: $    IS          is
294: $    integer     is_array(1)
295: $    PetscOffset i_is
296: $    int         ierr
297: $       call ISGetIndices(is,is_array,i_is,ierr)
298: $
299: $   Access first local entry in list
300: $      value = is_array(i_is + 1)
301: $
302: $      ...... other code
303: $       call ISRestoreIndices(is,is_array,i_is,ierr)

305:    See the Fortran chapter of the users manual and 
306:    petsc/src/is/examples/[tutorials,tests] for details.

308:    Level: intermediate

310: .seealso: ISGetIndices(), ISRestoreIndicesF90()
311: @*/
312: int ISRestoreIndices(IS is,int *ptr[])
313: {

319:   if (is->ops->restoreindices) {
320:     (*is->ops->restoreindices)(is,ptr);
321:   }
322:   return(0);
323: }

325: /*@C
326:    ISView - Displays an index set.

328:    Collective on IS

330:    Input Parameters:
331: +  is - the index set
332: -  viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.

334:    Level: intermediate

336: .seealso: PetscViewerASCIIOpen()
337: @*/
338: int ISView(IS is,PetscViewer viewer)
339: {

344:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(is->comm);
347: 
348:   (*is->ops->view)(is,viewer);
349:   return(0);
350: }

352: /*@
353:    ISSort - Sorts the indices of an index set.

355:    Collective on IS

357:    Input Parameters:
358: .  is - the index set

360:    Level: intermediate

362:    Concepts: index sets^sorting
363:    Concepts: sorting^index set

365: .seealso: ISSorted()
366: @*/
367: int ISSort(IS is)
368: {

373:   (*is->ops->sortindices)(is);
374:   return(0);
375: }

377: /*@C
378:    ISSorted - Checks the indices to determine whether they have been sorted.

380:    Collective on IS

382:    Input Parameter:
383: .  is - the index set

385:    Output Parameter:
386: .  flg - output flag, either PETSC_TRUE if the index set is sorted, 
387:          or PETSC_FALSE otherwise.

389:    Level: intermediate

391: .seealso: ISSort()
392: @*/
393: int ISSorted(IS is,PetscTruth *flg)
394: {

400:   (*is->ops->sorted)(is,flg);
401:   return(0);
402: }

404: /*@C
405:    ISDuplicate - Creates a duplicate copy of an index set.

407:    Collective on IS

409:    Input Parmeters:
410: .  is - the index set

412:    Output Parameters:
413: .  isnew - the copy of the index set

415:    Level: beginner

417:    Concepts: index sets^duplicating

419: .seealso: ISCreateGeneral()
420: @*/
421: int ISDuplicate(IS is,IS *newIS)
422: {

428:   (*is->ops->duplicate)(is,newIS);
429:   return(0);
430: }


433: /*MC
434:     ISGetIndicesF90 - Accesses the elements of an index set from Fortran90.
435:     The users should call ISRestoreIndicesF90() after having looked at the
436:     indices.  The user should NOT change the indices.

438:     Synopsis:
439:     ISGetIndicesF90(IS x,{Scalar, pointer :: xx_v(:)},integer ierr)

441:     Not collective

443:     Input Parameter:
444: .   x - index set

446:     Output Parameters:
447: +   xx_v - the Fortran90 pointer to the array
448: -   ierr - error code

450:     Example of Usage: 
451: .vb
452:     Scalar, pointer xx_v(:)
453:     ....
454:     call ISGetIndicesF90(x,xx_v,ierr)
455:     a = xx_v(3)
456:     call ISRestoreIndicesF90(x,xx_v,ierr)
457: .ve

459:     Notes:
460:     Not yet supported for all F90 compilers.

462:     Level: intermediate

464: .seealso:  ISRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices()

466:   Concepts: index sets^getting indices in f90
467:   Concepts: indices of index set in f90

469: M*/

471: /*MC
472:     ISRestoreIndicesF90 - Restores an index set to a usable state after
473:     a call to ISGetIndicesF90().

475:     Synopsis:
476:     ISRestoreIndicesF90(IS x,{Scalar, pointer :: xx_v(:)},integer ierr)

478:     Not collective

480:     Input Parameters:
481: .   x - index set
482: .   xx_v - the Fortran90 pointer to the array

484:     Output Parameter:
485: .   ierr - error code


488:     Example of Usage: 
489: .vb
490:     Scalar, pointer xx_v(:)
491:     ....
492:     call ISGetIndicesF90(x,xx_v,ierr)
493:     a = xx_v(3)
494:     call ISRestoreIndicesF90(x,xx_v,ierr)
495: .ve
496:    
497:     Notes:
498:     Not yet supported for all F90 compilers.

500:     Level: intermediate

502: .seealso:  ISGetIndicesF90(), ISGetIndices(), ISRestoreIndices()

504: M*/

506: /*MC
507:     ISBlockGetIndicesF90 - Accesses the elements of an index set from Fortran90.
508:     The users should call ISBlockRestoreIndicesF90() after having looked at the
509:     indices.  The user should NOT change the indices.

511:     Synopsis:
512:     ISBlockGetIndicesF90(IS x,{Scalar, pointer :: xx_v(:)},integer ierr)

514:     Not collective

516:     Input Parameter:
517: .   x - index set

519:     Output Parameters:
520: +   xx_v - the Fortran90 pointer to the array
521: -   ierr - error code
522:     Example of Usage: 
523: .vb
524:     Scalar, pointer xx_v(:)
525:     ....
526:     call ISBlockGetIndicesF90(x,xx_v,ierr)
527:     a = xx_v(3)
528:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
529: .ve

531:     Notes:
532:     Not yet supported for all F90 compilers

534:     Level: intermediate

536: .seealso:  ISBlockRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices(),
537:            ISRestoreIndices()

539:   Concepts: index sets^getting block indices in f90
540:   Concepts: indices of index set in f90
541:   Concepts: block^ indices of index set in f90

543: M*/

545: /*MC
546:     ISBlockRestoreIndicesF90 - Restores an index set to a usable state after
547:     a call to ISBlockGetIndicesF90().

549:     Synopsis:
550:     ISBlockRestoreIndicesF90(IS x,{Scalar, pointer :: xx_v(:)},integer ierr)

552:     Input Parameters:
553: +   x - index set
554: -   xx_v - the Fortran90 pointer to the array

556:     Output Parameter:
557: .   ierr - error code

559:     Example of Usage: 
560: .vb
561:     Scalar, pointer xx_v(:)
562:     ....
563:     call ISBlockGetIndicesF90(x,xx_v,ierr)
564:     a = xx_v(3)
565:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
566: .ve
567:    
568:     Notes:
569:     Not yet supported for all F90 compilers

571:     Level: intermediate

573: .seealso:  ISBlockGetIndicesF90(), ISGetIndices(), ISRestoreIndices(), ISRestoreIndicesF90()

575: M*/