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 private/isimpl.h
 7:  #include petscsys.h

  9: /* Logging support */
 10: PetscCookie  IS_COOKIE;

 14: /*@
 15:    ISIdentity - Determines whether index set is the identity mapping.

 17:    Collective on IS

 19:    Input Parmeters:
 20: .  is - the index set

 22:    Output Parameters:
 23: .  ident - PETSC_TRUE if an identity, else PETSC_FALSE

 25:    Level: intermediate

 27:    Concepts: identity mapping
 28:    Concepts: index sets^is identity

 30: .seealso: ISSetIdentity()
 31: @*/
 32: PetscErrorCode  ISIdentity(IS is,PetscTruth *ident)
 33: {

 39:   *ident = is->isidentity;
 40:   if (*ident) return(0);
 41:   if (is->ops->identity) {
 42:     (*is->ops->identity)(is,ident);
 43:   }
 44:   return(0);
 45: }

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

 52:    Collective on IS

 54:    Input Parmeters:
 55: .  is - the index set

 57:    Level: intermediate

 59:    Concepts: identity mapping
 60:    Concepts: index sets^is identity

 62: .seealso: ISIdentity()
 63: @*/
 64: PetscErrorCode  ISSetIdentity(IS is)
 65: {
 68:   is->isidentity = PETSC_TRUE;
 69:   return(0);
 70: }

 74: /*@
 75:    ISPermutation - PETSC_TRUE or PETSC_FALSE depending on whether the 
 76:    index set has been declared to be a permutation.

 78:    Collective on IS

 80:    Input Parmeters:
 81: .  is - the index set

 83:    Output Parameters:
 84: .  perm - PETSC_TRUE if a permutation, else PETSC_FALSE

 86:    Level: intermediate

 88:   Concepts: permutation
 89:   Concepts: index sets^is permutation

 91: .seealso: ISSetPermutation()
 92: @*/
 93: PetscErrorCode  ISPermutation(IS is,PetscTruth *perm)
 94: {
 98:   *perm = (PetscTruth) is->isperm;
 99:   return(0);
100: }

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

107:    Collective on IS

109:    Input Parmeters:
110: .  is - the index set

112:    Level: intermediate

114:   Concepts: permutation
115:   Concepts: index sets^permutation

117:    The debug version of the libraries (config/configure.py --with-debugging=1) checks if the 
118:   index set is actually a permutation. The optimized version just believes you.

120: .seealso: ISPermutation()
121: @*/
122: PetscErrorCode  ISSetPermutation(IS is)
123: {
126: #if defined(PETSC_USE_DEBUG)
127:   {
128:     PetscMPIInt    size;

131:     MPI_Comm_size(((PetscObject)is)->comm,&size);
132:     if (size == 1) {
133:       PetscInt       i,n,*idx;
134:       const PetscInt *iidx;
135: 
136:       ISGetSize(is,&n);
137:       PetscMalloc(n*sizeof(PetscInt),&idx);
138:       ISGetIndices(is,&iidx);
139:       PetscMemcpy(idx,iidx,n*sizeof(PetscInt));
140:       PetscSortInt(n,idx);
141:       for (i=0; i<n; i++) {
142:         if (idx[i] != i) SETERRQ(PETSC_ERR_ARG_WRONG,"Index set is not a permutation");
143:       }
144:       PetscFree(idx);
145:     }
146:   }
147: #endif
148:   is->isperm = PETSC_TRUE;
149:   return(0);
150: }

154: /*@
155:    ISDestroy - Destroys an index set.

157:    Collective on IS

159:    Input Parameters:
160: .  is - the index set

162:    Level: beginner

164: .seealso: ISCreateGeneral(), ISCreateStride(), ISCreateBlocked()
165: @*/
166: PetscErrorCode  ISDestroy(IS is)
167: {

172:   if (--((PetscObject)is)->refct > 0) return(0);
173:   (*is->ops->destroy)(is);
174:   return(0);
175: }

