Actual source code: general.c

petsc-3.11.4 2019-09-28
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:     ISAllGather(is,&istmp);
197:     ISSetPermutation(istmp);
198:     ISInvertPermutation(istmp,PETSC_DECIDE,&nistmp);
199:     ISDestroy(&istmp);
200:     /* get the part we need */
201:     if (nlocal == PETSC_DECIDE) nlocal = n;
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] = PetscMax(1,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:     PetscStackCallHDF5Return(dset_id,H5Dcreate2,(group, isname, inttype, filespace, H5P_DEFAULT, chunkspace, H5P_DEFAULT));
291:     PetscStackCallHDF5(H5Pclose,(chunkspace));
292:   } else {
293:     PetscStackCallHDF5Return(dset_id,H5Dopen2,(group, isname, H5P_DEFAULT));
294:     PetscStackCallHDF5(H5Dset_extent,(dset_id, dims));
295:   }
296:   PetscStackCallHDF5(H5Sclose,(filespace));

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

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

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

345:   ISGetIndices(is, &ind);
346:   PetscStackCallHDF5(H5Dwrite,(dset_id, inttype, memspace, filespace, plist_id, ind));
347:   PetscStackCallHDF5(H5Fflush,(file_id, H5F_SCOPE_GLOBAL));
348:   ISRestoreIndices(is, &ind);

350:   /* Close/release resources */
351:   PetscStackCallHDF5(H5Gclose,(group));
352:   PetscStackCallHDF5(H5Pclose,(plist_id));
353:   PetscStackCallHDF5(H5Sclose,(filespace));
354:   PetscStackCallHDF5(H5Sclose,(memspace));
355:   PetscStackCallHDF5(H5Dclose,(dset_id));
356:   PetscInfo1(is, "Wrote IS object with name %s\n", isname);
357:   return(0);
358: }
359: #endif

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

373:   /* ISGetLayout(is,&map); */
374:   PetscLayoutGetLocalSize(is->map, &n);
375:   PetscLayoutGetSize(is->map, &N);

377:   tr[0] = IS_FILE_CLASSID;
378:   tr[1] = N;

380:   PetscViewerBinaryGetSkipHeader(viewer,&skipHeader);
381:   if (!skipHeader) {
382:     PetscViewerBinaryWrite(viewer,tr,2,PETSC_INT,PETSC_FALSE);
383:   }

385:   PetscViewerBinaryGetUseMPIIO(viewer,&useMPIIO);
386: #if defined(PETSC_HAVE_MPIIO)
387:   if (useMPIIO) {
388:     MPI_File       mfdes;
389:     MPI_Offset     off;
390:     PetscMPIInt    lsize;
391:     PetscInt       rstart;
392:     const PetscInt *iarray;

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

408:   PetscViewerBinaryGetDescriptor(viewer,&fdes);
409:   MPI_Comm_rank(PetscObjectComm((PetscObject)is),&rank);
410:   MPI_Comm_size(PetscObjectComm((PetscObject)is),&size);

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

415:   PetscViewerFlowControlStart(viewer,&message_count,&flowcontrolcount);
416:   if (!rank) {
417:     PetscBinaryWrite(fdes,isa->idx,n,PETSC_INT,PETSC_FALSE);

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

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

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

456:     PetscObjectGetComm((PetscObject)viewer,&comm);
457:     MPI_Comm_rank(comm,&rank);
458:     MPI_Comm_size(comm,&size);

460:     PetscViewerGetFormat(viewer,&fmt);
461:     PetscViewerASCIIPushSynchronized(viewer);
462:     if (size > 1) {
463:       if (fmt == PETSC_VIEWER_ASCII_MATLAB) {
464:         const char* name;

466:         PetscObjectGetName((PetscObject)is,&name);
467:         PetscViewerASCIISynchronizedPrintf(viewer,"%s_%d = [...\n",name,rank);
468:         for (i=0; i<n; i++) {
469:           PetscViewerASCIISynchronizedPrintf(viewer,"%D\n",idx[i]+1);
470:         }
471:         PetscViewerASCIISynchronizedPrintf(viewer,"];\n");
472:       } else {
473:         PetscInt st = 0;

475:         if (fmt == PETSC_VIEWER_ASCII_INDEX) st = is->map->rstart;
476:         if (is->isperm) {
477:           PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Index set is permutation\n",rank);
478:         }
479:         PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Number of indices in set %D\n",rank,n);
480:         for (i=0; i<n; i++) {
481:           PetscViewerASCIISynchronizedPrintf(viewer,"[%d] %D %D\n",rank,i + st,idx[i]);
482:         }
483:       }
484:     } else {
485:       if (fmt == PETSC_VIEWER_ASCII_MATLAB) {
486:         const char* name;

488:         PetscObjectGetName((PetscObject)is,&name);
489:         PetscViewerASCIISynchronizedPrintf(viewer,"%s = [...\n",name);
490:         for (i=0; i<n; i++) {
491:           PetscViewerASCIISynchronizedPrintf(viewer,"%D\n",idx[i]+1);
492:         }
493:         PetscViewerASCIISynchronizedPrintf(viewer,"];\n");
494:       } else {
495:         PetscInt st = 0;

497:         if (fmt == PETSC_VIEWER_ASCII_INDEX) st = is->map->rstart;
498:         if (is->isperm) {
499:           PetscViewerASCIISynchronizedPrintf(viewer,"Index set is permutation\n");
500:         }
501:         PetscViewerASCIISynchronizedPrintf(viewer,"Number of indices in set %D\n",n);
502:         for (i=0; i<n; i++) {
503:           PetscViewerASCIISynchronizedPrintf(viewer,"%D %D\n",i + st,idx[i]);
504:         }
505:       }
506:     }
507:     PetscViewerFlush(viewer);
508:     PetscViewerASCIIPopSynchronized(viewer);
509:   } else if (isbinary) {
510:     ISView_General_Binary(is,viewer);
511:   } else if (ishdf5) {
512: #if defined(PETSC_HAVE_HDF5)
513:     ISView_General_HDF5(is,viewer);
514: #endif
515:   }
516:   return(0);
517: }

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

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

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

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

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

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

562: PetscErrorCode  ISToGeneral_General(IS is)
563: {
565:   return(0);
566: }

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

588: PETSC_INTERN PetscErrorCode ISSetUp_General(IS);

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

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

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

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

625:    Collective on MPI_Comm

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

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

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


643:    Level: beginner

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

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

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

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

664:    Collective on IS

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

672:    Level: beginner

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

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

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

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

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

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

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

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

720: PETSC_EXTERN PetscErrorCode ISCreate_General(IS is)
721: {
723:   IS_General     *sub;

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