Commits

Matt Knepley committed 5bfe0a2 Merge

Merged in jed/plex-scalable-partition-closure (pull request #11)

DMPlexCreatePartitionClosure: replace quadratic algorithm

Comments (0)

Files changed (5)

include/petscsys.h

 PETSC_EXTERN PetscErrorCode PetscSubcommSetType(PetscSubcomm,PetscSubcommType);
 PETSC_EXTERN PetscErrorCode PetscSubcommSetTypeGeneral(PetscSubcomm,PetscMPIInt,PetscMPIInt,PetscMPIInt);
 
+/*S
+   PetscSegBuffer - a segmented extendable buffer
+
+   Level: developer
+
+.seealso: PetscSegBufferCreate(), PetscSegBufferGet(), PetscSegBufferExtract(), PetscSegBufferDestroy()
+S*/
+typedef struct _n_PetscSegBuffer *PetscSegBuffer;
+PETSC_EXTERN PetscErrorCode PetscSegBufferCreate(PetscInt,PetscInt,PetscSegBuffer*);
+PETSC_EXTERN PetscErrorCode PetscSegBufferDestroy(PetscSegBuffer*);
+PETSC_EXTERN PetscErrorCode PetscSegBufferGet(PetscSegBuffer*,PetscInt,void*);
+PETSC_EXTERN PetscErrorCode PetscSegBufferExtractAlloc(PetscSegBuffer*,void*);
+PETSC_EXTERN PetscErrorCode PetscSegBufferExtractTo(PetscSegBuffer*,void*);
+PETSC_EXTERN PetscErrorCode PetscSegBufferExtractInPlace(PetscSegBuffer*,void*);
+
 /* Reset __FUNCT__ in case the user does not define it themselves */
 #undef __FUNCT__
 #define __FUNCT__ "User provided function"

src/dm/impls/plex/plex.c

 {
   /* const PetscInt  height = 0; */
   const PetscInt *partArray;
-  PetscInt       *allPoints, *partPoints = NULL;
-  PetscInt        rStart, rEnd, rank, maxPartSize = 0, newSize;
+  PetscInt       *allPoints, *packPoints;
+  PetscInt        rStart, rEnd, rank, pStart, pEnd, newSize;
   PetscErrorCode  ierr;
+  PetscBT         bt;
+  PetscSegBuffer  segpack,segpart;
 
   PetscFunctionBegin;
   ierr = PetscSectionGetChart(pointSection, &rStart, &rEnd);CHKERRQ(ierr);
   ierr = ISGetIndices(pointPartition, &partArray);CHKERRQ(ierr);
   ierr = PetscSectionCreate(PetscObjectComm((PetscObject)dm), section);CHKERRQ(ierr);
   ierr = PetscSectionSetChart(*section, rStart, rEnd);CHKERRQ(ierr);
+  ierr = DMPlexGetChart(dm,&pStart,&pEnd);CHKERRQ(ierr);
+  ierr = PetscBTCreate(pEnd-pStart,&bt);CHKERRQ(ierr);
+  ierr = PetscSegBufferCreate(sizeof(PetscInt),1000,&segpack);CHKERRQ(ierr);
+  ierr = PetscSegBufferCreate(sizeof(PetscInt),1000,&segpart);CHKERRQ(ierr);
   for (rank = rStart; rank < rEnd; ++rank) {
     PetscInt partSize = 0;
     PetscInt numPoints, offset, p;
 
       /* TODO Include support for height > 0 case */
       ierr = DMPlexGetTransitiveClosure(dm, point, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr);
-      /* Merge into existing points */
-      if (partSize+closureSize > maxPartSize) {
-        PetscInt *tmpPoints;
-
-        maxPartSize = PetscMax(partSize+closureSize, 2*maxPartSize);
-        ierr = PetscMalloc(maxPartSize * sizeof(PetscInt), &tmpPoints);CHKERRQ(ierr);
-        ierr = PetscMemcpy(tmpPoints, partPoints, partSize * sizeof(PetscInt));CHKERRQ(ierr);
-        ierr = PetscFree(partPoints);CHKERRQ(ierr);
-
-        partPoints = tmpPoints;
+      for (c=0; c<closureSize; c++) {
+        PetscInt cpoint = closure[c*2];
+        if (!PetscBTLookupSet(bt,cpoint-pStart)) {
+          PetscInt *pt;
+          partSize++;
+          ierr = PetscSegBufferGet(&segpart,1,&pt);CHKERRQ(ierr);
+          *pt = cpoint;
+        }
       }
-      for (c = 0; c < closureSize; ++c) partPoints[partSize+c] = closure[c*2];
-      partSize += closureSize;
-
-      ierr = PetscSortRemoveDupsInt(&partSize, partPoints);CHKERRQ(ierr);
       ierr = DMPlexRestoreTransitiveClosure(dm, point, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr);
     }
     ierr = PetscSectionSetDof(*section, rank, partSize);CHKERRQ(ierr);
+    ierr = PetscSegBufferGet(&segpack,partSize,&packPoints);CHKERRQ(ierr);
+    ierr = PetscSegBufferExtractTo(&segpart,packPoints);CHKERRQ(ierr);
+    ierr = PetscSortInt(partSize,packPoints);CHKERRQ(ierr);
+    for (p=0; p<partSize; p++) {ierr = PetscBTClear(bt,packPoints[p]-pStart);CHKERRQ(ierr);}
   }
+  ierr = PetscBTDestroy(&bt);CHKERRQ(ierr);
+  ierr = PetscSegBufferDestroy(&segpart);CHKERRQ(ierr);
+
   ierr = PetscSectionSetUp(*section);CHKERRQ(ierr);
   ierr = PetscSectionGetStorageSize(*section, &newSize);CHKERRQ(ierr);
   ierr = PetscMalloc(newSize * sizeof(PetscInt), &allPoints);CHKERRQ(ierr);
 
+  ierr = PetscSegBufferExtractInPlace(&segpack,&packPoints);CHKERRQ(ierr);
   for (rank = rStart; rank < rEnd; ++rank) {
-    PetscInt partSize = 0, newOffset;
-    PetscInt numPoints, offset, p;
-
-    ierr = PetscSectionGetDof(pointSection, rank, &numPoints);CHKERRQ(ierr);
-    ierr = PetscSectionGetOffset(pointSection, rank, &offset);CHKERRQ(ierr);
-    for (p = 0; p < numPoints; ++p) {
-      PetscInt  point   = partArray[offset+p], closureSize, c;
-      PetscInt *closure = NULL;
+    PetscInt numPoints, offset;
 
-      /* TODO Include support for height > 0 case */
-      ierr = DMPlexGetTransitiveClosure(dm, point, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr);
-      /* Merge into existing points */
-      for (c = 0; c < closureSize; ++c) partPoints[partSize+c] = closure[c*2];
-      partSize += closureSize;
-
-      ierr = PetscSortRemoveDupsInt(&partSize, partPoints);CHKERRQ(ierr);
-      ierr = DMPlexRestoreTransitiveClosure(dm, point, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr);
-    }
-    ierr = PetscSectionGetOffset(*section, rank, &newOffset);CHKERRQ(ierr);
-    ierr = PetscMemcpy(&allPoints[newOffset], partPoints, partSize * sizeof(PetscInt));CHKERRQ(ierr);
+    ierr = PetscSectionGetDof(*section, rank, &numPoints);CHKERRQ(ierr);
+    ierr = PetscSectionGetOffset(*section, rank, &offset);CHKERRQ(ierr);
+    ierr = PetscMemcpy(&allPoints[offset], packPoints, numPoints * sizeof(PetscInt));CHKERRQ(ierr);
+    packPoints += numPoints;
   }
+
+  ierr = PetscSegBufferDestroy(&segpack);CHKERRQ(ierr);
   ierr = ISRestoreIndices(pointPartition, &partArray);CHKERRQ(ierr);
-  ierr = PetscFree(partPoints);CHKERRQ(ierr);
   ierr = ISCreateGeneral(PetscObjectComm((PetscObject)dm), newSize, allPoints, PETSC_OWN_POINTER, partition);CHKERRQ(ierr);
   PetscFunctionReturn(0);
 }

src/sys/utils/makefile

 SOURCEC	  = arch.c fhost.c fuser.c memc.c mpiu.c psleep.c sortd.c sorti.c \
             str.c sortip.c pbarrier.c pdisplay.c ctable.c psplit.c \
             select.c mpimesg.c sseenabled.c mpitr.c  mpilong.c mathinf.c \
-            mpits.c
+            mpits.c segbuffer.c
 SOURCEF	  =
 SOURCEH	  = ../../../include/petscctable.h
 MANSEC	  = Sys

src/sys/utils/mpits.c

 }
 
 #if defined(PETSC_HAVE_MPI_IBARRIER)
