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