VRPTW with collects and deliveries

vrptw_collect_deliveries.py

#!/usr/bin/env python3
from ortools.constraint_solver import pywrapcp
from ortools.constraint_solver import routing_enums_pb2


def create_data_model(input_data=None):
    """Stores the data for the problem."""
    data = {}
    data['distance_matrix'] = [
        [
            0, 48349, 61297, 43122, 58087, 48400, 29370, 49387, 58984, 59945,
            57572, 43950, 28426, 28043, 28043, 41348, 43597
        ],
        [
            44926, 0, 11894, 8064, 27222, 10241, 28622, 6477, 29870, 41057,
            28459, 5578, 28392, 28009, 28009, 7609, 7214
        ],
        [
            58222, 12254, 0, 17514, 19163, 15709, 41917, 13627, 32970, 64286,
            19843, 16982, 41688, 41304, 41304, 18817, 16939
        ],
        [
            46900, 8066, 17196, 0, 15691, 6990, 30595, 5606, 21578, 32765,
            20166, 10333, 30365, 29982, 29982, 3476, 1387
        ],
        [
            56202, 18222, 19454, 15690, 0, 11355, 45304, 12531, 10896, 31239,
            7828, 21633, 45074, 44691, 44691, 18185, 15820
        ],
        [
            46296, 10250, 15501, 7127, 11404, 0, 36741, 4559, 16901, 28088,
            15490, 16479, 36511, 36128, 36128, 9621, 7847
        ],
        [
            35695, 30036, 72709, 31170, 73589, 37159, 0, 34590, 50431, 61618,
            49019, 25637, 4704, 5158, 5158, 28169, 30419
        ],
        [
            50241, 6733, 13590, 5639, 12542, 4552, 33936, 0, 23877, 35064,
            22466, 10145, 33706, 33323, 33323, 6943, 5065
        ],
        [
            56883, 29479, 33255, 21220, 10891, 16208, 49529, 21967, 0, 20883,
            5353, 29267, 49299, 48916, 48916, 22409, 22304
        ],
        [
            57692, 35694, 64740, 27435, 31235, 27714, 55743, 28181, 20885, 0,
            16859, 35481, 55513, 55130, 55130, 28624, 28519
        ],
        [
            54598, 27194, 20134, 18935, 7828, 13923, 47244, 19682, 5665, 16852,
            0, 26982, 47014, 46631, 46631, 20124, 20020
        ],
        [
            40881, 6290, 19144, 10267, 24957, 16256, 24576, 10606, 29527,
            40714, 28116, 0, 24346, 23963, 23963, 7266, 9516
        ],
        [
            34751, 29806, 42753, 30940, 72290, 36929, 4704, 34360, 50201,
            61388, 48789, 25407, 0, 1080, 1080, 27939, 30189
        ],
        [
            34368, 29423, 42370, 30557, 71907, 36546, 5158, 33977, 49817,
            61004, 48406, 25023, 1080, 0, 0, 27556, 29806
        ],
        [
            34368, 29423, 42370, 30557, 71907, 36546, 5158, 33977, 49817,
            61004, 48406, 25023, 1080, 0, 0, 27556, 29806
        ],
        [
            44443, 8089, 18213, 3302, 18078, 10088, 28139, 6624, 23396, 34583,
            21985, 7877, 27909, 27526, 27526, 0, 2452
        ],
        [
            46140, 7217, 16495, 1368, 16359, 8074, 29835, 4905, 22661, 33848,
            21250, 9573, 29606, 29222, 29222, 2842, 0
        ]
    ]
    data['time_matrix'] = [
        [
            0, 4081, 5512, 4029, 5343, 4521, 3639, 4487, 4877, 5696, 4927,
            3672, 3424, 3338, 3338, 3722, 4007
        ],
        [
            3884, 0, 1476, 1050, 1999, 1465, 3190, 887, 2555, 3725, 2604, 748,
            3057, 2970, 2970, 900, 960
        ],
        [
            5297, 1507, 0, 1841, 1116, 1854, 4603, 1434, 2353, 3056, 2292,
            2145, 4470, 4383, 4383, 1964, 1792
        ],
        [
            3823, 1051, 1896, 0, 1787, 1068, 3129, 786, 1847, 3017, 1897, 998,
            2996, 2910, 2910, 421, 231
        ],
        [
            5268, 2155, 1181, 1788, 0, 1501, 4711, 1465, 1410, 2867, 1243,
            2616, 4577, 4491, 4491, 2002, 1909
        ],
        [
            4414, 1467, 1896, 1068, 1487, 0, 3991, 777, 1899, 3069, 1949, 1859,
            3858, 3771, 3771, 1283, 1221
        ],
        [
            3587, 3276, 4081, 3171, 3941, 4033, 0, 3681, 4512, 5682, 4561,
            2868, 1058, 986, 986, 2857, 3142
        ],
        [
            4298, 924, 1522, 755, 1464, 773, 3605, 0, 2278, 3449, 2328, 1385,
            3471, 3385, 3385, 878, 706
        ],
        [
            4723, 2463, 2418, 1762, 1410, 1785, 4387, 2122, 0, 1578, 820, 2256,
            4254, 4168, 4168, 1679, 1914
        ],
        [
            5614, 3762, 3159, 3060, 2867, 3041, 5686, 3369, 1579, 0, 2076,
            3554, 5553, 5466, 5466, 2977, 3212
        ],
        [
            4871, 2611, 2358, 1910, 1243, 1933, 4535, 2270, 907, 2078, 0, 2404,
            4402, 4316, 4316, 1827, 2062
        ],
        [
            3290, 734, 2154, 873, 2454, 1735, 2596, 1335, 2214, 3384, 2264, 0,
            2463, 2377, 2377, 559, 844
        ],
        [
            3371, 3143, 4574, 3037, 4114, 3899, 1058, 3548, 4378, 5549, 4428,
            2734, 0, 243, 243, 2723, 3008
        ],
        [
            3285, 3056, 4487, 2951, 4028, 3813, 986, 3461, 4292, 5462, 4342,
            2648, 243, 0, 0, 2637, 2922
        ],
        [
            3285, 3056, 4487, 2951, 4028, 3813, 986, 3461, 4292, 5462, 4342,
            2648, 243, 0, 0, 2637, 2922
        ],
        [
            3570, 952, 1980, 420, 2072, 1382, 2876, 870, 1867, 3037, 1917, 744,
            2743, 2657, 2657, 0, 331
        ],
        [
            3795, 960, 1820, 228, 1912, 1220, 3101, 710, 1999, 3170, 2049, 970,
            2968, 2882, 2882, 375, 0
        ]
    ]
    data['num_vehicles'] = 6
    data['time_windows'] = [(0, 85800), (0, 86340), (0, 86340), (0, 86340),
                            (0, 86340), (0, 86340), (0, 86340), (28800, 82800),
                            (0, 86340), (0, 86340), (0, 86340), (0, 86340),
                            (0, 86340), (0, 86340), (0, 86340), (0, 86340),
                            (0, 86340)]
    data['depot'] = 0
    data['vehicle_capacities'] = [5, 5, 5, 5, 5, 5]
    data['deliveries'] = [
        0, # depot
        0,
        0,
        -1, # delivery
        -1, # delivery
        -1, # delivery
        -1, # delivery
        -1, # delivery
        0,
        -1, # delivery
        -1, # delivery
        0,
        -1, # delivery
        -1, # delivery
        -1, # delivery
        -1, # delivery
        -1 # delivery
    ]
    data['total_load'] = [
        0, # depot
        6, # impossible pickup
        1, # pickup
        -1, # delivery
        -1, # delivery
        -1, # delivery
        -1, # delivery
        -1, # delivery
        1, # pickup
        -1, # delivery
        -1, # delivery
        1, # pickup
        -1, # delivery
        -1, # delivery
        -1, # delivery
        -1, # delivery
        -1 # delivery
    ]
    data['driver_shift_times'] = [(0, 79200), (0, 79200), (0, 79200),
                                  (0, 79200), (0, 79200), (0, 79200)]

    assert len(data['distance_matrix']) == len(data['time_matrix'])
    assert len(data['time_matrix']) == len(data['time_windows'])
    assert len(data['distance_matrix']) == len(data['deliveries'])
    assert len(data['deliveries']) == len(data['total_load'])
    assert data['num_vehicles'] == len(data['driver_shift_times'])
    return data


