Actual source code: block.c

petsc-master 2019-12-08
Report Typos and Errors

  2: /*
  3:      Provides the functions for index sets (IS) defined by a list of integers.
  4:    These are for blocks of data, each block is indicated with a single integer.
  5: */
  6:  #include <petsc/private/isimpl.h>
  7:  #include <petscviewer.h>

  9: typedef struct {
 10:   PetscBool sorted;      /* are the blocks sorted? */
 11:   PetscBool allocated;   /* did we allocate the index array ourselves? */
 12:   PetscInt  *idx;
 13: } IS_Block;

 15: static PetscErrorCode ISDestroy_Block(IS is)
 16: {
 17:   IS_Block       *sub = (IS_Block*)is->data;

 21:   if (sub->allocated) {PetscFree(sub->idx);}
 22:   PetscObjectComposeFunction((PetscObject)is,"ISBlockSetIndices_C",NULL);
 23:   PetscObjectComposeFunction((PetscObject)is,"ISBlockGetIndices_C",NULL);
 24:   PetscObjectComposeFunction((PetscObject)is,"ISBlockRestoreIndices_C",NULL);
 25:   PetscObjectComposeFunction((PetscObject)is,"ISBlockGetSize_C",NULL);
 26:   PetscObjectComposeFunction((PetscObject)is,"ISBlockGetLocalSize_C",NULL);
 27:   PetscFree(is->data);
 28:   return(0);
 29: }

 31: static PetscErrorCode ISLocate_Block(IS is,PetscInt key,PetscInt *location)
 32: {
 33:   IS_Block       *sub = (IS_Block*)is->data;
 34:   PetscInt       numIdx, i, bs, bkey, mkey;
 35:   PetscBool      sorted;

 39:   PetscLayoutGetBlockSize(is->map,&bs);
 40:   PetscLayoutGetSize(is->map,&numIdx);
 41:   numIdx /= bs;
 42:   bkey    = key / bs;
 43:   mkey    = key % bs;
 44:   if (mkey < 0) {
 45:     bkey--;
 46:     mkey += bs;
 47:   }
 48:   ISGetInfo(is,IS_SORTED,IS_LOCAL,PETSC_TRUE,&sorted);
 49:   if (sorted) {
 50:     PetscFindInt(bkey,numIdx,sub->idx,location);
 51:   } else {
 52:     const PetscInt *idx = sub->idx;

 54:     *location = -1;
 55:     for (i = 0; i < numIdx; i++) {
 56:       if (idx[i] == bkey) {
 57:         *location = i;
 58:         break;
 59:       }
 60:     }
 61:   }
 62:   if (*location >= 0) {
 63:     *location = *location * bs + mkey;
 64:   }
 65:   return(0);
 66: }

 68: static PetscErrorCode ISGetIndices_Block(IS in,const PetscInt *idx[])
 69: {
 70:   IS_Block       *sub = (IS_Block*)in->data;
 72:   PetscInt       i,j,k,bs,n,*ii,*jj;

 75:   PetscLayoutGetBlockSize(in->map, &bs);
 76:   PetscLayoutGetLocalSize(in->map, &n);
 77:   n   /= bs;
 78:   if (bs == 1) *idx = sub->idx;
 79:   else {
 80:     PetscMalloc1(bs*n,&jj);
 81:     *idx = jj;
 82:     k    = 0;
 83:     ii   = sub->idx;
 84:     for (i=0; i<n; i++)
 85:       for (j=0; j<bs; j++)
 86:         jj[k++] = bs*ii[i] + j;
 87:   }
 88:   return(0);
 89: }

 91: static PetscErrorCode ISRestoreIndices_Block(IS is,const PetscInt *idx[])
 92: {
 93:   IS_Block       *sub = (IS_Block*)is->data;
 94:   PetscInt       bs;

 98:   PetscLayoutGetBlockSize(is->map, &bs);
 99:   if (bs != 1) {
100:     PetscFree(*(void**)idx);
101:   } else {
102:     if (*idx != sub->idx) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Must restore with value from ISGetIndices()");
103:   }
104:   return(0);
105: }