179: /*@
180:    ISInvertPermutation - Creates a new permutation that is the inverse of 
181:                          a given permutation.

183:    Collective on IS

185:    Input Parameter:
186: +  is - the index set
187: -  nlocal - number of indices on this processor in result (ignored for 1 proccessor) or
188:             use PETSC_DECIDE

190:    Output Parameter:
191: .  isout - the inverse permutation

193:    Level: intermediate

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

198:    Concepts: inverse permutation
199:    Concepts: permutation^inverse
200:    Concepts: index sets^inverting
201: @*/
202: PetscErrorCode  ISInvertPermutation(IS is,PetscInt nlocal,IS *isout)
203: {

209:   if (!is->isperm) SETERRQ(PETSC_ERR_ARG_WRONG,"Not a permutation, must call ISSetPermutation() on the IS first");
210:   (*is->ops->invertpermutation)(is,nlocal,isout);
211:   ISSetPermutation(*isout);
212:   return(0);
213: }

217: /*@
218:    ISGetSize - Returns the global length of an index set. 

220:    Not Collective

222:    Input Parameter:
223: .  is - the index set

225:    Output Parameter:
226: .  size - the global size

228:    Level: beginner

230:    Concepts: size^of index set
231:    Concepts: index sets^size

233: @*/
234: PetscErrorCode  ISGetSize(IS is,PetscInt *size)
235: {

241:   (*is->ops->getsize)(is,size);
242:   return(0);
243: }

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

250:    Not Collective

252:    Input Parameter:
253: .  is - the index set

255:    Output Parameter:
256: .  size - the local size

258:    Level: beginner

260:    Concepts: size^of index set
261:    Concepts: local size^of index set
262:    Concepts: index sets^local size
263:   
264: @*/
265: PetscErrorCode  ISGetLocalSize(IS is,PetscInt *size)
266: {

272:   (*is->ops->getlocalsize)(is,size);
273:   return(0);
274: }

278: /*@C
279:    ISGetIndices - Returns a pointer to the indices.  The user should call 
280:    ISRestoreIndices() after having looked at the indices.  The user should 
281:    NOT change the indices.

283:    Not Collective

285:    Input Parameter:
286: .  is - the index set

288:    Output Parameter:
289: .  ptr - the location to put the pointer to the indices

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:    Concepts: index sets^getting indices
311:    Concepts: indices of index set

313: .seealso: ISRestoreIndices(), ISGetIndicesF90()
314: @*/
315: PetscErrorCode  ISGetIndices(IS is,const PetscInt *ptr[])
316: {

322:   (*is->ops->getindices)(is,ptr);
323:   return(0);
324: }

328: /*@C
329:    ISRestoreIndices - Restores an index set to a usable state after a call 
330:                       to ISGetIndices().

332:    Not Collective

334:    Input Parameters:
335: +  is - the index set
336: -  ptr - the pointer obtained by ISGetIndices()

338:    Fortran Note:
339:    This routine is used differently from Fortran
340: $    IS          is
341: $    integer     is_array(1)
342: $    PetscOffset i_is
343: $    int         ierr
344: $       call ISGetIndices(is,is_array,i_is,ierr)
345: $
346: $   Access first local entry in list
347: $      value = is_array(i_is + 1)
348: $
349: $      ...... other code
350: $       call ISRestoreIndices(is,is_array,i_is,ierr)

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

355:    Level: intermediate

357: .seealso: ISGetIndices(), ISRestoreIndicesF90()
358: @*/
359: PetscErrorCode  ISRestoreIndices(IS is,const PetscInt *ptr[])
360: {

366:   if (is->ops->restoreindices) {
367:     (*is->ops->restoreindices)(is,ptr);
368:   }
369:   return(0);
370: }

374: /*@C
375:    ISView - Displays an index set.

377:    Collective on IS

379:    Input Parameters:
380: +  is - the index set
381: -  viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.

383:    Level: intermediate

385: .seealso: PetscViewerASCIIOpen()
386: @*/
387: PetscErrorCode  ISView(IS is,PetscViewer viewer)
388: {

393:   if (!viewer) {
394:     PetscViewerASCIIGetStdout(((PetscObject)is)->comm,&viewer);
395:   }
398: 
399:   (*is->ops->view)(is,viewer);
400:   return(0);
401: }

405: /*@
406:    ISSort - Sorts the indices of an index set.

408:    Collective on IS

410:    Input Parameters:
411: .  is - the index set

413:    Level: intermediate

415:    Concepts: index sets^sorting
416:    Concepts: sorting^index set

418: .seealso: ISSorted()
419: @*/
420: PetscErrorCode  ISSort(IS is)
421: {

426:   (*is->ops->sortindices)(is);
427:   return(0);
428: }

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

435:    Collective on IS

437:    Input Parameter:
438: .  is - the index set