def print_solution(data, manager, routing, solution):
    """Prints solution on console."""
    #objective
    print(f'Objective: {solution.ObjectiveValue()}')
    # Display dropped nodes.
    dropped_nodes = 'Dropped nodes:'
    for node in range(routing.Size()):
        if routing.IsStart(node) or routing.IsEnd(node):
            continue
        if solution.Value(routing.NextVar(node)) == node:
            dropped_nodes += f'{manager.IndexToNode(node)} '
    print(dropped_nodes)
    # Display routes
    total_distance = 0
    time_dimension = routing.GetDimensionOrDie('Time')
    capacity_dimension = routing.GetDimensionOrDie('Capacity')
    delivery_dimension = routing.GetDimensionOrDie('Delivery')
    for vehicle_id in range(data['num_vehicles']):
        index = routing.Start(vehicle_id)
        plan_output = f'Route for vehicle {vehicle_id}:\n'
        route_distance = 0
        while not routing.IsEnd(index):
            node_index = manager.IndexToNode(index)
            load_var = capacity_dimension.CumulVar(index)
            route_load = solution.Value(load_var)
            delivery_var = delivery_dimension.CumulVar(index)
            route_delivery = solution.Value(delivery_var)
            plan_output += f'{node_index} Load({route_load})Delivery({route_delivery}) -> '
            previous_index = index
            index = solution.Value(routing.NextVar(index))
            route_distance += routing.GetArcCostForVehicle(
                previous_index, index, vehicle_id)
        node_index = manager.IndexToNode(index)
        load_var = capacity_dimension.CumulVar(index)
        route_load = solution.Value(load_var)
        delivery_var = delivery_dimension.CumulVar(index)
        route_delivery = solution.Value(delivery_var)
        plan_output += f' {node_index} Load({route_load})Delivery({route_delivery})\n'
        plan_output += f'Distance of the route: {route_distance}m\n'
        print(plan_output)
        total_distance += route_distance
    print(f'Total distance of all routes: {total_distance}m')