107: static PetscErrorCode ISInvertPermutation_Block(IS is,PetscInt nlocal,IS *isout)
108: {
109:   IS_Block       *sub = (IS_Block*)is->data;
110:   PetscInt       i,*ii,bs,n,*idx = sub->idx;
111:   PetscMPIInt    size;

115:   MPI_Comm_size(PetscObjectComm((PetscObject)is),&size);
116:   PetscLayoutGetBlockSize(is->map, &bs);
117:   PetscLayoutGetLocalSize(is->map, &n);
118:   n   /= bs;
119:   if (size == 1) {
120:     PetscMalloc1(n,&ii);
121:     for (i=0; i<n; i++) ii[idx[i]] = i;
122:     ISCreateBlock(PETSC_COMM_SELF,bs,n,ii,PETSC_OWN_POINTER,isout);
123:     ISSetPermutation(*isout);
124:   } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No inversion written yet for block IS");
125:   return(0);
126: }

128: static PetscErrorCode ISView_Block(IS is, PetscViewer viewer)
129: {
130:   IS_Block       *sub = (IS_Block*)is->data;
132:   PetscInt       i,bs,n,*idx = sub->idx;
133:   PetscBool      iascii;

136:   PetscLayoutGetBlockSize(is->map, &bs);
137:   PetscLayoutGetLocalSize(is->map, &n);
138:   n   /= bs;
139:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
140:   if (iascii) {
141:     PetscViewerFormat fmt;

143:     PetscViewerGetFormat(viewer,&fmt);
144:     if (fmt == PETSC_VIEWER_ASCII_MATLAB) {
145:       IS             ist;
146:       const char     *name;
147:       const PetscInt *idx;
148:       PetscInt       n;

150:       PetscObjectGetName((PetscObject)is,&name);
151:       ISGetLocalSize(is,&n);
152:       ISGetIndices(is,&idx);
153:       ISCreateGeneral(PetscObjectComm((PetscObject)is),n,idx,PETSC_USE_POINTER,&ist);
154:       PetscObjectSetName((PetscObject)ist,name);
155:       ISView(ist,viewer);
156:       ISDestroy(&ist);
157:       ISRestoreIndices(is,&idx);
158:     } else {
159:       PetscBool isperm;

161:       ISGetInfo(is,IS_PERMUTATION,IS_GLOBAL,PETSC_FALSE,&isperm);
162:       if (isperm) {PetscViewerASCIIPrintf(viewer,"Block Index set is permutation\n");}
163:       PetscViewerASCIIPushSynchronized(viewer);
164:       PetscViewerASCIISynchronizedPrintf(viewer,"Block size %D\n",bs);
165:       PetscViewerASCIISynchronizedPrintf(viewer,"Number of block indices in set %D\n",n);
166:       PetscViewerASCIISynchronizedPrintf(viewer,"The first indices of each block are\n");
167:       for (i=0; i<n; i++) {
168:         PetscViewerASCIISynchronizedPrintf(viewer,"Block %D Index %D\n",i,idx[i]);
169:       }
170:       PetscViewerFlush(viewer);
171:       PetscViewerASCIIPopSynchronized(viewer);
172:     }
173:   }
174:   return(0);
175: }

177: static PetscErrorCode ISSort_Block(IS is)
178: {
179:   IS_Block       *sub = (IS_Block*)is->data;
180:   PetscInt       bs, n;

184:   PetscLayoutGetBlockSize(is->map, &bs);
185:   PetscLayoutGetLocalSize(is->map, &n);
186:   PetscSortInt(n/bs,sub->idx);
187:   return(0);
188: }

190: static PetscErrorCode ISSortRemoveDups_Block(IS is)
191: {
192:   IS_Block       *sub = (IS_Block*)is->data;
193:   PetscInt       bs, n, nb;
194:   PetscBool      sorted;

198:   PetscLayoutGetBlockSize(is->map, &bs);
199:   PetscLayoutGetLocalSize(is->map, &n);
200:   nb   = n/bs;
201:   ISGetInfo(is,IS_SORTED,IS_LOCAL,PETSC_TRUE,&sorted);
202:   if (sorted) {
203:     PetscSortedRemoveDupsInt(&nb,sub->idx);
204:   } else {
205:     PetscSortRemoveDupsInt(&nb,sub->idx);
206:   }
207:   PetscLayoutDestroy(&is->map);
208:   PetscLayoutCreateFromSizes(PetscObjectComm((PetscObject)is), nb*bs, PETSC_DECIDE, bs, &is->map);
209:   return(0);
210: }

212: static PetscErrorCode ISSorted_Block(IS is,PetscBool  *flg)
213: {

217:   ISGetInfo(is,IS_SORTED,IS_LOCAL,PETSC_TRUE,flg);
218:   return(0);
219: }

