Actual source code: data_bucket.h

  1: #pragma once

  3: #include <petsc/private/dmswarmimpl.h>

  5: #define DMSWARM_DATA_BUCKET_BUFFER_DEFAULT -1
  6: #define DMSWARM_DATAFIELD_POINT_ACCESS_GUARD

  8: /* Logging flag */
  9: #define DMSWARM_DATA_BUCKET_LOG

 11: typedef enum {
 12:   DATABUCKET_VIEW_STDOUT = 0,
 13:   DATABUCKET_VIEW_ASCII,
 14:   DATABUCKET_VIEW_BINARY,
 15:   DATABUCKET_VIEW_HDF5
 16: } DMSwarmDataBucketViewType;

 18: struct _p_DMSwarmDataField {
 19:   char         *registration_function;
 20:   PetscInt      L, bs;
 21:   PetscBool     active;
 22:   size_t        atomic_size;
 23:   char         *name; /* what are they called */
 24:   void         *data; /* the data - an array of structs */
 25:   PetscDataType petsc_type;
 26: };

 28: struct _p_DMSwarmDataBucket {
 29:   PetscInt          L;         /* number in use */
 30:   PetscInt          buffer;    /* memory buffer used for re-allocation */
 31:   PetscInt          allocated; /* number allocated, this will equal datafield->L */
 32:   PetscBool         finalised; /* DEPRECATED */
 33:   PetscInt          nfields;   /* how many fields of this type */
 34:   DMSwarmDataField *field;     /* the data */
 35: };

 37: #define DMSWARM_DATAFIELD_point_access(data, index, atomic_size)                (void *)((char *)(data) + (index) * (atomic_size))
 38: #define DMSWARM_DATAFIELD_point_access_offset(data, index, atomic_size, offset) (void *)((char *)(data) + (index) * (atomic_size) + (offset))

 40: PETSC_INTERN PetscErrorCode DMSwarmDataFieldStringInList(const char[], const PetscInt, const DMSwarmDataField[], PetscBool *);
 41: PETSC_INTERN PetscErrorCode DMSwarmDataFieldStringFindInList(const char[], const PetscInt, const DMSwarmDataField[], PetscInt *);

 43: PETSC_INTERN PetscErrorCode DMSwarmDataFieldCreate(const char[], const char[], const size_t, const PetscInt, DMSwarmDataField *);
 44: PETSC_INTERN PetscErrorCode DMSwarmDataFieldDestroy(DMSwarmDataField *);
 45: PETSC_INTERN PetscErrorCode DMSwarmDataBucketCreate(DMSwarmDataBucket *);
 46: PETSC_INTERN PetscErrorCode DMSwarmDataBucketDestroy(DMSwarmDataBucket *);
 47: PETSC_INTERN PetscErrorCode DMSwarmDataBucketQueryForActiveFields(DMSwarmDataBucket, PetscBool *);
 48: PETSC_INTERN PetscErrorCode DMSwarmDataBucketRegisterField(DMSwarmDataBucket, const char[], const char[], size_t, DMSwarmDataField *);

 50: PETSC_INTERN PetscErrorCode DMSwarmDataFieldGetNumEntries(DMSwarmDataField, PetscInt *);
 51: PETSC_INTERN PetscErrorCode DMSwarmDataFieldSetBlockSize(DMSwarmDataField, PetscInt);
 52: PETSC_INTERN PetscErrorCode DMSwarmDataFieldSetSize(DMSwarmDataField, const PetscInt);
 53: PETSC_INTERN PetscErrorCode DMSwarmDataFieldZeroBlock(DMSwarmDataField, const PetscInt, const PetscInt);
 54: PETSC_INTERN PetscErrorCode DMSwarmDataFieldGetAccess(const DMSwarmDataField);
 55: PETSC_INTERN PetscErrorCode DMSwarmDataFieldAccessPoint(const DMSwarmDataField, const PetscInt, void **);
 56: PETSC_INTERN PetscErrorCode DMSwarmDataFieldAccessPointOffset(const DMSwarmDataField, const size_t, const PetscInt, void **);
 57: PETSC_INTERN PetscErrorCode DMSwarmDataFieldRestoreAccess(DMSwarmDataField);
 58: PETSC_INTERN PetscErrorCode DMSwarmDataFieldVerifyAccess(const DMSwarmDataField, const size_t);
 59: PETSC_INTERN PetscErrorCode DMSwarmDataFieldGetAtomicSize(const DMSwarmDataField, size_t *);

 61: PETSC_INTERN PetscErrorCode DMSwarmDataFieldInsertPoint(const DMSwarmDataField, const PetscInt, const void *);
 62: PETSC_INTERN PetscErrorCode DMSwarmDataFieldCopyPoint(const PetscInt, const DMSwarmDataField, const PetscInt, const DMSwarmDataField);
 63: PETSC_INTERN PetscErrorCode DMSwarmDataFieldZeroPoint(const DMSwarmDataField, const PetscInt);

 65: PETSC_INTERN PetscErrorCode DMSwarmDataBucketFinalize(DMSwarmDataBucket);
 66: PETSC_INTERN PetscErrorCode DMSwarmDataBucketSetInitialSizes(DMSwarmDataBucket, const PetscInt, const PetscInt);
 67: PETSC_INTERN PetscErrorCode DMSwarmDataBucketSetSizes(DMSwarmDataBucket, const PetscInt, const PetscInt);
 68: PETSC_INTERN PetscErrorCode DMSwarmDataBucketGetSizes(DMSwarmDataBucket, PetscInt *, PetscInt *, PetscInt *);
 69: PETSC_INTERN PetscErrorCode DMSwarmDataBucketGetGlobalSizes(MPI_Comm, DMSwarmDataBucket, PetscInt *, PetscInt *, PetscInt *);
 70: PETSC_INTERN PetscErrorCode DMSwarmDataBucketGetDMSwarmDataFields(DMSwarmDataBucket, PetscInt *, DMSwarmDataField *[]);

 72: PETSC_INTERN PetscErrorCode DMSwarmDataBucketCopyPoint(const DMSwarmDataBucket, const PetscInt, const DMSwarmDataBucket, const PetscInt);
 73: PETSC_INTERN PetscErrorCode DMSwarmDataBucketCreateFromSubset(DMSwarmDataBucket, const PetscInt, const PetscInt[], DMSwarmDataBucket *);
 74: PETSC_INTERN PetscErrorCode DMSwarmDataBucketZeroPoint(const DMSwarmDataBucket, const PetscInt);

 76: PETSC_INTERN PetscErrorCode DMSwarmDataBucketView(MPI_Comm, DMSwarmDataBucket, const char[], DMSwarmDataBucketViewType);

 78: PETSC_INTERN PetscErrorCode DMSwarmDataBucketAddPoint(DMSwarmDataBucket);
 79: PETSC_INTERN PetscErrorCode DMSwarmDataBucketRemovePoint(DMSwarmDataBucket);
 80: PETSC_INTERN PetscErrorCode DMSwarmDataBucketRemovePointAtIndex(const DMSwarmDataBucket, const PetscInt);

 82: PETSC_INTERN PetscErrorCode DMSwarmDataBucketDuplicateFields(DMSwarmDataBucket, DMSwarmDataBucket *);
 83: PETSC_INTERN PetscErrorCode DMSwarmDataBucketInsertValues(DMSwarmDataBucket, DMSwarmDataBucket);

 85: /* helpers for parallel send/recv */
 86: PETSC_INTERN PetscErrorCode DMSwarmDataBucketCreatePackedArray(DMSwarmDataBucket, size_t *, void **);
 87: PETSC_INTERN PetscErrorCode DMSwarmDataBucketDestroyPackedArray(DMSwarmDataBucket, void **);
 88: PETSC_INTERN PetscErrorCode DMSwarmDataBucketFillPackedArray(DMSwarmDataBucket, const PetscInt, void *);
 89: PETSC_INTERN PetscErrorCode DMSwarmDataBucketInsertPackedArray(DMSwarmDataBucket, const PetscInt, void *);