1 static char help[] = "Example usage of extracting single cells with their associated fields from a swarm and putting it in a new swarm object\n"; 2 3 #include <petscdmplex.h> 4 #include <petscdmswarm.h> 5 #include <petscts.h> 6 7 typedef struct { 8 PetscInt particlesPerCell; /* The number of partices per cell */ 9 } AppCtx; 10 11 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 12 { 13 PetscErrorCode ierr; 14 15 PetscFunctionBeginUser; 16 options->particlesPerCell = 1; 17 18 ierr = PetscOptionsBegin(comm, "", "CellSwarm Options", "DMSWARM");CHKERRQ(ierr); 19 ierr = PetscOptionsInt("-particles_per_cell", "Number of particles per cell", "ex3.c", options->particlesPerCell, &options->particlesPerCell, NULL);CHKERRQ(ierr); 20 ierr = PetscOptionsEnd();CHKERRQ(ierr); 21 PetscFunctionReturn(0); 22 } 23 24 static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm, AppCtx *user) 25 { 26 PetscErrorCode ierr; 27 28 PetscFunctionBeginUser; 29 ierr = DMCreate(comm, dm);CHKERRQ(ierr); 30 ierr = DMSetType(*dm, DMPLEX);CHKERRQ(ierr); 31 ierr = DMSetFromOptions(*dm);CHKERRQ(ierr); 32 ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr); 33 PetscFunctionReturn(0); 34 } 35 36 static PetscErrorCode CreateParticles(DM dm, DM *sw, AppCtx *user) 37 { 38 PetscInt *cellid; 39 PetscInt dim, cStart, cEnd, c, Np = user->particlesPerCell, p; 40 PetscErrorCode ierr; 41 42 PetscFunctionBeginUser; 43 ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 44 ierr = DMCreate(PetscObjectComm((PetscObject) dm), sw);CHKERRQ(ierr); 45 ierr = DMSetType(*sw, DMSWARM);CHKERRQ(ierr); 46 ierr = DMSetDimension(*sw, dim);CHKERRQ(ierr); 47 ierr = DMSwarmSetType(*sw, DMSWARM_PIC);CHKERRQ(ierr); 48 ierr = DMSwarmSetCellDM(*sw, dm);CHKERRQ(ierr); 49 ierr = DMSwarmRegisterPetscDatatypeField(*sw, "kinematics", 2, PETSC_REAL);CHKERRQ(ierr); 50 ierr = DMSwarmFinalizeFieldRegister(*sw);CHKERRQ(ierr); 51 ierr = DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);CHKERRQ(ierr); 52 ierr = DMSwarmSetLocalSizes(*sw, (cEnd - cStart) * Np, 0);CHKERRQ(ierr); 53 ierr = DMSetFromOptions(*sw);CHKERRQ(ierr); 54 ierr = DMSwarmGetField(*sw, DMSwarmPICField_cellid, NULL, NULL, (void **) &cellid);CHKERRQ(ierr); 55 for (c = cStart; c < cEnd; ++c) { 56 for (p = 0; p < Np; ++p) { 57 const PetscInt n = c*Np + p; 58 59 cellid[n] = c; 60 } 61 } 62 ierr = DMSwarmRestoreField(*sw, DMSwarmPICField_cellid, NULL, NULL, (void **) &cellid);CHKERRQ(ierr); 63 ierr = PetscObjectSetName((PetscObject) *sw, "Particles");CHKERRQ(ierr); 64 ierr = DMViewFromOptions(*sw, NULL, "-sw_view");CHKERRQ(ierr); 65 PetscFunctionReturn(0); 66 } 67 68 int main(int argc,char **argv) 69 { 70 DM dm, sw, cellsw; /* Mesh and particle managers */ 71 MPI_Comm comm; 72 AppCtx user; 73 PetscErrorCode ierr; 74 75 ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr; 76 comm = PETSC_COMM_WORLD; 77 ierr = ProcessOptions(comm, &user);CHKERRQ(ierr); 78 ierr = CreateMesh(comm, &dm, &user);CHKERRQ(ierr); 79 ierr = CreateParticles(dm, &sw, &user);CHKERRQ(ierr); 80 ierr = DMSetApplicationContext(sw, &user);CHKERRQ(ierr); 81 ierr = DMCreate(comm, &cellsw);CHKERRQ(ierr); 82 ierr = DMSwarmGetCellSwarm(sw, 1, cellsw);CHKERRQ(ierr); 83 ierr = DMViewFromOptions(cellsw, NULL, "-subswarm_view");CHKERRQ(ierr); 84 ierr = DMSwarmRestoreCellSwarm(sw, 1, cellsw);CHKERRQ(ierr); 85 ierr = DMDestroy(&sw);CHKERRQ(ierr); 86 ierr = DMDestroy(&dm);CHKERRQ(ierr); 87 ierr = DMDestroy(&cellsw);CHKERRQ(ierr); 88 ierr = PetscFinalize(); 89 return ierr; 90 } 91 92 /*TEST 93 build: 94 requires: triangle !single !complex 95 test: 96 suffix: 1 97 args: -particles_per_cell 2 -dm_plex_box_faces 2,1 -dm_view -sw_view -subswarm_view 98 TEST*/ 99