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