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