Actual source code: general.c
petsc-3.5.4 2015-05-23
2: /*
3: Provides the functions for index sets (IS) defined by a list of integers.
4: */
5: #include <../src/vec/is/is/impls/general/general.h> /*I "petscis.h" I*/
6: #include <petscvec.h>
7: #include <petscviewer.h>
8: #include <petscviewerhdf5.h>
12: PetscErrorCode ISDuplicate_General(IS is,IS *newIS)
13: {
15: IS_General *sub = (IS_General*)is->data;
16: PetscInt n;
19: PetscLayoutGetLocalSize(is->map, &n);
20: ISCreateGeneral(PetscObjectComm((PetscObject) is), n, sub->idx, PETSC_COPY_VALUES, newIS);
21: return(0);
22: }
26: PetscErrorCode ISDestroy_General(IS is)
27: {
28: IS_General *is_general = (IS_General*)is->data;
32: if (is_general->allocated) {PetscFree(is_general->idx);}
33: PetscObjectComposeFunction((PetscObject)is,"ISGeneralSetIndices_C",0);
34: PetscFree(is->data);
35: return(0);
36: }
40: PetscErrorCode ISIdentity_General(IS is, PetscBool *ident)
41: {
42: IS_General *is_general = (IS_General*)is->data;
43: PetscInt i,n,*idx = is_general->idx;
47: PetscLayoutGetLocalSize(is->map, &n);
48: is->isidentity = PETSC_TRUE;
49: *ident = PETSC_TRUE;
50: for (i=0; i<n; i++) {
51: if (idx[i] != i) {
52: is->isidentity = PETSC_FALSE;
53: *ident = PETSC_FALSE;
54: break;
55: }
56: }
57: return(0);
58: }
62: static PetscErrorCode ISCopy_General(IS is,IS isy)
63: {
64: IS_General *is_general = (IS_General*)is->data,*isy_general = (IS_General*)isy->data;
65: PetscInt n, N, ny, Ny;
69: PetscLayoutGetLocalSize(is->map, &n);
70: PetscLayoutGetSize(is->map, &N);
71: PetscLayoutGetLocalSize(isy->map, &ny);
72: PetscLayoutGetSize(isy->map, &Ny);
73: if (n != ny || N != Ny) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Index sets incompatible");
74: isy_general->sorted = is_general->sorted;
75: PetscMemcpy(isy_general->idx,is_general->idx,n*sizeof(PetscInt));
76: return(0);
77: }
81: PetscErrorCode ISOnComm_General(IS is,MPI_Comm comm,PetscCopyMode mode,IS *newis)
82: {
84: IS_General *sub = (IS_General*)is->data;
85: PetscInt n;
88: if (mode == PETSC_OWN_POINTER) SETERRQ(comm,PETSC_ERR_ARG_WRONG,"Cannot use PETSC_OWN_POINTER");
89: PetscLayoutGetLocalSize(is->map, &n);
90: ISCreateGeneral(comm,n,sub->idx,mode,newis);
91: return(0);
92: }
96: static PetscErrorCode ISSetBlockSize_General(IS is,PetscInt bs)
97: {
98: #if defined(PETSC_USE_DEBUG)
99: IS_General *sub = (IS_General*)is->data;
100: PetscInt n;
101: #endif
105: PetscLayoutSetBlockSize(is->map, bs);
106: #if defined(PETSC_USE_DEBUG)
107: PetscLayoutGetLocalSize(is->map, &n);
108: {
109: PetscInt i,j;
110: for (i=0; i<n; i+=bs) {
111: for (j=0; j<bs; j++) {
112: if (sub->idx[i+j] != sub->idx[i]+j) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Index set does not have block structure, cannot set block size to %D",bs);
113: }
114: }
115: }
116: #endif
117: return(0);
118: }
122: static PetscErrorCode ISContiguousLocal_General(IS is,PetscInt gstart,PetscInt gend,PetscInt *start,PetscBool *contig)
123: {
124: IS_General *sub = (IS_General*)is->data;
125: PetscInt n,i,p;
129: *start = 0;
130: *contig = PETSC_TRUE;
131: PetscLayoutGetLocalSize(is->map, &n);
132: if (!n) return(0);
133: p = sub->idx[0];
134: if (p < gstart) goto nomatch;
135: *start = p - gstart;
136: if (n > gend-p) goto nomatch;
137: for (i=1; i<n; i++,p++) {
138: if (sub->idx[i] != p+1) goto nomatch;
139: }
140: return(0);
141: nomatch:
142: *start = -1;
143: *contig = PETSC_FALSE;
144: return(0);
145: }
149: PetscErrorCode ISGetIndices_General(IS in,const PetscInt *idx[])
150: {
151: IS_General *sub = (IS_General*)in->data;
154: *idx = sub->idx;
155: return(0);
156: }
160: PetscErrorCode ISRestoreIndices_General(IS in,const PetscInt *idx[])
161: {
162: IS_General *sub = (IS_General*)in->data;
165: if (*idx != sub->idx) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Must restore with value from ISGetIndices()");
166: return(0);
167: }
171: PetscErrorCode ISGetSize_General(IS is,PetscInt *size)
172: {
176: PetscLayoutGetSize(is->map, size);
177: return(0);
178: }
182: PetscErrorCode ISGetLocalSize_General(IS is,PetscInt *size)
183: {
187: PetscLayoutGetLocalSize(is->map, size);
188: return(0);
189: }
193: PetscErrorCode ISInvertPermutation_General(IS is,PetscInt nlocal,IS *isout)
194: {
195: IS_General *sub = (IS_General*)is->data;
196: PetscInt i,*ii,n,nstart;
197: const PetscInt *idx = sub->idx;
198: PetscMPIInt size;
199: IS istmp,nistmp;
203: PetscLayoutGetLocalSize(is->map, &n);
204: MPI_Comm_size(PetscObjectComm((PetscObject)is),&size);
205: if (size == 1) {
206: PetscMalloc1(n,&ii);
207: for (i=0; i<n; i++) ii[idx[i]] = i;
208: ISCreateGeneral(PETSC_COMM_SELF,n,ii,PETSC_OWN_POINTER,isout);
209: ISSetPermutation(*isout);
210: } else {
211: /* crude, nonscalable get entire IS on each processor */
212: if (nlocal == PETSC_DECIDE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Do not yet support nlocal of PETSC_DECIDE");
213: ISAllGather(is,&istmp);
214: ISSetPermutation(istmp);
215: ISInvertPermutation(istmp,PETSC_DECIDE,&nistmp);
216: ISDestroy(&istmp);
217: /* get the part we need */
218: MPI_Scan(&nlocal,&nstart,1,MPIU_INT,MPI_SUM,PetscObjectComm((PetscObject)is));
219: #if defined(PETSC_USE_DEBUG)
220: {
221: PetscInt N;
222: PetscMPIInt rank;
223: MPI_Comm_rank(PetscObjectComm((PetscObject)is),&rank);
224: PetscLayoutGetSize(is->map, &N);
225: if (rank == size-1) {
226: if (nstart != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Sum of nlocal lengths %d != total IS length %d",nstart,N);
227: }
228: }
229: #endif
230: nstart -= nlocal;
231: ISGetIndices(nistmp,&idx);
232: ISCreateGeneral(PetscObjectComm((PetscObject)is),nlocal,idx+nstart,PETSC_COPY_VALUES,isout);
233: ISRestoreIndices(nistmp,&idx);
234: ISDestroy(&nistmp);
235: }
236: return(0);
237: }
239: #if defined(PETSC_HAVE_HDF5)
242: PetscErrorCode ISView_General_HDF5(IS is, PetscViewer viewer)
243: {
244: hid_t filespace; /* file dataspace identifier */
245: hid_t chunkspace; /* chunk dataset property identifier */
246: hid_t plist_id; /* property list identifier */
247: hid_t dset_id; /* dataset identifier */
248: hid_t memspace; /* memory dataspace identifier */
249: hid_t inttype; /* int type (H5T_NATIVE_INT or H5T_NATIVE_LLONG) */
250: hid_t file_id, group;
251: herr_t status;
252: hsize_t dim, maxDims[3], dims[3], chunkDims[3], count[3],offset[3];
253: PetscInt bs, N, n, timestep, low;
254: const PetscInt *ind;
255: const char *isname;
256: PetscErrorCode ierr;
259: ISGetBlockSize(is,&bs);
260: PetscViewerHDF5OpenGroup(viewer, &file_id, &group);
261: PetscViewerHDF5GetTimestep(viewer, ×tep);
263: /* Create the dataspace for the dataset.
264: *
265: * dims - holds the current dimensions of the dataset
266: *
267: * maxDims - holds the maximum dimensions of the dataset (unlimited
268: * for the number of time steps with the current dimensions for the
269: * other dimensions; so only additional time steps can be added).
270: *
271: * chunkDims - holds the size of a single time step (required to
272: * permit extending dataset).
273: */
274: dim = 0;
275: if (timestep >= 0) {
276: dims[dim] = timestep+1;
277: maxDims[dim] = H5S_UNLIMITED;
278: chunkDims[dim] = 1;
279: ++dim;
280: }
281: ISGetSize(is, &N);
282: ISGetLocalSize(is, &n);
283: PetscHDF5IntCast(N/bs,dims + dim);
285: maxDims[dim] = dims[dim];
286: chunkDims[dim] = dims[dim];
287: ++dim;
288: if (bs >= 1) {
289: dims[dim] = bs;
290: maxDims[dim] = dims[dim];
291: chunkDims[dim] = dims[dim];
292: ++dim;
293: }
294: filespace = H5Screate_simple(dim, dims, maxDims);
295: if (filespace == -1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"Cannot H5Screate_simple()");
297: #if defined(PETSC_USE_64BIT_INDICES)
298: inttype = H5T_NATIVE_LLONG;
299: #else
300: inttype = H5T_NATIVE_INT;
301: #endif
303: /* Create the dataset with default properties and close filespace */
304: PetscObjectGetName((PetscObject) is, &isname);
305: if (!H5Lexists(group, isname, H5P_DEFAULT)) {
306: /* Create chunk */
307: chunkspace = H5Pcreate(H5P_DATASET_CREATE);
308: if (chunkspace == -1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "Cannot H5Pcreate()");
309: status = H5Pset_chunk(chunkspace, dim, chunkDims);CHKERRQ(status);
311: #if (H5_VERS_MAJOR * 10000 + H5_VERS_MINOR * 100 + H5_VERS_RELEASE >= 10800)
312: dset_id = H5Dcreate2(group, isname, inttype, filespace, H5P_DEFAULT, chunkspace, H5P_DEFAULT);
313: #else
314: dset_id = H5Dcreate(group, isname, inttype, filespace, H5P_DEFAULT);
315: #endif
316: if (dset_id == -1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "Cannot H5Dcreate2()");
317: status = H5Pclose(chunkspace);CHKERRQ(status);
318: } else {
319: dset_id = H5Dopen2(group, isname, H5P_DEFAULT);
320: status = H5Dset_extent(dset_id, dims);CHKERRQ(status);
321: }
322: status = H5Sclose(filespace);CHKERRQ(status);
324: /* Each process defines a dataset and writes it to the hyperslab in the file */
325: dim = 0;
326: if (timestep >= 0) {
327: count[dim] = 1;
328: ++dim;
329: }
330: PetscHDF5IntCast(n/bs,count + dim);
331: ++dim;
332: if (bs >= 1) {
333: count[dim] = bs;
334: ++dim;
335: }
336: if (n > 0) {
337: memspace = H5Screate_simple(dim, count, NULL);
338: if (memspace == -1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "Cannot H5Screate_simple()");
339: } else {
340: /* Can't create dataspace with zero for any dimension, so create null dataspace. */
341: memspace = H5Screate(H5S_NULL);
342: if (memspace == -1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "Cannot H5Screate()");
343: }
345: /* Select hyperslab in the file */
346: PetscLayoutGetRange(is->map, &low, NULL);
347: dim = 0;
348: if (timestep >= 0) {
349: offset[dim] = timestep;
350: ++dim;
351: }
352: PetscHDF5IntCast(low/bs,offset + dim);
353: ++dim;
354: if (bs >= 1) {
355: offset[dim] = 0;
356: ++dim;
357: }
358: if (n > 0) {
359: filespace = H5Dget_space(dset_id);
360: if (filespace == -1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "Cannot H5Dget_space()");
361: status = H5Sselect_hyperslab(filespace, H5S_SELECT_SET, offset, NULL, count, NULL);CHKERRQ(status);
362: } else {
363: /* Create null filespace to match null memspace. */
364: filespace = H5Screate(H5S_NULL);
365: if (filespace == -1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "Cannot H5Screate(H5S_NULL)");
366: }
368: /* Create property list for collective dataset write */
369: plist_id = H5Pcreate(H5P_DATASET_XFER);
370: if (plist_id == -1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "Cannot H5Pcreate()");
371: #if defined(PETSC_HAVE_H5PSET_FAPL_MPIO)
372: status = H5Pset_dxpl_mpio(plist_id, H5FD_MPIO_COLLECTIVE);CHKERRQ(status);
373: #endif
374: /* To write dataset independently use H5Pset_dxpl_mpio(plist_id, H5FD_MPIO_INDEPENDENT) */
376: ISGetIndices(is, &ind);
377: status = H5Dwrite(dset_id, inttype, memspace, filespace, plist_id, ind);CHKERRQ(status);
378: status = H5Fflush(file_id, H5F_SCOPE_GLOBAL);CHKERRQ(status);
379: ISGetIndices(is, &ind);
381: /* Close/release resources */
382: if (group != file_id) {status = H5Gclose(group);CHKERRQ(status);}
383: status = H5Pclose(plist_id);CHKERRQ(status);
384: status = H5Sclose(filespace);CHKERRQ(status);
385: status = H5Sclose(memspace);CHKERRQ(status);
386: status = H5Dclose(dset_id);CHKERRQ(status);
387: PetscInfo1(is, "Wrote IS object with name %s\n", isname);
388: return(0);
389: }
390: #endif
394: PetscErrorCode ISView_General_Binary(IS is,PetscViewer viewer)
395: {
397: IS_General *isa = (IS_General*) is->data;
398: PetscMPIInt rank,size,mesgsize,tag = ((PetscObject)viewer)->tag, mesglen;
399: PetscInt n,N,len,j,tr[2];
400: int fdes;
401: MPI_Status status;
402: PetscInt message_count,flowcontrolcount,*values;
405: PetscLayoutGetLocalSize(is->map, &n);
406: PetscLayoutGetSize(is->map, &N);
407: PetscViewerBinaryGetDescriptor(viewer,&fdes);
409: /* determine maximum message to arrive */
410: MPI_Comm_rank(PetscObjectComm((PetscObject)is),&rank);
411: MPI_Comm_size(PetscObjectComm((PetscObject)is),&size);
413: tr[0] = IS_FILE_CLASSID;
414: tr[1] = N;
415: PetscViewerBinaryWrite(viewer,tr,2,PETSC_INT,PETSC_FALSE);
416: MPI_Reduce(&n,&len,1,MPIU_INT,MPI_SUM,0,PetscObjectComm((PetscObject)is));
418: PetscViewerFlowControlStart(viewer,&message_count,&flowcontrolcount);
419: if (!rank) {
420: PetscBinaryWrite(fdes,isa->idx,n,PETSC_INT,PETSC_FALSE);
422: PetscMalloc1(len,&values);
423: PetscMPIIntCast(len,&mesgsize);
424: /* receive and save messages */
425: for (j=1; j<size; j++) {
426: PetscViewerFlowControlStepMaster(viewer,j,&message_count,flowcontrolcount);
427: MPI_Recv(values,mesgsize,MPIU_INT,j,tag,PetscObjectComm((PetscObject)is),&status);
428: MPI_Get_count(&status,MPIU_INT,&mesglen);
429: PetscBinaryWrite(fdes,values,(PetscInt)mesglen,PETSC_INT,PETSC_TRUE);
430: }
431: PetscViewerFlowControlEndMaster(viewer,&message_count);
432: PetscFree(values);
433: } else {
434: PetscViewerFlowControlStepWorker(viewer,rank,&message_count);
435: PetscMPIIntCast(n,&mesgsize);
436: MPI_Send(isa->idx,mesgsize,MPIU_INT,0,tag,PetscObjectComm((PetscObject)is));
437: PetscViewerFlowControlEndWorker(viewer,&message_count);
438: }
439: return(0);
440: }
444: PetscErrorCode ISView_General(IS is,PetscViewer viewer)
445: {
446: IS_General *sub = (IS_General*)is->data;
448: PetscInt i,n,*idx = sub->idx;
449: PetscBool iascii,isbinary,ishdf5;
452: PetscLayoutGetLocalSize(is->map, &n);
453: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
454: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
455: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
456: if (iascii) {
457: MPI_Comm comm;
458: PetscMPIInt rank,size;
460: PetscObjectGetComm((PetscObject)viewer,&comm);
461: MPI_Comm_rank(comm,&rank);
462: MPI_Comm_size(comm,&size);
464: PetscViewerASCIISynchronizedAllow(viewer,PETSC_TRUE);
465: if (size > 1) {
466: if (is->isperm) {
467: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Index set is permutation\n",rank);
468: }
469: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Number of indices in set %D\n",rank,n);
470: for (i=0; i<n; i++) {
471: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] %D %D\n",rank,i,idx[i]);
472: }
473: } else {
474: if (is->isperm) {
475: PetscViewerASCIISynchronizedPrintf(viewer,"Index set is permutation\n");
476: }
477: PetscViewerASCIISynchronizedPrintf(viewer,"Number of indices in set %D\n",n);
478: for (i=0; i<n; i++) {
479: PetscViewerASCIISynchronizedPrintf(viewer,"%D %D\n",i,idx[i]);
480: }
481: }
482: PetscViewerFlush(viewer);
483: PetscViewerASCIISynchronizedAllow(viewer,PETSC_FALSE);
484: } else if (isbinary) {
485: ISView_General_Binary(is,viewer);
486: } else if (ishdf5) {
487: #if defined(PETSC_HAVE_HDF5)
488: ISView_General_HDF5(is,viewer);
489: #endif
490: }
491: return(0);
492: }
496: PetscErrorCode ISSort_General(IS is)
497: {
498: IS_General *sub = (IS_General*)is->data;
499: PetscInt n;
503: if (sub->sorted) return(0);
504: PetscLayoutGetLocalSize(is->map, &n);
505: PetscSortInt(n,sub->idx);
506: sub->sorted = PETSC_TRUE;
507: return(0);
508: }
512: PetscErrorCode ISSortRemoveDups_General(IS is)
513: {
514: IS_General *sub = (IS_General*)is->data;
515: PetscInt n;
519: if (sub->sorted) return(0);
520: PetscLayoutGetLocalSize(is->map, &n);
521: PetscSortRemoveDupsInt(&n,sub->idx);
522: PetscLayoutSetLocalSize(is->map, n);
523: sub->sorted = PETSC_TRUE;
524: return(0);
525: }
529: PetscErrorCode ISSorted_General(IS is,PetscBool *flg)
530: {
531: IS_General *sub = (IS_General*)is->data;
534: *flg = sub->sorted;
535: return(0);
536: }
540: PetscErrorCode ISToGeneral_General(IS is)
541: {
543: return(0);
544: }
546: static struct _ISOps myops = { ISGetSize_General,
547: ISGetLocalSize_General,
548: ISGetIndices_General,
549: ISRestoreIndices_General,
550: ISInvertPermutation_General,
551: ISSort_General,
552: ISSortRemoveDups_General,
553: ISSorted_General,
554: ISDuplicate_General,
555: ISDestroy_General,
556: ISView_General,
557: ISLoad_Default,
558: ISIdentity_General,
559: ISCopy_General,
560: ISToGeneral_General,
561: ISOnComm_General,
562: ISSetBlockSize_General,
563: ISContiguousLocal_General};
567: PetscErrorCode ISCreateGeneral_Private(IS is)
568: {
570: IS_General *sub = (IS_General*)is->data;
571: const PetscInt *idx = sub->idx;
572: PetscBool sorted = PETSC_TRUE;
573: PetscInt n,i,min,max;
576: PetscLayoutGetLocalSize(is->map, &n);
577: PetscLayoutSetUp(is->map);
578: for (i=1; i<n; i++) {
579: if (idx[i] < idx[i-1]) {sorted = PETSC_FALSE; break;}
580: }
581: if (n) min = max = idx[0];
582: else min = max = 0;
583: for (i=1; i<n; i++) {
584: if (idx[i] < min) min = idx[i];
585: if (idx[i] > max) max = idx[i];
586: }
587: sub->sorted = sorted;
588: is->min = min;
589: is->max = max;
590: is->isperm = PETSC_FALSE;
591: is->isidentity = PETSC_FALSE;
592: ISViewFromOptions(is,NULL,"-is_view");
593: return(0);
594: }
598: /*@
599: ISCreateGeneral - Creates a data structure for an index set
600: containing a list of integers.
602: Collective on MPI_Comm
604: Input Parameters:
605: + comm - the MPI communicator
606: . n - the length of the index set
607: . idx - the list of integers
608: - mode - see PetscCopyMode for meaning of this flag.
610: Output Parameter:
611: . is - the new index set
613: Notes:
614: When the communicator is not MPI_COMM_SELF, the operations on IS are NOT
615: conceptually the same as MPI_Group operations. The IS are then
616: distributed sets of indices and thus certain operations on them are
617: collective.
620: Level: beginner
622: Concepts: index sets^creating
623: Concepts: IS^creating
625: .seealso: ISCreateStride(), ISCreateBlock(), ISAllGather()
626: @*/
627: PetscErrorCode ISCreateGeneral(MPI_Comm comm,PetscInt n,const PetscInt idx[],PetscCopyMode mode,IS *is)
628: {
632: ISCreate(comm,is);
633: ISSetType(*is,ISGENERAL);
634: ISGeneralSetIndices(*is,n,idx,mode);
635: return(0);
636: }
640: /*@
641: ISGeneralSetIndices - Sets the indices for an ISGENERAL index set
643: Collective on IS
645: Input Parameters:
646: + is - the index set
647: . n - the length of the index set
648: . idx - the list of integers
649: - mode - see PetscCopyMode for meaning of this flag.
651: Level: beginner
653: Concepts: index sets^creating
654: Concepts: IS^creating
656: .seealso: ISCreateGeneral(), ISCreateStride(), ISCreateBlock(), ISAllGather()
657: @*/
658: PetscErrorCode ISGeneralSetIndices(IS is,PetscInt n,const PetscInt idx[],PetscCopyMode mode)
659: {
663: PetscUseMethod(is,"ISGeneralSetIndices_C",(IS,PetscInt,const PetscInt[],PetscCopyMode),(is,n,idx,mode));
664: return(0);
665: }
669: PetscErrorCode ISGeneralSetIndices_General(IS is,PetscInt n,const PetscInt idx[],PetscCopyMode mode)
670: {
672: IS_General *sub = (IS_General*)is->data;
675: if (n < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"length < 0");
678: if (sub->allocated) {PetscFree(sub->idx);}
679: PetscLayoutSetLocalSize(is->map, n);
680: if (mode == PETSC_COPY_VALUES) {
681: PetscMalloc1(n,&sub->idx);
682: PetscLogObjectMemory((PetscObject)is,n*sizeof(PetscInt));
683: PetscMemcpy(sub->idx,idx,n*sizeof(PetscInt));
684: sub->allocated = PETSC_TRUE;
685: } else if (mode == PETSC_OWN_POINTER) {
686: sub->idx = (PetscInt*)idx;
687: sub->allocated = PETSC_TRUE;
688: } else {
689: sub->idx = (PetscInt*)idx;
690: sub->allocated = PETSC_FALSE;
691: }
692: ISCreateGeneral_Private(is);
693: return(0);
694: }
698: PETSC_EXTERN PetscErrorCode ISCreate_General(IS is)
699: {
701: IS_General *sub;
704: PetscMemcpy(is->ops,&myops,sizeof(myops));
705: PetscNewLog(is,&sub);
706: is->data = (void *) sub;
707: PetscObjectComposeFunction((PetscObject)is,"ISGeneralSetIndices_C",ISGeneralSetIndices_General);
708: return(0);
709: }