1*0f453b92SJed Brown #include <petscsys.h> 2*0f453b92SJed Brown 3*0f453b92SJed Brown /* Segmented (extendable) array implementation */ 4*0f453b92SJed Brown struct _n_PetscSegBuffer { 5*0f453b92SJed Brown PetscInt unitbytes; 6*0f453b92SJed Brown PetscInt alloc; 7*0f453b92SJed Brown PetscInt used; 8*0f453b92SJed Brown PetscInt tailused; 9*0f453b92SJed Brown PetscSegBuffer tail; 10*0f453b92SJed Brown union { /* Dummy types to ensure alignment */ 11*0f453b92SJed Brown PetscReal dummy_real; 12*0f453b92SJed Brown PetscInt dummy_int; 13*0f453b92SJed Brown char array[1]; 14*0f453b92SJed Brown } u; 15*0f453b92SJed Brown }; 16*0f453b92SJed Brown 17*0f453b92SJed Brown #undef __FUNCT__ 18*0f453b92SJed Brown #define __FUNCT__ "PetscSegBufferAlloc_Private" 19*0f453b92SJed Brown static PetscErrorCode PetscSegBufferAlloc_Private(PetscSegBuffer *seg,PetscInt count) 20*0f453b92SJed Brown { 21*0f453b92SJed Brown PetscErrorCode ierr; 22*0f453b92SJed Brown PetscSegBuffer newseg,s; 23*0f453b92SJed Brown PetscInt alloc; 24*0f453b92SJed Brown 25*0f453b92SJed Brown PetscFunctionBegin; 26*0f453b92SJed Brown s = *seg; 27*0f453b92SJed Brown /* Grow at least fast enough to hold next item, like Fibonacci otherwise (up to 1MB chunks) */ 28*0f453b92SJed Brown alloc = PetscMax(s->used+count,PetscMin(1000000/s->unitbytes+1,s->alloc+s->tailused)); 29*0f453b92SJed Brown ierr = PetscMalloc(offsetof(struct _n_PetscSegBuffer,u)+alloc*s->unitbytes,&newseg);CHKERRQ(ierr); 30*0f453b92SJed Brown ierr = PetscMemzero(newseg,offsetof(struct _n_PetscSegBuffer,u));CHKERRQ(ierr); 31*0f453b92SJed Brown 32*0f453b92SJed Brown newseg->unitbytes = s->unitbytes; 33*0f453b92SJed Brown newseg->tailused = s->used + s->tailused; 34*0f453b92SJed Brown newseg->tail = s; 35*0f453b92SJed Brown newseg->alloc = alloc; 36*0f453b92SJed Brown *seg = newseg; 37*0f453b92SJed Brown PetscFunctionReturn(0); 38*0f453b92SJed Brown } 39*0f453b92SJed Brown 40*0f453b92SJed Brown #undef __FUNCT__ 41*0f453b92SJed Brown #define __FUNCT__ "PetscSegBufferCreate" 42*0f453b92SJed Brown /*@C 43*0f453b92SJed Brown PetscSegBufferCreate - create segmented buffer 44*0f453b92SJed Brown 45*0f453b92SJed Brown Not Collective 46*0f453b92SJed Brown 47*0f453b92SJed Brown Input Arguments: 48*0f453b92SJed Brown + unitbytes - number of bytes that each entry will contain 49*0f453b92SJed Brown - expected - expected/typical number of entries 50*0f453b92SJed Brown 51*0f453b92SJed Brown Output Argument: 52*0f453b92SJed Brown . seg - segmented buffer object 53*0f453b92SJed Brown 54*0f453b92SJed Brown Level: developer 55*0f453b92SJed Brown 56*0f453b92SJed Brown .seealso: PetscSegBufferGet(), PetscSegBufferExtract(), PetscSegBufferDestroy() 57*0f453b92SJed Brown @*/ 58*0f453b92SJed Brown PetscErrorCode PetscSegBufferCreate(PetscInt unitbytes,PetscInt expected,PetscSegBuffer *seg) 59*0f453b92SJed Brown { 60*0f453b92SJed Brown PetscErrorCode ierr; 61*0f453b92SJed Brown 62*0f453b92SJed Brown PetscFunctionBegin; 63*0f453b92SJed Brown ierr = PetscMalloc(offsetof(struct _n_PetscSegBuffer,u)+expected*unitbytes,seg);CHKERRQ(ierr); 64*0f453b92SJed Brown ierr = PetscMemzero(*seg,offsetof(struct _n_PetscSegBuffer,u));CHKERRQ(ierr); 65*0f453b92SJed Brown 66*0f453b92SJed Brown (*seg)->unitbytes = unitbytes; 67*0f453b92SJed Brown (*seg)->alloc = expected; 68*0f453b92SJed Brown PetscFunctionReturn(0); 69*0f453b92SJed Brown } 70*0f453b92SJed Brown 71*0f453b92SJed Brown #undef __FUNCT__ 72*0f453b92SJed Brown #define __FUNCT__ "PetscSegBufferGet" 73*0f453b92SJed Brown /*@C 74*0f453b92SJed Brown PetscSegBufferGet - get new buffer space from a segmented buffer 75*0f453b92SJed Brown 76*0f453b92SJed Brown Not Collective 77*0f453b92SJed Brown 78*0f453b92SJed Brown Input Arguments: 79*0f453b92SJed Brown + seg - address of segmented buffer 80*0f453b92SJed Brown - count - number of entries needed 81*0f453b92SJed Brown 82*0f453b92SJed Brown Output Argument: 83*0f453b92SJed Brown . buf - address of new buffer for contiguous data 84*0f453b92SJed Brown 85*0f453b92SJed Brown Level: developer 86*0f453b92SJed Brown 87*0f453b92SJed Brown .seealso: PetscSegBufferCreate(), PetscSegBufferExtract(), PetscSegBufferDestroy() 88*0f453b92SJed Brown @*/ 89*0f453b92SJed Brown PetscErrorCode PetscSegBufferGet(PetscSegBuffer *seg,PetscInt count,void *buf) 90*0f453b92SJed Brown { 91*0f453b92SJed Brown PetscErrorCode ierr; 92*0f453b92SJed Brown PetscSegBuffer s; 93*0f453b92SJed Brown 94*0f453b92SJed Brown PetscFunctionBegin; 95*0f453b92SJed Brown s = *seg; 96*0f453b92SJed Brown if (PetscUnlikely(s->used + count > s->alloc)) {ierr = PetscSegBufferAlloc_Private(seg,count);CHKERRQ(ierr);} 97*0f453b92SJed Brown s = *seg; 98*0f453b92SJed Brown *(char**)buf = &s->u.array[s->used*s->unitbytes]; 99*0f453b92SJed Brown s->used += count; 100*0f453b92SJed Brown PetscFunctionReturn(0); 101*0f453b92SJed Brown } 102*0f453b92SJed Brown 103*0f453b92SJed Brown #undef __FUNCT__ 104*0f453b92SJed Brown #define __FUNCT__ "PetscSegBufferDestroy" 105*0f453b92SJed Brown /*@C 106*0f453b92SJed Brown PetscSegBufferDestroy - destroy segmented buffer 107*0f453b92SJed Brown 108*0f453b92SJed Brown Not Collective 109*0f453b92SJed Brown 110*0f453b92SJed Brown Input Arguments: 111*0f453b92SJed Brown . seg - address of segmented buffer object 112*0f453b92SJed Brown 113*0f453b92SJed Brown Level: developer 114*0f453b92SJed Brown 115*0f453b92SJed Brown .seealso: PetscSegBufferCreate() 116*0f453b92SJed Brown @*/ 117*0f453b92SJed Brown PetscErrorCode PetscSegBufferDestroy(PetscSegBuffer *seg) 118*0f453b92SJed Brown { 119*0f453b92SJed Brown PetscErrorCode ierr; 120*0f453b92SJed Brown PetscSegBuffer s; 121*0f453b92SJed Brown 122*0f453b92SJed Brown PetscFunctionBegin; 123*0f453b92SJed Brown for (s=*seg; s;) { 124*0f453b92SJed Brown PetscSegBuffer tail = s->tail; 125*0f453b92SJed Brown ierr = PetscFree(s);CHKERRQ(ierr); 126*0f453b92SJed Brown s = tail; 127*0f453b92SJed Brown } 128*0f453b92SJed Brown *seg = NULL; 129*0f453b92SJed Brown PetscFunctionReturn(0); 130*0f453b92SJed Brown } 131*0f453b92SJed Brown 132*0f453b92SJed Brown #undef __FUNCT__ 133*0f453b92SJed Brown #define __FUNCT__ "PetscSegBufferExtract" 134*0f453b92SJed Brown /*@C 135*0f453b92SJed Brown PetscSegBufferExtract - extract contiguous data and reset segmented buffer 136*0f453b92SJed Brown 137*0f453b92SJed Brown Not Collective 138*0f453b92SJed Brown 139*0f453b92SJed Brown Input Argument: 140*0f453b92SJed Brown . seg - segmented buffer 141*0f453b92SJed Brown 142*0f453b92SJed Brown Output Argument: 143*0f453b92SJed Brown . contiguous - address of new array containing contiguous data, caller frees with PetscFree() 144*0f453b92SJed Brown 145*0f453b92SJed Brown Level: developer 146*0f453b92SJed Brown 147*0f453b92SJed Brown .seealso: PetscSegBufferCreate(), PetscSegBufferGet(), PetscSegBufferDestroy() 148*0f453b92SJed Brown @*/ 149*0f453b92SJed Brown PetscErrorCode PetscSegBufferExtract(PetscSegBuffer *seg,void *contiguous) 150*0f453b92SJed Brown { 151*0f453b92SJed Brown PetscErrorCode ierr; 152*0f453b92SJed Brown PetscInt unitbytes; 153*0f453b92SJed Brown PetscSegBuffer s,t; 154*0f453b92SJed Brown char *contig,*ptr; 155*0f453b92SJed Brown 156*0f453b92SJed Brown PetscFunctionBegin; 157*0f453b92SJed Brown s = *seg; 158*0f453b92SJed Brown 159*0f453b92SJed Brown unitbytes = s->unitbytes; 160*0f453b92SJed Brown 161*0f453b92SJed Brown ierr = PetscMalloc((s->used+s->tailused)*unitbytes,&contig);CHKERRQ(ierr); 162*0f453b92SJed Brown ptr = contig + s->tailused*unitbytes; 163*0f453b92SJed Brown ierr = PetscMemcpy(ptr,s->u.array,s->used*unitbytes);CHKERRQ(ierr); 164*0f453b92SJed Brown for (t=s->tail; t;) { 165*0f453b92SJed Brown PetscSegBuffer tail = t->tail; 166*0f453b92SJed Brown ptr -= t->used*unitbytes; 167*0f453b92SJed Brown ierr = PetscMemcpy(ptr,t->u.array,t->used*unitbytes);CHKERRQ(ierr); 168*0f453b92SJed Brown ierr = PetscFree(t);CHKERRQ(ierr); 169*0f453b92SJed Brown t = tail; 170*0f453b92SJed Brown } 171*0f453b92SJed Brown if (ptr != contig) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Tail count does not match"); 172*0f453b92SJed Brown s->tailused = 0; 173*0f453b92SJed Brown s->tail = NULL; 174*0f453b92SJed Brown *(char**)contiguous = contig; 175*0f453b92SJed Brown PetscFunctionReturn(0); 176*0f453b92SJed Brown } 177