-/* Segmented (extendable) array implementation */
-typedef struct _SegArray *SegArray;
-struct _SegArray {
-  PetscInt unitbytes;
-  PetscInt alloc;
-  PetscInt used;
-  PetscInt tailused;
-  SegArray tail;
-  union {                       /* Dummy types to ensure alignment */
-    PetscReal dummy_real;
-    PetscInt  dummy_int;
-    char      array[1];
-  } u;
-};
-
-#undef __FUNCT__
-#define __FUNCT__ "SegArrayCreate"
-static PetscErrorCode SegArrayCreate(PetscInt unitbytes,PetscInt expected,SegArray *seg)
-{
-  PetscErrorCode ierr;
-
-  PetscFunctionBegin;
-  ierr = PetscMalloc(offsetof(struct _SegArray,u)+expected*unitbytes,seg);CHKERRQ(ierr);
-  ierr = PetscMemzero(*seg,offsetof(struct _SegArray,u));CHKERRQ(ierr);
-
-  (*seg)->unitbytes = unitbytes;
-  (*seg)->alloc     = expected;
-  PetscFunctionReturn(0);
-}
-
-#undef __FUNCT__
-#define __FUNCT__ "SegArrayAlloc_Private"
-static PetscErrorCode SegArrayAlloc_Private(SegArray *seg,PetscInt count)
-{
-  PetscErrorCode ierr;
-  SegArray       newseg,s;
-  PetscInt       alloc;
-
-  PetscFunctionBegin;
-  s = *seg;
-  /* Grow at least fast enough to hold next item, like Fibonacci otherwise (up to 1MB chunks) */
-  alloc = PetscMax(s->used+count,PetscMin(1000000/s->unitbytes+1,s->alloc+s->tailused));
-  ierr  = PetscMalloc(offsetof(struct _SegArray,u)+alloc*s->unitbytes,&newseg);CHKERRQ(ierr);
-  ierr  = PetscMemzero(newseg,offsetof(struct _SegArray,u));CHKERRQ(ierr);
-
-  newseg->unitbytes = s->unitbytes;
-  newseg->tailused  = s->used + s->tailused;
-  newseg->tail      = s;
-  newseg->alloc     = alloc;
-  *seg              = newseg;
-  PetscFunctionReturn(0);
-}
-
-#undef __FUNCT__
-#define __FUNCT__ "SegArrayGet"
-static PetscErrorCode SegArrayGet(SegArray *seg,PetscInt count,void *array)
-{
-  PetscErrorCode ierr;
-  SegArray       s;
-
-  PetscFunctionBegin;
-  s = *seg;
-  if (PetscUnlikely(s->used + count > s->alloc)) {ierr = SegArrayAlloc_Private(seg,count);CHKERRQ(ierr);}
-  s = *seg;
-  *(char**)array = &s->u.array[s->used*s->unitbytes];
-  s->used += count;
-  PetscFunctionReturn(0);
-}
-
-#undef __FUNCT__
-#define __FUNCT__ "SegArrayDestroy"
-static PetscErrorCode SegArrayDestroy(SegArray *seg)
-{
-  PetscErrorCode ierr;
-  SegArray       s;
-
-  PetscFunctionBegin;
-  for (s=*seg; s;) {
-    SegArray tail = s->tail;
-    ierr = PetscFree(s);CHKERRQ(ierr);
-    s = tail;
-  }
-  *seg = NULL;
-  PetscFunctionReturn(0);
-}
-
-#undef __FUNCT__
-#define __FUNCT__ "SegArrayExtract"
-/* Extracts contiguous data and resets segarray */
-static PetscErrorCode SegArrayExtract(SegArray *seg,void *contiguous)
-{
-  PetscErrorCode ierr;
-  PetscInt       unitbytes;
-  SegArray       s,t;
-  char           *contig,*ptr;
-
-  PetscFunctionBegin;
-  s = *seg;
-
-  unitbytes = s->unitbytes;
-
-  ierr = PetscMalloc((s->used+s->tailused)*unitbytes,&contig);CHKERRQ(ierr);
-  ptr  = contig + s->tailused*unitbytes;
-  ierr = PetscMemcpy(ptr,s->u.array,s->used*unitbytes);CHKERRQ(ierr);
-  for (t=s->tail; t;) {
-    SegArray tail = t->tail;
-    ptr -= t->used*unitbytes;
-    ierr = PetscMemcpy(ptr,t->u.array,t->used*unitbytes);CHKERRQ(ierr);
-    ierr = PetscFree(t);CHKERRQ(ierr);
-    t    = tail;
-  }
-  if (ptr != contig) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Tail count does not match");
-  s->tailused         = 0;
-  s->tail             = NULL;
-  *(char**)contiguous = contig;
-  PetscFunctionReturn(0);
-}
 
 #undef __FUNCT__
 #define __FUNCT__ "PetscCommBuildTwoSided_Ibarrier"
   PetscInt       i;
   char           *tdata;
   MPI_Request    *sendreqs,barrier;
