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 PETSCVEC_DLLEXPORT IS_COOKIE = 0;

 13: /*@C
 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 PETSCVEC_DLLEXPORT 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 PETSCVEC_DLLEXPORT ISSetIdentity(IS is)
 64: {
 67:   is->isidentity = PETSC_TRUE;
 68:   return(0);
 69: }

 73: /*@C
 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 PETSCVEC_DLLEXPORT 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 PETSCVEC_DLLEXPORT ISSetPermutation(IS is)
119: {
122:   is->isperm = PETSC_TRUE;
123:   return(0);
124: }

128: /*@C
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 PETSCVEC_DLLEXPORT ISDestroy(IS is)
141: {

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

153: /*@C
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 PETSCVEC_DLLEXPORT 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: }

189: #if defined(__cplusplus)
192: PetscErrorCode PETSCVEC_DLLEXPORT ISGetSizeNew(IS is,PetscInt *size)
193: {

199:   is->cops->getsize(size);
200:   return(0);
201: }
202: #endif

206: /*@
207:    ISGetSize - Returns the global length of an index set. 

209:    Not Collective

211:    Input Parameter:
212: .  is - the index set

214:    Output Parameter:
215: .  size - the global size

217:    Level: beginner

219:    Concepts: size^of index set
220:    Concepts: index sets^size

222: @*/
223: PetscErrorCode PETSCVEC_DLLEXPORT ISGetSize(IS is,PetscInt *size)
224: {

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

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

239:    Not Collective

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

244:    Output Parameter:
245: .  size - the local size

247:    Level: beginner

249:    Concepts: size^of index set
250:    Concepts: local size^of index set
251:    Concepts: index sets^local size
252:   
253: @*/
254: PetscErrorCode PETSCVEC_DLLEXPORT ISGetLocalSize(IS is,PetscInt *size)
255: {

261:   (*is->ops->getlocalsize)(is,size);
262:   return(0);
263: }

267: /*@C
268:    ISGetIndices - Returns a pointer to the indices.  The user should call 
269:    ISRestoreIndices() after having looked at the indices.  The user should 
270:    NOT change the indices.

272:    Not Collective

274:    Input Parameter:
275: .  is - the index set

277:    Output Parameter:
278: .  ptr - the location to put the pointer to the indices

280:    Fortran Note:
281:    This routine is used differently from Fortran
282: $    IS          is
283: $    integer     is_array(1)
284: $    PetscOffset i_is
285: $    int         ierr
286: $       call ISGetIndices(is,is_array,i_is,ierr)
287: $
288: $   Access first local entry in list
289: $      value = is_array(i_is + 1)
290: $
291: $      ...... other code
292: $       call ISRestoreIndices(is,is_array,i_is,ierr)

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

297:    Level: intermediate

299:    Concepts: index sets^getting indices
300:    Concepts: indices of index set

302: .seealso: ISRestoreIndices(), ISGetIndicesF90()
303: @*/
304: PetscErrorCode PETSCVEC_DLLEXPORT ISGetIndices(IS is,PetscInt *ptr[])
305: {

311:   (*is->ops->getindices)(is,ptr);
312:   return(0);
313: }

317: /*@C
318:    ISRestoreIndices - Restores an index set to a usable state after a call 
319:                       to ISGetIndices().

321:    Not Collective

323:    Input Parameters:
324: +  is - the index set
325: -  ptr - the pointer obtained by ISGetIndices()

327:    Fortran Note:
328:    This routine is used differently from Fortran
329: $    IS          is
330: $    integer     is_array(1)
331: $    PetscOffset i_is
332: $    int         ierr
333: $       call ISGetIndices(is,is_array,i_is,ierr)
334: $
335: $   Access first local entry in list
336: $      value = is_array(i_is + 1)
337: $
338: $      ...... other code
339: $       call ISRestoreIndices(is,is_array,i_is,ierr)

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

344:    Level: intermediate

346: .seealso: ISGetIndices(), ISRestoreIndicesF90()
347: @*/
348: PetscErrorCode PETSCVEC_DLLEXPORT ISRestoreIndices(IS is,PetscInt *ptr[])
349: {

355:   if (is->ops->restoreindices) {
356:     (*is->ops->restoreindices)(is,ptr);
357:   }
358:   return(0);
359: }

363: /*@C
364:    ISView - Displays an index set.

366:    Collective on IS

368:    Input Parameters:
369: +  is - the index set
370: -  viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.

372:    Level: intermediate

374: .seealso: PetscViewerASCIIOpen()
375: @*/
376: PetscErrorCode PETSCVEC_DLLEXPORT ISView(IS is,PetscViewer viewer)
377: {

382:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(is->comm);
385: 
386:   (*is->ops->view)(is,viewer);
387:   return(0);
388: }

392: /*@
393:    ISSort - Sorts the indices of an index set.

395:    Collective on IS

397:    Input Parameters:
398: .  is - the index set

400:    Level: intermediate

402:    Concepts: index sets^sorting
403:    Concepts: sorting^index set

405: .seealso: ISSorted()
406: @*/
407: PetscErrorCode PETSCVEC_DLLEXPORT ISSort(IS is)
408: {

413:   (*is->ops->sortindices)(is);
414:   return(0);
415: }

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

422:    Collective on IS

424:    Input Parameter:
425: .  is - the index set

427:    Output Parameter:
428: .  flg - output flag, either PETSC_TRUE if the index set is sorted, 
429:          or PETSC_FALSE otherwise.

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

435:    Level: intermediate

437: .seealso: ISSort()
438: @*/
439: PetscErrorCode PETSCVEC_DLLEXPORT ISSorted(IS is,PetscTruth *flg)
440: {

446:   (*is->ops->sorted)(is,flg);
447:   return(0);
448: }

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

455:    Collective on IS

457:    Input Parmeters:
458: .  is - the index set

460:    Output Parameters:
461: .  isnew - the copy of the index set

463:    Level: beginner

465:    Concepts: index sets^duplicating

467: .seealso: ISCreateGeneral()
468: @*/
469: PetscErrorCode PETSCVEC_DLLEXPORT ISDuplicate(IS is,IS *newIS)
470: {

476:   (*is->ops->duplicate)(is,newIS);
477:   return(0);
478: }


481: /*MC
482:     ISGetIndicesF90 - Accesses the elements of an index set from Fortran90.
483:     The users should call ISRestoreIndicesF90() after having looked at the
484:     indices.  The user should NOT change the indices.

486:     Synopsis:
487:     ISGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

489:     Not collective

491:     Input Parameter:
492: .   x - index set

494:     Output Parameters:
495: +   xx_v - the Fortran90 pointer to the array
496: -   ierr - error code

498:     Example of Usage: 
499: .vb
500:     PetscScalar, pointer xx_v(:)
501:     ....
502:     call ISGetIndicesF90(x,xx_v,ierr)
503:     a = xx_v(3)
504:     call ISRestoreIndicesF90(x,xx_v,ierr)
505: .ve

507:     Notes:
508:     Not yet supported for all F90 compilers.

510:     Level: intermediate

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

514:   Concepts: index sets^getting indices in f90
515:   Concepts: indices of index set in f90

517: M*/

519: /*MC
520:     ISRestoreIndicesF90 - Restores an index set to a usable state after
521:     a call to ISGetIndicesF90().

523:     Synopsis:
524:     ISRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

526:     Not collective

528:     Input Parameters:
529: .   x - index set
530: .   xx_v - the Fortran90 pointer to the array

532:     Output Parameter:
533: .   ierr - error code


536:     Example of Usage: 
537: .vb
538:     PetscScalar, pointer xx_v(:)
539:     ....
540:     call ISGetIndicesF90(x,xx_v,ierr)
541:     a = xx_v(3)
542:     call ISRestoreIndicesF90(x,xx_v,ierr)
543: .ve
544:    
545:     Notes:
546:     Not yet supported for all F90 compilers.

548:     Level: intermediate

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

552: M*/

554: /*MC
555:     ISBlockGetIndicesF90 - Accesses the elements of an index set from Fortran90.
556:     The users should call ISBlockRestoreIndicesF90() after having looked at the
557:     indices.  The user should NOT change the indices.

559:     Synopsis:
560:     ISBlockGetIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

562:     Not collective

564:     Input Parameter:
565: .   x - index set

567:     Output Parameters:
568: +   xx_v - the Fortran90 pointer to the array
569: -   ierr - error code
570:     Example of Usage: 
571: .vb
572:     PetscScalar, pointer xx_v(:)
573:     ....
574:     call ISBlockGetIndicesF90(x,xx_v,ierr)
575:     a = xx_v(3)
576:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
577: .ve

579:     Notes:
580:     Not yet supported for all F90 compilers

582:     Level: intermediate

584: .seealso:  ISBlockRestoreIndicesF90(), ISGetIndices(), ISRestoreIndices(),
585:            ISRestoreIndices()

587:   Concepts: index sets^getting block indices in f90
588:   Concepts: indices of index set in f90
589:   Concepts: block^ indices of index set in f90

591: M*/

593: /*MC
594:     ISBlockRestoreIndicesF90 - Restores an index set to a usable state after
595:     a call to ISBlockGetIndicesF90().

597:     Synopsis:
598:     ISBlockRestoreIndicesF90(IS x,{integer, pointer :: xx_v(:)},integer ierr)

600:     Input Parameters:
601: +   x - index set
602: -   xx_v - the Fortran90 pointer to the array

604:     Output Parameter:
605: .   ierr - error code

607:     Example of Usage: 
608: .vb
609:     PetscScalar, pointer xx_v(:)
610:     ....
611:     call ISBlockGetIndicesF90(x,xx_v,ierr)
612:     a = xx_v(3)
613:     call ISBlockRestoreIndicesF90(x,xx_v,ierr)
614: .ve
615:    
616:     Notes:
617:     Not yet supported for all F90 compilers

619:     Level: intermediate

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

623: M*/