221: static PetscErrorCode ISSortedLocal_Block(IS is,PetscBool *flg)
222: {
223:   IS_Block       *sub = (IS_Block*)is->data;
224:   PetscInt       n, bs, i, *idx;

228:   PetscLayoutGetLocalSize(is->map, &n);
229:   PetscLayoutGetBlockSize(is->map, &bs);
230:   n   /= bs;
231:   idx  = sub->idx;
232:   for (i = 1; i < n; i++) if (idx[i] < idx[i - 1]) break;
233:   if (i < n) *flg = PETSC_FALSE;
234:   else       *flg = PETSC_TRUE;
235:   return(0);
236: }

238: static PetscErrorCode ISUniqueLocal_Block(IS is,PetscBool *flg)
239: {
240:   IS_Block       *sub = (IS_Block*)is->data;
241:   PetscInt       n, bs, i, *idx, *idxcopy = NULL;
242:   PetscBool      sortedLocal;

246:   PetscLayoutGetLocalSize(is->map, &n);
247:   PetscLayoutGetBlockSize(is->map, &bs);
248:   n   /= bs;
249:   idx  = sub->idx;
250:   ISGetInfo(is,IS_SORTED,IS_LOCAL,PETSC_TRUE,&sortedLocal);
251:   if (!sortedLocal) {
252:     PetscMalloc1(n, &idxcopy);
253:     PetscArraycpy(idxcopy, idx, n);
254:     PetscSortInt(n, idxcopy);
255:     idx = idxcopy;
256:   }
257:   for (i = 1; i < n; i++) if (idx[i] == idx[i - 1]) break;
258:   if (i < n) *flg = PETSC_FALSE;
259:   else       *flg = PETSC_TRUE;
260:   PetscFree(idxcopy);
261:   return(0);
262: }

264: static PetscErrorCode ISPermutationLocal_Block(IS is,PetscBool *flg)
265: {
266:   IS_Block       *sub = (IS_Block*)is->data;
267:   PetscInt       n, bs, i, *idx, *idxcopy = NULL;
268:   PetscBool      sortedLocal;

272:   PetscLayoutGetLocalSize(is->map, &n);
273:   PetscLayoutGetBlockSize(is->map, &bs);
274:   n   /= bs;
275:   idx  = sub->idx;
276:   ISGetInfo(is,IS_SORTED,IS_LOCAL,PETSC_TRUE,&sortedLocal);
277:   if (!sortedLocal) {
278:     PetscMalloc1(n, &idxcopy);
279:     PetscArraycpy(idxcopy, idx, n);
280:     PetscSortInt(n, idxcopy);
281:     idx = idxcopy;
282:   }
283:   for (i = 0; i < n; i++) if (idx[i] != i) break;
284:   if (i < n) *flg = PETSC_FALSE;
285:   else       *flg = PETSC_TRUE;
286:   PetscFree(idxcopy);
287:   return(0);
288: }

290: static PetscErrorCode ISIntervalLocal_Block(IS is,PetscBool *flg)
291: {
292:   IS_Block       *sub = (IS_Block*)is->data;
293:   PetscInt       n, bs, i, *idx;

297:   PetscLayoutGetLocalSize(is->map, &n);
298:   PetscLayoutGetBlockSize(is->map, &bs);
299:   n   /= bs;
300:   idx  = sub->idx;
301:   for (i = 1; i < n; i++) if (idx[i] != idx[i - 1] + 1) break;
302:   if (i < n) *flg = PETSC_FALSE;
303:   else       *flg = PETSC_TRUE;
304:   return(0);
305: }

307: static PetscErrorCode ISDuplicate_Block(IS is,IS *newIS)
308: {
310:   IS_Block       *sub = (IS_Block*)is->data;
311:   PetscInt        bs, n;

314:   PetscLayoutGetBlockSize(is->map, &bs);
315:   PetscLayoutGetLocalSize(is->map, &n);
316:   n   /= bs;
317:   ISCreateBlock(PetscObjectComm((PetscObject)is),bs,n,sub->idx,PETSC_COPY_VALUES,newIS);
318:   return(0);
319: }