-  SegArray       segrank,segdata;
+  PetscSegBuffer segrank,segdata;
 
   PetscFunctionBegin;
   ierr  = PetscCommGetNewTag(comm,&tag);CHKERRQ(ierr);
   for (i=0; i<nto; i++) {
     ierr = MPI_Issend((void*)(tdata+count*unitbytes*i),count,dtype,toranks[i],tag,comm,sendreqs+i);CHKERRQ(ierr);
   }
-  ierr = SegArrayCreate(sizeof(PetscMPIInt),4,&segrank);CHKERRQ(ierr);
-  ierr = SegArrayCreate(unitbytes,4*count,&segdata);CHKERRQ(ierr);
+  ierr = PetscSegBufferCreate(sizeof(PetscMPIInt),4,&segrank);CHKERRQ(ierr);
+  ierr = PetscSegBufferCreate(unitbytes,4*count,&segdata);CHKERRQ(ierr);
 
   nrecvs  = 0;
   barrier = MPI_REQUEST_NULL;
     if (flag) {                 /* incoming message */
       PetscMPIInt *recvrank;
       void        *buf;
-      ierr      = SegArrayGet(&segrank,1,&recvrank);CHKERRQ(ierr);
-      ierr      = SegArrayGet(&segdata,count,&buf);CHKERRQ(ierr);
+      ierr      = PetscSegBufferGet(&segrank,1,&recvrank);CHKERRQ(ierr);
+      ierr      = PetscSegBufferGet(&segdata,count,&buf);CHKERRQ(ierr);
       *recvrank = status.MPI_SOURCE;
       ierr      = MPI_Recv(buf,count,dtype,status.MPI_SOURCE,tag,comm,MPI_STATUS_IGNORE);CHKERRQ(ierr);
       nrecvs++;
     }
   }
   *nfrom = nrecvs;
