1 2 #include "data_bucket.h" 3 4 /* string helpers */ 5 #undef __FUNCT__ 6 #define __FUNCT__ "StringInList" 7 PetscErrorCode StringInList(const char name[],const PetscInt N,const DataField gfield[],PetscBool *val) 8 { 9 PetscInt i; 10 11 *val = PETSC_FALSE; 12 for (i=0; i<N; i++) { 13 if (strcmp( name, gfield[i]->name) == 0 ) { 14 *val = PETSC_TRUE; 15 PetscFunctionReturn(0); 16 } 17 } 18 PetscFunctionReturn(0); 19 } 20 21 #undef __FUNCT__ 22 #define __FUNCT__ "StringFindInList" 23 PetscErrorCode StringFindInList(const char name[],const PetscInt N,const DataField gfield[],PetscInt *index) 24 { 25 PetscInt i; 26 27 *index = -1; 28 for (i=0; i<N; i++) { 29 if (strcmp( name, gfield[i]->name ) == 0) { 30 *index = i; 31 PetscFunctionReturn(0); 32 } 33 } 34 PetscFunctionReturn(0); 35 } 36 37 #undef __FUNCT__ 38 #define __FUNCT__ "DataFieldCreate" 39 PetscErrorCode DataFieldCreate(const char registeration_function[],const char name[],const size_t size,const PetscInt L,DataField *DF) 40 { 41 DataField df; 42 43 df = malloc( sizeof(struct _p_DataField) ); 44 memset( df, 0, sizeof(struct _p_DataField) ); 45 46 47 asprintf( &df->registeration_function, "%s", registeration_function ); 48 asprintf( &df->name, "%s", name ); 49 df->atomic_size = size; 50 df->L = L; 51 df->bs = 1; 52 53 df->data = malloc( size * L ); /* allocate something so we don't have to reallocate */ 54 memset( df->data, 0, size * L ); 55 56 *DF = df; 57 PetscFunctionReturn(0); 58 } 59 60 #undef __FUNCT__ 61 #define __FUNCT__ "DataFieldDestroy" 62 PetscErrorCode DataFieldDestroy(DataField *DF) 63 { 64 DataField df = *DF; 65 66 free(df->registeration_function); 67 free(df->name); 68 free(df->data); 69 free(df); 70 71 *DF = NULL; 72 PetscFunctionReturn(0); 73 } 74 75 /* data bucket */ 76 #undef __FUNCT__ 77 #define __FUNCT__ "DataBucketCreate" 78 PetscErrorCode DataBucketCreate(DataBucket *DB) 79 { 80 DataBucket db; 81 82 83 db = malloc( sizeof(struct _p_DataBucket) ); 84 memset( db, 0, sizeof(struct _p_DataBucket) ); 85 86 db->finalised = PETSC_FALSE; 87 88 /* create empty spaces for fields */ 89 db->L = -1; 90 db->buffer = 1; 91 db->allocated = 1; 92 93 db->nfields = 0; 94 db->field = malloc(sizeof(DataField)); 95 96 *DB = db; 97 PetscFunctionReturn(0); 98 } 99 100 #undef __FUNCT__ 101 #define __FUNCT__ "DataBucketDestroy" 102 PetscErrorCode DataBucketDestroy(DataBucket *DB) 103 { 104 DataBucket db = *DB; 105 PetscInt f; 106 PetscErrorCode ierr; 107 108 /* release fields */ 109 for (f=0; f<db->nfields; f++) { 110 ierr = DataFieldDestroy(&db->field[f]);CHKERRQ(ierr); 111 } 112 113 /* this will catch the initially allocated objects in the event that no fields are registered */ 114 if (db->field != NULL) { 115 free(db->field); 116 } 117 118 free(db); 119 120 *DB = NULL; 121 PetscFunctionReturn(0); 122 } 123 124 #undef __FUNCT__ 125 #define __FUNCT__ "DataBucketQueryForActiveFields" 126 PetscErrorCode DataBucketQueryForActiveFields(DataBucket db,PetscBool *any_active_fields) 127 { 128 PetscInt f; 129 *any_active_fields = PETSC_FALSE; 130 for (f=0; f<db->nfields; f++) { 131 if (db->field[f]->active) { 132 *any_active_fields = PETSC_TRUE; 133 PetscFunctionReturn(0); 134 } 135 } 136 PetscFunctionReturn(0); 137 } 138 139 #undef __FUNCT__ 140 #define __FUNCT__ "DataBucketRegisterField" 141 PetscErrorCode DataBucketRegisterField( 142 DataBucket db, 143 const char registeration_function[], 144 const char field_name[], 145 size_t atomic_size, DataField *_gfield) 146 { 147 PetscBool val; 148 DataField *field,fp; 149 PetscErrorCode ierr; 150 151 /* check we haven't finalised the registration of fields */ 152 /* 153 if(db->finalised==PETSC_TRUE) { 154 printf("ERROR: DataBucketFinalize() has been called. Cannot register more fields\n"); 155 ERROR(); 156 } 157 */ 158 159 /* check for repeated name */ 160 ierr = StringInList( field_name, db->nfields, (const DataField*)db->field, &val );CHKERRQ(ierr); 161 if (val == PETSC_TRUE) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Field %s already exists. Cannot add same field twice",field_name); 162 163 /* create new space for data */ 164 field = realloc( db->field, sizeof(DataField)*(db->nfields+1)); 165 db->field = field; 166 167 /* add field */ 168 ierr = DataFieldCreate( registeration_function, field_name, atomic_size, db->allocated, &fp );CHKERRQ(ierr); 169 db->field[ db->nfields ] = fp; 170 171 db->nfields++; 172 173 if (_gfield != NULL) { 174 *_gfield = fp; 175 } 176 PetscFunctionReturn(0); 177 } 178 179 /* 180 #define DataBucketRegisterField(db,name,size,k) {\ 181 char *location;\ 182 asprintf(&location,"Registered by %s() at line %d within file %s", __FUNCTION__, __LINE__, __FILE__);\ 183 _DataBucketRegisterField( (db), location, (name), (size), (k) );\ 184 free(location);\ 185 } 186 */ 187 188 #undef __FUNCT__ 189 #define __FUNCT__ "DataBucketGetDataFieldByName" 190 PetscErrorCode DataBucketGetDataFieldByName(DataBucket db,const char name[],DataField *gfield) 191 { 192 PetscErrorCode ierr; 193 PetscInt idx; 194 PetscBool found; 195 196 ierr = StringInList(name,db->nfields,(const DataField*)db->field,&found);CHKERRQ(ierr); 197 if (!found) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Cannot find DataField with name %s",name); 198 199 ierr = StringFindInList(name,db->nfields,(const DataField*)db->field,&idx);CHKERRQ(ierr); 200 201 *gfield = db->field[idx]; 202 PetscFunctionReturn(0); 203 } 204 205 #undef __FUNCT__ 206 #define __FUNCT__ "DataBucketQueryDataFieldByName" 207 PetscErrorCode DataBucketQueryDataFieldByName(DataBucket db,const char name[],PetscBool *found) 208 { 209 PetscErrorCode ierr; 210 *found = PETSC_FALSE; 211 ierr = StringInList(name,db->nfields,(const DataField*)db->field,found);CHKERRQ(ierr); 212 PetscFunctionReturn(0); 213 } 214 215 #undef __FUNCT__ 216 #define __FUNCT__ "DataBucketFinalize" 217 PetscErrorCode DataBucketFinalize(DataBucket db) 218 { 219 db->finalised = PETSC_TRUE; 220 PetscFunctionReturn(0); 221 } 222 223 #undef __FUNCT__ 224 #define __FUNCT__ "DataFieldGetNumEntries" 225 PetscErrorCode DataFieldGetNumEntries(DataField df,PetscInt *sum) 226 { 227 *sum = df->L; 228 PetscFunctionReturn(0); 229 } 230 231 #undef __FUNCT__ 232 #define __FUNCT__ "DataFieldSetBlockSize" 233 PetscErrorCode DataFieldSetBlockSize(DataField df,PetscInt blocksize) 234 { 235 df->bs = blocksize; 236 PetscFunctionReturn(0); 237 } 238 239 #undef __FUNCT__ 240 #define __FUNCT__ "DataFieldSetSize" 241 PetscErrorCode DataFieldSetSize(DataField df,const PetscInt new_L) 242 { 243 void *tmp_data; 244 245 if (new_L <= 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Cannot set size of DataField to be <= 0"); 246 247 if (new_L == df->L) PetscFunctionReturn(0); 248 249 if (new_L > df->L) { 250 251 tmp_data = realloc( df->data, df->atomic_size * (new_L) ); 252 df->data = tmp_data; 253 254 /* init new contents */ 255 memset( ( ((char*)df->data)+df->L*df->atomic_size), 0, (new_L-df->L)*df->atomic_size ); 256 257 } else { 258 /* reallocate pointer list, add +1 in case new_L = 0 */ 259 tmp_data = realloc( df->data, df->atomic_size * (new_L+1) ); 260 df->data = tmp_data; 261 } 262 263 df->L = new_L; 264 PetscFunctionReturn(0); 265 } 266 267 #undef __FUNCT__ 268 #define __FUNCT__ "DataFieldZeroBlock" 269 PetscErrorCode DataFieldZeroBlock(DataField df,const PetscInt start,const PetscInt end) 270 { 271 if (start > end) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_USER,"Cannot zero a block of entries if start(%D) > end(%D)",start,end); 272 273 if (start < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Cannot zero a block of entries if start(%D) < 0",start); 274 275 if (end > df->L) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_USER,"Cannot zero a block of entries if end(%D) >= array size(%D)",end,df->L); 276 277 memset( ( ((char*)df->data)+start*df->atomic_size), 0, (end-start)*df->atomic_size ); 278 PetscFunctionReturn(0); 279 } 280 281 /* 282 A negative buffer value will simply be ignored and the old buffer value will be used. 283 */ 284 #undef __FUNCT__ 285 #define __FUNCT__ "DataBucketSetSizes" 286 PetscErrorCode DataBucketSetSizes(DataBucket db,const PetscInt L,const PetscInt buffer) 287 { 288 PetscInt current_allocated,new_used,new_unused,new_buffer,new_allocated,f; 289 PetscBool any_active_fields; 290 PetscErrorCode ierr; 291 292 if (db->finalised == PETSC_FALSE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"You must call DataBucketFinalize() before DataBucketSetSizes()"); 293 294 ierr = DataBucketQueryForActiveFields(db,&any_active_fields);CHKERRQ(ierr); 295 if (any_active_fields) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Cannot safely re-size as at least one DataField is currently being accessed"); 296 297 current_allocated = db->allocated; 298 299 new_used = L; 300 new_unused = current_allocated - new_used; 301 new_buffer = db->buffer; 302 if (buffer >= 0) { /* update the buffer value */ 303 new_buffer = buffer; 304 } 305 new_allocated = new_used + new_buffer; 306 307 /* action */ 308 if (new_allocated > current_allocated) { 309 /* increase size to new_used + new_buffer */ 310 for (f=0; f<db->nfields; f++) { 311 ierr = DataFieldSetSize( db->field[f], new_allocated );CHKERRQ(ierr); 312 } 313 314 db->L = new_used; 315 db->buffer = new_buffer; 316 db->allocated = new_used + new_buffer; 317 } 318 else { 319 if (new_unused > 2 * new_buffer) { 320 321 /* shrink array to new_used + new_buffer */ 322 for (f=0; f<db->nfields; f++) { 323 ierr = DataFieldSetSize( db->field[f], new_allocated );CHKERRQ(ierr); 324 } 325 326 db->L = new_used; 327 db->buffer = new_buffer; 328 db->allocated = new_used + new_buffer; 329 } 330 else { 331 db->L = new_used; 332 db->buffer = new_buffer; 333 } 334 } 335 336 /* zero all entries from db->L to db->allocated */ 337 for (f=0; f<db->nfields; f++) { 338 DataField field = db->field[f]; 339 ierr = DataFieldZeroBlock(field, db->L,db->allocated);CHKERRQ(ierr); 340 } 341 PetscFunctionReturn(0); 342 } 343 344 #undef __FUNCT__ 345 #define __FUNCT__ "DataBucketSetInitialSizes" 346 PetscErrorCode DataBucketSetInitialSizes(DataBucket db,const PetscInt L,const PetscInt buffer) 347 { 348 PetscInt f; 349 PetscErrorCode ierr; 350 351 ierr = DataBucketSetSizes(db,L,buffer);CHKERRQ(ierr); 352 353 for (f=0; f<db->nfields; f++) { 354 DataField field = db->field[f]; 355 ierr = DataFieldZeroBlock(field,0,db->allocated);CHKERRQ(ierr); 356 } 357 PetscFunctionReturn(0); 358 } 359 360 #undef __FUNCT__ 361 #define __FUNCT__ "DataBucketGetSizes" 362 PetscErrorCode DataBucketGetSizes(DataBucket db,PetscInt *L,PetscInt *buffer,PetscInt *allocated) 363 { 364 if (L) { *L = db->L; } 365 if (buffer) { *buffer = db->buffer; } 366 if (allocated) { *allocated = db->allocated; } 367 PetscFunctionReturn(0); 368 } 369 370 #undef __FUNCT__ 371 #define __FUNCT__ "DataBucketGetGlobalSizes" 372 PetscErrorCode DataBucketGetGlobalSizes(MPI_Comm comm,DataBucket db,PetscInt *L,PetscInt *buffer,PetscInt *allocated) 373 { 374 PetscInt _L,_buffer,_allocated; 375 PetscInt ierr; 376 377 _L = db->L; 378 _buffer = db->buffer; 379 _allocated = db->allocated; 380 381 if (L) { ierr = MPI_Allreduce(&_L,L,1,MPIU_INT,MPI_SUM,comm);CHKERRQ(ierr); } 382 if (buffer) { ierr = MPI_Allreduce(&_buffer,buffer,1,MPIU_INT,MPI_SUM,comm);CHKERRQ(ierr); } 383 if (allocated) { ierr = MPI_Allreduce(&_allocated,allocated,1,MPIU_INT,MPI_SUM,comm);CHKERRQ(ierr); } 384 PetscFunctionReturn(0); 385 } 386 387 #undef __FUNCT__ 388 #define __FUNCT__ "DataBucketGetDataFields" 389 PetscErrorCode DataBucketGetDataFields(DataBucket db,PetscInt *L,DataField *fields[]) 390 { 391 if (L) { *L = db->nfields; } 392 if (fields) { *fields = db->field; } 393 PetscFunctionReturn(0); 394 } 395 396 #undef __FUNCT__ 397 #define __FUNCT__ "DataFieldGetAccess" 398 PetscErrorCode DataFieldGetAccess(const DataField gfield) 399 { 400 if (gfield->active) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Field \"%s\" is already active. You must call DataFieldRestoreAccess()",gfield->name); 401 402 gfield->active = PETSC_TRUE; 403 PetscFunctionReturn(0); 404 } 405 406 #undef __FUNCT__ 407 #define __FUNCT__ "DataFieldAccessPoint" 408 PetscErrorCode DataFieldAccessPoint(const DataField gfield,const PetscInt pid,void **ctx_p) 409 { 410 #ifdef DATAFIELD_POINT_ACCESS_GUARD 411 /* debug mode */ 412 /* check point is valid */ 413 if (pid < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"index must be >= 0"); 414 if (pid >= gfield->L) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"index must be < %D",gfield->L); 415 416 if (gfield->active == PETSC_FALSE) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Field \"%s\" is not active. You must call DataFieldGetAccess() before point data can be retrivied",gfield->name); 417 #endif 418 419 //*ctx_p = (void*)( ((char*)gfield->data) + pid * gfield->atomic_size); 420 *ctx_p = __DATATFIELD_point_access(gfield->data,pid,gfield->atomic_size); 421 PetscFunctionReturn(0); 422 } 423 424 #undef __FUNCT__ 425 #define __FUNCT__ "DataFieldAccessPointOffset" 426 PetscErrorCode DataFieldAccessPointOffset(const DataField gfield,const size_t offset,const PetscInt pid,void **ctx_p) 427 { 428 #ifdef DATAFIELD_POINT_ACCESS_GUARD 429 /* debug mode */ 430 431 /* check point is valid */ 432 /* if( offset < 0 ) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"offset must be >= 0");*//* Note compiler realizes this can never happen with an unsigned PetscInt */ 433 if (offset >= gfield->atomic_size) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"offset must be < %zu",gfield->atomic_size); 434 435 /* check point is valid */ 436 if (pid < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"index must be >= 0"); 437 if (pid >= gfield->L) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"index must be < %D",gfield->L); 438 439 if (gfield->active == PETSC_FALSE) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Field \"%s\" is not active. You must call DataFieldGetAccess() before point data can be retrivied",gfield->name); 440 #endif 441 442 *ctx_p = __DATATFIELD_point_access_offset(gfield->data,pid,gfield->atomic_size,offset); 443 PetscFunctionReturn(0); 444 } 445 446 #undef __FUNCT__ 447 #define __FUNCT__ "DataFieldRestoreAccess" 448 PetscErrorCode DataFieldRestoreAccess(DataField gfield) 449 { 450 if (gfield->active == PETSC_FALSE) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Field \"%s\" is not active. You must call DataFieldGetAccess()", gfield->name ); 451 452 gfield->active = PETSC_FALSE; 453 PetscFunctionReturn(0); 454 } 455 456 #undef __FUNCT__ 457 #define __FUNCT__ "DataFieldVerifyAccess" 458 PetscErrorCode DataFieldVerifyAccess(const DataField gfield,const size_t size) 459 { 460 #ifdef DATAFIELD_POINT_ACCESS_GUARD 461 if (gfield->atomic_size != size) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_USER,"Field \"%s\" must be mapped to %zu bytes, your intended structure is %zu bytes in length.", 462 gfield->name, gfield->atomic_size, size ); 463 #endif 464 PetscFunctionReturn(0); 465 } 466 467 #undef __FUNCT__ 468 #define __FUNCT__ "DataFieldGetAtomicSize" 469 PetscErrorCode DataFieldGetAtomicSize(const DataField gfield,size_t *size) 470 { 471 if (size) { *size = gfield->atomic_size; } 472 PetscFunctionReturn(0); 473 } 474 475 #undef __FUNCT__ 476 #define __FUNCT__ "DataFieldGetEntries" 477 PetscErrorCode DataFieldGetEntries(const DataField gfield,void **data) 478 { 479 if (data) { 480 *data = gfield->data; 481 } 482 PetscFunctionReturn(0); 483 } 484 485 #undef __FUNCT__ 486 #define __FUNCT__ "DataFieldRestoreEntries" 487 PetscErrorCode DataFieldRestoreEntries(const DataField gfield,void **data) 488 { 489 if (data) { 490 *data = NULL; 491 } 492 PetscFunctionReturn(0); 493 } 494 495 /* y = x */ 496 #undef __FUNCT__ 497 #define __FUNCT__ "DataBucketCopyPoint" 498 PetscErrorCode DataBucketCopyPoint(const DataBucket xb,const PetscInt pid_x, 499 const DataBucket yb,const PetscInt pid_y) 500 { 501 PetscInt f; 502 PetscErrorCode ierr; 503 504 for (f=0; f<xb->nfields; f++) { 505 void *dest; 506 void *src; 507 508 ierr = DataFieldGetAccess( xb->field[f] );CHKERRQ(ierr); 509 if (xb != yb) { ierr = DataFieldGetAccess( yb->field[f] );CHKERRQ(ierr); } 510 511 ierr = DataFieldAccessPoint( xb->field[f],pid_x, &src );CHKERRQ(ierr); 512 ierr = DataFieldAccessPoint( yb->field[f],pid_y, &dest );CHKERRQ(ierr); 513 514 memcpy( dest, src, xb->field[f]->atomic_size ); 515 516 ierr = DataFieldRestoreAccess( xb->field[f] );CHKERRQ(ierr); 517 if (xb != yb) { ierr = DataFieldRestoreAccess( yb->field[f] );CHKERRQ(ierr); } 518 } 519 PetscFunctionReturn(0); 520 } 521 522 #undef __FUNCT__ 523 #define __FUNCT__ "DataBucketCreateFromSubset" 524 PetscErrorCode DataBucketCreateFromSubset(DataBucket DBIn,const PetscInt N,const PetscInt list[],DataBucket *DB) 525 { 526 PetscInt nfields; 527 DataField *fields; 528 DataBucketCreate(DB); 529 PetscInt f,L,buffer,allocated,p; 530 PetscErrorCode ierr; 531 532 /* copy contents of DBIn */ 533 ierr = DataBucketGetDataFields(DBIn,&nfields,&fields);CHKERRQ(ierr); 534 ierr = DataBucketGetSizes(DBIn,&L,&buffer,&allocated);CHKERRQ(ierr); 535 536 for (f=0; f<nfields; f++) { 537 ierr = DataBucketRegisterField(*DB,"DataBucketCreateFromSubset",fields[f]->name,fields[f]->atomic_size,NULL);CHKERRQ(ierr); 538 } 539 ierr = DataBucketFinalize(*DB);CHKERRQ(ierr); 540 541 ierr = DataBucketSetSizes(*DB,L,buffer);CHKERRQ(ierr); 542 543 /* now copy the desired guys from DBIn => DB */ 544 for (p=0; p<N; p++) { 545 ierr = DataBucketCopyPoint(DBIn,list[p], *DB,p);CHKERRQ(ierr); 546 } 547 PetscFunctionReturn(0); 548 } 549 550 // insert into an exisitng location 551 #undef __FUNCT__ 552 #define __FUNCT__ "DataFieldInsertPoint" 553 PetscErrorCode DataFieldInsertPoint(const DataField field,const PetscInt index,const void *ctx) 554 { 555 556 #ifdef DATAFIELD_POINT_ACCESS_GUARD 557 /* check point is valid */ 558 if (index < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"index must be >= 0"); 559 if (index >= field->L) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"index must be < %D",field->L); 560 #endif 561 562 // memcpy( (void*)((char*)field->data + index*field->atomic_size), ctx, field->atomic_size ); 563 memcpy( __DATATFIELD_point_access(field->data,index,field->atomic_size), ctx, field->atomic_size ); 564 PetscFunctionReturn(0); 565 } 566 567 // remove data at index - replace with last point 568 #undef __FUNCT__ 569 #define __FUNCT__ "DataBucketRemovePointAtIndex" 570 PetscErrorCode DataBucketRemovePointAtIndex(const DataBucket db,const PetscInt index) 571 { 572 PetscInt f; 573 PetscBool any_active_fields; 574 PetscErrorCode ierr; 575 576 #ifdef DATAFIELD_POINT_ACCESS_GUARD 577 /* check point is valid */ 578 if (index < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"index must be >= 0"); 579 if (index >= db->allocated) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"index must be < %D",db->L+db->buffer); 580 #endif 581 582 ierr = DataBucketQueryForActiveFields(db,&any_active_fields);CHKERRQ(ierr); 583 if (any_active_fields) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Cannot safely remove point as at least one DataField is currently being accessed"); 584 585 if (index >= db->L) { /* this point is not in the list - no need to error, but I will anyway */ 586 SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_USER,"You should not be trying to remove point at index=%D since it's < db->L = %D", index, db->L ); 587 } 588 589 if (index != db->L-1) { /* not last point in list */ 590 for (f=0; f<db->nfields; f++) { 591 DataField field = db->field[f]; 592 593 /* copy then remove */ 594 ierr = DataFieldCopyPoint( db->L-1,field, index,field );CHKERRQ(ierr); 595 596 //DataFieldZeroPoint(field,index); 597 } 598 } 599 600 /* decrement size */ 601 /* this will zero out an crap at the end of the list */ 602 ierr = DataBucketRemovePoint(db);CHKERRQ(ierr); 603 PetscFunctionReturn(0); 604 } 605 606 /* copy x into y */ 607 #undef __FUNCT__ 608 #define __FUNCT__ "DataFieldCopyPoint" 609 PetscErrorCode DataFieldCopyPoint(const PetscInt pid_x,const DataField field_x, 610 const PetscInt pid_y,const DataField field_y ) 611 { 612 613 #ifdef DATAFIELD_POINT_ACCESS_GUARD 614 /* check point is valid */ 615 if (pid_x < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"(IN) index must be >= 0"); 616 if (pid_x >= field_x->L) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"(IN) index must be < %D",field_x->L); 617 618 if (pid_y < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"(OUT) index must be >= 0"); 619 if (pid_y >= field_y->L) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"(OUT) index must be < %D",field_y->L); 620 621 if( field_y->atomic_size != field_x->atomic_size ) { 622 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"atomic size must match"); 623 } 624 #endif 625 /* 626 memcpy( (void*)((char*)field_y->data + pid_y*field_y->atomic_size), 627 (void*)((char*)field_x->data + pid_x*field_x->atomic_size), 628 field_x->atomic_size ); 629 */ 630 memcpy( __DATATFIELD_point_access(field_y->data,pid_y,field_y->atomic_size), 631 __DATATFIELD_point_access(field_x->data,pid_x,field_x->atomic_size), 632 field_y->atomic_size ); 633 PetscFunctionReturn(0); 634 } 635 636 637 // zero only the datafield at this point 638 #undef __FUNCT__ 639 #define __FUNCT__ "DataFieldZeroPoint" 640 PetscErrorCode DataFieldZeroPoint(const DataField field,const PetscInt index) 641 { 642 #ifdef DATAFIELD_POINT_ACCESS_GUARD 643 /* check point is valid */ 644 if (index < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"index must be >= 0"); 645 if (index >= field->L) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"index must be < %D",field->L); 646 #endif 647 648 // memset( (void*)((char*)field->data + index*field->atomic_size), 0, field->atomic_size ); 649 memset( __DATATFIELD_point_access(field->data,index,field->atomic_size), 0, field->atomic_size ); 650 PetscFunctionReturn(0); 651 } 652 653 // zero ALL data for this point 654 #undef __FUNCT__ 655 #define __FUNCT__ "DataBucketZeroPoint" 656 PetscErrorCode DataBucketZeroPoint(const DataBucket db,const PetscInt index) 657 { 658 PetscInt f; 659 PetscErrorCode ierr; 660 661 /* check point is valid */ 662 if (index < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"index must be >= 0"); 663 if (index >= db->allocated) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"index must be < %D",db->allocated); 664 665 for (f=0; f<db->nfields; f++) { 666 DataField field = db->field[f]; 667 668 ierr = DataFieldZeroPoint(field,index);CHKERRQ(ierr); 669 } 670 PetscFunctionReturn(0); 671 } 672 673 /* increment */ 674 #undef __FUNCT__ 675 #define __FUNCT__ "DataBucketAddPoint" 676 PetscErrorCode DataBucketAddPoint(DataBucket db) 677 { 678 PetscErrorCode ierr; 679 680 ierr = DataBucketSetSizes(db,db->L+1,-1);CHKERRQ(ierr); 681 PetscFunctionReturn(0); 682 } 683 684 /* decrement */ 685 #undef __FUNCT__ 686 #define __FUNCT__ "DataBucketRemovePoint" 687 PetscErrorCode DataBucketRemovePoint(DataBucket db) 688 { 689 PetscErrorCode ierr; 690 691 ierr = DataBucketSetSizes(db,db->L-1,-1);CHKERRQ(ierr); 692 PetscFunctionReturn(0); 693 } 694 695 #undef __FUNCT__ 696 #define __FUNCT__ "_DataFieldViewBinary" 697 PetscErrorCode _DataFieldViewBinary(DataField field,FILE *fp) 698 { 699 fprintf(fp,"<DataField>\n"); 700 fprintf(fp,"%d\n", field->L); 701 fprintf(fp,"%zu\n",field->atomic_size); 702 fprintf(fp,"%s\n", field->registeration_function); 703 fprintf(fp,"%s\n", field->name); 704 705 fwrite(field->data, field->atomic_size, field->L, fp); 706 /* 707 printf(" ** wrote %zu bytes for DataField \"%s\" \n", field->atomic_size * field->L, field->name ); 708 */ 709 fprintf(fp,"\n</DataField>\n"); 710 PetscFunctionReturn(0); 711 } 712 713 #undef __FUNCT__ 714 #define __FUNCT__ "_DataBucketRegisterFieldFromFile" 715 PetscErrorCode _DataBucketRegisterFieldFromFile(FILE *fp,DataBucket db) 716 { 717 PetscBool val; 718 DataField *field; 719 720 DataField gfield; 721 char dummy[100]; 722 char registeration_function[5000]; 723 char field_name[5000]; 724 PetscInt L; 725 size_t atomic_size,strL; 726 PetscErrorCode ierr; 727 728 /* check we haven't finalised the registration of fields */ 729 /* 730 if(db->finalised==PETSC_TRUE) { 731 printf("ERROR: DataBucketFinalize() has been called. Cannot register more fields\n"); 732 ERROR(); 733 } 734 */ 735 736 737 /* read file contents */ 738 fgets(dummy,99,fp); //printf("read(header): %s", dummy ); 739 740 fscanf( fp, "%d\n",&L); //printf("read(L): %d\n", L); 741 742 fscanf( fp, "%zu\n",&atomic_size); //printf("read(size): %zu\n",atomic_size); 743 744 fgets(registeration_function,4999,fp); //printf("read(reg func): %s", registeration_function ); 745 strL = strlen(registeration_function); 746 if (strL > 1) { 747 registeration_function[strL-1] = 0; 748 } 749 750 fgets(field_name,4999,fp); //printf("read(name): %s", field_name ); 751 strL = strlen(field_name); 752 if (strL > 1) { 753 field_name[strL-1] = 0; 754 } 755 756 #ifdef DATA_BUCKET_LOG 757 PetscPrintf(PETSC_COMM_SELF," ** read L=%D; atomic_size=%zu; reg_func=\"%s\"; name=\"%s\" \n", L,atomic_size,registeration_function,field_name); 758 #endif 759 760 761 /* check for repeated name */ 762 ierr = StringInList( field_name, db->nfields, (const DataField*)db->field, &val );CHKERRQ(ierr); 763 if (val == PETSC_TRUE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Cannot add same field twice"); 764 765 /* create new space for data */ 766 field = realloc( db->field, sizeof(DataField)*(db->nfields+1)); 767 db->field = field; 768 769 /* add field */ 770 ierr = DataFieldCreate( registeration_function, field_name, atomic_size, L, &gfield );CHKERRQ(ierr); 771 772 /* copy contents of file */ 773 fread(gfield->data, gfield->atomic_size, gfield->L, fp); 774 #ifdef DATA_BUCKET_LOG 775 PetscPrintf(PETSC_COMM_SELF," ** read %zu bytes for DataField \"%s\" \n", gfield->atomic_size * gfield->L, field_name ); 776 #endif 777 /* finish reading meta data */ 778 fgets(dummy,99,fp); //printf("read(header): %s", dummy ); 779 fgets(dummy,99,fp); //printf("read(header): %s", dummy ); 780 781 db->field[ db->nfields ] = gfield; 782 783 db->nfields++; 784 PetscFunctionReturn(0); 785 } 786 787 #undef __FUNCT__ 788 #define __FUNCT__ "_DataBucketViewAscii_HeaderWrite_v00" 789 PetscErrorCode _DataBucketViewAscii_HeaderWrite_v00(FILE *fp) 790 { 791 fprintf(fp,"<DataBucketHeader>\n"); 792 fprintf(fp,"type=DataBucket\n"); 793 fprintf(fp,"format=ascii\n"); 794 fprintf(fp,"version=0.0\n"); 795 fprintf(fp,"options=\n"); 796 fprintf(fp,"</DataBucketHeader>\n"); 797 PetscFunctionReturn(0); 798 } 799 800 #undef __FUNCT__ 801 #define __FUNCT__ "_DataBucketViewAscii_HeaderRead_v00" 802 PetscErrorCode _DataBucketViewAscii_HeaderRead_v00(FILE *fp) 803 { 804 char dummy[100]; 805 size_t strL; 806 807 // header open 808 fgets(dummy,99,fp); //printf("read(header): %s", dummy ); 809 810 // type 811 fgets(dummy,99,fp); //printf("read(header): %s", dummy ); 812 strL = strlen(dummy); 813 if (strL > 1) { dummy[strL-1] = 0; } 814 if (strcmp(dummy,"type=DataBucket") != 0) { 815 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Data file doesn't contain a DataBucket type"); 816 } 817 818 // format 819 fgets(dummy,99,fp); //printf("read(header): %s", dummy ); 820 821 // version 822 fgets(dummy,99,fp); //printf("read(header): %s", dummy ); 823 strL = strlen(dummy); 824 if (strL > 1) { dummy[strL-1] = 0; } 825 if (strcmp(dummy,"version=0.0") != 0) { 826 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"DataBucket file must be parsed with version=0.0 : You tried %s", dummy); 827 } 828 829 // options 830 fgets(dummy,99,fp); //printf("read(header): %s", dummy ); 831 // header close 832 fgets(dummy,99,fp); //printf("read(header): %s", dummy ); 833 PetscFunctionReturn(0); 834 } 835 836 #undef __FUNCT__ 837 #define __FUNCT__ "_DataBucketLoadFromFileBinary_SEQ" 838 PetscErrorCode _DataBucketLoadFromFileBinary_SEQ(const char filename[],DataBucket *_db) 839 { 840 DataBucket db; 841 FILE *fp; 842 PetscInt L,buffer,f,nfields; 843 PetscErrorCode ierr; 844 845 846 #ifdef DATA_BUCKET_LOG 847 PetscPrintf(PETSC_COMM_SELF,"** DataBucketLoadFromFile **\n"); 848 #endif 849 850 /* open file */ 851 fp = fopen(filename,"rb"); 852 if (fp == NULL) { 853 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_FILE_OPEN,"Cannot open file with name %s", filename); 854 } 855 856 /* read header */ 857 ierr = _DataBucketViewAscii_HeaderRead_v00(fp);CHKERRQ(ierr); 858 859 fscanf(fp,"%d\n%d\n%d\n",&L,&buffer,&nfields); 860 861 ierr = DataBucketCreate(&db);CHKERRQ(ierr); 862 863 for (f=0; f<nfields; f++) { 864 ierr = _DataBucketRegisterFieldFromFile(fp,db);CHKERRQ(ierr); 865 } 866 fclose(fp); 867 868 ierr = DataBucketFinalize(db);CHKERRQ(ierr); 869 870 /* 871 DataBucketSetSizes(db,L,buffer); 872 */ 873 db->L = L; 874 db->buffer = buffer; 875 db->allocated = L + buffer; 876 877 *_db = db; 878 PetscFunctionReturn(0); 879 } 880 881 #undef __FUNCT__ 882 #define __FUNCT__ "DataBucketLoadFromFile" 883 PetscErrorCode DataBucketLoadFromFile(MPI_Comm comm,const char filename[],DataBucketViewType type,DataBucket *db) 884 { 885 PetscMPIInt nproc,rank; 886 PetscErrorCode ierr; 887 888 ierr = MPI_Comm_size(comm,&nproc);CHKERRQ(ierr); 889 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 890 891 #ifdef DATA_BUCKET_LOG 892 PetscPrintf(PETSC_COMM_SELF,"** DataBucketLoadFromFile **\n"); 893 #endif 894 if (type == DATABUCKET_VIEW_STDOUT) { 895 896 } else if (type == DATABUCKET_VIEW_ASCII) { 897 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Cannot be implemented as we don't know the underlying particle data structure"); 898 } else if (type == DATABUCKET_VIEW_BINARY) { 899 if (nproc == 1) { 900 ierr = _DataBucketLoadFromFileBinary_SEQ(filename,db);CHKERRQ(ierr); 901 } else { 902 char *name; 903 904 asprintf(&name,"%s_p%1.5d",filename, rank ); 905 ierr = _DataBucketLoadFromFileBinary_SEQ(name,db);CHKERRQ(ierr); 906 free(name); 907 } 908 } else { 909 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Unknown viewer requested"); 910 } 911 PetscFunctionReturn(0); 912 } 913 914 #undef __FUNCT__ 915 #define __FUNCT__ "_DataBucketViewBinary" 916 PetscErrorCode _DataBucketViewBinary(DataBucket db,const char filename[]) 917 { 918 FILE *fp = NULL; 919 PetscInt f; 920 PetscErrorCode ierr; 921 922 fp = fopen(filename,"wb"); 923 if (fp == NULL) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Cannot open file with name %s", filename); 924 925 /* db header */ 926 ierr =_DataBucketViewAscii_HeaderWrite_v00(fp);CHKERRQ(ierr); 927 928 /* meta-data */ 929 fprintf(fp,"%d\n%d\n%d\n", db->L,db->buffer,db->nfields); 930 931 for (f=0; f<db->nfields; f++) { 932 /* load datafields */ 933 ierr = _DataFieldViewBinary(db->field[f],fp);CHKERRQ(ierr); 934 } 935 936 fclose(fp); 937 PetscFunctionReturn(0); 938 } 939 940 #undef __FUNCT__ 941 #define __FUNCT__ "DataBucketView_SEQ" 942 PetscErrorCode DataBucketView_SEQ(DataBucket db,const char filename[],DataBucketViewType type) 943 { 944 PetscErrorCode ierr; 945 946 switch (type) { 947 case DATABUCKET_VIEW_STDOUT: 948 { 949 PetscInt f; 950 double memory_usage_total = 0.0; 951 952 PetscPrintf(PETSC_COMM_SELF,"DataBucketView(SEQ): (\"%s\")\n",filename); 953 PetscPrintf(PETSC_COMM_SELF," L = %D \n", db->L ); 954 PetscPrintf(PETSC_COMM_SELF," buffer = %D \n", db->buffer ); 955 PetscPrintf(PETSC_COMM_SELF," allocated = %D \n", db->allocated ); 956 957 PetscPrintf(PETSC_COMM_SELF," nfields registered = %D \n", db->nfields ); 958 for (f=0; f<db->nfields; f++) { 959 double memory_usage_f = (double)(db->field[f]->atomic_size * db->allocated) * 1.0e-6; 960 961 PetscPrintf(PETSC_COMM_SELF," [%3D]: field name ==>> %30s : Mem. usage = %1.2e (MB) \n", f, db->field[f]->name, memory_usage_f ); 962 PetscPrintf(PETSC_COMM_SELF," blocksize = %D \n", db->field[f]->bs ); 963 if (db->field[f]->bs != 1) { 964 PetscPrintf(PETSC_COMM_SELF," atomic size = %zu [full block, bs=%D]\n", db->field[f]->atomic_size,db->field[f]->bs ); 965 PetscPrintf(PETSC_COMM_SELF," atomic size/item = %zu \n", db->field[f]->atomic_size/db->field[f]->bs ); 966 } else { 967 PetscPrintf(PETSC_COMM_SELF," atomic size = %zu \n", db->field[f]->atomic_size ); 968 } 969 memory_usage_total += memory_usage_f; 970 } 971 PetscPrintf(PETSC_COMM_SELF," Total mem. usage = %1.2e (MB) \n", memory_usage_total ); 972 } 973 break; 974 975 case DATABUCKET_VIEW_ASCII: 976 { 977 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Cannot be implemented as we don't know the underlying particle data structure"); 978 } 979 break; 980 981 case DATABUCKET_VIEW_BINARY: 982 { 983 ierr = _DataBucketViewBinary(db,filename);CHKERRQ(ierr); 984 } 985 break; 986 987 case DATABUCKET_VIEW_HDF5: 988 { 989 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No HDF5 support"); 990 } 991 break; 992 993 default: 994 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Unknown viewer method requested"); 995 break; 996 } 997 PetscFunctionReturn(0); 998 } 999 1000 #undef __FUNCT__ 1001 #define __FUNCT__ "DataBucketView_MPI" 1002 PetscErrorCode DataBucketView_MPI(MPI_Comm comm,DataBucket db,const char filename[],DataBucketViewType type) 1003 { 1004 PetscErrorCode ierr; 1005 1006 switch (type) { 1007 case DATABUCKET_VIEW_STDOUT: 1008 { 1009 PetscInt f; 1010 PetscInt L,buffer,allocated; 1011 double memory_usage_total,memory_usage_total_local = 0.0; 1012 PetscMPIInt rank; 1013 1014 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 1015 1016 ierr = DataBucketGetGlobalSizes(comm,db,&L,&buffer,&allocated);CHKERRQ(ierr); 1017 1018 for (f=0; f<db->nfields; f++) { 1019 double memory_usage_f = (double)(db->field[f]->atomic_size * db->allocated) * 1.0e-6; 1020 1021 memory_usage_total_local += memory_usage_f; 1022 } 1023 ierr = MPI_Allreduce(&memory_usage_total_local,&memory_usage_total,1,MPI_DOUBLE,MPI_SUM,comm);CHKERRQ(ierr); 1024 1025 if (rank == 0) { 1026 PetscPrintf(comm,"DataBucketView(MPI): (\"%s\")\n",filename); 1027 PetscPrintf(comm," L = %D \n", L ); 1028 PetscPrintf(comm," buffer (max) = %D \n", buffer ); 1029 PetscPrintf(comm," allocated = %D \n", allocated ); 1030 1031 PetscPrintf(comm," nfields registered = %D \n", db->nfields ); 1032 for (f=0; f<db->nfields; f++) { 1033 double memory_usage_f = (double)(db->field[f]->atomic_size * db->allocated) * 1.0e-6; 1034 1035 PetscPrintf(PETSC_COMM_SELF," [%3D]: field name ==>> %30s : Mem. usage = %1.2e (MB) : rank0\n", f, db->field[f]->name, memory_usage_f ); 1036 } 1037 1038 PetscPrintf(PETSC_COMM_SELF," Total mem. usage = %1.2e (MB) : collective\n", memory_usage_total ); 1039 } 1040 1041 } 1042 break; 1043 1044 case DATABUCKET_VIEW_ASCII: 1045 { 1046 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Cannot be implemented as we don't know the underlying data structure"); 1047 } 1048 break; 1049 1050 case DATABUCKET_VIEW_BINARY: 1051 { 1052 char *name; 1053 PetscMPIInt rank; 1054 1055 /* create correct extension */ 1056 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 1057 asprintf(&name,"%s_p%1.5d",filename, rank ); 1058 1059 ierr = _DataBucketViewBinary(db,name);CHKERRQ(ierr); 1060 1061 free(name); 1062 } 1063 break; 1064 1065 case DATABUCKET_VIEW_HDF5: 1066 { 1067 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No support for HDF5"); 1068 } 1069 break; 1070 1071 default: 1072 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Unknown viewer method requested"); 1073 break; 1074 } 1075 PetscFunctionReturn(0); 1076 } 1077 1078 #undef __FUNCT__ 1079 #define __FUNCT__ "DataBucketView" 1080 PetscErrorCode DataBucketView(MPI_Comm comm,DataBucket db,const char filename[],DataBucketViewType type) 1081 { 1082 PetscMPIInt nproc; 1083 PetscErrorCode ierr; 1084 1085 ierr = MPI_Comm_size(comm,&nproc);CHKERRQ(ierr); 1086 if (nproc == 1) { 1087 ierr =DataBucketView_SEQ(db,filename,type);CHKERRQ(ierr); 1088 } else { 1089 ierr = DataBucketView_MPI(comm,db,filename,type);CHKERRQ(ierr); 1090 } 1091 PetscFunctionReturn(0); 1092 } 1093 1094 #undef __FUNCT__ 1095 #define __FUNCT__ "DataBucketDuplicateFields" 1096 PetscErrorCode DataBucketDuplicateFields(DataBucket dbA,DataBucket *dbB) 1097 { 1098 DataBucket db2; 1099 PetscInt f; 1100 PetscErrorCode ierr; 1101 1102 ierr = DataBucketCreate(&db2);CHKERRQ(ierr); 1103 1104 /* copy contents from dbA into db2 */ 1105 for (f=0; f<dbA->nfields; f++) { 1106 DataField field; 1107 size_t atomic_size; 1108 char *name; 1109 1110 field = dbA->field[f]; 1111 1112 atomic_size = field->atomic_size; 1113 name = field->name; 1114 1115 ierr = DataBucketRegisterField(db2,"DataBucketDuplicateFields",name,atomic_size,NULL);CHKERRQ(ierr); 1116 } 1117 ierr = DataBucketFinalize(db2);CHKERRQ(ierr); 1118 ierr = DataBucketSetInitialSizes(db2,0,1000);CHKERRQ(ierr); 1119 1120 /* set pointer */ 1121 *dbB = db2; 1122 PetscFunctionReturn(0); 1123 } 1124 1125 /* 1126 Insert points from db2 into db1 1127 db1 <<== db2 1128 */ 1129 #undef __FUNCT__ 1130 #define __FUNCT__ "DataBucketInsertValues" 1131 PetscErrorCode DataBucketInsertValues(DataBucket db1,DataBucket db2) 1132 { 1133 PetscInt n_mp_points1,n_mp_points2; 1134 PetscInt n_mp_points1_new,p; 1135 PetscErrorCode ierr; 1136 1137 ierr = DataBucketGetSizes(db1,&n_mp_points1,0,0);CHKERRQ(ierr); 1138 ierr = DataBucketGetSizes(db2,&n_mp_points2,0,0);CHKERRQ(ierr); 1139 1140 n_mp_points1_new = n_mp_points1 + n_mp_points2; 1141 ierr = DataBucketSetSizes(db1,n_mp_points1_new,-1);CHKERRQ(ierr); 1142 1143 for (p=0; p<n_mp_points2; p++) { 1144 // db1 <<== db2 // 1145 ierr =DataBucketCopyPoint( db2,p, db1,(n_mp_points1 + p) );CHKERRQ(ierr); 1146 } 1147 PetscFunctionReturn(0); 1148 } 1149 1150 /* helpers for parallel send/recv */ 1151 #undef __FUNCT__ 1152 #define __FUNCT__ "DataBucketCreatePackedArray" 1153 PetscErrorCode DataBucketCreatePackedArray(DataBucket db,size_t *bytes,void **buf) 1154 { 1155 PetscInt f; 1156 size_t sizeof_marker_contents; 1157 void *buffer; 1158 1159 sizeof_marker_contents = 0; 1160 for (f=0; f<db->nfields; f++) { 1161 DataField df = db->field[f]; 1162 1163 sizeof_marker_contents += df->atomic_size; 1164 } 1165 1166 buffer = malloc(sizeof_marker_contents); 1167 memset(buffer,0,sizeof_marker_contents); 1168 1169 if (bytes) { *bytes = sizeof_marker_contents; } 1170 if (buf) { *buf = buffer; } 1171 PetscFunctionReturn(0); 1172 } 1173 1174 #undef __FUNCT__ 1175 #define __FUNCT__ "DataBucketDestroyPackedArray" 1176 PetscErrorCode DataBucketDestroyPackedArray(DataBucket db,void **buf) 1177 { 1178 if (buf) { 1179 free(*buf); 1180 *buf = NULL; 1181 } 1182 PetscFunctionReturn(0); 1183 } 1184 1185 #undef __FUNCT__ 1186 #define __FUNCT__ "DataBucketFillPackedArray" 1187 PetscErrorCode DataBucketFillPackedArray(DataBucket db,const PetscInt index,void *buf) 1188 { 1189 PetscInt f; 1190 void *data,*data_p; 1191 size_t asize,offset; 1192 1193 offset = 0; 1194 for (f=0; f<db->nfields; f++) { 1195 DataField df = db->field[f]; 1196 1197 asize = df->atomic_size; 1198 1199 data = (void*)( df->data ); 1200 data_p = (void*)( (char*)data + index*asize ); 1201 1202 memcpy( (void*)((char*)buf + offset), data_p, asize); 1203 offset = offset + asize; 1204 } 1205 PetscFunctionReturn(0); 1206 } 1207 1208 #undef __FUNCT__ 1209 #define __FUNCT__ "DataBucketInsertPackedArray" 1210 PetscErrorCode DataBucketInsertPackedArray(DataBucket db,const PetscInt idx,void *data) 1211 { 1212 PetscInt f; 1213 void *data_p; 1214 size_t offset; 1215 PetscErrorCode ierr; 1216 1217 offset = 0; 1218 for (f=0; f<db->nfields; f++) { 1219 DataField df = db->field[f]; 1220 1221 data_p = (void*)( (char*)data + offset ); 1222 1223 ierr = DataFieldInsertPoint(df, idx, (void*)data_p );CHKERRQ(ierr); 1224 offset = offset + df->atomic_size; 1225 } 1226 PetscFunctionReturn(0); 1227 } 1228