OR-Tools VRP Solver Under-performing

Question

I am trying to test or-tools routing solver to solve the basic TSP problem, but I have not been able to get it properly working. I have a distance matrix and a bunch of greedy solutions already generated before the problem is sent over to the routing solver. As an example, I have set up a problem using the sample python code from this website shared at the bottom.

In the example, I have 10 cities with an asymmetric distance matrix. The best greedy solution (closest-first starting from different cities) is stored as an initial solution in data. I have two functions: solve_from_initial_route() and solve_from_scratch() that solve the same problem with or without the information of an initial solution and also yield the same result. The solver exhibits some surprising behavior here:

• The function solve_from_initial_route() yields the initial greedy solution as the final solution and exits immediately (4-5 ms) without any attempt to solve and generate runtime logs (even though search logging is enabled).
• The function solve_from_scratch() also yields the same greedy solution as the final solution and does produce runtime logs showing that it has evaluated a lot of options. But the funny thing is that no matter how long I run the solver for, the solution is always the same. The solver is somehow not acting smartly and is always evaluating worse options. On the other hand, a genetic algorithm run on the same problem produces a much better solution than the greedy initial solution in under 1 second!
• I have also tried all other local search options and for randomized algorithms like simulated annealing, you would expect to see some variation in the result over different time limits, but that doesn't happen and the solution is always the same as the initial greedy solution!
• Also, the routing solver goes on even though it should technically have evaluated all combinations. With 10 cities there are a total of 10! = 3628800 sequences to be evaluated in total, but with a higher limit than that the solver keeps running till it has hit the search limit and not the actual problem limit.

It is possible that I have not set all options correctly or am missing something in my code. I would appreciate any help in making the solver work as expected.

Thanks!

!pip install ortools
from __future__ import print_function
from ortools.constraint_solver import pywrapcp
from ortools.constraint_solver import routing_enums_pb2

def create_data_model():
    """Stores the data for the problem."""
    data = {}
    data['distance_matrix'] = [
        [0, 227543, 133934, 200896, 106495, 163222, 75896, 139494, 46460, 102942],
        [135873, 0, 15673, 174874, 80474, 197318, 109993, 232377, 139343, 46665],
        [229482, 15673, 0, 88692, 183092, 125214, 214714, 153718, 247723, 140274],
        [108503, 174151, 80542, 0, 15674, 169948, 82622, 205007, 111973, 49550],
        [195308, 94193, 167348, 21174, 0, 105716, 169428, 134221, 198779, 136356],
        [77835, 203602, 109992, 176954, 82554, 0, 15660, 174340, 81306, 79000],
        [172835, 119784, 213500, 94785, 189185, 21172, 0, 96019, 190024, 174000],
        [48413, 232967, 139358, 206320, 111919, 168647, 81321, 0, 15662, 108366],
        [141422, 153773, 247490, 128774, 204928, 101504, 174329, 15662, 0, 201374],
        [104492, 139205, 45595, 143494, 49093, 165938, 78612, 200997, 107963, 0]
    ]
    data['initial_routes'] = [
        [8, 7, 6, 5, 4, 3, 2, 1]
    ]
    data['num_vehicles'] = 1
    data['start_idx'] = [0]
    data['end_idx'] = [9]
    return data

def print_solution(data, manager, routing, solution):
    """Prints solution on console."""
    max_route_distance = 0
    for vehicle_id in range(data['num_vehicles']):
        index = routing.Start(vehicle_id)
        plan_output = 'Route for vehicle {}:\n'.format(vehicle_id)
        route_distance = 0
        while not routing.IsEnd(index):
            plan_output += ' {} -> '.format(manager.IndexToNode(index))
            previous_index = index
            index = solution.Value(routing.NextVar(index))
            route_distance += routing.GetArcCostForVehicle(
                previous_index, index, vehicle_id)
        plan_output += '{}\n'.format(manager.IndexToNode(index))
        plan_output += 'Distance of the route: {}m\n'.format(route_distance)
        print(plan_output)
        max_route_distance = max(route_distance, max_route_distance)
    print('Maximum of the route distances: {}m'.format(max_route_distance))

