Actual source code: general.c

petsc-3.9.4 2018-09-11
Report Typos and Errors

  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>
  6:  #include <petscvec.h>
  7:  #include <petscviewer.h>
  8: #include <petscviewerhdf5.h>

 10: static PetscErrorCode ISDuplicate_General(IS is,IS *newIS)
 11: {
 13:   IS_General     *sub = (IS_General*)is->data;
 14:   PetscInt       n;

 17:   PetscLayoutGetLocalSize(is->map, &n);
 18:   ISCreateGeneral(PetscObjectComm((PetscObject) is), n, sub->idx, PETSC_COPY_VALUES, newIS);
 19:   return(0);
 20: }

 22: static PetscErrorCode ISDestroy_General(IS is)
 23: {
 24:   IS_General     *is_general = (IS_General*)is->data;

 28:   if (is_general->allocated) {PetscFree(is_general->idx);}
 29:   PetscObjectComposeFunction((PetscObject)is,"ISGeneralSetIndices_C",NULL);
 30:   PetscFree(is->data);
 31:   return(0);
 32: }

 34: static PetscErrorCode ISIdentity_General(IS is, PetscBool *ident)
 35: {
 36:   IS_General *is_general = (IS_General*)is->data;
 37:   PetscInt   i,n,*idx = is_general->idx;

 41:   PetscLayoutGetLocalSize(is->map, &n);
 42:   is->isidentity = PETSC_TRUE;
 43:   *ident         = PETSC_TRUE;
 44:   for (i=0; i<n; i++) {
 45:     if (idx[i] != i) {
 46:       is->isidentity = PETSC_FALSE;
 47:       *ident         = PETSC_FALSE;
 48:       break;
 49:     }
 50:   }
 51:   return(0);
 52: }

 54: static PetscErrorCode ISCopy_General(IS is,IS isy)
 55: {
 56:   IS_General     *is_general = (IS_General*)is->data,*isy_general = (IS_General*)isy->data;
 57:   PetscInt       n, N, ny, Ny;

 61:   PetscLayoutGetLocalSize(is->map, &n);
 62:   PetscLayoutGetSize(is->map, &N);
 63:   PetscLayoutGetLocalSize(isy->map, &ny);
 64:   PetscLayoutGetSize(isy->map, &Ny);
 65:   if (n != ny || N != Ny) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Index sets incompatible");
 66:   isy_general->sorted = is_general->sorted;
 67:   PetscMemcpy(isy_general->idx,is_general->idx,n*sizeof(PetscInt));
 68:   return(0);
 69: }

 71: static PetscErrorCode ISOnComm_General(IS is,MPI_Comm comm,PetscCopyMode mode,IS *newis)
 72: {
 74:   IS_General     *sub = (IS_General*)is->data;
 75:   PetscInt       n;

 78:   if (mode == PETSC_OWN_POINTER) SETERRQ(comm,PETSC_ERR_ARG_WRONG,"Cannot use PETSC_OWN_POINTER");
 79:   PetscLayoutGetLocalSize(is->map, &n);
 80:   ISCreateGeneral(comm,n,sub->idx,mode,newis);
 81:   return(0);
 82: }

 84: static PetscErrorCode ISSetBlockSize_General(IS is,PetscInt bs)
 85: {

 89:   PetscLayoutSetBlockSize(is->map, bs);
 90:   return(0);
 91: }

 93: static PetscErrorCode ISContiguousLocal_General(IS is,PetscInt gstart,PetscInt gend,PetscInt *start,PetscBool *contig)
 94: {
 95:   IS_General *sub = (IS_General*)is->data;
 96:   PetscInt   n,i,p;

100:   *start  = 0;
101:   *contig = PETSC_TRUE;
102:   PetscLayoutGetLocalSize(is->map, &n);
103:   if (!n) return(0);
104:   p = sub->idx[0];
105:   if (p < gstart) goto nomatch;
106:   *start = p - gstart;
107:   if (n > gend-p) goto nomatch;
108:   for (i=1; i<n; i++,p++) {
109:     if (sub->idx[i] != p+1) goto nomatch;
110:   }
111:   return(0);
112: nomatch:
113:   *start  = -1;
114:   *contig = PETSC_FALSE;
115:   return(0);
116: }

