Actual source code: petscdmswarm.h

  1: #ifndef PETSCDMSWARM_H
  2: #define PETSCDMSWARM_H

  4: #include <petscdm.h>
  5: #include <petscdt.h>

  7: /* SUBMANSEC = DMSwarm */

  9: /*E
 10:    DMSwarmType - Defines the type of `DMSWARM`

 12:    Values:
 13: +  `DMSWARM_BASIC` - defines N entries of varied data-types which the user may register.
 14: -  `DMSWARM_PIC` - suitable for particle-in-cell methods. Configured as `DMSWARM_PIC`, the swarm will be aware of, another `DM` which serves as the background mesh.
 15:    Fields specific to particle-in-cell methods are registered by default. These include spatial coordinates, a unique identifier, a cell index and an index for
 16:    the owning rank. The background mesh will (by default) define the spatial decomposition of the points defined in the swarm. `DMSWARM_PIC` provides support
 17:    for particle-in-cell operations such as defining initial point coordinates, communicating particles between sub-domains, projecting particle data fields on to the mesh.

 19:    Level: beginner

 21: .seealso: `DMSWARM`, `DMSwarmSetType()`
 22: E*/
 23: typedef enum {
 24:   DMSWARM_BASIC = 0,
 25:   DMSWARM_PIC
 26: } DMSwarmType;

 28: typedef enum {
 29:   DMSWARM_MIGRATE_BASIC = 0,
 30:   DMSWARM_MIGRATE_DMCELLNSCATTER,
 31:   DMSWARM_MIGRATE_DMCELLEXACT,
 32:   DMSWARM_MIGRATE_USER
 33: } DMSwarmMigrateType;

 35: typedef enum {
 36:   DMSWARM_COLLECT_BASIC = 0,
 37:   DMSWARM_COLLECT_DMDABOUNDINGBOX,
 38:   DMSWARM_COLLECT_GENERAL,
 39:   DMSWARM_COLLECT_USER
 40: } DMSwarmCollectType;

 42: /*E
 43:    DMSwarmPICLayoutType - Defines the method used to define particle coordinates within each cell. The layouts are constructured using the reference cell geometry

 45:    Values:
 46: +  `DMSWARMPIC_LAYOUT_REGULAR` - defines points on a regular ijk mesh. When using `DMSWARMPIC_LAYOUT_REGULAR`, the fill_param defines the number of points in each spatial direction.
 47: .  `DMSWARMPIC_LAYOUT_GAUSS` -  defines points using an npoint Gauss-Legendre tensor product quadrature rule.
 48:    When using `DMSWARMPIC_LAYOUT_GAUSS`, the fill_param defines the number of quadrature points in each spatial direction.

 50: -  `DMSWARMPIC_LAYOUT_SUBDIVISION` - defines points on the centroid of a sub-divided reference cell.
 51:    When using `DMSWARMPIC_LAYOUT_SUBDIVISION`, the fill_param defines the number times the reference cell is sub-divided.

 53:    Level: beginner

 55: .seealso: `DMSWARM`, `DM`, `DMSwarmInsertPointsUsingCellDM()`
 56: E*/
 57: typedef enum {
 58:   DMSWARMPIC_LAYOUT_REGULAR = 0,
 59:   DMSWARMPIC_LAYOUT_GAUSS,
 60:   DMSWARMPIC_LAYOUT_SUBDIVISION
 61: } DMSwarmPICLayoutType;

 63: PETSC_EXTERN const char *DMSwarmTypeNames[];
 64: PETSC_EXTERN const char *DMSwarmMigrateTypeNames[];
 65: PETSC_EXTERN const char *DMSwarmCollectTypeNames[];

 67: PETSC_EXTERN const char DMSwarmField_pid[];
 68: PETSC_EXTERN const char DMSwarmField_rank[];
 69: PETSC_EXTERN const char DMSwarmPICField_coor[];
 70: PETSC_EXTERN const char DMSwarmPICField_cellid[];

 72: PETSC_EXTERN PetscErrorCode DMSwarmCreateGlobalVectorFromField(DM, const char[], Vec *);
 73: PETSC_EXTERN PetscErrorCode DMSwarmDestroyGlobalVectorFromField(DM, const char[], Vec *);
 74: PETSC_EXTERN PetscErrorCode DMSwarmCreateLocalVectorFromField(DM, const char[], Vec *);
 75: PETSC_EXTERN PetscErrorCode DMSwarmDestroyLocalVectorFromField(DM, const char[], Vec *);

 77: PETSC_EXTERN PetscErrorCode DMSwarmInitializeFieldRegister(DM);
 78: PETSC_EXTERN PetscErrorCode DMSwarmFinalizeFieldRegister(DM);
 79: PETSC_EXTERN PetscErrorCode DMSwarmSetLocalSizes(DM, PetscInt, PetscInt);
 80: PETSC_EXTERN PetscErrorCode DMSwarmRegisterPetscDatatypeField(DM, const char[], PetscInt, PetscDataType);
 81: PETSC_EXTERN PetscErrorCode DMSwarmRegisterUserStructField(DM, const char[], size_t);
 82: PETSC_EXTERN PetscErrorCode DMSwarmRegisterUserDatatypeField(DM, const char[], size_t, PetscInt);
 83: PETSC_EXTERN PetscErrorCode DMSwarmGetField(DM, const char[], PetscInt *, PetscDataType *, void **);
 84: PETSC_EXTERN PetscErrorCode DMSwarmRestoreField(DM, const char[], PetscInt *, PetscDataType *, void **);

 86: PETSC_EXTERN PetscErrorCode DMSwarmVectorDefineField(DM, const char[]);

 88: PETSC_EXTERN PetscErrorCode DMSwarmAddPoint(DM);
 89: PETSC_EXTERN PetscErrorCode DMSwarmAddNPoints(DM, PetscInt);
 90: PETSC_EXTERN PetscErrorCode DMSwarmRemovePoint(DM);
 91: PETSC_EXTERN PetscErrorCode DMSwarmRemovePointAtIndex(DM, PetscInt);
 92: PETSC_EXTERN PetscErrorCode DMSwarmCopyPoint(DM dm, PetscInt, PetscInt);

 94: PETSC_EXTERN PetscErrorCode DMSwarmGetLocalSize(DM, PetscInt *);
 95: PETSC_EXTERN PetscErrorCode DMSwarmGetSize(DM, PetscInt *);
 96: PETSC_EXTERN PetscErrorCode DMSwarmGetMigrateType(DM, DMSwarmMigrateType *);
 97: PETSC_EXTERN PetscErrorCode DMSwarmSetMigrateType(DM, DMSwarmMigrateType);
 98: PETSC_EXTERN PetscErrorCode DMSwarmMigrate(DM, PetscBool);