321: static PetscErrorCode ISCopy_Block(IS is,IS isy)
322: {
323:   IS_Block       *is_block = (IS_Block*)is->data,*isy_block = (IS_Block*)isy->data;
324:   PetscInt       bs, n, N, bsy, ny, Ny;

328:   PetscLayoutGetBlockSize(is->map, &bs);
329:   PetscLayoutGetLocalSize(is->map, &n);
330:   PetscLayoutGetSize(is->map, &N);
331:   PetscLayoutGetBlockSize(isy->map, &bsy);
332:   PetscLayoutGetLocalSize(isy->map, &ny);
333:   PetscLayoutGetSize(isy->map, &Ny);
334:   if (n != ny || N != Ny || bs != bsy) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Index sets incompatible");
335:   PetscArraycpy(isy_block->idx,is_block->idx,(n/bs));
336:   return(0);
337: }

339: static PetscErrorCode ISOnComm_Block(IS is,MPI_Comm comm,PetscCopyMode mode,IS *newis)
340: {
342:   IS_Block       *sub = (IS_Block*)is->data;
343:   PetscInt       bs, n;

346:   if (mode == PETSC_OWN_POINTER) SETERRQ(comm,PETSC_ERR_ARG_WRONG,"Cannot use PETSC_OWN_POINTER");
347:   PetscLayoutGetBlockSize(is->map, &bs);
348:   PetscLayoutGetLocalSize(is->map, &n);
349:   ISCreateBlock(comm,bs,n/bs,sub->idx,mode,newis);
350:   return(0);
351: }

353: static PetscErrorCode ISSetBlockSize_Block(IS is,PetscInt bs)
354: {

358:   PetscLayoutSetBlockSize(is->map, bs);
359:   return(0);
360: }

362: static PetscErrorCode ISToGeneral_Block(IS inis)
363: {
364:   IS_Block       *sub   = (IS_Block*)inis->data;
365:   PetscInt       bs,n;
366:   const PetscInt *idx;

370:   ISGetBlockSize(inis,&bs);
371:   ISGetLocalSize(inis,&n);
372:   ISGetIndices(inis,&idx);
373:   if (bs == 1) {
374:     PetscCopyMode mode = sub->allocated ? PETSC_OWN_POINTER : PETSC_USE_POINTER;
375:     sub->allocated = PETSC_FALSE; /* prevent deallocation when changing the subtype*/
376:     ISSetType(inis,ISGENERAL);
377:     ISGeneralSetIndices(inis,n,idx,mode);
378:   } else {
379:     ISSetType(inis,ISGENERAL);
380:     ISGeneralSetIndices(inis,n,idx,PETSC_OWN_POINTER);
381:   }
382:   return(0);
383: }


386: static struct _ISOps myops = { ISGetIndices_Block,
387:                                ISRestoreIndices_Block,
388:                                ISInvertPermutation_Block,
389:                                ISSort_Block,
390:                                ISSortRemoveDups_Block,
391:                                ISSorted_Block,
392:                                ISDuplicate_Block,
393:                                ISDestroy_Block,
394:                                ISView_Block,
395:                                ISLoad_Default,
396:                                ISCopy_Block,
397:                                ISToGeneral_Block,
398:                                ISOnComm_Block,
399:                                ISSetBlockSize_Block,
400:                                0,
401:                                ISLocate_Block,
402:                                /* we can have specialized local routines for determining properties,
403:                                 * but unless the block size is the same on each process (which is not guaranteed at
404:                                 * the moment), then trying to do something specialized for global properties is too
405:                                 * complicated */
406:                                ISSortedLocal_Block,
407:                                NULL,
408:                                ISUniqueLocal_Block,
409:                                NULL,
410:                                ISPermutationLocal_Block,
411:                                NULL,
412:                                ISIntervalLocal_Block,
413:                                NULL};

