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 dim; /* The topological mesh dimension */ 9 PetscBool simplex; /* Flag for simplices or tensor cells */ 10 char filename[PETSC_MAX_PATH_LEN]; /* Name of the mesh filename if any */ 11 PetscInt particlesPerCell; /* The number of partices per cell */ 12 } AppCtx; 13 14 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 15 { 16 PetscErrorCode ierr; 17 18 PetscFunctionBeginUser; 19 options->dim = 2; 20 options->simplex = PETSC_TRUE; 21 options->particlesPerCell = 1; 22 ierr = PetscStrcpy(options->filename, "");CHKERRQ(ierr); 23 ierr = PetscOptionsBegin(comm, "", "CellSwarm Options", "DMSWARM");CHKERRQ(ierr); 24 ierr = PetscOptionsInt("-dim", "The topological mesh dimension", "ex3.c", options->dim, &options->dim, NULL);CHKERRQ(ierr); 25 ierr = PetscOptionsBool("-simplex", "The flag for simplices or tensor cells", "ex3.c", options->simplex, &options->simplex, NULL);CHKERRQ(ierr); 26 ierr = PetscOptionsString("-mesh", "Name of the mesh filename if any", "ex3.c", options->filename, options->filename, PETSC_MAX_PATH_LEN, NULL);CHKERRQ(ierr); 27 ierr = PetscOptionsInt("-particles_per_cell", "Number of particles per cell", "ex3.c", options->particlesPerCell, &options->particlesPerCell, NULL);CHKERRQ(ierr); 28 ierr = PetscOptionsEnd();CHKERRQ(ierr); 29 PetscFunctionReturn(0); 30 } 31 32 static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm, AppCtx *user) 33 { 34 PetscBool flg; 35 PetscErrorCode ierr; 36 37 PetscFunctionBeginUser; 38 ierr = PetscStrcmp(user->filename, "", &flg);CHKERRQ(ierr); 39 if (flg) { 40 ierr = DMPlexCreateBoxMesh(comm, user->dim, user->simplex, NULL, NULL, NULL, NULL, PETSC_TRUE, dm);CHKERRQ(ierr); 41 } else { 42 ierr = DMPlexCreateFromFile(comm, user->filename, PETSC_TRUE, dm);CHKERRQ(ierr); 43 ierr = DMGetDimension(*dm, &user->dim);CHKERRQ(ierr); 44 } 45 { 46 DM distributedMesh = NULL; 47 48 ierr = DMPlexDistribute(*dm, 0, NULL, &distributedMesh);CHKERRQ(ierr); 49 if (distributedMesh) { 50 ierr = DMDestroy(dm);CHKERRQ(ierr); 51 *dm = distributedMesh; 52 } 53 } 54 ierr = DMLocalizeCoordinates(*dm);CHKERRQ(ierr); /* needed for periodic */ 55 ierr = DMSetFromOptions(*dm);CHKERRQ(ierr); 56 ierr = PetscObjectSetName((PetscObject) *dm, "Mesh");CHKERRQ(ierr); 57 ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr); 58 PetscFunctionReturn(0); 59 } 60 61 static PetscErrorCode CreateParticles(DM dm, DM *sw, AppCtx *user) 62 { 63 PetscInt *cellid; 64 PetscInt dim, cStart, cEnd, c, Np = user->particlesPerCell, p; 65 PetscErrorCode ierr; 66 67 PetscFunctionBeginUser; 68 ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 69 ierr = DMCreate(PetscObjectComm((PetscObject) dm), sw);CHKERRQ(ierr); 70 ierr = DMSetType(*sw, DMSWARM);CHKERRQ(ierr); 71 ierr = DMSetDimension(*sw, dim);CHKERRQ(ierr); 72 ierr = DMSwarmSetType(*sw, DMSWARM_PIC);CHKERRQ(ierr); 73 ierr = DMSwarmSetCellDM(*sw, dm);CHKERRQ(ierr); 74 ierr = DMSwarmRegisterPetscDatatypeField(*sw, "kinematics", 2, PETSC_REAL);CHKERRQ(ierr); 75 ierr = DMSwarmFinalizeFieldRegister(*sw);CHKERRQ(ierr); 76 ierr = DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);CHKERRQ(ierr); 77 ierr = DMSwarmSetLocalSizes(*sw, (cEnd - cStart) * Np, 0);CHKERRQ(ierr); 78 ierr = DMSetFromOptions(*sw);CHKERRQ(ierr); 79 ierr = DMSwarmGetField(*sw, DMSwarmPICField_cellid, NULL, NULL, (void **) &cellid);CHKERRQ(ierr); 80 for (c = cStart; c < cEnd; ++c) { 81 for (p = 0; p < Np; ++p) { 82 const PetscInt n = c*Np + p; 83 84 cellid[n] = c; 85 } 86 } 87 ierr = DMSwarmRestoreField(*sw, DMSwarmPICField_cellid, NULL, NULL, (void **) &cellid);CHKERRQ(ierr); 88 ierr = PetscObjectSetName((PetscObject) *sw, "Particles");CHKERRQ(ierr); 89 ierr = DMViewFromOptions(*sw, NULL, "-sw_view");CHKERRQ(ierr); 90 PetscFunctionReturn(0); 91 } 92 93 int main(int argc,char **argv) 94 { 95 DM dm, sw, cellsw; /* Mesh and particle managers */ 96 MPI_Comm comm; 97 AppCtx user; 98 PetscErrorCode ierr; 99 100 ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr; 101 comm = PETSC_COMM_WORLD; 102 ierr = ProcessOptions(comm, &user);CHKERRQ(ierr); 103 ierr = CreateMesh(comm, &dm, &user);CHKERRQ(ierr); 104 ierr = CreateParticles(dm, &sw, &user);CHKERRQ(ierr); 105 ierr = DMSetApplicationContext(sw, &user);CHKERRQ(ierr); 106 ierr = DMCreate(comm, &cellsw);CHKERRQ(ierr); 107 ierr = DMSwarmGetCellSwarm(sw, 1, cellsw);CHKERRQ(ierr); 108 ierr = DMViewFromOptions(cellsw, NULL, "-subswarm_view");CHKERRQ(ierr); 109 ierr = DMSwarmRestoreCellSwarm(sw, 1, cellsw);CHKERRQ(ierr); 110 ierr = DMDestroy(&sw);CHKERRQ(ierr); 111 ierr = DMDestroy(&dm);CHKERRQ(ierr); 112 ierr = DMDestroy(&cellsw); 113 ierr = PetscFinalize(); 114 return ierr; 115 } 116 117 /*TEST 118 build: 119 requires: triangle !single !complex 120 test: 121 suffix: 1 122 args: -particles_per_cell 2 -dm_plex_box_faces 2,1 -dm_view -sw_view -subswarm_view 123 TEST*/ 124