-  ierr   = SegArrayExtract(&segrank,fromranks);CHKERRQ(ierr);
-  ierr   = SegArrayDestroy(&segrank);CHKERRQ(ierr);
-  ierr   = SegArrayExtract(&segdata,fromdata);CHKERRQ(ierr);
-  ierr   = SegArrayDestroy(&segdata);CHKERRQ(ierr);
+  ierr   = PetscSegBufferExtractAlloc(&segrank,fromranks);CHKERRQ(ierr);
+  ierr   = PetscSegBufferDestroy(&segrank);CHKERRQ(ierr);
+  ierr   = PetscSegBufferExtractAlloc(&segdata,fromdata);CHKERRQ(ierr);
+  ierr   = PetscSegBufferDestroy(&segdata);CHKERRQ(ierr);
   PetscFunctionReturn(0);
 }
 #endif

src/sys/utils/segbuffer.c

+#include <petscsys.h>
+
+/* Segmented (extendable) array implementation */
+struct _n_PetscSegBuffer {
+  PetscInt unitbytes;
+  PetscInt alloc;
+  PetscInt used;
+  PetscInt tailused;
+  PetscSegBuffer tail;
+  union {                       /* Dummy types to ensure alignment */
+    PetscReal dummy_real;
+    PetscInt  dummy_int;
+    char      array[1];
+  } u;
+};
+
+#undef __FUNCT__
+#define __FUNCT__ "PetscSegBufferAlloc_Private"
+static PetscErrorCode PetscSegBufferAlloc_Private(PetscSegBuffer *seg,PetscInt count)
+{
+  PetscErrorCode ierr;
+  PetscSegBuffer newseg,s;
+  PetscInt       alloc;
+
+  PetscFunctionBegin;
+  s = *seg;
+  /* Grow at least fast enough to hold next item, like Fibonacci otherwise (up to 1MB chunks) */
+  alloc = PetscMax(s->used+count,PetscMin(1000000/s->unitbytes+1,s->alloc+s->tailused));
+  ierr  = PetscMalloc(offsetof(struct _n_PetscSegBuffer,u)+alloc*s->unitbytes,&newseg);CHKERRQ(ierr);
+  ierr  = PetscMemzero(newseg,offsetof(struct _n_PetscSegBuffer,u));CHKERRQ(ierr);
+
+  newseg->unitbytes = s->unitbytes;
+  newseg->tailused  = s->used + s->tailused;
+  newseg->tail      = s;
+  newseg->alloc     = alloc;
+  *seg              = newseg;
+  PetscFunctionReturn(0);
+}
+
+#undef __FUNCT__
+#define __FUNCT__ "PetscSegBufferCreate"
+/*@C
+   PetscSegBufferCreate - create segmented buffer
+
+   Not Collective
+
+   Input Arguments:
++  unitbytes - number of bytes that each entry will contain
+-  expected - expected/typical number of entries
+
+   Output Argument:
+.  seg - segmented buffer object
+
+   Level: developer
+
+.seealso: PetscSegBufferGet(), PetscSegBufferExtractAlloc(), PetscSegBufferExtractTo(), PetscSegBufferExtractInPlace(), PetscSegBufferDestroy()
+@*/
+PetscErrorCode PetscSegBufferCreate(PetscInt unitbytes,PetscInt expected,PetscSegBuffer *seg)
+{
+  PetscErrorCode ierr;
+
+  PetscFunctionBegin;
+  ierr = PetscMalloc(offsetof(struct _n_PetscSegBuffer,u)+expected*unitbytes,seg);CHKERRQ(ierr);
+  ierr = PetscMemzero(*seg,offsetof(struct _n_PetscSegBuffer,u));CHKERRQ(ierr);
+
+  (*seg)->unitbytes = unitbytes;
+  (*seg)->alloc     = expected;
+  PetscFunctionReturn(0);
+}
+
+#undef __FUNCT__
+#define __FUNCT__ "PetscSegBufferGet"
+/*@C
+   PetscSegBufferGet - get new buffer space from a segmented buffer
+
+   Not Collective
+
+   Input Arguments:
++  seg - address of segmented buffer
+-  count - number of entries needed
+
+   Output Argument:
+.  buf - address of new buffer for contiguous data
+
+   Level: developer
+
+.seealso: PetscSegBufferCreate(), PetscSegBufferExtractAlloc(), PetscSegBufferExtractTo(), PetscSegBufferExtractInPlace(), PetscSegBufferDestroy()
+@*/
+PetscErrorCode PetscSegBufferGet(PetscSegBuffer *seg,PetscInt count,void *buf)
+{
+  PetscErrorCode ierr;
+  PetscSegBuffer s;
+
+  PetscFunctionBegin;
+  s = *seg;
+  if (PetscUnlikely(s->used + count > s->alloc)) {ierr = PetscSegBufferAlloc_Private(seg,count);CHKERRQ(ierr);}
+  s = *seg;
+  *(char**)buf = &s->u.array[s->used*s->unitbytes];
+  s->used += count;
+  PetscFunctionReturn(0);
+}
+
+#undef __FUNCT__
+#define __FUNCT__ "PetscSegBufferDestroy"
+/*@C
+   PetscSegBufferDestroy - destroy segmented buffer
+
+   Not Collective
+
+   Input Arguments:
+.  seg - address of segmented buffer object
+
+   Level: developer
+
+.seealso: PetscSegBufferCreate()
+@*/
+PetscErrorCode PetscSegBufferDestroy(PetscSegBuffer *seg)
+{
+  PetscErrorCode ierr;
+  PetscSegBuffer s;
+
+  PetscFunctionBegin;
+  for (s=*seg; s;) {
+    PetscSegBuffer tail = s->tail;
+    ierr = PetscFree(s);CHKERRQ(ierr);
+    s = tail;
+  }
+  *seg = NULL;
+  PetscFunctionReturn(0);
+}
+
+#undef __FUNCT__
+#define __FUNCT__ "PetscSegBufferExtractTo"
+/*@C
+   PetscSegBufferExtractTo - extract contiguous data to provided buffer and reset segmented buffer
+
+   Not Collective
+
+   Input Argument:
++  seg - segmented buffer
+-  contig - allocated buffer to hold contiguous data
+
+   Level: developer
+
+.seealso: PetscSegBufferCreate(), PetscSegBufferGet(), PetscSegBufferDestroy(), PetscSegBufferExtractAlloc(), PetscSegBufferExtractInPlace()
+@*/
+PetscErrorCode PetscSegBufferExtractTo(PetscSegBuffer *seg,void *contig)
+{
+  PetscErrorCode ierr;
+  PetscInt       unitbytes;
+  PetscSegBuffer s,t;
+  char           *ptr;
+
+  PetscFunctionBegin;
+  s = *seg;
+
+  unitbytes = s->unitbytes;
+
+  ptr  = ((char*)contig) + s->tailused*unitbytes;
+  ierr = PetscMemcpy(ptr,s->u.array,s->used*unitbytes);CHKERRQ(ierr);
+  for (t=s->tail; t;) {
+    PetscSegBuffer tail = t->tail;
+    ptr -= t->used*unitbytes;
+    ierr = PetscMemcpy(ptr,t->u.array,t->used*unitbytes);CHKERRQ(ierr);
+    ierr = PetscFree(t);CHKERRQ(ierr);
+    t    = tail;
+  }
+  if (ptr != contig) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Tail count does not match");
+  s->used             = 0;
+  s->tailused         = 0;
+  s->tail             = NULL;
+  PetscFunctionReturn(0);
+}
+
+#undef __FUNCT__
+#define __FUNCT__ "PetscSegBufferExtractAlloc"
+/*@C
+   PetscSegBufferExtractAlloc - extract contiguous data to new allocation and reset segmented buffer
+
+   Not Collective
+
+   Input Argument:
+.  seg - segmented buffer
+
+   Output Argument:
+.  contiguous - address of new array containing contiguous data, caller frees with PetscFree()
+
+   Level: developer
+
+   Developer Notes: 'seg' argument is a pointer so that implementation could reallocate, though this is not currently done
+
+.seealso: PetscSegBufferCreate(), PetscSegBufferGet(), PetscSegBufferDestroy(), PetscSegBufferExtractTo(), PetscSegBufferExtractInPlace()
+@*/
+PetscErrorCode PetscSegBufferExtractAlloc(PetscSegBuffer *seg,void *contiguous)
+{
+  PetscErrorCode ierr;
+  PetscSegBuffer s;
+  void           *contig;
+
+  PetscFunctionBegin;
+  s = *seg;
+
+  ierr = PetscMalloc((s->used+s->tailused)*s->unitbytes,&contig);CHKERRQ(ierr);
+  ierr = PetscSegBufferExtractTo(seg,contig);CHKERRQ(ierr);
+  *(void**)contiguous = contig;
+  PetscFunctionReturn(0);
+}
+
+#undef __FUNCT__
+#define __FUNCT__ "PetscSegBufferExtractInPlace"
+/*@C
+   PetscSegBufferExtractInPlace - extract in-place contiguous representation of data and reset segmented buffer for reuse
+
+   Collective
+
+   Input Arguments:
+.  seg - segmented buffer object
+
+   Output Arguments:
+.  contig - address of pointer to contiguous memory
+
+   Level: developer
+
+.seealso: PetscSegBufferExtractAlloc(), PetscSegBufferExtractTo()
+@*/
+PetscErrorCode PetscSegBufferExtractInPlace(PetscSegBuffer *seg,void *contig)
+{
+  PetscErrorCode ierr;
+
+  PetscFunctionBegin;
+  if (!(*seg)->tail) {
+    *(char**)contig = (*seg)->u.array;
+  } else {
+    PetscSegBuffer s = *seg,newseg;
+
+    ierr = PetscSegBufferCreate(s->unitbytes,s->used+s->tailused,&newseg);CHKERRQ(ierr);
+    ierr = PetscSegBufferExtractTo(seg,newseg->u.array);CHKERRQ(ierr);
+    ierr = PetscSegBufferDestroy(seg);CHKERRQ(ierr);
+    *seg = newseg;
+    *(void**)contig = newseg->u.array;
+  }
+  PetscFunctionReturn(0);
+}