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