def solve_from_initial_route():
    """Solve the CVRP problem."""
    # Instantiate the data problem.
    data = create_data_model()

    # Create the routing index manager.
    manager = pywrapcp.RoutingIndexManager(len(data['distance_matrix']),
                                           data['num_vehicles'], data['start_idx'],
                                           data['end_idx'])

    # Create Routing Model.
    routing = pywrapcp.RoutingModel(manager)


    # Create and register a transit callback.
    def distance_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)
        return data['distance_matrix'][from_node][to_node]

    transit_callback_index = routing.RegisterTransitCallback(distance_callback)

    # Define cost of each arc.
    routing.SetArcCostEvaluatorOfAllVehicles(transit_callback_index)

    initial_solution = routing.ReadAssignmentFromRoutes(data['initial_routes'],True)
    print('Initial solution:')
    print_solution(data, manager, routing, initial_solution)

    # Set default search parameters.
    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.seconds = 2
    search_parameters.lns_time_limit.seconds = 1
    search_parameters.solution_limit = 15000
    search_parameters.log_search = True

    # Solve the problem.
    solution = routing.SolveFromAssignmentWithParameters(initial_solution, search_parameters)

    # Print solution on console.
    if solution:
        print('Solution after search:')
        print_solution(data, manager, routing, solution)

def solve_from_scratch():
    """Solve the CVRP problem."""
    # Instantiate the data problem.
    data = create_data_model()

    # Create the routing index manager
    manager = pywrapcp.RoutingIndexManager(len(data['distance_matrix']),
                                           data['num_vehicles'], data['start_idx'],
                                           data['end_idx'])

    # Create Routing Model.
    routing = pywrapcp.RoutingModel(manager)

    # Create and register a transit callback.
    def distance_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)
        return data['distance_matrix'][from_node][to_node]

    transit_callback_index = routing.RegisterTransitCallback(distance_callback)

    # Define cost of each arc.
    routing.SetArcCostEvaluatorOfAllVehicles(transit_callback_index)

    # Set default search parameters.
    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.seconds = 2
    search_parameters.lns_time_limit.seconds = 1
    search_parameters.solution_limit = 150000
    search_parameters.log_search = True

    # Solve the problem.
    solution = routing.SolveWithParameters(search_parameters)

    # Print solution on console.
    if solution:
        print('Solution after search:')
        print_solution(data, manager, routing, solution)

if __name__ == '__main__':
    #solve_from_initial_route()
    solve_from_scratch()

Answer 1

If I run your code, it prints

Solution after search:
Route for vehicle 0:
 0 ->  8 ->  7 ->  6 ->  5 ->  4 ->  3 ->  2 ->  1 -> 9
Distance of the route: 411223m

This is the optimal path from 0 to 9 using all the nodes (I checked that using the tsp_sat code). And it is found in less than 1 ms.

Now, on the log part,

Solution #5265 (783010, objective minimum = 411223, objective maximum = 1103173, time = 1998 ms, branches = 26227, failures = 14386, depth = 33, OrOpt<3>, neighbors = 351904, filtered neighbors = 5265, accepted neighbors = 5265, memory used = 35.63 MB, limit = 99%)

GLS penalizes the cost function, so the actual value 783010 is not a true distance, but a penalized distance.

Now, solve_from_initial_route() hits a known bug

Here is the correct solve code

# Create the routing index manager.
manager = pywrapcp.RoutingIndexManager(len(data['distance_matrix']),
                                       data['num_vehicles'], 
                                       data['start_idx'],
                                       data['end_idx'])

# Create Routing Model.
routing = pywrapcp.RoutingModel(manager)


# Create and register a transit callback.
def distance_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)
    return data['distance_matrix'][from_node][to_node]

transit_callback_index = routing.RegisterTransitCallback(distance_callback)

# Define cost of each arc.
routing.SetArcCostEvaluatorOfAllVehicles(transit_callback_index)

# Set default search parameters.
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.seconds = 1
search_parameters.lns_time_limit.seconds = 1
search_parameters.solution_limit = 15000
search_parameters.log_search = True

routing.CloseModelWithParameters(search_parameters)

initial_solution = routing.ReadAssignmentFromRoutes(data['initial_routes'],
                                                    True)
print('Initial solution:')
print_solution(data, manager, routing, initial_solution)


# Solve the problem.
solution =  routing.SolveFromAssignmentWithParameters(initial_solution, search_parameters)

This initializes the parameters correctly, and the search finds the optimal solution.