118: static PetscErrorCode ISLocate_General(IS is,PetscInt key,PetscInt *location)
119: {
120:   IS_General     *sub = (IS_General*)is->data;
121:   PetscInt       numIdx, i;

125:   PetscLayoutGetLocalSize(is->map,&numIdx);
126:   if (sub->sorted) { PetscFindInt(key,numIdx,sub->idx,location);}
127:   else {
128:     const PetscInt *idx = sub->idx;

130:     *location = -1;
131:     for (i = 0; i < numIdx; i++) {
132:       if (idx[i] == key) {
133:         *location = i;
134:         return(0);
135:       }
136:     }
137:   }
138:   return(0);
139: }

141: static PetscErrorCode ISGetIndices_General(IS in,const PetscInt *idx[])
142: {
143:   IS_General *sub = (IS_General*)in->data;

146:   *idx = sub->idx;
147:   return(0);
148: }

150: static PetscErrorCode ISRestoreIndices_General(IS in,const PetscInt *idx[])
151: {
152:   IS_General *sub = (IS_General*)in->data;

155:   if (*idx != sub->idx) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Must restore with value from ISGetIndices()");
156:   return(0);
157: }

159: static PetscErrorCode ISGetSize_General(IS is,PetscInt *size)
160: {

164:   PetscLayoutGetSize(is->map, size);
165:   return(0);
166: }

168: static PetscErrorCode ISGetLocalSize_General(IS is,PetscInt *size)
169: {

173:   PetscLayoutGetLocalSize(is->map, size);
174:   return(0);
175: }

