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