def add_capacity_constraints(routing, manager, data):
    """Adds capacity constraint"""
    def delivery_callback(from_index):
        from_node = manager.IndexToNode(from_index)
        return data['deliveries'][from_node]

    delivery_callback_index = routing.RegisterUnaryTransitCallback(
        delivery_callback)
    delivery = 'Delivery'
    routing.AddDimensionWithVehicleCapacity(
        delivery_callback_index,
        0,  # null capacity slack
        data['vehicle_capacities'],  # vehicle maximum capacities
        False,  # don't start cumul to zero
        delivery)

    def demand_callback(from_index):
        from_node = manager.IndexToNode(from_index)
        return data['total_load'][from_node]

    demand_evaluator_index = routing.RegisterUnaryTransitCallback(
        demand_callback)
    capacity = 'Capacity'
    routing.AddDimensionWithVehicleCapacity(
        demand_evaluator_index,
        0,  # null capacity slack
        data['vehicle_capacities'],  # vehicle maximum capacities
        False,  # start cumul to zero
        capacity)

    capacity_dimension = routing.GetDimensionOrDie(capacity)
    delivery_dimension = routing.GetDimensionOrDie(delivery)
    for i in range(data['num_vehicles']):
        index = routing.Start(i)
        #capacity_dimension.SetCumulVarSoftLowerBound(
        #    index=i,
        #    lower_bound=data['vehicle_capacities'][i],
        #    coefficient=100_000_000)
        routing.solver().Add(
            delivery_dimension.CumulVar(index) == capacity_dimension.CumulVar(index))

        end_index = routing.End(i)
        routing.AddVariableMinimizedByFinalizer(
            delivery_dimension.CumulVar(end_index))
        routing.AddVariableMinimizedByFinalizer(
            capacity_dimension.CumulVar(end_index))

    # capacity_dimension.SetGlobalSpanCostCoefficient(max(data['vehicle_capacities']))


