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 src/vec/is/isimpl.h

  8: /* Logging support */
  9: PetscCookie  IS_COOKIE = 0;

 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: .seealso: ISPermutation()
117: @*/
118: PetscErrorCode  ISSetPermutation(IS is)
119: {
122:   is->isperm = PETSC_TRUE;
123:   return(0);
124: }

128: /*@
129:    ISDestroy - Destroys an index set.

131:    Collective on IS

133:    Input Parameters:
134: .  is - the index set

136:    Level: beginner

138: .seealso: ISCreateGeneral(), ISCreateStride(), ISCreateBlocked()
139: @*/
140: PetscErrorCode  ISDestroy(IS is)
141: {

146:   if (--is->refct > 0) return(0);
147:   (*is->ops->destroy)(is);
148:   return(0);
149: }

153: /*@
154:    ISInvertPermutation - Creates a new permutation that is the inverse of 
155:                          a given permutation.

157:    Collective on IS

159:    Input Parameter:
160: +  is - the index set
161: -  nlocal - number of indices on this processor in result (ignored for 1 proccessor) or
162:             use PETSC_DECIDE

164:    Output Parameter:
165: .  isout - the inverse permutation

167:    Level: intermediate

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

172:    Concepts: inverse permutation
173:    Concepts: permutation^inverse
174:    Concepts: index sets^inverting
175: @*/
176: PetscErrorCode  ISInvertPermutation(IS is,PetscInt nlocal,IS *isout)
177: {

183:   if (!is->isperm) SETERRQ(PETSC_ERR_ARG_WRONG,"not a permutation");
184:   (*is->ops->invertpermutation)(is,nlocal,isout);
185:   ISSetPermutation(*isout);
186:   return(0);
187: }

191: /*@
192:    ISGetSize - Returns the global length of an index set. 

194:    Not Collective

196:    Input Parameter:
197: .  is - the index set

199:    Output Parameter:
200: .  size - the global size

202:    Level: beginner

204:    Concepts: size^of index set
205:    Concepts: index sets^size

207: @*/
208: PetscErrorCode  ISGetSize(IS is,PetscInt *size)
209: {

215:   (*is->ops->getsize)(is,size);
216:   return(0);
217: }

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

224:    Not Collective

226:    Input Parameter:
227: .  is - the index set

229:    Output Parameter:
230: .  size - the local size

232:    Level: beginner

234:    Concepts: size^of index set
235:    Concepts: local size^of index set
236:    Concepts: index sets^local size
237:   
238: @*/
239: PetscErrorCode  ISGetLocalSize(IS is,PetscInt *size)
240: {

246:   (*is->ops->getlocalsize)(is,size);
247:   return(0);
248: }

252: /*@C
253:    ISGetIndices - Returns a pointer to the indices.  The user should call 
254:    ISRestoreIndices() after having looked at the indices.  The user should 
255:    NOT change the indices.

257:    Not Collective

259:    Input Parameter:
260: .  is - the index set

262:    Output Parameter:
263: .  ptr - the location to put the pointer to the indices

265:    Fortran Note:
266:    This routine is used differently from Fortran
267: $    IS          is
268: $    integer     is_array(1)
269: $    PetscOffset i_is
270: $    int         ierr
271: $       call ISGetIndices(is,is_array,i_is,ierr)
272: $
273: $   Access first local entry in list
274: $      value = is_array(i_is + 1)
275: $
276: $      ...... other code
277: $       call ISRestoreIndices(is,is_array,i_is,ierr)

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

282:    Level: intermediate

284:    Concepts: index sets^getting indices
285:    Concepts: indices of index set

287: .seealso: ISRestoreIndices(), ISGetIndicesF90()
288: @*/
289: PetscErrorCode  ISGetIndices(IS is,PetscInt *ptr[])
290: {

296:   (*is->ops->getindices)(is,ptr);
297:   return(0);
298: }

302: /*@C
303:    ISRestoreIndices - Restores an index set to a usable state after a call 
304:                       to ISGetIndices().

306:    Not Collective

308:    Input Parameters:
309: +  is - the index set
310: -  ptr - the pointer obtained by ISGetIndices()

312:    Fortran Note:
313:    This routine is used differently from Fortran
314: $    IS          is
315: $    integer     is_array(1)
316: $    PetscOffset i_is
317: $    int         ierr
318: $       call ISGetIndices(is,is_array,i_is,ierr)
319: $
320: $   Access first local entry in list
321: $      value = is_array(i_is + 1)
322: $
323: $      ...... other code
324: $       call ISRestoreIndices(is,is_array,i_is,ierr)

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

329:    Level: intermediate

331: .seealso: ISGetIndices(), ISRestoreIndicesF90()
332: @*/
333: PetscErrorCode  ISRestoreIndices(IS is,PetscInt *ptr[])
334: {

340:   if (is->ops->restoreindices) {
341:     (*is->ops->restoreindices)(is,ptr);
342:   }
343:   return(0);
344: }

348: /*@C
349:    ISView - Displays an index set.

351:    Collective on IS

353:    Input Parameters:
354: +  is - the index set
355: -  viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.

357:    Level: intermediate

359: .seealso: PetscViewerASCIIOpen()
360: @*/
361: PetscErrorCode  ISView(IS is,PetscViewer viewer)
362: {

367:   if (!viewer) {
368:     PetscViewerASCIIGetStdout(is->comm,&viewer);
369:   }
372: 
373:   (*is->ops->view)(is,viewer);
374:   return(0);
375: }

