Solve TSP without crossing through the object

I have a grid of points which are forming a cube. The cube looks like this:

Now I want to create a path through all the points, so I had the following method which applies TSP and optimizes the path with 2-Opt:

def tsp_2opt(points):
    # Return input when only one point is given
    if len(points) <= 1:
        return points

    # Create a distance matrix between all points
    dist_matrix = np.sqrt(np.sum((points[:, np.newaxis] - points) ** 2, axis=2))

    # Create an initial tour using the nearest neighbor algorithm
    n = len(points)
    unvisited = set(range(n))
    curr_point = 0  # start with the first point
    tour = [curr_point]
    unvisited.remove(curr_point)
    while unvisited:
        next_point = min(unvisited, key=lambda x: dist_matrix[curr_point, x])
        if len(unvisited) == 1:
            tour.append(next_point)
            break
        tour.append(next_point)
        unvisited.remove(next_point)
        curr_point = next_point

    # Use 2-Opt algorithm to improve the tour
    improved = True
    while improved:
        improved = False
        for i in range(n-2):
            for j in range(i+2, n-1):
                if dist_matrix[tour[i], tour[j]] + dist_matrix[tour[i+1], tour[j+1]] < dist_matrix[tour[i], tour[i+1]] + dist_matrix[tour[j], tour[j+1]]:
                    tour[i+1:j+1] = reversed(tour[i+1:j+1])
                    improved = True

    return points[np.array(tour)]

This resulted in:

This is almost what I want, except that this breaks one rule that my model should prevent. The path is not allowed to go through the object, so the last part it crosses from middle to middle, because it has no other option to reach that last point. I tried setting the distance matrix to infinity for points on a different surface, so that it would not try to go there, but in my case it has no other option.
These are my data points:

[[  0. 100. 100.]
 [  0. 100.  50.]
 [  0. 100.   0.]
 [ 50. 100. 100.]
 [ 50. 100.  50.]
 [ 50. 100.   0.]
 [100. 100. 100.]
 [100. 100.  50.]
 [100. 100.   0.]
 [100.   0.   0.]
 [100.  50.   0.]
 [100. 100.   0.]
 [100.   0.  50.]
 [100.  50.  50.]
 [100. 100.  50.]
 [100.   0. 100.]
 [100.  50. 100.]
 [100. 100. 100.]
 [  0.   0.   0.]
 [  0.  50.   0.]
 [  0. 100.   0.]
 [  0.   0.  50.]
 [  0.  50.  50.]
 [  0. 100.  50.]
 [  0.   0. 100.]
 [  0.  50. 100.]
 [  0. 100. 100.]
 [  0.   0. 100.]
 [  0.  50. 100.]
 [  0. 100. 100.]
 [ 50.   0. 100.]
 [ 50.  50. 100.]
 [ 50. 100. 100.]
 [100.   0. 100.]
 [100.  50. 100.]
 [100. 100. 100.]
 [  0.   0.   0.]
 [  0.  50.   0.]
 [  0. 100.   0.]
 [ 50.   0.   0.]
 [ 50.  50.   0.]
 [ 50. 100.   0.]
 [100.   0.   0.]
 [100.  50.   0.]
 [100. 100.   0.]
 [  0.   0. 100.]
 [  0.   0.  50.]
 [  0.   0.   0.]
 [ 50.   0. 100.]
 [ 50.   0.  50.]
 [ 50.   0.   0.]
 [100.   0. 100.]
 [100.   0.  50.]
 [100.   0.   0.]]

First I created a check that it is only allowed to move one axis at the time, so that these cross paths would not occur, but the problem is that that solution only works for cubes. I need to make this work for other shapes as well.

Because you do not constrain the distance from the current point to the third point. You should have used 3D search, but you used 2D, so there is a constraint missing.

The constraint from the first point to the second point is that the distance is equal to 50, and the constraint from the first point to the third point is that the distance is greater than 50 and less than 100.

One observation: as currently coded, your algorithm (both the first part where you build the initial path and the second part where you try to optimize that path) are perfectly fine with producing paths with edges going "through" the object.

Your goal seems to be to have a final path that lies on the surface of your shape. One way to achieve that could be to only deal with paths already lying on the surface (both initial path and optimization).

(note that this constraint may prevent exploring some possible optimizations at the end, I haven't looked into that).

for example:

• either add a second connected matrix, or an adjacency list for each node, or set dist to "-1" or to "None" for non forbidden paths,
• on top of checking the distance between nodes, check that going to that neighbor is allowed (e.g: your condition to create the initial path will be a little bit more complex than the one-liner min(unvisited, key=lambda x: dist_matrix[curr_point, x]), and the if condition of your optimizing loop will check a few more things)

In order to generate the initial path, you will also need to change your algorithm :

greedily selecting one neighbor to your current node may lead you to dead ends (e.g: current in the middle of a face, with all surrounding nodes explored, but some other nodes still unvisited on the opposite face)

You may try :

• backtracking
• try to include unvisited nodes by expanding a path locally (e.g: replace a sgement ...a--b... with a path ...a--e--f--b...)