415: /*@
416:    ISBlockSetIndices - The indices are relative to entries, not blocks.

418:    Collective on IS

420:    Input Parameters:
421: +  is - the index set
422: .  bs - number of elements in each block, one for each block and count of block not indices
423: .   n - the length of the index set (the number of blocks)
424: .  idx - the list of integers, these are by block, not by location
425: -  mode - see PetscCopyMode, only PETSC_COPY_VALUES and PETSC_OWN_POINTER are supported


428:    Notes:
429:    When the communicator is not MPI_COMM_SELF, the operations on the
430:    index sets, IS, are NOT conceptually the same as MPI_Group operations.
431:    The index sets are then distributed sets of indices and thus certain operations
432:    on them are collective.

434:    Example:
435:    If you wish to index the values {0,1,4,5}, then use
436:    a block size of 2 and idx of {0,2}.

438:    Level: beginner

440: .seealso: ISCreateStride(), ISCreateGeneral(), ISAllGather()
441: @*/
442: PetscErrorCode  ISBlockSetIndices(IS is,PetscInt bs,PetscInt n,const PetscInt idx[],PetscCopyMode mode)
443: {

447:   ISClearInfoCache(is,PETSC_FALSE);
448:   PetscUseMethod(is,"ISBlockSetIndices_C",(IS,PetscInt,PetscInt,const PetscInt[],PetscCopyMode),(is,bs,n,idx,mode));
449:   return(0);
450: }

452: static PetscErrorCode  ISBlockSetIndices_Block(IS is,PetscInt bs,PetscInt n,const PetscInt idx[],PetscCopyMode mode)
453: {
455:   PetscInt       i,min,max;
456:   IS_Block       *sub = (IS_Block*)is->data;
457:   PetscLayout    map;

460:   if (bs < 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"block size < 1");
461:   if (n < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"length < 0");

464:   PetscLayoutCreateFromSizes(PetscObjectComm((PetscObject)is),n*bs,is->map->N,bs,&map);
465:   PetscLayoutDestroy(&is->map);
466:   is->map = map;

468:   if (sub->allocated) {PetscFree(sub->idx);}
469:   if (mode == PETSC_COPY_VALUES) {
470:     PetscMalloc1(n,&sub->idx);
471:     PetscLogObjectMemory((PetscObject)is,n*sizeof(PetscInt));
472:     PetscArraycpy(sub->idx,idx,n);
473:     sub->allocated = PETSC_TRUE;
474:   } else if (mode == PETSC_OWN_POINTER) {
475:     sub->idx = (PetscInt*) idx;
476:     PetscLogObjectMemory((PetscObject)is,n*sizeof(PetscInt));
477:     sub->allocated = PETSC_TRUE;
478:   } else if (mode == PETSC_USE_POINTER) {
479:     sub->idx = (PetscInt*) idx;
480:     sub->allocated = PETSC_FALSE;
481:   }

483:   if (n) {
484:     min = max = idx[0];
485:     for (i=1; i<n; i++) {
486:       if (idx[i] < min) min = idx[i];
487:       if (idx[i] > max) max = idx[i];
488:     }
489:     is->min = bs*min;
490:     is->max = bs*max+bs-1;
491:   } else {
492:     is->min = PETSC_MAX_INT;
493:     is->max = PETSC_MIN_INT;
494:   }
495:   return(0);
496: }

498: /*@
499:    ISCreateBlock - Creates a data structure for an index set containing
500:    a list of integers. The indices are relative to entries, not blocks.

502:    Collective

504:    Input Parameters:
505: +  comm - the MPI communicator
506: .  bs - number of elements in each block
507: .  n - the length of the index set (the number of blocks)
508: .  idx - the list of integers, one for each block and count of block not indices
509: -  mode - see PetscCopyMode, only PETSC_COPY_VALUES and PETSC_OWN_POINTER are supported in this routine

511:    Output Parameter:
512: .  is - the new index set

514:    Notes:
515:    When the communicator is not MPI_COMM_SELF, the operations on the
516:    index sets, IS, are NOT conceptually the same as MPI_Group operations.
517:    The index sets are then distributed sets of indices and thus certain operations
518:    on them are collective.

520:    Example:
521:    If you wish to index the values {0,1,6,7}, then use
522:    a block size of 2 and idx of {0,3}.

524:    Level: beginner

526: .seealso: ISCreateStride(), ISCreateGeneral(), ISAllGather()
527: @*/
528: PetscErrorCode  ISCreateBlock(MPI_Comm comm,PetscInt bs,PetscInt n,const PetscInt idx[],PetscCopyMode mode,IS *is)
529: {

534:   if (bs < 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"block size < 1");
535:   if (n < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"length < 0");

538:   ISCreate(comm,is);
539:   ISSetType(*is,ISBLOCK);
540:   ISBlockSetIndices(*is,bs,n,idx,mode);
541:   return(0);
542: }

544: static PetscErrorCode  ISBlockGetIndices_Block(IS is,const PetscInt *idx[])
545: {
546:   IS_Block *sub = (IS_Block*)is->data;

549:   *idx = sub->idx;
550:   return(0);
551: }

