I am working with a mesh of triangles (in 3D, although I doubt it makes a difference). The mesh is given as list of lists, each list containing the indices of the three vertices of a triangle in said mesh.

For example:

0 2 6
0 12 15
2 66 10234

I want to find the neighbours of each point, i.e., to have an array where the first row contains the neighbours of point 0, the second the neighbours of point 1 and so on.

How do I make it more efficient?

points is a list of lists containing the point coordinates, while triangles is a list of lists containing the triangles in the mesh:

from numpy import ma

pointsNb = points.shape[0]

trianglesNb = triangles.shape[0]

neighbours = [[] for i in range(pointsNb)] # for each point in points, make a list of its neighbours

# find neighbours
for i in range(0, trianglesNb):
    triangle = ma.array(triangles[i][:], mask=False)
    for j, point in enumerate(triangles[i]):
        triangle[j] = ma.masked
        for possible_neighbour in triangle[~triangle.mask]:
            if possible_neighbour not in neighbours[point]:
        triangle[j] = point
  • \$\begingroup\$ Please add a definition of ma. \$\endgroup\$
    – vnp
    Commented Aug 19, 2014 at 17:15
  • \$\begingroup\$ ma is the masked property for arrays in numpy \$\endgroup\$
    – John
    Commented Aug 19, 2014 at 22:07

1 Answer 1

  1. There's no documentation. What is this code supposed to do? How am I supposed to call it?

  2. There's no functional decomposition. Code is easier to understand and test if it's organized into functions with well-defined inputs and outputs.

  3. There's no inherent order to the neighbours of a point, so the natural representation of the neighbours is a set, not a list. Also, if you used a set, then you'd have a more efficient membership test, and you could update it with new points without having to check whether they are already present.

So I would write:

from collections import defaultdict

def find_neighbours(triangles):
    """Given an array of triangles (triples of points), return a graph in
    adjacency representation: that is, a map from point to the set of
    its neighbours.

    graph = defaultdict(set)
    for x, y, z in triangles:
        graph[x].update((y, z))
        graph[y].update((x, z))
        graph[z].update((x, y))
    return graph

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.