Shortest Hamilton path

1,128 views
Skip to first unread message

Riyas M. J

unread,
Jan 24, 2022, 4:52:23 AM1/24/22
to or-tools-discuss
Dear all,

I would like to find the shortest path that connects all the points once. Theoretically, it is called the shortest hamiltonian path (SHP). The SHP is similar to the Traveling Salesperson Problem (TSP). In TSP, the endpoint should be equal to the starting point (the salesman gets back to his office). Whereas, in SHP, the starting and endpoint can be different. 

I would like to perform the SHP for around 1000 points in Python programming. Even though I am able to follow the python example of TSP (TPS using python), I am finding a hard time performing SHP in python. 

Kindly help me in this regard. I have tried to achieve this using various tools including the networkx package, Delaunay triangulation using scipy, and brute force. This is part of my project and I am running out of time after all the failures. So please help me achieve SHP using or-tools in python. 
Thank you.

Best regards
MJ

Frederic Didier

unread,
Jan 24, 2022, 5:17:26 AM1/24/22
to or-tools...@googlegroups.com
You can transform a SHP problem into a TSP one simply by adding one extra node connected to all other nodes with arcs of weight zero.
You can then reuse the TSP example. 
That said, 1k nodes starts to be big, especially if you feed the fully connected graph.

--
You received this message because you are subscribed to the Google Groups "or-tools-discuss" group.
To unsubscribe from this group and stop receiving emails from it, send an email to or-tools-discu...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/or-tools-discuss/f4aa7104-97fa-4381-bc61-2239ae888b9an%40googlegroups.com.

Riyas MJ

unread,
Jan 24, 2022, 5:42:24 AM1/24/22
to or-tools-discuss
Thank you @Frederic Didier for the response.

1000 points in a fully connected graph is heavy indeed. However, this is a part of a bigger project and it can handle the computation. After exploring various other options, I found this as the best for my need. 

Did you mean by adding another point in the distance matrix and specifying zero as the distance to each other points? I am giving an example of what I am trying to say. 

## If the below is the original distance-matrix
distance_matrix = {
0: {0: 0, 1: 125, 2: 50, 3: 90, 4: 28, 5: 257},
 1: {0: 125, 1: 0, 2: 75, 3: 38, 4: 154, 5: 134},
 2: {0: 50, 1: 75, 2: 0, 3: 40, 4: 78, 5: 207},
 3: {0: 90, 1: 38, 2: 40, 3: 0, 4: 118, 5: 166},
 4: {0: 28, 1: 154, 2: 78, 3: 118, 4: 0, 5: 285},
 5: {0: 257, 1: 134, 2: 207, 3: 166, 4: 285, 5: 0}}

## Add an additional point to and assign zero distance with all other points such as
distance_matrix = {
0: {0: 0, 1: 125, 2: 50, 3: 90, 4: 28, 5: 257, 6: 0},
 1: {0: 125, 1: 0, 2: 75, 3: 38, 4: 154, 5: 134, 6: 0},
 2: {0: 50, 1: 75, 2: 0, 3: 40, 4: 78, 5: 207, 6: 0},
 3: {0: 90, 1: 38, 2: 40, 3: 0, 4: 118, 5: 166, 6: 0},
 4: {0: 28, 1: 154, 2: 78, 3: 118, 4: 0, 5: 285, 6: 0},
 5: {0: 257, 1: 134, 2: 207, 3: 166, 4: 285, 5: 0, 6:0},
 6: {0: 0, 1: 0, 2: 0, 3: 0, 4: 0, 5: 0, 6: 0}}

Thank you
MJ

Mizux Seiha

unread,
Jan 24, 2022, 6:13:52 AM1/24/22
to or-tools-discuss
yes and don't forget to set 6 as your depot when instantiating the RoutingIndexManager

note: an other way is to perform a check in your transit callback if you don't want/can touch the matrix
```python
def my_transit_callback(from_index, to_index):
"""Returns the distance between the two nodes."""

# Convert from routing variable Index to distance matrix NodeIndex.

from_node = manager.IndexToNode(from_index)

to_node = manager.IndexToNode(to_index)
if from_node == 6 or to_node == 6: return 0
return data['distance_matrix'][from_node][to_node] # you can keep you matrix 6x6 aka idx range [0,5]x[0,5]

```

Riyas MJ

unread,
Jan 24, 2022, 8:48:51 AM1/24/22
to or-tools-discuss
Thank you for the confirmation @mizu. 
However, it didn't result well. The result I am getting is a TSP without the final joining with the initial point. It is a path but this is not the shortest Hamiltonian path. The result shows the inherent nature of TSP to come back to its origin and that results in a longer journey than the shortest path covering all the points. A sample image is attached. From the image, it is clear that the initial and final points are 1 node difference (in every trial).

I am also attaching the actual working code of the project if any reference is required. 
The project actually divides a large area into several small ones and makes 5 random points on each of them. Thereafter, it is trying to find the shortest path touching every node.
Thank you.
tsp_fail.png
working_code.py

Mizux Seiha

