1 2 #ifndef __DATA_BUCKET_H__ 3 #define __DATA_BUCKET_H__ 4 5 6 #include <petsc.h> 7 #include <petsc/private/dmswarmimpl.h> /*I "petscdmswarm.h" I*/ 8 9 10 #define DEFAULT -32654789 11 12 #define DATAFIELD_POINT_ACCESS_GUARD 13 14 /* Logging flag */ 15 #define DATA_BUCKET_LOG 16 17 18 typedef enum { DATABUCKET_VIEW_STDOUT=0, DATABUCKET_VIEW_ASCII, DATABUCKET_VIEW_BINARY, DATABUCKET_VIEW_HDF5 } DataBucketViewType; 19 20 21 22 struct _p_DataField { 23 char *registeration_function; 24 PetscInt L; 25 PetscBool active; 26 size_t atomic_size; 27 char *name; /* what are they called */ 28 void *data; /* the data - an array of structs */ 29 PetscDataType petsc_type; 30 }; 31 32 struct _p_DataBucket { 33 PetscInt L; /* number in use */ 34 PetscInt buffer; /* memory buffer used for re-allocation */ 35 PetscInt allocated; /* number allocated, this will equal datafield->L */ 36 PetscBool finalised; /* DEPRECIATED */ 37 PetscInt nfields; /* how many fields of this type */ 38 DataField *field; /* the data */ 39 }; 40 41 #define ERROR() {\ 42 printf("ERROR: %s() from line %d in %s !!\n", __FUNCTION__, __LINE__, __FILE__);\ 43 exit(EXIT_FAILURE);\ 44 } 45 46 #define __DATATFIELD_point_access(data,index,atomic_size) (void*)((char*)(data) + (index)*(atomic_size)) 47 #define __DATATFIELD_point_access_offset(data,index,atomic_size,offset) (void*)((char*)(data) + (index)*(atomic_size) + (offset)) 48 49 50 51 void StringInList( const char name[], const PetscInt N, const DataField gfield[], PetscBool *val ); 52 void StringFindInList( const char name[], const PetscInt N, const DataField gfield[], PetscInt *index ); 53 54 void DataFieldCreate( const char registeration_function[], const char name[], const size_t size, const PetscInt L, DataField *DF ); 55 void DataFieldDestroy( DataField *DF ); 56 void DataBucketCreate( DataBucket *DB ); 57 void DataBucketDestroy( DataBucket *DB ); 58 void _DataBucketRegisterField( 59 DataBucket db, 60 const char registeration_function[], 61 const char field_name[], 62 size_t atomic_size, DataField *_gfield ); 63 64 65 #define DataBucketRegisterField(db,name,size,k) {\ 66 char *location;\ 67 asprintf(&location,"Registered by %s() at line %d within file %s", __FUNCTION__, __LINE__, __FILE__);\ 68 _DataBucketRegisterField( (db), location, (name), (size), (k) );\ 69 free(location);\ 70 } 71 72 void DataFieldGetNumEntries(DataField df, PetscInt *sum); 73 void DataFieldSetSize( DataField df, const PetscInt new_L ); 74 void DataFieldZeroBlock( DataField df, const PetscInt start, const PetscInt end ); 75 void DataFieldGetAccess( const DataField gfield ); 76 void DataFieldAccessPoint( const DataField gfield, const PetscInt pid, void **ctx_p ); 77 void DataFieldAccessPointOffset( const DataField gfield, const size_t offset, const PetscInt pid, void **ctx_p ); 78 void DataFieldRestoreAccess( DataField gfield ); 79 void DataFieldVerifyAccess( const DataField gfield, const size_t size); 80 void DataFieldGetAtomicSize(const DataField gfield,size_t *size); 81 82 void DataFieldGetEntries(const DataField gfield,void **data); 83 void DataFieldRestoreEntries(const DataField gfield,void **data); 84 85 void DataFieldInsertPoint( const DataField field, const PetscInt index, const void *ctx ); 86 void DataFieldCopyPoint( const PetscInt pid_x, const DataField field_x, 87 const PetscInt pid_y, const DataField field_y ); 88 void DataFieldZeroPoint( const DataField field, const PetscInt index ); 89 90 void DataBucketGetDataFieldByName(DataBucket db,const char name[],DataField *gfield); 91 void DataBucketQueryDataFieldByName(DataBucket db,const char name[],PetscBool *found); 92 void DataBucketFinalize(DataBucket db); 93 void DataBucketSetInitialSizes( DataBucket db, const PetscInt L, const PetscInt buffer ); 94 void DataBucketSetSizes( DataBucket db, const PetscInt L, const PetscInt buffer ); 95 void DataBucketGetSizes( DataBucket db, PetscInt *L, PetscInt *buffer, PetscInt *allocated ); 96 void DataBucketGetGlobalSizes(MPI_Comm comm, DataBucket db, PetscInt *L, PetscInt *buffer, PetscInt *allocated ); 97 void DataBucketGetDataFields( DataBucket db, PetscInt *L, DataField *fields[] ); 98 99 void DataBucketCopyPoint( const DataBucket xb, const PetscInt pid_x, 100 const DataBucket yb, const PetscInt pid_y ); 101 void DataBucketCreateFromSubset( DataBucket DBIn, const PetscInt N, const PetscInt list[], DataBucket *DB ); 102 void DataBucketZeroPoint( const DataBucket db, const PetscInt index ); 103 104 //void DataBucketLoadFromFile(const char filename[], DataBucketViewType type, DataBucket *db); 105 void DataBucketLoadFromFile(MPI_Comm comm,const char filename[], DataBucketViewType type, DataBucket *db); 106 //void DataBucketView(DataBucket db,const char filename[],DataBucketViewType type); 107 void DataBucketView(MPI_Comm comm,DataBucket db,const char filename[],DataBucketViewType type); 108 109 void DataBucketAddPoint( DataBucket db ); 110 void DataBucketRemovePoint( DataBucket db ); 111 void DataBucketRemovePointAtIndex( const DataBucket db, const PetscInt index ); 112 113 void DataBucketDuplicateFields(DataBucket dbA,DataBucket *dbB); 114 void DataBucketInsertValues(DataBucket db1,DataBucket db2); 115 116 /* helpers for parallel send/recv */ 117 void DataBucketCreatePackedArray(DataBucket db,size_t *bytes,void **buf); 118 void DataBucketDestroyPackedArray(DataBucket db,void **buf); 119 void DataBucketFillPackedArray(DataBucket db,const PetscInt index,void *buf); 120 void DataBucketInsertPackedArray(DataBucket db,const PetscInt idx,void *data); 121 122 123 #endif 124 125