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,bs; 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 __DATATFIELD_point_access(data,index,atomic_size) (void*)((char*)(data) + (index)*(atomic_size)) 40 #define __DATATFIELD_point_access_offset(data,index,atomic_size,offset) (void*)((char*)(data) + (index)*(atomic_size) + (offset)) 41 42 43 44 PetscErrorCode StringInList(const char name[],const PetscInt N,const DataField gfield[],PetscBool *val); 45 PetscErrorCode StringFindInList(const char name[],const PetscInt N,const DataField gfield[],PetscInt *index); 46 47 PetscErrorCode DataFieldCreate(const char registeration_function[],const char name[],const size_t size,const PetscInt L,DataField *DF); 48 PetscErrorCode DataFieldDestroy(DataField *DF); 49 PetscErrorCode DataBucketCreate(DataBucket *DB); 50 PetscErrorCode DataBucketDestroy(DataBucket *DB); 51 PetscErrorCode DataBucketQueryForActiveFields(DataBucket db,PetscBool *any_active_fields); 52 PetscErrorCode DataBucketRegisterField( 53 DataBucket db, 54 const char registeration_function[], 55 const char field_name[], 56 size_t atomic_size,DataField *_gfield ); 57 58 59 PetscErrorCode DataFieldGetNumEntries(DataField df,PetscInt *sum); 60 PetscErrorCode DataFieldSetBlockSize(DataField df,PetscInt blocksize); 61 PetscErrorCode DataFieldSetSize(DataField df,const PetscInt new_L); 62 PetscErrorCode DataFieldZeroBlock(DataField df,const PetscInt start,const PetscInt end); 63 PetscErrorCode DataFieldGetAccess(const DataField gfield); 64 PetscErrorCode DataFieldAccessPoint(const DataField gfield,const PetscInt pid,void **ctx_p); 65 PetscErrorCode DataFieldAccessPointOffset(const DataField gfield,const size_t offset,const PetscInt pid,void **ctx_p); 66 PetscErrorCode DataFieldRestoreAccess(DataField gfield); 67 PetscErrorCode DataFieldVerifyAccess(const DataField gfield,const size_t size); 68 PetscErrorCode DataFieldGetAtomicSize(const DataField gfield,size_t *size); 69 70 PetscErrorCode DataFieldGetEntries(const DataField gfield,void **data); 71 PetscErrorCode DataFieldRestoreEntries(const DataField gfield,void **data); 72 73 PetscErrorCode DataFieldInsertPoint(const DataField field,const PetscInt index,const void *ctx); 74 PetscErrorCode DataFieldCopyPoint(const PetscInt pid_x,const DataField field_x,const PetscInt pid_y,const DataField field_y); 75 PetscErrorCode DataFieldZeroPoint(const DataField field,const PetscInt index); 76 77 PetscErrorCode DataBucketGetDataFieldByName(DataBucket db,const char name[],DataField *gfield); 78 PetscErrorCode DataBucketQueryDataFieldByName(DataBucket db,const char name[],PetscBool *found); 79 PetscErrorCode DataBucketFinalize(DataBucket db); 80 PetscErrorCode DataBucketSetInitialSizes(DataBucket db,const PetscInt L,const PetscInt buffer); 81 PetscErrorCode DataBucketSetSizes(DataBucket db,const PetscInt L,const PetscInt buffer); 82 PetscErrorCode DataBucketGetSizes(DataBucket db,PetscInt *L,PetscInt *buffer,PetscInt *allocated); 83 PetscErrorCode DataBucketGetGlobalSizes(MPI_Comm comm,DataBucket db,PetscInt *L,PetscInt *buffer,PetscInt *allocated); 84 PetscErrorCode DataBucketGetDataFields(DataBucket db,PetscInt *L,DataField *fields[]); 85 86 PetscErrorCode DataBucketCopyPoint(const DataBucket xb,const PetscInt pid_x,const DataBucket yb,const PetscInt pid_y); 87 PetscErrorCode DataBucketCreateFromSubset(DataBucket DBIn,const PetscInt N,const PetscInt list[],DataBucket *DB); 88 PetscErrorCode DataBucketZeroPoint(const DataBucket db,const PetscInt index); 89 90 PetscErrorCode DataBucketView(MPI_Comm comm,DataBucket db,const char filename[],DataBucketViewType type); 91 92 PetscErrorCode DataBucketAddPoint(DataBucket db); 93 PetscErrorCode DataBucketRemovePoint(DataBucket db); 94 PetscErrorCode DataBucketRemovePointAtIndex(const DataBucket db,const PetscInt index); 95 96 PetscErrorCode DataBucketDuplicateFields(DataBucket dbA,DataBucket *dbB); 97 PetscErrorCode DataBucketInsertValues(DataBucket db1,DataBucket db2); 98 99 /* helpers for parallel send/recv */ 100 PetscErrorCode DataBucketCreatePackedArray(DataBucket db,size_t *bytes,void **buf); 101 PetscErrorCode DataBucketDestroyPackedArray(DataBucket db,void **buf); 102 PetscErrorCode DataBucketFillPackedArray(DataBucket db,const PetscInt index,void *buf); 103 PetscErrorCode DataBucketInsertPackedArray(DataBucket db,const PetscInt idx,void *data); 104 105 106 #endif 107 108