def add_time_window_constraints(routing, manager, data):
    # _, service_time = get_service_time(input_data)
    def time_callback(from_index, to_index):

        from_node = manager.IndexToNode(from_index)
        to_node = manager.IndexToNode(to_index)
        return data['time_matrix'][from_node][
            to_node]  #+ int(service_time[from_node]*60)

    transit_callback_index_time = routing.RegisterTransitCallback(
        time_callback)

    time = 'Time'
    routing.AddDimension(
            transit_callback_index_time,
            86400,
            86400,
            False,
            time)

    time_dimension = routing.GetDimensionOrDie(time)

    for location_idx, time_window in enumerate(data['time_windows']):
        if location_idx == 0:
            continue
        index = manager.NodeToIndex(location_idx)
        time_dimension.CumulVar(index).SetRange(time_window[0], time_window[1])
        routing.AddToAssignment(time_dimension.SlackVar(index))
        # Add time window constraints for each vehicle start node.
    for vehicle_id in range(data['num_vehicles']):
        index = routing.Start(vehicle_id)
        time_dimension.CumulVar(index).SetRange(
            data['driver_shift_times'][vehicle_id][0],
            data['driver_shift_times'][vehicle_id][1])

        routing.AddToAssignment(time_dimension.SlackVar(index))
        # routing.AddToAssignment(time_dimension.SlackVar(routing.End(vehicle_id)))

    # Instantiate route start and end times to produce feasible times.
    for i in range(data['num_vehicles']):
        routing.AddVariableMinimizedByFinalizer(
            time_dimension.CumulVar(routing.Start(i)))
        routing.AddVariableMinimizedByFinalizer(
            time_dimension.CumulVar(routing.End(i)))

    routing.SetArcCostEvaluatorOfAllVehicles(transit_callback_index_time)
    # time_dimension.SetGlobalSpanCostCoefficient(1)


def add_distance_callback(routing, manager, data):
    def distance_callback(from_index, to_index):
        from_node = manager.IndexToNode(from_index)
        to_node = manager.IndexToNode(to_index)
        return data['distance_matrix'][from_node][to_node]

    transit_callback_index_dist = routing.RegisterTransitCallback(
        distance_callback)
    routing.SetArcCostEvaluatorOfAllVehicles(transit_callback_index_dist)

    dimension_name = 'Distance'
    routing.AddDimension(transit_callback_index_dist, 0, 250000, True,
                         dimension_name)
    distance_dimension = routing.GetDimensionOrDie(dimension_name)
    # distance_dimension.SetGlobalSpanCostCoefficient(40)


def Run_Routing():
    data = create_data_model()
    manager = pywrapcp.RoutingIndexManager(
        len(data['distance_matrix']),
        data['num_vehicles'],
        data['depot'],
    )
    routing = pywrapcp.RoutingModel(manager)

    add_distance_callback(routing, manager, data)
    add_capacity_constraints(routing, manager, data)
    add_time_window_constraints(routing, manager, data)

    penalty = 100_000_000
    for node in range(1, len(data['distance_matrix'])):
        routing.AddDisjunction([manager.NodeToIndex(node)], penalty)

    search_parameters = pywrapcp.DefaultRoutingSearchParameters()

    search_parameters.first_solution_strategy = (
        routing_enums_pb2.FirstSolutionStrategy.PATH_CHEAPEST_ARC)

    solution = routing.SolveWithParameters(search_parameters)

    if solution:
        return print_solution(data, manager, routing, solution)
    else:
        print(solution)


Run_Routing()

zz_readme.md

Potential output:

%./abhik.py
Objective: 100034154
Dropped nodes:1 
Route for vehicle 0:
0 Load(0)Delivery(0) ->  0 Load(0)Delivery(0)
Distance of the route: 0m

Route for vehicle 1:
0 Load(0)Delivery(0) ->  0 Load(0)Delivery(0)
Distance of the route: 0m

Route for vehicle 2:
0 Load(0)Delivery(0) ->  0 Load(0)Delivery(0)
Distance of the route: 0m

Route for vehicle 3:
0 Load(5)Delivery(5) -> 9 Load(5)Delivery(5) -> 8 Load(4)Delivery(4) -> 10 Load(5)Delivery(4) -> 4 Load(4)Delivery(3) -> 2 Load(3)Delivery(2) -> 5 Load(4)Delivery(2) -> 7 Load(3)Delivery(1) -> 11 Load(2)Delivery(0) ->  0 Load(3)Delivery(0)
Distance of the route: 17825m

Route for vehicle 4:
0 Load(3)Delivery(3) -> 15 Load(3)Delivery(3) -> 16 Load(2)Delivery(2) -> 3 Load(1)Delivery(1) ->  0 Load(0)Delivery(0)
Distance of the route: 8104m

Route for vehicle 5:
0 Load(4)Delivery(4) -> 6 Load(4)Delivery(4) -> 12 Load(3)Delivery(3) -> 13 Load(2)Delivery(2) -> 14 Load(1)Delivery(1) ->  0 Load(0)Delivery(0)
Distance of the route: 8225m

Total distance of all routes: 34154m

If instead of units/slots we use those those 2 dimensions of "deliveries" and "total load" both with volume and weight, is it possible for the solver to take both under consideration in order to return a feasible solution?