Quickly finding the number of neighbors that cells has

[cococastano]

I am trying to write a cellular automata routine where only the states of cells at the outer surface of my unstructured grid are changed with each time step. I saw a post by @banesullivan that I thought may work if I could determine whether cells had less than 6 neighbors. Unfortunately this implementation would take too long given the number of time steps I need to take. There must be some way to query the number of neighbors a cell has, i.e. the number of faces a cell shares. I have looked around a lot for this and am fairly new to pyvista. I would really appreciate some guidance! Its possible I am just not familiar with the syntax

I tried using the GetCellNeighbors method but the neighbor list is being returned empty:

cell_idx = 1000  # try for neighbors of cell 1000
cell = grid.GetCell(cell_idx)
pt_ids = cell.GetPointIds()  # get point ids of cell of interest
neighbors = pv.vtk.vtkIdList()  # make vtk list to write neighbors to (?)
grid.GetCellNeighbors(cell_idx,pt_ids,neighbors)
print(neighbors)

Thank you, Nic

[banesullivan]

This example might help on how to properly use GetCellNeighbors: https://lorensen.github.io/VTKExamples/site/Cxx/PolyData/CellPointNeighbors/

We don't implement any ways to do this out of the box with PyVista. The closest thing I can think of is to follow what I outline here: https://github.com/pyvista/pyvista-support/issues/96#issuecomment-571864471 . There are some enumerations (python for-loops) in there, so performance may be an issue for you... also that code is only valid for PolyData.

---

Here is an example that uses our extract_points method - this performs a selection on the mesh for all cells containing the give point IDs. Simply feed it the point IDs of the cell of interest and you'll get that cell's neighbors. I believe this should have good performance because we are only running the vtkExtractSelection filter under the hood.

import pyvista as pv
from pyvista import examples
import numpy as np

def get_neighbor_cell_ids(grid, cell_idx):
    """Helper to get neighbor cell IDs."""
    cell = grid.GetCell(cell_idx)
    pids = pv.vtk_id_list_to_array(cell.GetPointIds())
    neighbors = set(grid.extract_points(pids)["vtkOriginalCellIds"])
    neighbors.discard(cell_idx)
    return np.array(list(neighbors))

grid = examples.load_hexbeam()
cell_idx = 10
neighbors = get_neighbor_cell_ids(grid, cell_idx)

# And plot to show
p = pv.Plotter(notebook=0)
p.add_mesh(grid.extract_all_edges(), color='k', label="whole mesh")
p.add_mesh(grid.extract_cells(neighbors), color=True, opacity=0.5, label="neighbors")
p.add_mesh(grid.extract_cells(cell_idx), color='pink', opacity=0.75, label="the cell")
p.add_legend()
p.show()