unread,
Jan 25, 2022, 3:14:27 AM1/25/22
to or-tools-discuss
AFAIK you are using only the initial search
```py
    # Setting first solution heuristic.
    search_parameters = pywrapcp.DefaultRoutingSearchParameters()
    search_parameters.first_solution_strategy = (
        routing_enums_pb2.FirstSolutionStrategy.PATH_CHEAPEST_ARC)

    # Solve the problem.
    solution = routing.SolveWithParameters(search_parameters)
```
If  I were you I would try to also enable the GLS few seconds
```python
    search_parameters.local_search_metaheuristic = (
        routing_enums_pb2.LocalSearchMetaheuristic.GUIDED_LOCAL_SEARCH)
    search_parameters.time_limit.FromSeconds(30)
```

Mizux Seiha

unread,
Jan 25, 2022, 4:56:39 AM1/25/22
to or-tools-discuss
Figure_1.png
After few test it seems you had few other mistakes:
* depot should be the last node of index len(PointsAll)
* my_route should be route[0][1:-1]
* the code provided overwrite the distance_matrix you hacked few line above -> erasing the dummy node in the distance matrix

my diff:

```diff
%diff -u --color working_code.py plop.py
--- working_code.py        2022-01-25 10:55:21.805758120 +0100
+++ plop.py        2022-01-25 10:49:10.149326981 +0100
@@ -91,7 +91,8 @@
 ip_poly = 'POLYGON((7511352.944174381 7112117.215453198,7511653.217271762 7108850.245242323,7515833.018061557 7109162.529626478,7515556.765888274 7112537.598119422,7511352.944174381 7112117.215453198))'
 
 ip_shape = shapely.wkt.loads(ip_poly)
-polygonAll = grids(ip_shape, grid_size)
+ip_grid_size = 800
+polygonAll = grids(ip_shape, ip_grid_size)
 finalize = False
 len(polygonAll)
 
@@ -125,16 +126,6 @@
     return distances
 
 
-distance_matrix = compute_euclidean_distance_matrix(data['locations'])
-
-## Adding the dummy point at the end with zero distance to every point
-dict_dummy_pt = {}
-for num in range(len(pointsAll) + 1):
-    dict_dummy_pt[num] = 0
-distance_matrix[len(pointsAll)] = dict_dummy_pt
-for num in range(len(pointsAll) + 1):
-    distance_matrix[num][len(pointsAll)] = 0
-
 ## o-r Tools
 from ortools.constraint_solver import routing_enums_pb2
 from ortools.constraint_solver import pywrapcp
@@ -144,7 +135,7 @@
 def create_data_model():
     """Stores the data for the problem."""
     data['num_vehicles'] = 1
-    data['depot'] = 0
+    data['depot'] = len(pointsAll)
     return data
 
 
@@ -169,8 +160,21 @@
     """Entry point of the program."""
     # Instantiate the data problem.
     data = create_data_model()
+    print(f'depot: {data["depot"]}')
     distance_matrix = compute_euclidean_distance_matrix(data['locations'])
 
+    ## Adding the dummy point at the end with zero distance to every point
+    print(f'before distance matrix len: {len(distance_matrix)}')
+    dict_dummy_pt = {}
+    for num in range(len(pointsAll) + 1):
+        dict_dummy_pt[num] = 0
+    distance_matrix[len(pointsAll)] = dict_dummy_pt
+    for num in range(len(pointsAll) + 1):
+        distance_matrix[num][len(pointsAll)] = 0
+
+    print(f'after distance matrix len: {len(distance_matrix)}')
+    print(f'pointsAll len: {len(pointsAll)}')
+
     # Create the routing index manager.
     manager = pywrapcp.RoutingIndexManager(len(distance_matrix),
                                            data['num_vehicles'], data['depot'])
@@ -194,7 +198,10 @@

     search_parameters = pywrapcp.DefaultRoutingSearchParameters()
     search_parameters.first_solution_strategy = (
         routing_enums_pb2.FirstSolutionStrategy.PATH_CHEAPEST_ARC)
-
+    search_parameters.local_search_metaheuristic = (
+        routing_enums_pb2.LocalSearchMetaheuristic.GUIDED_LOCAL_SEARCH)
+    search_parameters.time_limit.FromSeconds(10)
+    search_parameters.log_search = True

     # Solve the problem.
     solution = routing.SolveWithParameters(search_parameters)
 
@@ -222,7 +229,7 @@
     #     x =0
     # # print(routes)
     global my_route
-    my_route = routes[0][:-1]
+    my_route = routes[0][1:-1]
 
     # print(my_route)
```

makansij

unread,
Jun 12, 2023, 5:18:13 PM6/12/23
to or-tools-discuss
Hi,

I tried to examine this example to answer my own question posted here: https://groups.google.com/g/or-tools-discuss/c/5VXTnoumQlQ

But I don't see how you obtained "routes" variable in your code above.  Can you tell me how to obtain "routes"?  Or where it is defined/declared?  

Thanks.

Reply all
Reply to author
Forward
0 new messages