Actual source code: block.c
petsc-3.14.6 2021-03-30
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: if (n) {
81: PetscMalloc1(bs*n,&jj);
82: *idx = jj;
83: k = 0;
84: ii = sub->idx;
85: for (i=0; i<n; i++)
86: for (j=0; j<bs; j++)
87: jj[k++] = bs*ii[i] + j;
88: } else {
89: /* do not malloc for zero size because F90Array1dCreate() inside ISRestoreArrayF90() does not keep array when zero length array */
90: *idx = NULL;
91: }
92: }
93: return(0);
94: }
96: static PetscErrorCode ISRestoreIndices_Block(IS is,const PetscInt *idx[])
97: {
98: IS_Block *sub = (IS_Block*)is->data;
99: PetscInt bs;
103: PetscLayoutGetBlockSize(is->map, &bs);
104: if (bs != 1) {
105: PetscFree(*(void**)idx);
106: } else {
107: /* F90Array1dCreate() inside ISRestoreArrayF90() does not keep array when zero length array */
108: if (is->map->n > 0 && *idx != sub->idx) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Must restore with value from ISGetIndices()");
109: }
110: return(0);
111: }
113: static PetscErrorCode ISInvertPermutation_Block(IS is,PetscInt nlocal,IS *isout)
114: {
115: IS_Block *sub = (IS_Block*)is->data;
116: PetscInt i,*ii,bs,n,*idx = sub->idx;
117: PetscMPIInt size;
121: MPI_Comm_size(PetscObjectComm((PetscObject)is),&size);
122: PetscLayoutGetBlockSize(is->map, &bs);
123: PetscLayoutGetLocalSize(is->map, &n);
124: n /= bs;
125: if (size == 1) {
126: PetscMalloc1(n,&ii);
127: for (i=0; i<n; i++) ii[idx[i]] = i;
128: ISCreateBlock(PETSC_COMM_SELF,bs,n,ii,PETSC_OWN_POINTER,isout);
129: ISSetPermutation(*isout);
130: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No inversion written yet for block IS");
131: return(0);
132: }
134: static PetscErrorCode ISView_Block(IS is, PetscViewer viewer)
135: {
136: IS_Block *sub = (IS_Block*)is->data;
138: PetscInt i,bs,n,*idx = sub->idx;
139: PetscBool iascii,ibinary;
142: PetscLayoutGetBlockSize(is->map, &bs);
143: PetscLayoutGetLocalSize(is->map, &n);
144: n /= bs;
145: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
146: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&ibinary);
147: if (iascii) {
148: PetscViewerFormat fmt;
150: PetscViewerGetFormat(viewer,&fmt);
151: if (fmt == PETSC_VIEWER_ASCII_MATLAB) {
152: IS ist;
153: const char *name;
154: const PetscInt *idx;
155: PetscInt n;
157: PetscObjectGetName((PetscObject)is,&name);
158: ISGetLocalSize(is,&n);
159: ISGetIndices(is,&idx);
160: ISCreateGeneral(PetscObjectComm((PetscObject)is),n,idx,PETSC_USE_POINTER,&ist);
161: PetscObjectSetName((PetscObject)ist,name);
162: ISView(ist,viewer);
163: ISDestroy(&ist);
164: ISRestoreIndices(is,&idx);
165: } else {
166: PetscBool isperm;
168: ISGetInfo(is,IS_PERMUTATION,IS_GLOBAL,PETSC_FALSE,&isperm);
169: if (isperm) {PetscViewerASCIIPrintf(viewer,"Block Index set is permutation\n");}
170: PetscViewerASCIIPushSynchronized(viewer);
171: PetscViewerASCIISynchronizedPrintf(viewer,"Block size %D\n",bs);
172: PetscViewerASCIISynchronizedPrintf(viewer,"Number of block indices in set %D\n",n);
173: PetscViewerASCIISynchronizedPrintf(viewer,"The first indices of each block are\n");
174: for (i=0; i<n; i++) {
175: PetscViewerASCIISynchronizedPrintf(viewer,"Block %D Index %D\n",i,idx[i]);
176: }
177: PetscViewerFlush(viewer);
178: PetscViewerASCIIPopSynchronized(viewer);
179: }
180: } else if (ibinary) {
181: ISView_Binary(is,viewer);
182: }
183: return(0);
184: }
186: static PetscErrorCode ISSort_Block(IS is)
187: {
188: IS_Block *sub = (IS_Block*)is->data;
189: PetscInt bs, n;
193: PetscLayoutGetBlockSize(is->map, &bs);
194: PetscLayoutGetLocalSize(is->map, &n);
195: PetscIntSortSemiOrdered(n/bs,sub->idx);
196: return(0);
197: }
199: static PetscErrorCode ISSortRemoveDups_Block(IS is)
200: {
201: IS_Block *sub = (IS_Block*)is->data;
202: PetscInt bs, n, nb;
203: PetscBool sorted;
207: PetscLayoutGetBlockSize(is->map, &bs);
208: PetscLayoutGetLocalSize(is->map, &n);
209: nb = n/bs;
210: ISGetInfo(is,IS_SORTED,IS_LOCAL,PETSC_TRUE,&sorted);
211: if (sorted) {
212: PetscSortedRemoveDupsInt(&nb,sub->idx);
213: } else {
214: PetscSortRemoveDupsInt(&nb,sub->idx);
215: }
216: PetscLayoutDestroy(&is->map);
217: PetscLayoutCreateFromSizes(PetscObjectComm((PetscObject)is), nb*bs, PETSC_DECIDE, bs, &is->map);
218: return(0);
219: }
221: static PetscErrorCode ISSorted_Block(IS is,PetscBool *flg)
222: {
226: ISGetInfo(is,IS_SORTED,IS_LOCAL,PETSC_TRUE,flg);
227: return(0);
228: }
230: static PetscErrorCode ISSortedLocal_Block(IS is,PetscBool *flg)
231: {
232: IS_Block *sub = (IS_Block*)is->data;
233: PetscInt n, bs, i, *idx;
237: PetscLayoutGetLocalSize(is->map, &n);
238: PetscLayoutGetBlockSize(is->map, &bs);
239: n /= bs;
240: idx = sub->idx;
241: for (i = 1; i < n; i++) if (idx[i] < idx[i - 1]) break;
242: if (i < n) *flg = PETSC_FALSE;
243: else *flg = PETSC_TRUE;
244: return(0);
245: }
247: static PetscErrorCode ISUniqueLocal_Block(IS is,PetscBool *flg)
248: {
249: IS_Block *sub = (IS_Block*)is->data;
250: PetscInt n, bs, i, *idx, *idxcopy = NULL;
251: PetscBool sortedLocal;
255: PetscLayoutGetLocalSize(is->map, &n);
256: PetscLayoutGetBlockSize(is->map, &bs);
257: n /= bs;
258: idx = sub->idx;
259: ISGetInfo(is,IS_SORTED,IS_LOCAL,PETSC_TRUE,&sortedLocal);
260: if (!sortedLocal) {
261: PetscMalloc1(n, &idxcopy);
262: PetscArraycpy(idxcopy, idx, n);
263: PetscIntSortSemiOrdered(n, idxcopy);
264: idx = idxcopy;
265: }
266: for (i = 1; i < n; i++) if (idx[i] == idx[i - 1]) break;
267: if (i < n) *flg = PETSC_FALSE;
268: else *flg = PETSC_TRUE;
269: PetscFree(idxcopy);
270: return(0);
271: }
273: static PetscErrorCode ISPermutationLocal_Block(IS is,PetscBool *flg)
274: {
275: IS_Block *sub = (IS_Block*)is->data;
276: PetscInt n, bs, i, *idx, *idxcopy = NULL;
277: PetscBool sortedLocal;
281: PetscLayoutGetLocalSize(is->map, &n);
282: PetscLayoutGetBlockSize(is->map, &bs);
283: n /= bs;
284: idx = sub->idx;
285: ISGetInfo(is,IS_SORTED,IS_LOCAL,PETSC_TRUE,&sortedLocal);
286: if (!sortedLocal) {
287: PetscMalloc1(n, &idxcopy);
288: PetscArraycpy(idxcopy, idx, n);
289: PetscIntSortSemiOrdered(n, idxcopy);
290: idx = idxcopy;
291: }
292: for (i = 0; i < n; i++) if (idx[i] != i) break;
293: if (i < n) *flg = PETSC_FALSE;
294: else *flg = PETSC_TRUE;
295: PetscFree(idxcopy);
296: return(0);
297: }
299: static PetscErrorCode ISIntervalLocal_Block(IS is,PetscBool *flg)
300: {
301: IS_Block *sub = (IS_Block*)is->data;
302: PetscInt n, bs, i, *idx;
306: PetscLayoutGetLocalSize(is->map, &n);
307: PetscLayoutGetBlockSize(is->map, &bs);
308: n /= bs;
309: idx = sub->idx;
310: for (i = 1; i < n; i++) if (idx[i] != idx[i - 1] + 1) break;
311: if (i < n) *flg = PETSC_FALSE;
312: else *flg = PETSC_TRUE;
313: return(0);
314: }
316: static PetscErrorCode ISDuplicate_Block(IS is,IS *newIS)
317: {
319: IS_Block *sub = (IS_Block*)is->data;
320: PetscInt bs, n;
323: PetscLayoutGetBlockSize(is->map, &bs);
324: PetscLayoutGetLocalSize(is->map, &n);
325: n /= bs;
326: ISCreateBlock(PetscObjectComm((PetscObject)is),bs,n,sub->idx,PETSC_COPY_VALUES,newIS);
327: return(0);
328: }
330: static PetscErrorCode ISCopy_Block(IS is,IS isy)
331: {
332: IS_Block *is_block = (IS_Block*)is->data,*isy_block = (IS_Block*)isy->data;
333: PetscInt bs, n, N, bsy, ny, Ny;
337: PetscLayoutGetBlockSize(is->map, &bs);
338: PetscLayoutGetLocalSize(is->map, &n);
339: PetscLayoutGetSize(is->map, &N);
340: PetscLayoutGetBlockSize(isy->map, &bsy);
341: PetscLayoutGetLocalSize(isy->map, &ny);
342: PetscLayoutGetSize(isy->map, &Ny);
343: if (n != ny || N != Ny || bs != bsy) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Index sets incompatible");
344: PetscArraycpy(isy_block->idx,is_block->idx,(n/bs));
345: return(0);
346: }
348: static PetscErrorCode ISOnComm_Block(IS is,MPI_Comm comm,PetscCopyMode mode,IS *newis)
349: {
351: IS_Block *sub = (IS_Block*)is->data;
352: PetscInt bs, n;
355: if (mode == PETSC_OWN_POINTER) SETERRQ(comm,PETSC_ERR_ARG_WRONG,"Cannot use PETSC_OWN_POINTER");
356: PetscLayoutGetBlockSize(is->map, &bs);
357: PetscLayoutGetLocalSize(is->map, &n);
358: ISCreateBlock(comm,bs,n/bs,sub->idx,mode,newis);
359: return(0);
360: }
362: static PetscErrorCode ISSetBlockSize_Block(IS is,PetscInt bs)
363: {
367: if (is->map->bs > 0 && bs != is->map->bs) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Cannot change blocksize %D (to %D) if ISType is ISBLOCK",is->map->bs,bs);
368: PetscLayoutSetBlockSize(is->map, bs);
369: return(0);
370: }
372: static PetscErrorCode ISToGeneral_Block(IS inis)
373: {
374: IS_Block *sub = (IS_Block*)inis->data;
375: PetscInt bs,n;
376: const PetscInt *idx;
380: ISGetBlockSize(inis,&bs);
381: ISGetLocalSize(inis,&n);
382: ISGetIndices(inis,&idx);
383: if (bs == 1) {
384: PetscCopyMode mode = sub->allocated ? PETSC_OWN_POINTER : PETSC_USE_POINTER;
385: sub->allocated = PETSC_FALSE; /* prevent deallocation when changing the subtype*/
386: ISSetType(inis,ISGENERAL);
387: ISGeneralSetIndices(inis,n,idx,mode);
388: } else {
389: ISSetType(inis,ISGENERAL);
390: ISGeneralSetIndices(inis,n,idx,PETSC_OWN_POINTER);
391: }
392: return(0);
393: }
396: static struct _ISOps myops = { ISGetIndices_Block,
397: ISRestoreIndices_Block,
398: ISInvertPermutation_Block,
399: ISSort_Block,
400: ISSortRemoveDups_Block,
401: ISSorted_Block,
402: ISDuplicate_Block,
403: ISDestroy_Block,
404: ISView_Block,
405: ISLoad_Default,
406: ISCopy_Block,
407: ISToGeneral_Block,
408: ISOnComm_Block,
409: ISSetBlockSize_Block,
410: NULL,
411: ISLocate_Block,
412: /* we can have specialized local routines for determining properties,
413: * but unless the block size is the same on each process (which is not guaranteed at
414: * the moment), then trying to do something specialized for global properties is too
415: * complicated */
416: ISSortedLocal_Block,
417: NULL,
418: ISUniqueLocal_Block,
419: NULL,
420: ISPermutationLocal_Block,
421: NULL,
422: ISIntervalLocal_Block,
423: NULL};
425: /*@
426: ISBlockSetIndices - Set integers representing blocks of indices in an index set.
428: Collective on IS
430: Input Parameters:
431: + is - the index set
432: . bs - number of elements in each block
433: . n - the length of the index set (the number of blocks)
434: . idx - the list of integers, one for each block, the integers contain the index of the first index of each block divided by the block size
435: - mode - see PetscCopyMode, only PETSC_COPY_VALUES and PETSC_OWN_POINTER are supported
438: Notes:
439: When the communicator is not MPI_COMM_SELF, the operations on the
440: index sets, IS, are NOT conceptually the same as MPI_Group operations.
441: The index sets are then distributed sets of indices and thus certain operations
442: on them are collective.
444: Example:
445: If you wish to index the values {0,1,4,5}, then use
446: a block size of 2 and idx of {0,2}.
448: Level: beginner
450: .seealso: ISCreateStride(), ISCreateGeneral(), ISAllGather(), ISCreateBlock(), ISBLOCK, ISGeneralSetIndices()
451: @*/
452: PetscErrorCode ISBlockSetIndices(IS is,PetscInt bs,PetscInt n,const PetscInt idx[],PetscCopyMode mode)
453: {
457: ISClearInfoCache(is,PETSC_FALSE);
458: PetscUseMethod(is,"ISBlockSetIndices_C",(IS,PetscInt,PetscInt,const PetscInt[],PetscCopyMode),(is,bs,n,idx,mode));
459: return(0);
460: }
462: static PetscErrorCode ISBlockSetIndices_Block(IS is,PetscInt bs,PetscInt n,const PetscInt idx[],PetscCopyMode mode)
463: {
465: PetscInt i,min,max;
466: IS_Block *sub = (IS_Block*)is->data;
467: PetscLayout map;
470: if (bs < 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"block size < 1");
471: if (n < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"length < 0");
474: PetscLayoutCreateFromSizes(PetscObjectComm((PetscObject)is),n*bs,is->map->N,bs,&map);
475: PetscLayoutDestroy(&is->map);
476: is->map = map;
478: if (sub->allocated) {PetscFree(sub->idx);}
479: if (mode == PETSC_COPY_VALUES) {
480: PetscMalloc1(n,&sub->idx);
481: PetscLogObjectMemory((PetscObject)is,n*sizeof(PetscInt));
482: PetscArraycpy(sub->idx,idx,n);
483: sub->allocated = PETSC_TRUE;
484: } else if (mode == PETSC_OWN_POINTER) {
485: sub->idx = (PetscInt*) idx;
486: PetscLogObjectMemory((PetscObject)is,n*sizeof(PetscInt));
487: sub->allocated = PETSC_TRUE;
488: } else if (mode == PETSC_USE_POINTER) {
489: sub->idx = (PetscInt*) idx;
490: sub->allocated = PETSC_FALSE;
491: }
493: if (n) {
494: min = max = idx[0];
495: for (i=1; i<n; i++) {
496: if (idx[i] < min) min = idx[i];
497: if (idx[i] > max) max = idx[i];
498: }
499: is->min = bs*min;
500: is->max = bs*max+bs-1;
501: } else {
502: is->min = PETSC_MAX_INT;
503: is->max = PETSC_MIN_INT;
504: }
505: return(0);
506: }
508: /*@
509: ISCreateBlock - Creates a data structure for an index set containing
510: a list of integers. Each integer represents a fixed block size set of indices.
512: Collective
514: Input Parameters:
515: + comm - the MPI communicator
516: . bs - number of elements in each block
517: . n - the length of the index set (the number of blocks)
518: . idx - the list of integers, one for each block, the integers contain the index of the first entry of each block divided by the block size
519: - mode - see PetscCopyMode, only PETSC_COPY_VALUES and PETSC_OWN_POINTER are supported in this routine
521: Output Parameter:
522: . is - the new index set
524: Notes:
525: When the communicator is not MPI_COMM_SELF, the operations on the
526: index sets, IS, are NOT conceptually the same as MPI_Group operations.
527: The index sets are then distributed sets of indices and thus certain operations
528: on them are collective.
530: Example:
531: If you wish to index the values {0,1,6,7}, then use
532: a block size of 2 and idx of {0,3}.
534: Level: beginner
536: .seealso: ISCreateStride(), ISCreateGeneral(), ISAllGather(), ISBlockSetIndices(), ISBLOCK, ISGENERAL
537: @*/
538: PetscErrorCode ISCreateBlock(MPI_Comm comm,PetscInt bs,PetscInt n,const PetscInt idx[],PetscCopyMode mode,IS *is)
539: {
544: if (bs < 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"block size < 1");
545: if (n < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"length < 0");
548: ISCreate(comm,is);
549: ISSetType(*is,ISBLOCK);
550: ISBlockSetIndices(*is,bs,n,idx,mode);
551: return(0);
552: }
554: static PetscErrorCode ISBlockGetIndices_Block(IS is,const PetscInt *idx[])
555: {
556: IS_Block *sub = (IS_Block*)is->data;
559: *idx = sub->idx;
560: return(0);
561: }
563: static PetscErrorCode ISBlockRestoreIndices_Block(IS is,const PetscInt *idx[])
564: {
566: return(0);
567: }
569: /*@C
570: ISBlockGetIndices - Gets the indices associated with each block.
572: Not Collective
574: Input Parameter:
575: . is - the index set
577: Output Parameter:
578: . idx - the integer indices, one for each block and count of block not indices
580: Level: intermediate
582: .seealso: ISGetIndices(), ISBlockRestoreIndices(), ISBLOCK, ISBlockSetIndices(), ISCreateBlock()
583: @*/
584: PetscErrorCode ISBlockGetIndices(IS is,const PetscInt *idx[])
585: {
589: PetscUseMethod(is,"ISBlockGetIndices_C",(IS,const PetscInt*[]),(is,idx));
590: return(0);
591: }
593: /*@C
594: ISBlockRestoreIndices - Restores the indices associated with each block.
596: Not Collective
598: Input Parameter:
599: . is - the index set
601: Output Parameter:
602: . idx - the integer indices
604: Level: intermediate
606: .seealso: ISRestoreIndices(), ISBlockGetIndices()
607: @*/
608: PetscErrorCode ISBlockRestoreIndices(IS is,const PetscInt *idx[])
609: {
613: PetscUseMethod(is,"ISBlockRestoreIndices_C",(IS,const PetscInt*[]),(is,idx));
614: return(0);
615: }
617: /*@
618: ISBlockGetLocalSize - Returns the local number of blocks in the index set.
620: Not Collective
622: Input Parameter:
623: . is - the index set
625: Output Parameter:
626: . size - the local number of blocks
628: Level: intermediate
630: .seealso: ISGetBlockSize(), ISBlockGetSize(), ISGetSize(), ISCreateBlock(), ISBLOCK
631: @*/
632: PetscErrorCode ISBlockGetLocalSize(IS is,PetscInt *size)
633: {
637: PetscUseMethod(is,"ISBlockGetLocalSize_C",(IS,PetscInt*),(is,size));
638: return(0);
639: }
641: static PetscErrorCode ISBlockGetLocalSize_Block(IS is,PetscInt *size)
642: {
643: PetscInt bs, n;
647: PetscLayoutGetBlockSize(is->map, &bs);
648: PetscLayoutGetLocalSize(is->map, &n);
649: *size = n/bs;
650: return(0);
651: }
653: /*@
654: ISBlockGetSize - Returns the global number of blocks in the index set.
656: Not Collective
658: Input Parameter:
659: . is - the index set
661: Output Parameter:
662: . size - the global number of blocks
664: Level: intermediate
666: .seealso: ISGetBlockSize(), ISBlockGetLocalSize(), ISGetSize(), ISCreateBlock(), ISBLOCK
667: @*/
668: PetscErrorCode ISBlockGetSize(IS is,PetscInt *size)
669: {
673: PetscUseMethod(is,"ISBlockGetSize_C",(IS,PetscInt*),(is,size));
674: return(0);
675: }
677: static PetscErrorCode ISBlockGetSize_Block(IS is,PetscInt *size)
678: {
679: PetscInt bs, N;
683: PetscLayoutGetBlockSize(is->map, &bs);
684: PetscLayoutGetSize(is->map, &N);
685: *size = N/bs;
686: return(0);
687: }
689: PETSC_EXTERN PetscErrorCode ISCreate_Block(IS is)
690: {
692: IS_Block *sub;
695: PetscNewLog(is,&sub);
696: is->data = (void *) sub;
697: PetscMemcpy(is->ops,&myops,sizeof(myops));
698: PetscObjectComposeFunction((PetscObject)is,"ISBlockSetIndices_C",ISBlockSetIndices_Block);
699: PetscObjectComposeFunction((PetscObject)is,"ISBlockGetIndices_C",ISBlockGetIndices_Block);
700: PetscObjectComposeFunction((PetscObject)is,"ISBlockRestoreIndices_C",ISBlockRestoreIndices_Block);
701: PetscObjectComposeFunction((PetscObject)is,"ISBlockGetSize_C",ISBlockGetSize_Block);
702: PetscObjectComposeFunction((PetscObject)is,"ISBlockGetLocalSize_C",ISBlockGetLocalSize_Block);
703: return(0);
704: }