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

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

 15:    Collective on IS

 17:    Input Parmeters:
 18: .  is - the index set

 20:    Output Parameters:
 21: .  ident - PETSC_TRUE if an identity, else PETSC_FALSE

 23:    Level: intermediate

 25:    Concepts: identity mapping
 26:    Concepts: index sets^is identity

 28: .seealso: ISSetIdentity()
 29: @*/
 30: int ISIdentity(IS is,PetscTruth *ident)
 31: {

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

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

 50:    Collective on IS

 52:    Input Parmeters:
 53: .  is - the index set

 55:    Level: intermediate

 57:    Concepts: identity mapping
 58:    Concepts: index sets^is identity

 60: .seealso: ISIdentity()
 61: @*/
 62: int ISSetIdentity(IS is)
 63: {
 66:   is->isidentity = PETSC_TRUE;
 67:   return(0);
 68: }

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

 76:    Collective on IS

 78:    Input Parmeters:
 79: .  is - the index set

 81:    Output Parameters:
 82: .  perm - PETSC_TRUE if a permutation, else PETSC_FALSE

 84:    Level: intermediate

 86:   Concepts: permutation
 87:   Concepts: index sets^is permutation

 89: .seealso: ISSetPermutation()
 90: @*/
 91: int ISPermutation(IS is,PetscTruth *perm)
 92: {
 96:   *perm = (PetscTruth) is->isperm;
 97:   return(0);
 98: }

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

105:    Collective on IS

107:    Input Parmeters:
108: .  is - the index set

110:    Level: intermediate

112:   Concepts: permutation
113:   Concepts: index sets^permutation

115: .seealso: ISPermutation()
116: @*/
117: int ISSetPermutation(IS is)
118: {
121:   is->isperm = PETSC_TRUE;
122:   return(0);
123: }

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

130:    Collective on IS

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

135:    Level: beginner

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

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

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

150:   (*is->ops->destroy)(is);
151:   return(0);
152: }

156: /*@C
157:    ISInvertPermutation - Creates a new permutation that is the inverse of 
158:                          a given permutation.

160:    Collective on IS

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

167:    Output Parameter:
168: .  isout - the inverse permutation

170:    Level: intermediate

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

175:    Concepts: inverse permutation
176:    Concepts: permutation^inverse
177:    Concepts: index sets^inverting
178: @*/
179: int ISInvertPermutation(IS is,int nlocal,IS *isout)
180: {

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

192: #if defined(__cplusplus)
195: int ISGetSizeNew(IS is,int *size)
196: {

202:   is->cops->getsize(size);
203:   return(0);
204: }
205: #endif

209: /*@
210:    ISGetSize - Returns the global length of an index set. 

212:    Not Collective

214:    Input Parameter:
215: .  is - the index set

217:    Output Parameter:
218: .  size - the global size

220:    Level: beginner

222:    Concepts: size^of index set
223:    Concepts: index sets^size

225: @*/
226: int ISGetSize(IS is,int *size)
227: {

233:   (*is->ops->getsize)(is,size);
234:   return(0);
235: }

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

242:    Not Collective

244:    Input Parameter:
245: .  is - the index set

247:    Output Parameter:
248: .  size - the local size

250:    Level: beginner

252:    Concepts: size^of index set
253:    Concepts: local size^of index set
254:    Concepts: index sets^local size
255:   
256: @*/
257: int ISGetLocalSize(IS is,int *size)
258: {

264:   (*is->ops->getlocalsize)(is,size);
265:   return(0);
266: }

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

275:    Not Collective

277:    Input Parameter:
278: .  is - the index set

280:    Output Parameter:
281: .  ptr - the location to put the pointer to the indices

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

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

300:    Level: intermediate

302:    Concepts: index sets^getting indices
303:    Concepts: indices of index set

305: .seealso: ISRestoreIndices(), ISGetIndicesF90()
306: @*/
307: int ISGetIndices(IS is,int *ptr[])
308: {

314:   (*is->ops->getindices)(is,ptr);
315:   return(0);
316: }

320: /*@C
321:    ISRestoreIndices - Restores an index set to a usable state after a call 
322:                       to ISGetIndices().

324:    Not Collective

326:    Input Parameters:
327: +  is - the index set
328: -  ptr - the pointer obtained by ISGetIndices()

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

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

347:    Level: intermediate

349: .seealso: ISGetIndices(), ISRestoreIndicesF90()
350: @*/
351: int ISRestoreIndices(IS is,int *ptr[])
352: {

358:   if (is->ops->restoreindices) {
359:     (*is->ops->restoreindices)(is,ptr);
360:   }
361:   return(0);
362: }

366: /*@C
367:    ISView - Displays an index set.

369:    Collective on IS

371:    Input Parameters:
372: +  is - the index set
373: -  viewer - viewer used to display the set, for example PETSC_VIEWER_STDOUT_SELF.

375:    Level: intermediate

377: .seealso: PetscViewerASCIIOpen()
378: @*/
379: int ISView(IS is,PetscViewer viewer)
380: {

385:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(is->comm);
388: 
389:   (*is->ops->view)(is,viewer);
390:   return(0);
391: }

395: /*@
396:    ISSort - Sorts the indices of an index set.

398:    Collective on IS

400:    Input Parameters:
401: .  is - the index set

403:    Level: intermediate

405:    Concepts: index sets^sorting
406:    Concepts: sorting^index set

408: .seealso: ISSorted()
409: @*/
410: int ISSort(IS is)
411: {

416:   (*is->ops->sortindices)(is);
417:   return(0);
418: }

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

425:    Collective on IS

427:    Input Parameter:
428: .  is - the index set

430:    Output Parameter:
431: .  flg - output flag, either PETSC_TRUE if the index set is sorted, 
432:          or PETSC_FALSE otherwise.

434:    Level: intermediate

436: .seealso: ISSort()
437: @*/
438: int ISSorted(IS is,PetscTruth *flg)
439: {

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

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

454:    Collective on IS

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

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

462:    Level: beginner

464:    Concepts: index sets^duplicating

466: .seealso: ISCreateGeneral()
467: @*/
468: int ISDuplicate(IS is,IS *newIS)
469: {

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


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

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

488:     Not collective

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

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

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

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

509:     Level: intermediate

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

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

516: M*/

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

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

525:     Not collective

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

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


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

547:     Level: intermediate

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

551: M*/

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

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

561:     Not collective

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

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

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

581:     Level: intermediate

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

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

590: M*/

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

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

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

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

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

618:     Level: intermediate

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

622: M*/