100: PETSC_EXTERN PetscErrorCode DMSwarmCollectViewCreate(DM);
101: PETSC_EXTERN PetscErrorCode DMSwarmCollectViewDestroy(DM);
102: PETSC_EXTERN PetscErrorCode DMSwarmSetCellDM(DM, DM);
103: PETSC_EXTERN PetscErrorCode DMSwarmGetCellDM(DM, DM *);

105: PETSC_EXTERN PetscErrorCode DMSwarmSetType(DM, DMSwarmType);

107: PETSC_EXTERN PetscErrorCode DMSwarmSetPointsUniformCoordinates(DM, PetscReal *, PetscReal *, PetscInt *, InsertMode);
108: PETSC_EXTERN PetscErrorCode DMSwarmSetPointCoordinates(DM, PetscInt, PetscReal *, PetscBool, InsertMode);
109: PETSC_EXTERN PetscErrorCode DMSwarmInsertPointsUsingCellDM(DM, DMSwarmPICLayoutType, PetscInt);
110: PETSC_EXTERN PetscErrorCode DMSwarmSetPointCoordinatesCellwise(DM, PetscInt, PetscReal *);
111: PETSC_EXTERN PetscErrorCode DMSwarmSetPointCoordinatesRandom(DM, PetscInt);
112: PETSC_EXTERN PetscErrorCode DMSwarmViewFieldsXDMF(DM, const char *, PetscInt, const char **);
113: PETSC_EXTERN PetscErrorCode DMSwarmViewXDMF(DM, const char *);

115: PETSC_EXTERN PetscErrorCode DMSwarmSortGetAccess(DM);
116: PETSC_EXTERN PetscErrorCode DMSwarmSortRestoreAccess(DM);
117: PETSC_EXTERN PetscErrorCode DMSwarmSortGetPointsPerCell(DM, PetscInt, PetscInt *, PetscInt **);
118: PETSC_EXTERN PetscErrorCode DMSwarmSortGetNumberOfPointsPerCell(DM, PetscInt, PetscInt *);
119: PETSC_EXTERN PetscErrorCode DMSwarmSortGetIsValid(DM, PetscBool *);
120: PETSC_EXTERN PetscErrorCode DMSwarmSortGetSizes(DM, PetscInt *, PetscInt *);

122: PETSC_EXTERN PetscErrorCode DMSwarmProjectFields(DM, PetscInt, const char **, Vec **, PetscBool);
123: PETSC_EXTERN PetscErrorCode DMSwarmCreateMassMatrixSquare(DM, DM, Mat *);

125: PETSC_EXTERN PetscErrorCode DMSwarmGetCellSwarm(DM, PetscInt, DM);
126: PETSC_EXTERN PetscErrorCode DMSwarmRestoreCellSwarm(DM, PetscInt, DM);
127: PETSC_EXTERN PetscErrorCode DMSwarmGetNumSpecies(DM, PetscInt *);
128: PETSC_EXTERN PetscErrorCode DMSwarmSetNumSpecies(DM, PetscInt);
129: PETSC_EXTERN PetscErrorCode DMSwarmGetCoordinateFunction(DM, PetscErrorCode (**)(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar[], void *));
130: PETSC_EXTERN PetscErrorCode DMSwarmSetCoordinateFunction(DM, PetscErrorCode (*)(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar[], void *));
131: PETSC_EXTERN PetscErrorCode DMSwarmGetVelocityFunction(DM, PetscErrorCode (**)(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar[], void *));
132: PETSC_EXTERN PetscErrorCode DMSwarmSetVelocityFunction(DM, PetscErrorCode (*)(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar[], void *));
133: PETSC_EXTERN PetscErrorCode DMSwarmComputeLocalSize(DM, PetscInt, PetscProbFunc);
134: PETSC_EXTERN PetscErrorCode DMSwarmComputeLocalSizeFromOptions(DM);
135: PETSC_EXTERN PetscErrorCode DMSwarmInitializeCoordinates(DM);
136: PETSC_EXTERN PetscErrorCode DMSwarmInitializeVelocities(DM, PetscProbFunc, const PetscReal[]);
137: PETSC_EXTERN PetscErrorCode DMSwarmInitializeVelocitiesFromOptions(DM, const PetscReal[]);

139: #endif