177: static PetscErrorCode ISInvertPermutation_General(IS is,PetscInt nlocal,IS *isout)
178: {
179:   IS_General     *sub = (IS_General*)is->data;
180:   PetscInt       i,*ii,n,nstart;
181:   const PetscInt *idx = sub->idx;
182:   PetscMPIInt    size;
183:   IS             istmp,nistmp;

187:   PetscLayoutGetLocalSize(is->map, &n);
188:   MPI_Comm_size(PetscObjectComm((PetscObject)is),&size);
189:   if (size == 1) {
190:     PetscMalloc1(n,&ii);
191:     for (i=0; i<n; i++) ii[idx[i]] = i;
192:     ISCreateGeneral(PETSC_COMM_SELF,n,ii,PETSC_OWN_POINTER,isout);
193:     ISSetPermutation(*isout);
194:   } else {
195:     /* crude, nonscalable get entire IS on each processor */
196:     if (nlocal == PETSC_DECIDE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Do not yet support nlocal of PETSC_DECIDE");
197:     ISAllGather(is,&istmp);
198:     ISSetPermutation(istmp);
199:     ISInvertPermutation(istmp,PETSC_DECIDE,&nistmp);
200:     ISDestroy(&istmp);
201:     /* get the part we need */
202:     MPI_Scan(&nlocal,&nstart,1,MPIU_INT,MPI_SUM,PetscObjectComm((PetscObject)is));
203: #if defined(PETSC_USE_DEBUG)
204:     {
205:       PetscInt    N;
206:       PetscMPIInt rank;
207:       MPI_Comm_rank(PetscObjectComm((PetscObject)is),&rank);
208:       PetscLayoutGetSize(is->map, &N);
209:       if (rank == size-1) {
210:         if (nstart != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Sum of nlocal lengths %d != total IS length %d",nstart,N);
211:       }
212:     }
213: #endif
214:     nstart -= nlocal;
215:     ISGetIndices(nistmp,&idx);
216:     ISCreateGeneral(PetscObjectComm((PetscObject)is),nlocal,idx+nstart,PETSC_COPY_VALUES,isout);
217:     ISRestoreIndices(nistmp,&idx);
218:     ISDestroy(&nistmp);
219:   }
220:   return(0);
221: }

223: #if defined(PETSC_HAVE_HDF5)
224: static PetscErrorCode ISView_General_HDF5(IS is, PetscViewer viewer)
225: {
226:   hid_t           filespace;  /* file dataspace identifier */
227:   hid_t           chunkspace; /* chunk dataset property identifier */
228:   hid_t           plist_id;   /* property list identifier */
229:   hid_t           dset_id;    /* dataset identifier */
230:   hid_t           memspace;   /* memory dataspace identifier */
231:   hid_t           inttype;    /* int type (H5T_NATIVE_INT or H5T_NATIVE_LLONG) */
232:   hid_t           file_id, group;
233:   hsize_t         dim, maxDims[3], dims[3], chunkDims[3], count[3],offset[3];
234:   PetscInt        bs, N, n, timestep, low;
235:   const PetscInt *ind;
236:   const char     *isname;
237:   PetscErrorCode  ierr;

240:   ISGetBlockSize(is,&bs);
241:   PetscViewerHDF5OpenGroup(viewer, &file_id, &group);
242:   PetscViewerHDF5GetTimestep(viewer, &timestep);

244:   /* Create the dataspace for the dataset.
245:    *
246:    * dims - holds the current dimensions of the dataset
247:    *
248:    * maxDims - holds the maximum dimensions of the dataset (unlimited
249:    * for the number of time steps with the current dimensions for the
250:    * other dimensions; so only additional time steps can be added).
251:    *
252:    * chunkDims - holds the size of a single time step (required to
253:    * permit extending dataset).
254:    */
255:   dim = 0;
256:   if (timestep >= 0) {
257:     dims[dim]      = timestep+1;
258:     maxDims[dim]   = H5S_UNLIMITED;
259:     chunkDims[dim] = 1;
260:     ++dim;
261:   }
262:   ISGetSize(is, &N);
263:   ISGetLocalSize(is, &n);
264:   PetscHDF5IntCast(N/bs,dims + dim);

266:   maxDims[dim]   = dims[dim];
267:   chunkDims[dim] = dims[dim];
268:   ++dim;
269:   if (bs >= 1) {
270:     dims[dim]      = bs;
271:     maxDims[dim]   = dims[dim];
272:     chunkDims[dim] = dims[dim];
273:     ++dim;
274:   }
275:   PetscStackCallHDF5Return(filespace,H5Screate_simple,(dim, dims, maxDims));

277: #if defined(PETSC_USE_64BIT_INDICES)
278:   inttype = H5T_NATIVE_LLONG;
279: #else
280:   inttype = H5T_NATIVE_INT;
281: #endif

283:   /* Create the dataset with default properties and close filespace */
284:   PetscObjectGetName((PetscObject) is, &isname);
285:   if (!H5Lexists(group, isname, H5P_DEFAULT)) {
286:     /* Create chunk */
287:     PetscStackCallHDF5Return(chunkspace,H5Pcreate,(H5P_DATASET_CREATE));
288:     PetscStackCallHDF5(H5Pset_chunk,(chunkspace, dim, chunkDims));

290: #if (H5_VERS_MAJOR * 10000 + H5_VERS_MINOR * 100 + H5_VERS_RELEASE >= 10800)
291:     PetscStackCallHDF5Return(dset_id,H5Dcreate2,(group, isname, inttype, filespace, H5P_DEFAULT, chunkspace, H5P_DEFAULT));
292: #else
293:     PetscStackCallHDF5Return(dset_id,H5Dcreate,(group, isname, inttype, filespace, H5P_DEFAULT));
294: #endif
295:     PetscStackCallHDF5(H5Pclose,(chunkspace));
296:   } else {
297:     PetscStackCallHDF5Return(dset_id,H5Dopen2,(group, isname, H5P_DEFAULT));
298:     PetscStackCallHDF5(H5Dset_extent,(dset_id, dims));
299:   }
300:   PetscStackCallHDF5(H5Sclose,(filespace));

302:   /* Each process defines a dataset and writes it to the hyperslab in the file */
303:   dim = 0;
304:   if (timestep >= 0) {
305:     count[dim] = 1;
306:     ++dim;
307:   }
308:   PetscHDF5IntCast(n/bs,count + dim);
309:   ++dim;
310:   if (bs >= 1) {
311:     count[dim] = bs;
312:     ++dim;
313:   }
314:   if (n > 0) {
315:     PetscStackCallHDF5Return(memspace,H5Screate_simple,(dim, count, NULL));
316:   } else {
317:     /* Can't create dataspace with zero for any dimension, so create null dataspace. */
318:     PetscStackCallHDF5Return(memspace,H5Screate,(H5S_NULL));
319:   }

321:   /* Select hyperslab in the file */
322:   PetscLayoutGetRange(is->map, &low, NULL);
323:   dim  = 0;
324:   if (timestep >= 0) {
325:     offset[dim] = timestep;
326:     ++dim;
327:   }
328:   PetscHDF5IntCast(low/bs,offset + dim);
329:   ++dim;
330:   if (bs >= 1) {
331:     offset[dim] = 0;
332:     ++dim;
333:   }
334:   if (n > 0) {
335:     PetscStackCallHDF5Return(filespace,H5Dget_space,(dset_id));
336:     PetscStackCallHDF5(H5Sselect_hyperslab,(filespace, H5S_SELECT_SET, offset, NULL, count, NULL));
337:   } else {
338:     /* Create null filespace to match null memspace. */
339:     PetscStackCallHDF5Return(filespace,H5Screate,(H5S_NULL));
340:   }

342:   /* Create property list for collective dataset write */
343:   PetscStackCallHDF5Return(plist_id,H5Pcreate,(H5P_DATASET_XFER));
344: #if defined(PETSC_HAVE_H5PSET_FAPL_MPIO)
345:   PetscStackCallHDF5(H5Pset_dxpl_mpio,(plist_id, H5FD_MPIO_COLLECTIVE));
346: #endif
347:   /* To write dataset independently use H5Pset_dxpl_mpio(plist_id, H5FD_MPIO_INDEPENDENT) */

349:   ISGetIndices(is, &ind);
350:   PetscStackCallHDF5(H5Dwrite,(dset_id, inttype, memspace, filespace, plist_id, ind));
351:   PetscStackCallHDF5(H5Fflush,(file_id, H5F_SCOPE_GLOBAL));
352:   ISRestoreIndices(is, &ind);

354:   /* Close/release resources */
355:   if (group != file_id) PetscStackCallHDF5(H5Gclose,(group));
356:   PetscStackCallHDF5(H5Pclose,(plist_id));
357:   PetscStackCallHDF5(H5Sclose,(filespace));
358:   PetscStackCallHDF5(H5Sclose,(memspace));
359:   PetscStackCallHDF5(H5Dclose,(dset_id));
360:   PetscInfo1(is, "Wrote IS object with name %s\n", isname);
361:   return(0);
362: }
363: #endif

365: static PetscErrorCode ISView_General_Binary(IS is,PetscViewer viewer)
366: {
368:   PetscBool      skipHeader,useMPIIO;
369:   IS_General     *isa = (IS_General*) is->data;
370:   PetscMPIInt    rank,size,mesgsize,tag = ((PetscObject)viewer)->tag, mesglen;
371:   PetscInt       n,N,len,j,tr[2];
372:   int            fdes;
373:   MPI_Status     status;
374:   PetscInt       message_count,flowcontrolcount,*values;

377:   /* ISGetLayout(is,&map); */
378:   PetscLayoutGetLocalSize(is->map, &n);
379:   PetscLayoutGetSize(is->map, &N);

381:   tr[0] = IS_FILE_CLASSID;
382:   tr[1] = N;

384:   PetscViewerBinaryGetSkipHeader(viewer,&skipHeader);
385:   if (!skipHeader) {
386:     PetscViewerBinaryWrite(viewer,tr,2,PETSC_INT,PETSC_FALSE);
387:   }

389:   PetscViewerBinaryGetUseMPIIO(viewer,&useMPIIO);
390: #if defined(PETSC_HAVE_MPIIO)
391:   if (useMPIIO) {
392:     MPI_File       mfdes;
393:     MPI_Offset     off;
394:     PetscMPIInt    lsize;
395:     PetscInt       rstart;
396:     const PetscInt *iarray;

398:     PetscMPIIntCast(n,&lsize);
399:     PetscViewerBinaryGetMPIIODescriptor(viewer,&mfdes);
400:     PetscViewerBinaryGetMPIIOOffset(viewer,&off);
401:     PetscLayoutGetRange(is->map,&rstart,NULL);
402:     off += rstart*(MPI_Offset)sizeof(PetscInt); /* off is MPI_Offset, not PetscMPIInt */
403:     MPI_File_set_view(mfdes,off,MPIU_INT,MPIU_INT,(char*)"native",MPI_INFO_NULL);
404:     ISGetIndices(is,&iarray);
405:     MPIU_File_write_all(mfdes,(void*)iarray,lsize,MPIU_INT,MPI_STATUS_IGNORE);
406:     ISRestoreIndices(is,&iarray);
407:     PetscViewerBinaryAddMPIIOOffset(viewer,N*(MPI_Offset)sizeof(PetscInt));
408:     return(0);
409:   }
410: #endif

412:   PetscViewerBinaryGetDescriptor(viewer,&fdes);
413:   MPI_Comm_rank(PetscObjectComm((PetscObject)is),&rank);
414:   MPI_Comm_size(PetscObjectComm((PetscObject)is),&size);

416:   /* determine maximum message to arrive */
417:   MPI_Reduce(&n,&len,1,MPIU_INT,MPI_MAX,0,PetscObjectComm((PetscObject)is));

419:   PetscViewerFlowControlStart(viewer,&message_count,&flowcontrolcount);
420:   if (!rank) {
421:     PetscBinaryWrite(fdes,isa->idx,n,PETSC_INT,PETSC_FALSE);

423:     PetscMalloc1(len,&values);
424:     PetscMPIIntCast(len,&mesgsize);
425:     /* receive and save messages */
426:     for (j=1; j<size; j++) {
427:       PetscViewerFlowControlStepMaster(viewer,j,&message_count,flowcontrolcount);
428:       MPI_Recv(values,mesgsize,MPIU_INT,j,tag,PetscObjectComm((PetscObject)is),&status);
429:       MPI_Get_count(&status,MPIU_INT,&mesglen);
430:       PetscBinaryWrite(fdes,values,(PetscInt)mesglen,PETSC_INT,PETSC_TRUE);
431:     }
432:     PetscViewerFlowControlEndMaster(viewer,&message_count);
433:     PetscFree(values);
434:   } else {
435:     PetscViewerFlowControlStepWorker(viewer,rank,&message_count);
436:     PetscMPIIntCast(n,&mesgsize);
437:     MPI_Send(isa->idx,mesgsize,MPIU_INT,0,tag,PetscObjectComm((PetscObject)is));
438:     PetscViewerFlowControlEndWorker(viewer,&message_count);
439:   }
440:   return(0);
441: }

443: static PetscErrorCode ISView_General(IS is,PetscViewer viewer)
444: {
445:   IS_General     *sub = (IS_General*)is->data;
447:   PetscInt       i,n,*idx = sub->idx;
448:   PetscBool      iascii,isbinary,ishdf5,matl;

451:   PetscLayoutGetLocalSize(is->map, &n);
452:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
453:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
454:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
455:   if (iascii) {
456:     MPI_Comm          comm;
457:     PetscMPIInt       rank,size;
458:     PetscViewerFormat fmt;

460:     PetscObjectGetComm((PetscObject)viewer,&comm);
461:     MPI_Comm_rank(comm,&rank);
462:     MPI_Comm_size(comm,&size);

464:     PetscViewerGetFormat(viewer,&fmt);
465:     matl = (PetscBool)(fmt == PETSC_VIEWER_ASCII_MATLAB);
466:     PetscViewerASCIIPushSynchronized(viewer);
467:     if (size > 1) {
468:       if (matl) {
469:         const char* name;

471:         PetscObjectGetName((PetscObject)is,&name);
472:         PetscViewerASCIISynchronizedPrintf(viewer,"%s_%d = [...\n",name,rank);
473:         for (i=0; i<n; i++) {
474:           PetscViewerASCIISynchronizedPrintf(viewer,"%D\n",idx[i]+1);
475:         }
476:         PetscViewerASCIISynchronizedPrintf(viewer,"];\n");
477:       } else {
478:         if (is->isperm) {
479:           PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Index set is permutation\n",rank);
480:         }
481:         PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Number of indices in set %D\n",rank,n);
482:         for (i=0; i<n; i++) {
483:           PetscViewerASCIISynchronizedPrintf(viewer,"[%d] %D %D\n",rank,i,idx[i]);
484:         }
485:       }
486:     } else {
487:       if (matl) {
488:         const char* name;

490:         PetscObjectGetName((PetscObject)is,&name);
491:         PetscViewerASCIISynchronizedPrintf(viewer,"%s = [...\n",name);
492:         for (i=0; i<n; i++) {
493:           PetscViewerASCIISynchronizedPrintf(viewer,"%D\n",idx[i]+1);
494:         }
495:         PetscViewerASCIISynchronizedPrintf(viewer,"];\n");
496:       } else {
497:         if (is->isperm) {
498:           PetscViewerASCIISynchronizedPrintf(viewer,"Index set is permutation\n");
499:         }
500:         PetscViewerASCIISynchronizedPrintf(viewer,"Number of indices in set %D\n",n);
501:         for (i=0; i<n; i++) {
502:           PetscViewerASCIISynchronizedPrintf(viewer,"%D %D\n",i,idx[i]);
503:         }
504:       }
505:     }
506:     PetscViewerFlush(viewer);
507:     PetscViewerASCIIPopSynchronized(viewer);
508:   } else if (isbinary) {
509:     ISView_General_Binary(is,viewer);
510:   } else if (ishdf5) {
511: #if defined(PETSC_HAVE_HDF5)
512:     ISView_General_HDF5(is,viewer);
513: #endif
514:   }
515:   return(0);
516: }

518: static PetscErrorCode ISSort_General(IS is)
519: {
520:   IS_General     *sub = (IS_General*)is->data;
521:   PetscInt       n;

525:   if (sub->sorted) return(0);
526:   PetscLayoutGetLocalSize(is->map, &n);
527:   PetscSortInt(n,sub->idx);
528:   sub->sorted = PETSC_TRUE;
529:   return(0);
530: }

532: static PetscErrorCode ISSortRemoveDups_General(IS is)
533: {
534:   IS_General     *sub = (IS_General*)is->data;
535:   PetscInt       n;

539:   PetscLayoutGetLocalSize(is->map, &n);
540:   if (sub->sorted) {
541:     PetscSortedRemoveDupsInt(&n,sub->idx);
542:   } else {
543:     PetscSortRemoveDupsInt(&n,sub->idx);
544:   }
545:   PetscLayoutSetLocalSize(is->map, n);
546:   PetscLayoutSetSize(is->map, PETSC_DECIDE);
547:   PetscLayoutSetUp(is->map);
548:   sub->sorted = PETSC_TRUE;
549:   return(0);
550: }

552: static PetscErrorCode ISSorted_General(IS is,PetscBool  *flg)
553: {
554:   IS_General *sub = (IS_General*)is->data;

557:   *flg = sub->sorted;
558:   return(0);
559: }

561: PetscErrorCode  ISToGeneral_General(IS is)
562: {
564:   return(0);
565: }

567: static struct _ISOps myops = { ISGetSize_General,
568:                                ISGetLocalSize_General,
569:                                ISGetIndices_General,
570:                                ISRestoreIndices_General,
571:                                ISInvertPermutation_General,
572:                                ISSort_General,
573:                                ISSortRemoveDups_General,
574:                                ISSorted_General,
575:                                ISDuplicate_General,
576:                                ISDestroy_General,
577:                                ISView_General,
578:                                ISLoad_Default,
579:                                ISIdentity_General,
580:                                ISCopy_General,
581:                                ISToGeneral_General,
582:                                ISOnComm_General,
583:                                ISSetBlockSize_General,
584:                                ISContiguousLocal_General,
585:                                ISLocate_General};

587: PETSC_INTERN PetscErrorCode ISSetUp_General(IS);

589: PetscErrorCode ISSetUp_General(IS is)
590: {
592:   IS_General     *sub = (IS_General*)is->data;
593:   const PetscInt *idx = sub->idx;
594:   PetscInt       n,i,min,max;

597:   PetscLayoutGetLocalSize(is->map, &n);

599:   sub->sorted = PETSC_TRUE;
600:   for (i=1; i<n; i++) {
601:     if (idx[i] < idx[i-1]) {sub->sorted = PETSC_FALSE; break;}
602:   }
603:   if (n) {
604:     min = max = idx[0];
605:     for (i=1; i<n; i++) {
606:       if (idx[i] < min) min = idx[i];
607:       if (idx[i] > max) max = idx[i];
608:     }
609:     is->min = min;
610:     is->max = max;
611:   } else {
612:     is->min = PETSC_MAX_INT;
613:     is->max = PETSC_MIN_INT;
614:   }
615:   is->isperm     = PETSC_FALSE;
616:   is->isidentity = PETSC_FALSE;
617:   return(0);
618: }

620: /*@
621:    ISCreateGeneral - Creates a data structure for an index set
622:    containing a list of integers.

624:    Collective on MPI_Comm

626:    Input Parameters:
627: +  comm - the MPI communicator
628: .  n - the length of the index set
629: .  idx - the list of integers
630: -  mode - PETSC_COPY_VALUES, PETSC_OWN_POINTER, or PETSC_USE_POINTER; see PetscCopyMode for meaning of this flag.

632:    Output Parameter:
633: .  is - the new index set

635:    Notes:
636:    When the communicator is not MPI_COMM_SELF, the operations on IS are NOT
637:    conceptually the same as MPI_Group operations. The IS are then
638:    distributed sets of indices and thus certain operations on them are
639:    collective.


642:    Level: beginner

644:   Concepts: index sets^creating
645:   Concepts: IS^creating

647: .seealso: ISCreateStride(), ISCreateBlock(), ISAllGather(), PETSC_COPY_VALUES, PETSC_OWN_POINTER, PETSC_USE_POINTER, PetscCopyMode
648: @*/
649: PetscErrorCode  ISCreateGeneral(MPI_Comm comm,PetscInt n,const PetscInt idx[],PetscCopyMode mode,IS *is)
650: {

654:   ISCreate(comm,is);
655:   ISSetType(*is,ISGENERAL);
656:   ISGeneralSetIndices(*is,n,idx,mode);
657:   return(0);
658: }

660: /*@
661:    ISGeneralSetIndices - Sets the indices for an ISGENERAL index set

663:    Collective on IS

665:    Input Parameters:
666: +  is - the index set
667: .  n - the length of the index set
668: .  idx - the list of integers
669: -  mode - see PetscCopyMode for meaning of this flag.

671:    Level: beginner

673:   Concepts: index sets^creating
674:   Concepts: IS^creating

676: .seealso: ISCreateGeneral(), ISCreateStride(), ISCreateBlock(), ISAllGather()
677: @*/
678: PetscErrorCode  ISGeneralSetIndices(IS is,PetscInt n,const PetscInt idx[],PetscCopyMode mode)
679: {

683:   PetscUseMethod(is,"ISGeneralSetIndices_C",(IS,PetscInt,const PetscInt[],PetscCopyMode),(is,n,idx,mode));
684:   return(0);
685: }

687: PetscErrorCode  ISGeneralSetIndices_General(IS is,PetscInt n,const PetscInt idx[],PetscCopyMode mode)
688: {
690:   IS_General     *sub = (IS_General*)is->data;

693:   if (n < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"length < 0");

696:   PetscLayoutSetLocalSize(is->map,n);
697:   PetscLayoutSetUp(is->map);

699:   if (sub->allocated) {PetscFree(sub->idx);}
700:   if (mode == PETSC_COPY_VALUES) {
701:     PetscMalloc1(n,&sub->idx);
702:     PetscLogObjectMemory((PetscObject)is,n*sizeof(PetscInt));
703:     PetscMemcpy(sub->idx,idx,n*sizeof(PetscInt));
704:     sub->allocated = PETSC_TRUE;
705:   } else if (mode == PETSC_OWN_POINTER) {
706:     sub->idx = (PetscInt*)idx;
707:     PetscLogObjectMemory((PetscObject)is,n*sizeof(PetscInt));
708:     sub->allocated = PETSC_TRUE;
709:   } else {
710:     sub->idx = (PetscInt*)idx;
711:     sub->allocated = PETSC_FALSE;
712:   }

714:   ISSetUp_General(is);
715:   ISViewFromOptions(is,NULL,"-is_view");
716:   return(0);
717: }

719: PETSC_EXTERN PetscErrorCode ISCreate_General(IS is)
720: {
722:   IS_General     *sub;

725:   PetscNewLog(is,&sub);
726:   is->data = (void *) sub;
727:   PetscMemcpy(is->ops,&myops,sizeof(myops));
728:   PetscObjectComposeFunction((PetscObject)is,"ISGeneralSetIndices_C",ISGeneralSetIndices_General);
729:   return(0);
730: }