379: /*@
380:    ISSort - Sorts the indices of an index set.

382:    Collective on IS

384:    Input Parameters:
385: .  is - the index set

387:    Level: intermediate

389:    Concepts: index sets^sorting
390:    Concepts: sorting^index set

392: .seealso: ISSorted()
393: @*/
394: PetscErrorCode  ISSort(IS is)
395: {

400:   (*is->ops->sortindices)(is);
401:   return(0);
402: }

406: /*@
407:    ISSorted - Checks the indices to determine whether they have been sorted.

409:    Collective on IS

411:    Input Parameter:
412: .  is - the index set

414:    Output Parameter:
415: .  flg - output flag, either PETSC_TRUE if the index set is sorted, 
416:          or PETSC_FALSE otherwise.

418:    Notes: For parallel IS objects this only indicates if the local part of the IS
419:           is sorted. So some processors may return PETSC_TRUE while others may 
420:           return PETSC_FALSE.

422:    Level: intermediate

424: .seealso: ISSort()
425: @*/
426: PetscErrorCode  ISSorted(IS is,PetscTruth *flg)
427: {

433:   (*is->ops->sorted)(is,flg);
434:   return(0);
435: }

439: /*@
440:    ISDuplicate - Creates a duplicate copy of an index set.

442:    Collective on IS

444:    Input Parmeters:
445: .  is - the index set

447:    Output Parameters:
448: .  isnew - the copy of the index set

450:    Level: beginner

452:    Concepts: index sets^duplicating

454: .seealso: ISCreateGeneral()
455: @*/
456: PetscErrorCode  ISDuplicate(IS is,IS *newIS)
457: {

463:   (*is->ops->duplicate)(is,newIS);
464:   return(0);
465: }


468: /*MC
469:     ISGetIndicesF90 - Accesses the elements of an index set from Fortran90.
470:     The users should call ISRestoreIndicesF90() after having looked at the
471:     indices.  The user should NOT change the indices.

473:     Synopsis:
474:     ISGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

476:     Not collective

478:     Input Parameter:
479: .   x - index set

481:     Output Parameters:
482: +   xx_v - the Fortran90 pointer to the array
483: -   ierr - error code

485:     Example of Usage: 
486: .vb
487:     PetscScalar, pointer xx_v(:)
488:     ....
489:     call ISGetIndicesF90(x,xx_v,ierr)
490:     a = xx_v(3)
491:     call ISRestoreIndicesF90(x,xx_v,ierr)
492: .ve

494:     Notes:
495:     Not yet supported for all F90 compilers.

497:     Level: intermediate

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

501:   Concepts: index sets^getting indices in f90
502:   Concepts: indices of index set in f90

504: M*/

506: /*MC
507:     ISRestoreIndicesF90 - Restores an index set to a usable state after
508:     a call to ISGetIndicesF90().

510:     Synopsis:
511:     ISRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

513:     Not collective

515:     Input Parameters:
516: .   x - index set
517: .   xx_v - the Fortran90 pointer to the array

519:     Output Parameter:
520: .   ierr - error code


523:     Example of Usage: 
524: .vb
525:     PetscScalar, pointer xx_v(:)
526:     ....
527:     call ISGetIndicesF90(x,xx_v,ierr)
528:     a = xx_v(3)
529:     call ISRestoreIndicesF90(x,xx_v,ierr)
530: .ve
531:    
532:     Notes:
533:     Not yet supported for all F90 compilers.

535:     Level: intermediate

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

539: M*/

541: /*MC
542:     ISBlockGetIndicesF90 - Accesses the elements of an index set from Fortran90.
543:     The users should call ISBlockRestoreIndicesF90() after having looked at the
544:     indices.  The user should NOT change the indices.

546:     Synopsis:
547:     ISBlockGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

549:     Not collective

551:     Input Parameter:
552: .   x - index set

554:     Output Parameters:
555: +   xx_v - the Fortran90 pointer to the array
556: -   ierr - error code
557:     Example of Usage: 
558: .vb
559:     PetscScalar, pointer xx_v(:)
560:     ....
561:     call ISBlockGetIndicesF90(x,xx_v,ierr)
562:     a = xx_v(3)
563:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
564: .ve

566:     Notes:
567:     Not yet supported for all F90 compilers

569:     Level: intermediate

571: .seealso:  ISBlockRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices(),
572:            ISRestoreIndices()

574:   Concepts: index sets^getting block indices in f90
575:   Concepts: indices of index set in f90
576:   Concepts: block^ indices of index set in f90

578: M*/

580: /*MC
581:     ISBlockRestoreIndicesF90 - Restores an index set to a usable state after
582:     a call to ISBlockGetIndicesF90().

584:     Synopsis:
585:     ISBlockRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

587:     Input Parameters:
588: +   x - index set
589: -   xx_v - the Fortran90 pointer to the array

591:     Output Parameter:
592: .   ierr - error code

594:     Example of Usage: 
595: .vb
596:     PetscScalar, pointer xx_v(:)
597:     ....
598:     call ISBlockGetIndicesF90(x,xx_v,ierr)
599:     a = xx_v(3)
600:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
601: .ve
602:    
603:     Notes:
604:     Not yet supported for all F90 compilers

606:     Level: intermediate

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

610: M*/