1 2 #define PETSCDM_DLL 3 #include <petsc/private/dmswarmimpl.h> /*I "petscdmswarm.h" I*/ 4 #include <petscsf.h> 5 #include <petscdmda.h> 6 #include <petscdmplex.h> 7 8 /* 9 Error chceking macto to ensure the swarm type is correct and that a cell DM has been set 10 */ 11 #define DMSWARMPICVALID(dm) \ 12 { \ 13 DM_Swarm *_swarm = (DM_Swarm*)(dm)->data; \ 14 if (_swarm->swarm_type != DMSWARM_PIC) SETERRQ(PetscObjectComm((PetscObject)(dm)),PETSC_ERR_SUP,"Only valid for DMSwarm-PIC. You must call DMSwarmSetType(dm,DMSWARM_PIC)"); \ 15 else \ 16 if (!_swarm->dmcell) SETERRQ(PetscObjectComm((PetscObject)(dm)),PETSC_ERR_SUP,"Only valid for DMSwarmPIC if the cell DM is set. You must call DMSwarmSetCellDM(dm,celldm)"); \ 17 } 18 19 /* Coordinate insertition/addition API */ 20 /*@C 21 DMSwarmSetPointsUniformCoordinates - Set point coordinates in a DMSwarm on a regular (ijk) grid 22 23 Collective on DM 24 25 Input parameters: 26 + dm - the DMSwarm 27 . min - minimum coordinate values in the x, y, z directions (array of length dim) 28 . max - maximum coordinate values in the x, y, z directions (array of length dim) 29 . npoints - number of points in each spatial direction (array of length dim) 30 - mode - indicates whether to append points to the swarm (ADD_VALUES), or over-ride existing points (INSERT_VALUES) 31 32 Level: beginner 33 34 Notes: 35 When using mode = INSERT_VALUES, this method will reset the number of particles in the DMSwarm 36 to be npoints[0]*npoints[1] (2D) or npoints[0]*npoints[1]*npoints[2] (3D). When using mode = ADD_VALUES, 37 new points will be appended to any already existing in the DMSwarm 38 39 .seealso: DMSwarmSetType(), DMSwarmSetCellDM(), DMSwarmType 40 @*/ 41 PETSC_EXTERN PetscErrorCode DMSwarmSetPointsUniformCoordinates(DM dm,PetscReal min[],PetscReal max[],PetscInt npoints[],InsertMode mode) 42 { 43 PetscErrorCode ierr; 44 PetscReal gmin[] = {PETSC_MAX_REAL ,PETSC_MAX_REAL, PETSC_MAX_REAL}; 45 PetscReal gmax[] = {PETSC_MIN_REAL, PETSC_MIN_REAL, PETSC_MIN_REAL}; 46 PetscInt i,j,k,N,bs,b,n_estimate,n_curr,n_new_est,p,n_found; 47 Vec coorlocal; 48 const PetscScalar *_coor; 49 DM celldm; 50 PetscReal dx[3]; 51 Vec pos; 52 PetscScalar *_pos; 53 PetscReal *swarm_coor; 54 PetscInt *swarm_cellid; 55 PetscSF sfcell = NULL; 56 const PetscSFNode *LA_sfcell; 57 58 PetscFunctionBegin; 59 DMSWARMPICVALID(dm); 60 ierr = DMSwarmGetCellDM(dm,&celldm);CHKERRQ(ierr); 61 ierr = DMGetCoordinatesLocal(celldm,&coorlocal);CHKERRQ(ierr); 62 ierr = VecGetSize(coorlocal,&N);CHKERRQ(ierr); 63 ierr = VecGetBlockSize(coorlocal,&bs);CHKERRQ(ierr); 64 N = N / bs; 65 ierr = VecGetArrayRead(coorlocal,&_coor);CHKERRQ(ierr); 66 for (i=0; i<N; i++) { 67 for (b=0; b<bs; b++) { 68 gmin[b] = PetscMin(gmin[b],_coor[bs*i+b]); 69 gmax[b] = PetscMax(gmax[b],_coor[bs*i+b]); 70 } 71 } 72 ierr = VecRestoreArrayRead(coorlocal,&_coor);CHKERRQ(ierr); 73 74 for (b=0; b<bs; b++) { 75 dx[b] = (max[b] - min[b])/((PetscReal)(npoints[b]-1)); 76 } 77 78 /* determine number of points living in the bounding box */ 79 n_estimate = 0; 80 if (bs == 2) { npoints[2] = 1; } 81 for (k=0; k<npoints[2]; k++) { 82 for (j=0; j<npoints[1]; j++) { 83 for (i=0; i<npoints[0]; i++) { 84 PetscReal xp[] = {0.0,0.0,0.0}; 85 PetscInt ijk[3]; 86 PetscBool point_inside = PETSC_TRUE; 87 88 ijk[0] = i; 89 ijk[1] = j; 90 ijk[2] = k; 91 for (b=0; b<bs; b++) { 92 xp[b] = min[b] + ijk[b] * dx[b]; 93 } 94 for (b=0; b<bs; b++) { 95 if (xp[b] < gmin[b]) { point_inside = PETSC_FALSE; } 96 if (xp[b] > gmax[b]) { point_inside = PETSC_FALSE; } 97 } 98 if (point_inside) { n_estimate++; } 99 } 100 } 101 } 102 103 /* create candidate list */ 104 ierr = VecCreate(PetscObjectComm((PetscObject)dm),&pos);CHKERRQ(ierr); 105 ierr = VecSetSizes(pos,bs*n_estimate,PETSC_DECIDE);CHKERRQ(ierr); 106 ierr = VecSetBlockSize(pos,bs);CHKERRQ(ierr); 107 ierr = VecSetFromOptions(pos);CHKERRQ(ierr); 108 ierr = VecGetArray(pos,&_pos);CHKERRQ(ierr); 109 110 n_estimate = 0; 111 for (k=0; k<npoints[2]; k++) { 112 for (j=0; j<npoints[1]; j++) { 113 for (i=0; i<npoints[0]; i++) { 114 PetscReal xp[] = {0.0,0.0,0.0}; 115 PetscInt ijk[3]; 116 PetscBool point_inside = PETSC_TRUE; 117 118 ijk[0] = i; 119 ijk[1] = j; 120 ijk[2] = k; 121 for (b=0; b<bs; b++) { 122 xp[b] = min[b] + ijk[b] * dx[b]; 123 } 124 for (b=0; b<bs; b++) { 125 if (xp[b] < gmin[b]) { point_inside = PETSC_FALSE; } 126 if (xp[b] > gmax[b]) { point_inside = PETSC_FALSE; } 127 } 128 if (point_inside) { 129 for (b=0; b<bs; b++) { 130 _pos[bs*n_estimate+b] = xp[b]; 131 } 132 n_estimate++; 133 } 134 } 135 } 136 } 137 ierr = VecRestoreArray(pos,&_pos);CHKERRQ(ierr); 138 139 /* locate points */ 140 ierr = DMLocatePoints(celldm,pos,DM_POINTLOCATION_NONE,&sfcell);CHKERRQ(ierr); 141 142 ierr = PetscSFGetGraph(sfcell, NULL, NULL, NULL, &LA_sfcell);CHKERRQ(ierr); 143 n_found = 0; 144 for (p=0; p<n_estimate; p++) { 145 if (LA_sfcell[p].index != DMLOCATEPOINT_POINT_NOT_FOUND) { 146 n_found++; 147 } 148 } 149 150 /* adjust size */ 151 if (mode == ADD_VALUES) { 152 ierr = DMSwarmGetLocalSize(dm,&n_curr);CHKERRQ(ierr); 153 n_new_est = n_curr + n_found; 154 ierr = DMSwarmSetLocalSizes(dm,n_new_est,-1);CHKERRQ(ierr); 155 } 156 if (mode == INSERT_VALUES) { 157 n_curr = 0; 158 n_new_est = n_found; 159 ierr = DMSwarmSetLocalSizes(dm,n_new_est,-1);CHKERRQ(ierr); 160 } 161 162 /* initialize new coords, cell owners, pid */ 163 ierr = VecGetArrayRead(pos,&_coor);CHKERRQ(ierr); 164 ierr = DMSwarmGetField(dm,DMSwarmPICField_coor,NULL,NULL,(void**)&swarm_coor);CHKERRQ(ierr); 165 ierr = DMSwarmGetField(dm,DMSwarmPICField_cellid,NULL,NULL,(void**)&swarm_cellid);CHKERRQ(ierr); 166 n_found = 0; 167 for (p=0; p<n_estimate; p++) { 168 if (LA_sfcell[p].index != DMLOCATEPOINT_POINT_NOT_FOUND) { 169 for (b=0; b<bs; b++) { 170 swarm_coor[bs*(n_curr + n_found) + b] = _coor[bs*p+b]; 171 } 172 swarm_cellid[n_curr + n_found] = LA_sfcell[p].index; 173 n_found++; 174 } 175 } 176 ierr = DMSwarmRestoreField(dm,DMSwarmPICField_cellid,NULL,NULL,(void**)&swarm_cellid);CHKERRQ(ierr); 177 ierr = DMSwarmRestoreField(dm,DMSwarmPICField_coor,NULL,NULL,(void**)&swarm_coor);CHKERRQ(ierr); 178 ierr = VecRestoreArrayRead(pos,&_coor);CHKERRQ(ierr); 179 180 ierr = PetscSFDestroy(&sfcell);CHKERRQ(ierr); 181 ierr = VecDestroy(&pos);CHKERRQ(ierr); 182 183 PetscFunctionReturn(0); 184 } 185 186 /*@C 187 DMSwarmSetPointCoordinates - Set point coordinates in a DMSwarm from a user defined list 188 189 Collective on DM 190 191 Input parameters: 192 + dm - the DMSwarm 193 . npoints - the number of points to insert 194 . coor - the coordinate values 195 . redundant - if set to PETSC_TRUE, it is assumed that npoints and coor[] are only valid on rank 0 and should be broadcast to other ranks 196 - mode - indicates whether to append points to the swarm (ADD_VALUES), or over-ride existing points (INSERT_VALUES) 197 198 Level: beginner 199 200 Notes: 201 If the user has specified redundant = PETSC_FALSE, the cell DM will attempt to locate the coordinates provided by coor[] within 202 its sub-domain. If they any values within coor[] are not located in the sub-domain, they will be ignored and will not get 203 added to the DMSwarm. 204 205 .seealso: DMSwarmSetType(), DMSwarmSetCellDM(), DMSwarmType, DMSwarmSetPointsUniformCoordinates() 206 @*/ 207 PETSC_EXTERN PetscErrorCode DMSwarmSetPointCoordinates(DM dm,PetscInt npoints,PetscReal coor[],PetscBool redundant,InsertMode mode) 208 { 209 PetscErrorCode ierr; 210 PetscReal gmin[] = {PETSC_MAX_REAL ,PETSC_MAX_REAL, PETSC_MAX_REAL}; 211 PetscReal gmax[] = {PETSC_MIN_REAL, PETSC_MIN_REAL, PETSC_MIN_REAL}; 212 PetscInt i,N,bs,b,n_estimate,n_curr,n_new_est,p,n_found; 213 Vec coorlocal; 214 const PetscScalar *_coor; 215 DM celldm; 216 Vec pos; 217 PetscScalar *_pos; 218 PetscReal *swarm_coor; 219 PetscInt *swarm_cellid; 220 PetscSF sfcell = NULL; 221 const PetscSFNode *LA_sfcell; 222 PetscReal *my_coor; 223 PetscInt my_npoints; 224 PetscMPIInt rank; 225 MPI_Comm comm; 226 227 PetscFunctionBegin; 228 DMSWARMPICVALID(dm); 229 ierr = PetscObjectGetComm((PetscObject)dm,&comm);CHKERRQ(ierr); 230 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 231 232 ierr = DMSwarmGetCellDM(dm,&celldm);CHKERRQ(ierr); 233 ierr = DMGetCoordinatesLocal(celldm,&coorlocal);CHKERRQ(ierr); 234 ierr = VecGetSize(coorlocal,&N);CHKERRQ(ierr); 235 ierr = VecGetBlockSize(coorlocal,&bs);CHKERRQ(ierr); 236 N = N / bs; 237 ierr = VecGetArrayRead(coorlocal,&_coor);CHKERRQ(ierr); 238 for (i=0; i<N; i++) { 239 for (b=0; b<bs; b++) { 240 gmin[b] = PetscMin(gmin[b],_coor[bs*i+b]); 241 gmax[b] = PetscMax(gmax[b],_coor[bs*i+b]); 242 } 243 } 244 ierr = VecRestoreArrayRead(coorlocal,&_coor);CHKERRQ(ierr); 245 246 /* broadcast points from rank 0 if requested */ 247 if (redundant) { 248 my_npoints = npoints; 249 ierr = MPI_Bcast(&my_npoints,1,MPIU_INT,0,comm);CHKERRQ(ierr); 250 251 if (rank > 0) { /* allocate space */ 252 ierr = PetscMalloc1(my_npoints,&my_coor);CHKERRQ(ierr); 253 } else { 254 my_coor = coor; 255 } 256 ierr = MPI_Bcast(my_coor,bs*my_npoints,MPIU_REAL,0,comm);CHKERRQ(ierr); 257 } else { 258 my_npoints = npoints; 259 my_coor = coor; 260 } 261 262 /* determine the number of points living in the bounding box */ 263 n_estimate = 0; 264 for (i=0; i<my_npoints; i++) { 265 PetscBool point_inside = PETSC_TRUE; 266 267 for (b=0; b<bs; b++) { 268 if (my_coor[bs*i+b] < gmin[b]) { point_inside = PETSC_FALSE; } 269 if (my_coor[bs*i+b] > gmax[b]) { point_inside = PETSC_FALSE; } 270 } 271 if (point_inside) { n_estimate++; } 272 } 273 274 /* create candidate list */ 275 ierr = VecCreate(PetscObjectComm((PetscObject)dm),&pos);CHKERRQ(ierr); 276 ierr = VecSetSizes(pos,bs*n_estimate,PETSC_DECIDE);CHKERRQ(ierr); 277 ierr = VecSetBlockSize(pos,bs);CHKERRQ(ierr); 278 ierr = VecSetFromOptions(pos);CHKERRQ(ierr); 279 ierr = VecGetArray(pos,&_pos);CHKERRQ(ierr); 280 281 n_estimate = 0; 282 for (i=0; i<my_npoints; i++) { 283 PetscBool point_inside = PETSC_TRUE; 284 285 for (b=0; b<bs; b++) { 286 if (my_coor[bs*i+b] < gmin[b]) { point_inside = PETSC_FALSE; } 287 if (my_coor[bs*i+b] > gmax[b]) { point_inside = PETSC_FALSE; } 288 } 289 if (point_inside) { 290 for (b=0; b<bs; b++) { 291 _pos[bs*n_estimate+b] = my_coor[bs*i+b]; 292 } 293 n_estimate++; 294 } 295 } 296 ierr = VecRestoreArray(pos,&_pos);CHKERRQ(ierr); 297 298 /* locate points */ 299 ierr = DMLocatePoints(celldm,pos,DM_POINTLOCATION_NONE,&sfcell);CHKERRQ(ierr); 300 301 ierr = PetscSFGetGraph(sfcell, NULL, NULL, NULL, &LA_sfcell);CHKERRQ(ierr); 302 n_found = 0; 303 for (p=0; p<n_estimate; p++) { 304 if (LA_sfcell[p].index != DMLOCATEPOINT_POINT_NOT_FOUND) { 305 n_found++; 306 } 307 } 308 309 /* adjust size */ 310 if (mode == ADD_VALUES) { 311 ierr = DMSwarmGetLocalSize(dm,&n_curr);CHKERRQ(ierr); 312 n_new_est = n_curr + n_found; 313 ierr = DMSwarmSetLocalSizes(dm,n_new_est,-1);CHKERRQ(ierr); 314 } 315 if (mode == INSERT_VALUES) { 316 n_curr = 0; 317 n_new_est = n_found; 318 ierr = DMSwarmSetLocalSizes(dm,n_new_est,-1);CHKERRQ(ierr); 319 } 320 321 /* initialize new coords, cell owners, pid */ 322 ierr = VecGetArrayRead(pos,&_coor);CHKERRQ(ierr); 323 ierr = DMSwarmGetField(dm,DMSwarmPICField_coor,NULL,NULL,(void**)&swarm_coor);CHKERRQ(ierr); 324 ierr = DMSwarmGetField(dm,DMSwarmPICField_cellid,NULL,NULL,(void**)&swarm_cellid);CHKERRQ(ierr); 325 n_found = 0; 326 for (p=0; p<n_estimate; p++) { 327 if (LA_sfcell[p].index != DMLOCATEPOINT_POINT_NOT_FOUND) { 328 for (b=0; b<bs; b++) { 329 swarm_coor[bs*(n_curr + n_found) + b] = _coor[bs*p+b]; 330 } 331 swarm_cellid[n_curr + n_found] = LA_sfcell[p].index; 332 n_found++; 333 } 334 } 335 ierr = DMSwarmRestoreField(dm,DMSwarmPICField_cellid,NULL,NULL,(void**)&swarm_cellid);CHKERRQ(ierr); 336 ierr = DMSwarmRestoreField(dm,DMSwarmPICField_coor,NULL,NULL,(void**)&swarm_coor);CHKERRQ(ierr); 337 ierr = VecRestoreArrayRead(pos,&_coor);CHKERRQ(ierr); 338 339 if (redundant) { 340 if (rank > 0) { 341 ierr = PetscFree(my_coor);CHKERRQ(ierr); 342 } 343 } 344 ierr = PetscSFDestroy(&sfcell);CHKERRQ(ierr); 345 ierr = VecDestroy(&pos);CHKERRQ(ierr); 346 347 PetscFunctionReturn(0); 348 } 349 350 extern PetscErrorCode private_DMSwarmInsertPointsUsingCellDM_DA(DM,DM,DMSwarmPICLayoutType,PetscInt); 351 extern PetscErrorCode private_DMSwarmInsertPointsUsingCellDM_PLEX(DM,DM,DMSwarmPICLayoutType,PetscInt); 352 353 /*@C 354 DMSwarmInsertPointsUsingCellDM - Insert point coordinates within each cell 355 356 Not collective 357 358 Input parameters: 359 + dm - the DMSwarm 360 . layout_type - method used to fill each cell with the cell DM 361 - fill_param - parameter controlling how many points per cell are added (the meaning of this parameter is dependent on the layout type) 362 363 Level: beginner 364 365 Notes: 366 The insert method will reset any previous defined points within the DMSwarm 367 368 .seealso: DMSwarmPICLayoutType, DMSwarmSetType(), DMSwarmSetCellDM(), DMSwarmType 369 @*/ 370 PETSC_EXTERN PetscErrorCode DMSwarmInsertPointsUsingCellDM(DM dm,DMSwarmPICLayoutType layout_type,PetscInt fill_param) 371 { 372 PetscErrorCode ierr; 373 DM celldm; 374 PetscBool isDA,isPLEX; 375 376 PetscFunctionBegin; 377 DMSWARMPICVALID(dm); 378 ierr = DMSwarmGetCellDM(dm,&celldm);CHKERRQ(ierr); 379 ierr = PetscObjectTypeCompare((PetscObject)celldm,DMDA,&isDA);CHKERRQ(ierr); 380 ierr = PetscObjectTypeCompare((PetscObject)celldm,DMPLEX,&isPLEX);CHKERRQ(ierr); 381 if (isDA) { 382 ierr = private_DMSwarmInsertPointsUsingCellDM_DA(dm,celldm,layout_type,fill_param);CHKERRQ(ierr); 383 } else if (isPLEX) { 384 ierr = private_DMSwarmInsertPointsUsingCellDM_PLEX(dm,celldm,layout_type,fill_param);CHKERRQ(ierr); 385 } else SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Only supported for cell DMs of type DMDA and DMPLEX"); 386 387 PetscFunctionReturn(0); 388 } 389 390 /* 391 PETSC_EXTERN PetscErrorCode DMSwarmAddPointCoordinatesCellWise(DM dm,PetscInt cell,PetscInt npoints,PetscReal xi[],PetscBool proximity_initialization) 392 { 393 PetscFunctionBegin; 394 PetscFunctionReturn(0); 395 } 396 */ 397 398 /* Field projection API */ 399 /* 400 PETSC_EXTERN PetscErrorCode DMSwarmProjectFields(DM dm,PetscInt project_type,PetscInt nfields,const char *fieldnames[],Vec *fields) 401 { 402 PetscFunctionBegin; 403 PetscFunctionReturn(0); 404 } 405 */ 406 407 /*@C 408 DMSwarmCreatePointPerCellCount - Count the number of points within all cells in the cell DM 409 410 Not collective 411 412 Input parameter: 413 . dm - the DMSwarm 414 415 Output parameters: 416 + ncells - the number of cells in the cell DM (optional argument, pass NULL to ignore) 417 - count - array of length ncells containing the number of points per cell 418 419 Level: beginner 420 421 Notes: 422 The array count is allocated internally and must be free'd by the user. 423 424 .seealso: DMSwarmSetType(), DMSwarmSetCellDM(), DMSwarmType 425 @*/ 426 PETSC_EXTERN PetscErrorCode DMSwarmCreatePointPerCellCount(DM dm,PetscInt *ncells,PetscInt **count) 427 { 428 PetscErrorCode ierr; 429 PetscBool isvalid; 430 PetscInt nel; 431 PetscInt *sum; 432 433 PetscFunctionBegin; 434 ierr = DMSwarmSortGetIsValid(dm,&isvalid);CHKERRQ(ierr); 435 nel = 0; 436 if (isvalid) { 437 PetscInt e; 438 439 ierr = DMSwarmSortGetSizes(dm,&nel,NULL);CHKERRQ(ierr); 440 441 ierr = PetscMalloc1(nel,&sum);CHKERRQ(ierr); 442 for (e=0; e<nel; e++) { 443 ierr = DMSwarmSortGetNumberOfPointsPerCell(dm,e,&sum[e]);CHKERRQ(ierr); 444 } 445 } else { 446 DM celldm; 447 PetscBool isda,isplex,isshell; 448 PetscInt p,npoints; 449 PetscInt *swarm_cellid; 450 451 /* get the number of cells */ 452 ierr = DMSwarmGetCellDM(dm,&celldm);CHKERRQ(ierr); 453 ierr = PetscObjectTypeCompare((PetscObject)celldm,DMDA,&isda);CHKERRQ(ierr); 454 ierr = PetscObjectTypeCompare((PetscObject)celldm,DMPLEX,&isplex);CHKERRQ(ierr); 455 ierr = PetscObjectTypeCompare((PetscObject)celldm,DMSHELL,&isshell);CHKERRQ(ierr); 456 if (isda) { 457 PetscInt _nel,_npe; 458 const PetscInt *_element; 459 460 ierr = DMDAGetElements(celldm,&_nel,&_npe,&_element);CHKERRQ(ierr); 461 nel = _nel; 462 ierr = DMDARestoreElements(celldm,&_nel,&_npe,&_element);CHKERRQ(ierr); 463 } else if (isplex) { 464 PetscInt ps,pe; 465 466 ierr = DMPlexGetHeightStratum(celldm,0,&ps,&pe);CHKERRQ(ierr); 467 nel = pe - ps; 468 } else if (isshell) { 469 PetscErrorCode (*method_DMShellGetNumberOfCells)(DM,PetscInt*); 470 471 ierr = PetscObjectQueryFunction((PetscObject)celldm,"DMGetNumberOfCells_C",&method_DMShellGetNumberOfCells);CHKERRQ(ierr); 472 if (method_DMShellGetNumberOfCells) { 473 ierr = method_DMShellGetNumberOfCells(celldm,&nel);CHKERRQ(ierr); 474 } else SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Cannot determine the number of cells for the DMSHELL object. User must provide a method via PetscObjectComposeFunction( (PetscObject)shelldm, \"DMGetNumberOfCells_C\", your_function_to_compute_number_of_cells );"); 475 } else SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Cannot determine the number of cells for a DM not of type DA, PLEX or SHELL"); 476 477 ierr = PetscMalloc1(nel,&sum);CHKERRQ(ierr); 478 ierr = PetscMemzero(sum,sizeof(PetscInt)*nel);CHKERRQ(ierr); 479 ierr = DMSwarmGetLocalSize(dm,&npoints);CHKERRQ(ierr); 480 ierr = DMSwarmGetField(dm,DMSwarmPICField_cellid,NULL,NULL,(void**)&swarm_cellid);CHKERRQ(ierr); 481 for (p=0; p<npoints; p++) { 482 if (swarm_cellid[p] != DMLOCATEPOINT_POINT_NOT_FOUND) { 483 sum[ swarm_cellid[p] ]++; 484 } 485 } 486 ierr = DMSwarmRestoreField(dm,DMSwarmPICField_cellid,NULL,NULL,(void**)&swarm_cellid);CHKERRQ(ierr); 487 } 488 if (ncells) { *ncells = nel; } 489 *count = sum; 490 PetscFunctionReturn(0); 491 } 492