440:    Output Parameter:
441: .  flg - output flag, either PETSC_TRUE if the index set is sorted, 
442:          or PETSC_FALSE otherwise.

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

448:    Level: intermediate

450: .seealso: ISSort()
451: @*/
452: PetscErrorCode  ISSorted(IS is,PetscTruth *flg)
453: {

459:   (*is->ops->sorted)(is,flg);
460:   return(0);
461: }

465: /*@
466:    ISDuplicate - Creates a duplicate copy of an index set.

468:    Collective on IS

470:    Input Parmeters:
471: .  is - the index set

473:    Output Parameters:
474: .  isnew - the copy of the index set

476:    Level: beginner

478:    Concepts: index sets^duplicating

480: .seealso: ISCreateGeneral()
481: @*/
482: PetscErrorCode  ISDuplicate(IS is,IS *newIS)
483: {

489:   (*is->ops->duplicate)(is,newIS);
490:   return(0);
491: }


494: /*MC
495:     ISGetIndicesF90 - Accesses the elements of an index set from Fortran90.
496:     The users should call ISRestoreIndicesF90() after having looked at the
497:     indices.  The user should NOT change the indices.

499:     Synopsis:
500:     ISGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

502:     Not collective

504:     Input Parameter:
505: .   x - index set

507:     Output Parameters:
508: +   xx_v - the Fortran90 pointer to the array
509: -   ierr - error code

511:     Example of Usage: 
512: .vb
513:     PetscScalar, pointer xx_v(:)
514:     ....
515:     call ISGetIndicesF90(x,xx_v,ierr)
516:     a = xx_v(3)
517:     call ISRestoreIndicesF90(x,xx_v,ierr)
518: .ve

520:     Notes:
521:     Not yet supported for all F90 compilers.

523:     Level: intermediate

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

527:   Concepts: index sets^getting indices in f90
528:   Concepts: indices of index set in f90

530: M*/

532: /*MC
533:     ISRestoreIndicesF90 - Restores an index set to a usable state after
534:     a call to ISGetIndicesF90().

536:     Synopsis:
537:     ISRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

539:     Not collective

541:     Input Parameters:
542: .   x - index set
543: .   xx_v - the Fortran90 pointer to the array

545:     Output Parameter:
546: .   ierr - error code


549:     Example of Usage: 
550: .vb
551:     PetscScalar, pointer xx_v(:)
552:     ....
553:     call ISGetIndicesF90(x,xx_v,ierr)
554:     a = xx_v(3)
555:     call ISRestoreIndicesF90(x,xx_v,ierr)
556: .ve
557:    
558:     Notes:
559:     Not yet supported for all F90 compilers.

561:     Level: intermediate

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

565: M*/

567: /*MC
568:     ISBlockGetIndicesF90 - Accesses the elements of an index set from Fortran90.
569:     The users should call ISBlockRestoreIndicesF90() after having looked at the
570:     indices.  The user should NOT change the indices.

572:     Synopsis:
573:     ISBlockGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

575:     Not collective

577:     Input Parameter:
578: .   x - index set

580:     Output Parameters:
581: +   xx_v - the Fortran90 pointer to the array
582: -   ierr - error code
583:     Example of Usage: 
584: .vb
585:     PetscScalar, pointer xx_v(:)
586:     ....
587:     call ISBlockGetIndicesF90(x,xx_v,ierr)
588:     a = xx_v(3)
589:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
590: .ve

592:     Notes:
593:     Not yet supported for all F90 compilers

595:     Level: intermediate

597: .seealso:  ISBlockRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices(),
598:            ISRestoreIndices()

600:   Concepts: index sets^getting block indices in f90
601:   Concepts: indices of index set in f90
602:   Concepts: block^ indices of index set in f90

604: M*/

606: /*MC
607:     ISBlockRestoreIndicesF90 - Restores an index set to a usable state after
608:     a call to ISBlockGetIndicesF90().

610:     Synopsis:
611:     ISBlockRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

613:     Input Parameters:
614: +   x - index set
615: -   xx_v - the Fortran90 pointer to the array

617:     Output Parameter:
618: .   ierr - error code

620:     Example of Usage: 
621: .vb
622:     PetscScalar, pointer xx_v(:)
623:     ....
624:     call ISBlockGetIndicesF90(x,xx_v,ierr)
625:     a = xx_v(3)
626:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
627: .ve
628:    
629:     Notes:
630:     Not yet supported for all F90 compilers

632:     Level: intermediate

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

636: M*/