553: static PetscErrorCode  ISBlockRestoreIndices_Block(IS is,const PetscInt *idx[])
554: {
556:   return(0);
557: }

559: /*@C
560:    ISBlockGetIndices - Gets the indices associated with each block.

562:    Not Collective

564:    Input Parameter:
565: .  is - the index set

567:    Output Parameter:
568: .  idx - the integer indices, one for each block and count of block not indices

570:    Level: intermediate

572: .seealso: ISGetIndices(), ISBlockRestoreIndices()
573: @*/
574: PetscErrorCode  ISBlockGetIndices(IS is,const PetscInt *idx[])
575: {

579:   PetscUseMethod(is,"ISBlockGetIndices_C",(IS,const PetscInt*[]),(is,idx));
580:   return(0);
581: }

583: /*@C
584:    ISBlockRestoreIndices - Restores the indices associated with each block.

586:    Not Collective

588:    Input Parameter:
589: .  is - the index set

591:    Output Parameter:
592: .  idx - the integer indices

594:    Level: intermediate

596: .seealso: ISRestoreIndices(), ISBlockGetIndices()
597: @*/
598: PetscErrorCode  ISBlockRestoreIndices(IS is,const PetscInt *idx[])
599: {

603:   PetscUseMethod(is,"ISBlockRestoreIndices_C",(IS,const PetscInt*[]),(is,idx));
604:   return(0);
605: }

607: /*@
608:    ISBlockGetLocalSize - Returns the local number of blocks in the index set.

610:    Not Collective

612:    Input Parameter:
613: .  is - the index set

615:    Output Parameter:
616: .  size - the local number of blocks

618:    Level: intermediate


621: .seealso: ISGetBlockSize(), ISBlockGetSize(), ISGetSize(), ISCreateBlock()
622: @*/
623: PetscErrorCode  ISBlockGetLocalSize(IS is,PetscInt *size)
624: {

628:   PetscUseMethod(is,"ISBlockGetLocalSize_C",(IS,PetscInt*),(is,size));
629:   return(0);
630: }

632: static PetscErrorCode  ISBlockGetLocalSize_Block(IS is,PetscInt *size)
633: {
634:   PetscInt       bs, n;

638:   PetscLayoutGetBlockSize(is->map, &bs);
639:   PetscLayoutGetLocalSize(is->map, &n);
640:   *size = n/bs;
641:   return(0);
642: }

644: /*@
645:    ISBlockGetSize - Returns the global number of blocks in the index set.

647:    Not Collective

649:    Input Parameter:
650: .  is - the index set

652:    Output Parameter:
653: .  size - the global number of blocks

655:    Level: intermediate


658: .seealso: ISGetBlockSize(), ISBlockGetLocalSize(), ISGetSize(), ISCreateBlock()
659: @*/
660: PetscErrorCode  ISBlockGetSize(IS is,PetscInt *size)
661: {

665:   PetscUseMethod(is,"ISBlockGetSize_C",(IS,PetscInt*),(is,size));
666:   return(0);
667: }

669: static PetscErrorCode  ISBlockGetSize_Block(IS is,PetscInt *size)
670: {
671:   PetscInt       bs, N;

675:   PetscLayoutGetBlockSize(is->map, &bs);
676:   PetscLayoutGetSize(is->map, &N);
677:   *size = N/bs;
678:   return(0);
679: }

681: PETSC_EXTERN PetscErrorCode ISCreate_Block(IS is)
682: {
684:   IS_Block       *sub;

687:   PetscNewLog(is,&sub);
688:   is->data = (void *) sub;
689:   PetscMemcpy(is->ops,&myops,sizeof(myops));
690:   PetscObjectComposeFunction((PetscObject)is,"ISBlockSetIndices_C",ISBlockSetIndices_Block);
691:   PetscObjectComposeFunction((PetscObject)is,"ISBlockGetIndices_C",ISBlockGetIndices_Block);
692:   PetscObjectComposeFunction((PetscObject)is,"ISBlockRestoreIndices_C",ISBlockRestoreIndices_Block);
693:   PetscObjectComposeFunction((PetscObject)is,"ISBlockGetSize_C",ISBlockGetSize_Block);
694:   PetscObjectComposeFunction((PetscObject)is,"ISBlockGetLocalSize_C",ISBlockGetLocalSize_Block);
695:   return(0);
696: }