1 #pragma once 2 3 #include <petscvec.h> /*I "petscvec.h" I*/ 4 #include <petscmat.h> /*I "petscmat.h" I*/ 5 #include <petscdmswarm.h> /*I "petscdmswarm.h" I*/ 6 #include <petsc/private/dmimpl.h> 7 8 PETSC_EXTERN PetscBool SwarmProjcite; 9 PETSC_EXTERN const char SwarmProjCitation[]; 10 11 PETSC_EXTERN PetscLogEvent DMSWARM_Migrate; 12 PETSC_EXTERN PetscLogEvent DMSWARM_SetSizes; 13 PETSC_EXTERN PetscLogEvent DMSWARM_AddPoints; 14 PETSC_EXTERN PetscLogEvent DMSWARM_RemovePoints; 15 PETSC_EXTERN PetscLogEvent DMSWARM_Sort; 16 PETSC_EXTERN PetscLogEvent DMSWARM_DataExchangerTopologySetup; 17 PETSC_EXTERN PetscLogEvent DMSWARM_DataExchangerBegin; 18 PETSC_EXTERN PetscLogEvent DMSWARM_DataExchangerEnd; 19 PETSC_EXTERN PetscLogEvent DMSWARM_DataExchangerSendCount; 20 PETSC_EXTERN PetscLogEvent DMSWARM_DataExchangerPack; 21 22 typedef struct _p_CellDMInfo *CellDMInfo; 23 struct _p_CellDMInfo { 24 DM dm; // The cell DM 25 PetscInt Nf; // The number of DM fields 26 char **dmFields; // Swarm fields defining this DM 27 char *coordField; // Swarm field for coordinates on this DM 28 CellDMInfo next; // Next struct down in the stack 29 }; 30 31 /* 32 Error checking to ensure the swarm type is correct and that a cell DM has been set 33 */ 34 #define DMSWARMPICVALID(obj) \ 35 do { \ 36 DM_Swarm *_swarm = (DM_Swarm *)(obj)->data; \ 37 PetscCheck(_swarm->swarm_type == DMSWARM_PIC, PetscObjectComm((PetscObject)(dm)), PETSC_ERR_SUP, "Valid only for DMSwarm-PIC. You must call DMSwarmSetType(dm,DMSWARM_PIC)"); \ 38 PetscCheck(_swarm->cellinfo && _swarm->cellinfo[0].dm, PetscObjectComm((PetscObject)(dm)), PETSC_ERR_SUP, "Valid only for DMSwarmPIC if the cell DM is set. You must call DMSwarmSetCellDM() or DMSwarmPushCellDM()"); \ 39 } while (0) 40 41 typedef struct { 42 DMSwarmDataBucket db; 43 PetscInt refct; 44 PetscBool field_registration_initialized; 45 PetscBool field_registration_finalized; 46 /* DMSwarmProjectMethod *swarm_project;*/ /* swarm, geometry, result */ 47 48 /* PetscInt overlap; */ 49 /* PetscErrorCode (*update_overlap)(void); */ 50 51 const char **vec_field_names; 52 PetscInt vec_field_num; 53 PetscInt vec_field_bs, vec_field_nlocal; 54 const char *coord_name; 55 56 PetscBool issetup; 57 DMSwarmType swarm_type; 58 DMSwarmMigrateType migrate_type; 59 DMSwarmCollectType collect_type; 60 61 CellDMInfo cellinfo; // The stack of cell DMs 62 63 PetscBool migrate_error_on_missing_point; 64 65 PetscBool collect_view_active; 66 PetscInt collect_view_reset_nlocal; 67 DMSwarmSort sort_context; 68 69 /* Support for PIC */ 70 PetscInt Ns; /* The number of particle species */ 71 72 PetscSimplePointFn *coordFunc; /* Function to set particle coordinates */ 73 PetscSimplePointFn *velFunc; /* Function to set particle velocities */ 74 75 /* Debugging */ 76 PetscInt printCoords; 77 PetscInt printWeights; 78 } DM_Swarm; 79 80 typedef struct { 81 PetscInt point_index; 82 PetscInt cell_index; 83 } SwarmPoint; 84 85 struct _p_DMSwarmSort { 86 PetscBool isvalid; 87 PetscInt ncells, npoints; 88 PetscInt *pcell_offsets; 89 SwarmPoint *list; 90 }; 91 92 PETSC_INTERN PetscErrorCode DMSwarmMigrate_Push_Basic(DM, PetscBool); 93 PETSC_INTERN PetscErrorCode DMSwarmMigrate_CellDMScatter(DM, PetscBool); 94 PETSC_INTERN PetscErrorCode DMSwarmMigrate_CellDMExact(DM, PetscBool); 95 96 PETSC_EXTERN PetscErrorCode DMSwarmReplace_Internal(DM, DM *); 97