1 static const char help[] = "Test initialization and migration with swarm."; 2 3 #include <petscdmplex.h> 4 #include <petscdmswarm.h> 5 6 typedef struct { 7 PetscReal L[3]; /* Dimensions of the mesh bounding box */ 8 } AppCtx; 9 10 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 11 { 12 PetscFunctionBeginUser; 13 PetscOptionsBegin(comm, "", "Swarm configuration options", "DMSWARM"); 14 PetscOptionsEnd(); 15 PetscFunctionReturn(PETSC_SUCCESS); 16 } 17 18 static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm, AppCtx *user) 19 { 20 PetscReal low[3], high[3]; 21 PetscInt cdim, d; 22 23 PetscFunctionBeginUser; 24 PetscCall(DMCreate(comm, dm)); 25 PetscCall(DMSetType(*dm, DMPLEX)); 26 PetscCall(DMSetFromOptions(*dm)); 27 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 28 PetscCall(DMGetCoordinateDim(*dm, &cdim)); 29 PetscCall(DMGetBoundingBox(*dm, low, high)); 30 PetscPrintf(PETSC_COMM_WORLD, "dim: %" PetscInt_FMT "\n", cdim); 31 for (d = 0; d < cdim; ++d) user->L[d] = high[d] - low[d]; 32 PetscFunctionReturn(PETSC_SUCCESS); 33 } 34 35 /* 36 This function initializes all particles on rank 0. 37 They are sent to other ranks to test migration across non nearest neighbors 38 */ 39 static PetscErrorCode CreateSwarm(DM dm, DM *sw, AppCtx *user) 40 { 41 PetscInt particleInitSize = 10; 42 PetscReal *coords, upper[3], lower[3]; 43 PetscInt *cellid, rank, size, Np, d, dim; 44 MPI_Comm comm; 45 46 PetscFunctionBegin; 47 comm = PETSC_COMM_WORLD; 48 MPI_Comm_rank(comm, &rank); 49 MPI_Comm_size(comm, &size); 50 PetscCall(DMGetBoundingBox(dm, lower, upper)); 51 PetscCall(DMCreate(PETSC_COMM_WORLD, sw)); 52 PetscCall(DMGetDimension(dm, &dim)); 53 PetscCall(DMSetType(*sw, DMSWARM)); 54 PetscCall(DMSetDimension(*sw, dim)); 55 PetscCall(DMSwarmSetType(*sw, DMSWARM_PIC)); 56 PetscCall(DMSwarmSetCellDM(*sw, dm)); 57 PetscCall(DMSwarmFinalizeFieldRegister(*sw)); 58 PetscCall(DMSwarmSetLocalSizes(*sw, rank == 0 ? particleInitSize : 0, 0)); 59 PetscCall(DMSetFromOptions(*sw)); 60 PetscCall(DMSwarmGetField(*sw, DMSwarmPICField_coor, NULL, NULL, (void **)&coords)); 61 PetscCall(DMSwarmGetField(*sw, DMSwarmPICField_cellid, NULL, NULL, (void **)&cellid)); 62 PetscCall(DMSwarmGetLocalSize(*sw, &Np)); 63 for (PetscInt p = 0; p < Np; ++p) { 64 for (PetscInt d = 0; d < dim; ++d) { coords[p * dim + d] = 0.5 * (upper[d] - lower[d]); } 65 coords[p * dim + 1] = (upper[1] - lower[1]) / particleInitSize * p + lower[1]; 66 cellid[p] = 0; 67 } 68 PetscCall(DMSwarmRestoreField(*sw, DMSwarmPICField_coor, NULL, NULL, (void **)&coords)); 69 PetscCall(DMSwarmRestoreField(*sw, DMSwarmPICField_cellid, NULL, NULL, (void **)&cellid)); 70 PetscCall(DMViewFromOptions(*sw, NULL, "-sw_view")); 71 PetscFunctionReturn(PETSC_SUCCESS); 72 } 73 74 /* 75 Configure the swarm on rank 0 with all particles 76 located there, then migrate where they need to be 77 */ 78 static PetscErrorCode CheckMigrate(DM sw) 79 { 80 Vec preMigrate, postMigrate, tmp; 81 PetscInt preSize, postSize; 82 PetscReal prenorm, postnorm; 83 84 PetscFunctionBeginUser; 85 PetscCall(DMSwarmCreateGlobalVectorFromField(sw, DMSwarmPICField_coor, &tmp)); 86 PetscCall(VecDuplicate(tmp, &preMigrate)); 87 PetscCall(VecCopy(tmp, preMigrate)); 88 PetscCall(DMSwarmDestroyGlobalVectorFromField(sw, DMSwarmPICField_coor, &tmp)); 89 PetscCall(DMSwarmMigrate(sw, PETSC_TRUE)); 90 PetscCall(DMSwarmCreateGlobalVectorFromField(sw, DMSwarmPICField_coor, &tmp)); 91 PetscCall(VecDuplicate(tmp, &postMigrate)); 92 PetscCall(VecCopy(tmp, postMigrate)); 93 PetscCall(DMSwarmDestroyGlobalVectorFromField(sw, DMSwarmPICField_coor, &tmp)); 94 PetscCall(VecGetSize(preMigrate, &preSize)); 95 PetscCall(VecGetSize(postMigrate, &postSize)); 96 PetscCheck(preSize == postSize, PETSC_COMM_WORLD, PETSC_ERR_ARG_SIZ, "Particles either lost or duplicated. Pre migrate global size %" PetscInt_FMT " != Post migrate size %" PetscInt_FMT "", preSize, postSize); 97 PetscCall(VecNorm(preMigrate, NORM_2, &prenorm)); 98 PetscCall(VecNorm(postMigrate, NORM_2, &postnorm)); 99 PetscCheck(PetscAbsReal(prenorm - postnorm) < PETSC_MACHINE_EPSILON, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Particle coordinates corrupted in migrate with abs(norm(pre) - norm(post)) = %.16g", PetscAbsReal(prenorm - postnorm)); 100 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Migrate check passes.\n")); 101 PetscCall(VecDestroy(&preMigrate)); 102 PetscCall(VecDestroy(&postMigrate)); 103 PetscFunctionReturn(PETSC_SUCCESS); 104 } 105 106 /* 107 Checks added points persist on migrate 108 */ 109 static PetscErrorCode CheckPointInsertion(DM sw) 110 { 111 PetscInt rank, size, Np_pre, Np_post; 112 MPI_Comm comm; 113 114 PetscFunctionBeginUser; 115 comm = PETSC_COMM_WORLD; 116 MPI_Comm_rank(comm, &rank); 117 MPI_Comm_size(comm, &size); 118 PetscCall(PetscPrintf(comm, "Basic point insertion check...\n")); 119 PetscCall(DMSwarmGetSize(sw, &Np_pre)); 120 if (rank == 0) PetscCall(DMSwarmAddPoint(sw)); 121 PetscCall(DMSwarmGetSize(sw, &Np_post)); 122 PetscCheck(Np_post == (Np_pre + 1), PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Particle insertion failed. Global size pre insertion: %" PetscInt_FMT " global size post insertion: %" PetscInt_FMT, Np_pre, Np_post); 123 PetscCall(CheckMigrate(sw)); 124 PetscCall(PetscPrintf(comm, "Basic point insertion check passes.\n")); 125 PetscFunctionReturn(PETSC_SUCCESS); 126 } 127 128 /* 129 Checks tie breaking works properly when a particle 130 is located at a shared boundary. The higher rank should 131 recieve the particle while the lower rank deletes it. 132 133 TODO: Currently only works for 2 procs. 134 */ 135 static PetscErrorCode CheckPointInsertion_Boundary(DM sw) 136 { 137 PetscInt rank, size, Np_loc_pre, Np_loc_post, dim; 138 PetscReal lbox_low[3], lbox_high[3], gbox_low[3], gbox_high[3]; 139 MPI_Comm comm; 140 DM cdm; 141 142 PetscFunctionBeginUser; 143 comm = PETSC_COMM_WORLD; 144 MPI_Comm_rank(comm, &rank); 145 MPI_Comm_size(comm, &size); 146 PetscCall(PetscPrintf(comm, "Rank boundary point insertion check...\n")); 147 PetscCall(DMSwarmGetCellDM(sw, &cdm)); 148 PetscCall(DMGetDimension(cdm, &dim)); 149 PetscCall(DMGetBoundingBox(cdm, gbox_low, gbox_high)); 150 if (rank == 0) { 151 PetscReal *coords; 152 PetscInt adjacentdim, Np; 153 154 PetscCall(DMGetLocalBoundingBox(cdm, lbox_low, lbox_high)); 155 // find a face that belongs to the neighbor. 156 for (PetscInt d = 0; d < dim; ++d) { 157 if ((gbox_high[d] - lbox_high[d]) != 0.) adjacentdim = d; 158 } 159 PetscCall(DMSwarmAddPoint(sw)); 160 PetscCall(DMSwarmGetLocalSize(sw, &Np)); 161 PetscCall(DMSwarmGetField(sw, DMSwarmPICField_coor, NULL, NULL, (void **)&coords)); 162 for (PetscInt d = 0; d < dim; ++d) coords[(Np - 1) * dim + d] = 0.5 * (lbox_high[d] - lbox_low[d]); 163 coords[(Np - 1) * dim + adjacentdim] = lbox_high[adjacentdim]; 164 PetscCall(DMSwarmRestoreField(sw, DMSwarmPICField_coor, NULL, NULL, (void **)&coords)); 165 } 166 PetscCall(DMSwarmGetLocalSize(sw, &Np_loc_pre)); 167 PetscCall(CheckMigrate(sw)); 168 PetscCall(DMSwarmGetLocalSize(sw, &Np_loc_post)); 169 if (rank == 0) PetscCheck(Np_loc_pre == (Np_loc_post + 1), PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Migration tie breaking failed on rank %" PetscInt_FMT ". Particle on boundary not sent.", rank); 170 if (rank == 1) PetscCheck(Np_loc_pre == (Np_loc_post - 1), PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Migration tie breaking failed on rank %" PetscInt_FMT ". Particle on boundary not recieved.", rank); 171 PetscCall(PetscPrintf(comm, "Rank boundary point insertion check passes.\n")); 172 PetscFunctionReturn(PETSC_SUCCESS); 173 } 174 175 int main(int argc, char **argv) 176 { 177 DM dm, sw; 178 PetscInt rank, size; 179 MPI_Comm comm; 180 AppCtx user; 181 182 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 183 comm = PETSC_COMM_WORLD; 184 PetscCall(ProcessOptions(comm, &user)); 185 PetscCall(CreateMesh(comm, &dm, &user)); 186 PetscCall(CreateSwarm(dm, &sw, &user)); 187 PetscCall(CheckMigrate(sw)); 188 PetscCall(CheckPointInsertion(sw)); 189 PetscCall(CheckPointInsertion_Boundary(sw)); 190 PetscCall(DMDestroy(&sw)); 191 PetscCall(DMDestroy(&dm)); 192 PetscFinalize(); 193 return 0; 194 } 195 196 /*TEST 197 # Swarm does not handle complex or quad 198 build: 199 requires: !complex double 200 201 test: 202 suffix: swarm_migrate_hash 203 requires: ctetgen 204 nsize: 2 205 args: -dm_plex_simplex 0 -dm_distribute_overlap 1 -dm_plex_hash_location true -dm_plex_box_faces 10,10,10 -dm_plex_box_lower 0.,0.,0. -dm_plex_box_upper 1.,1.,10. -dm_plex_box_bd none,none,none -dm_plex_dim 3 206 filter: grep -v marker | grep -v atomic | grep -v usage 207 test: 208 suffix: swarm_migrate_scan 209 requires: ctetgen 210 nsize: 2 211 args: -dm_plex_simplex 0 -dm_distribute_overlap 1 -dm_plex_hash_location false -dm_plex_box_faces 10,10,10 -dm_plex_box_lower 0.,0.,0. -dm_plex_box_upper 1.,1.,10. -dm_plex_box_bd none,none,none -dm_plex_dim 3 212 filter: grep -v marker | grep -v atomic | grep -v usage 213 TEST*/ 214