Chained Pickup and Delivery Problem 2024-05-02

A question came up a last fall on the OR-Tools forum about how to enforce a chained pickup and delivery problem. The standard pickup and delivery problem is to pick something up from one place, and then deliver that thing to some other place. This question asked about an interesting and fairly common variation of the standard PDP: a chained sequence of pickups and deliveries.

The original question was titled: “VRPPD with Pickup-Delivery-Pickup-Delivery chain (passenger journey with interchanges)”. The full forum exchange can be found here.

The question describes a routing problem for passengers where pickups and deliveries are sequenced. My earliest consulting project with OR-Tools faced this exact problem, with pickups and deliveries of patients to and from medical clinics. The patients were picked up at their homes, dropped off at the clinic for their appointment, then picked up from the clinic and dropped off back at their homes after the appointment was complete.

In my case, there were known time windows to help with the sequencing. The patients had appointment times, and the appointments had known durations. The clinic delivery time had to be before the appointment, and the clinic pickup time had to be after the appointment end.

In the question posed on the forum, the time windows are not known in advance; all that is known is the order of the trip segments. The example given was a trip from. A to B, and then from B to C. The pickup at A has to be followed by the delivery at B, which is in turn followed by the pickup at B and the final delivery at C. Similar to the normal pickup and delivery problem, all nodes must be active, or all nodes must be inactive. That is, if the pickup at A does not happen, then obviously the delivery at B cannot happen, but it is also true that the entire second leg from B to C also cannot happen.

There is no simple API call for this case, so I thought it would be interesting to figure out how to do it, especially given that I had done something very similar years ago.

I created a dummy distance matrix for 17 nodes, leaving node 0 as the depot node. Then I chose to define the pickup and delivery sequences as a list of nodes. So for example, the list [1, 6, 2, 10] means pickup at 1, deliver at 6; pickup at 6, deliver at 2; pickup at 2, deliver at 10, which completes the sequence.

The next issue the came to mind was that OR-Tools cannot have duplicated visits to a single node. So at the very least the program will have to create dummy nodes for each of the interior nodes in a sequence. Continuing the prior example, the sequence [1, 6, 2, 10] will have to have a dummy node for the pickup at 6 and the pickup at 2, since the actual nodes 6 and 2 are used for the deliveries. More generally, since there is no limit on the numbers of sequences, each sequence needs dummy nodes for each of its nodes, with duplicated dummies for the interior nodes.

To handle this, I set up some simple mapping dictionaries between the dummy nodes and the original nodes, and then iterated through all of the sequences in the problem set to generate the dummies and map them to the reals.

In code, that looks like:

Example python for creating dummy nodes to avoid duplicate nodes

    node_sequences: typing.List[typing.List[int]] = [
        [1, 6, 2, 10],
        [1, 6, 1],
        [1, 10],
        [4, 3, 5, 9],
    ]
    data["original_node_sequences"] = node_sequences
    dummy_node_sequences: typing.List[typing.List[int]] = []
    dropoff_pickup_alternatives: typing.Dict[int, int] = {}
    dummy_to_real: typing.Dict[int, int] = {0: 0}
    real_to_dummy: typing.Dict[int, typing.List[int]] = {}

    num_nodes = len(data["distance_matrix"])

    # create the dummy node sequence, that pairs up from/to and creates no duplicates
    for sequence in node_sequences:
        dummy_sequence: typing.List[int] = []
        first = sequence[0]
        last = sequence[-1]

        for node in sequence:
            dummy_node = node
            if node in real_to_dummy:
                dummy_node = num_nodes
                num_nodes += 1
            else:
                real_to_dummy[node] = []
            real_to_dummy[node].append(dummy_node)
            assert dummy_node not in dummy_to_real
            dummy_to_real[dummy_node] = node
            dummy_sequence.append(dummy_node)
            if node not in [first, last]:
                # internal node to a sequence, will be a dropoff then a pickup
                dummy_pickup = num_nodes
                num_nodes += 1
                dropoff_pickup_alternatives[dummy_node] = dummy_pickup
                dummy_to_real[dummy_pickup] = node

        dummy_node_sequences.append(dummy_sequence)

The obvious consequence is that the more and longer sequences that are fed into the problem, the more difficult the problem will be to solve because the number of nodes will increase. A problem like [1, 2], [3, 4] will be easier to solve than a problem like [1, 2, 3, 4] because the former needs no dummy nodes, and the later will needs at least 2. The larger the problem, the more difficult it is to solve.

I’m somebody who likes to see a problem in order to get in there and fix it. So following along with that logic, I am first going to present what happens if you just load up the sequences as defined into a standard (go look at the examples in the repo) formulation of the pickup and delivery problem.

This test case uses a 17-node distance matrix I generated randomly, four vehicles, and the six required trip sequences shown in the figure:

Test case trip sequences

[1, 6, 2, 10],
[1, 6, 1],
[1, 10],
[4, 3, 5, 9],
[7, 8, 15, 11],
[13, 12, 16, 14]

After making sure that the solver can drop any node by putting in a very large disjunction value for each, the solver ran and quickly produced the following set of routes:

Test case results

Serving sequence [1, 6, 2, 10]
   1 pickup time 6 -> 6 deliver time 14 on vehicle 2
   6 pickup time 14 -> 2 deliver time 17 on vehicle 2
   2 pickup time 17 -> 10 deliver time 21 on vehicle 2
Serving sequence [1, 6, 1]
   1 pickup time 6 -> 6 deliver time 14 on vehicle 2
   6 pickup time 6 -> 1 deliver time 14 on vehicle 3
Serving sequence [1, 10]
   1 pickup time 6 -> 10 deliver time 21 on vehicle 2
Serving sequence [4, 3, 5, 9]
   4 pickup time 16 -> 3 deliver time 17 on vehicle 3
   3 pickup time 17 -> 5 deliver time 24 on vehicle 3
   5 pickup time 12 -> 9 deliver time 25 on vehicle 2
Serving sequence [7, 8, 15, 11]
   7 pickup time 2 -> 8 deliver time 6 on vehicle 0
   8 pickup time 6 -> 15 deliver time 18 on vehicle 0
   15 pickup time 18 -> 11 deliver time 21 on vehicle 0
Serving sequence [13, 12, 16, 14]
   13 pickup time 13 -> 12 deliver time 22 on vehicle 0
   12 pickup time 4 -> 16 deliver time 12 on vehicle 1
   16 pickup time 12 -> 14 deliver time 14 on vehicle 1
There are 0 dropped nodes.  They are: []

The first trip chain result is okay. There isn’t any time to do anything at the various destinations, but vehicle 2 serves each of them in order.

The second trip sequence shows the problem, however. The passenger at 1 joins the passenger from the first trip sequence (also at 1) on vehicle 2, and is delivered to node 6 at time 14. However, according to the trip plan, that same passenger needs to be picked up at node 6 at time 6 in vehicle 3! This is impossible, of course, and needs to be corrected by setting up the constraints properly.

The other sequences show similar problems. They tend to work okay if the sequence is served by the same vehicle, but when the vehicles change, the timing gets messed up.

So the question is clear: how to order the sequences of pickups and deliveries such that any vehicle can handle any leg of any sequence.

The insight here is that there needs to be some quantity that is uniform across all vehicles. Typically this is the time dimension. Unlike a distance traveled or an elapsed time dimension, absolute time is usually defined by allowing a very large or unlimited slack value, and carefully accounting for travel time and any service time at the nodes inside of the dimension callbacks. If necessary, the starting value of time can be set not some non-zero value, so that vehicles can start at different times (for example if different vehicles have different operating hours).

I’ve seen this come up multiple times on the forum, so it is worth belaboring a bit. A node count dimension, a distance dimension, and an elapsed time dimension are all specific to a single vehicle.

Consider a node count dimension. When a vehicle visits its first node, the node count increments to 1. When the vehicle visits the second node, it increments to 2. These values are completely independent of any other vehicle’s node count.

Now consider elapsed time. A vehicle starts at the depot at time zero, then travels to the first node, picks up something, etc. The dimension tracks the time difference between when the vehicle starts, and when it performs different actions. It does not represent time for any other vehicle, because it starts from zero when the vehicle starts moving.

Elapsed time is very close to absolute time. The only difference is that with absolute time, the starting time for the vehicle is not fixed at zero, but rather allowed to vary depending on when the vehicle starts moving. So if a vehicle starts at 8am and time is tracked in minutes, then the start of an absolute time dimension would be 480.

If all vehicles start at exactly the same time, and if none of the pickups and deliveries have any fixed, absolute time windows (such as pickup after 10AM or deliver by 1PM), then there is no practical difference between elapsed time per vehicle and absolute time. It is as if all of the vehicles synchronize their watches at exactly 8AM.

That is how I am tracking absolute time here. Because I don’t really have real pickups and dropoffs, I don’t care about when those happen. So I am letting time start at 0, at which time all of the vehicles leave the depot.

This gives me the time dimension:

A dimension to track absolute time in minutes, with all vehicles starting at time zero.

# Add Time Window constraint.
# [START time_constraint]
time_evaluator_index = routing.RegisterTransitCallback(
    partial(create_time_evaluator(data), manager)
)
dimension_name = "Time"  # in minutes
routing.AddDimension(
    time_evaluator_index,
    480,  # one day is 8 hours, with infinite slack
    480,  # vehicle maximum time is one day
    True,  # start cumul to zero
    dimension_name,
)
time_dimension = routing.GetDimensionOrDie(dimension_name)

With this time dimension, I have something that makes sense to compare across different vehicles. If vehicle A visits node 1 at time 10, that is the same quantity for vehicle B, because they synchronized their clocks when they left the depot. With this, it is now possible to compare quantities between different vehicles. That is, it is possible to require that dropoff events happen earlier in “time” than the next pickup event for the person. This code is shown below, and it carefully links up the dropoff events from one sequence with the pickup and dropoff events of the next sequence.

How to maintain ordering of nodes served by different vehicles

time_dimension = routing.GetDimensionOrDie("Time")
# Add sequential pdp constraints each location in each pdp chain
solver = routing.solver()

for sequence in data["pickups_deliveries"]:
    prior_index = None
    prior_node = None
    for next_node in sequence:
        next_index = manager.NodeToIndex(next_node)
        if prior_index is None:
            prior_index = next_index
            prior_node = next_node
            continue
        # prior index has to happen before next index in Time, which
        # is common across all vehicles
        solver.Add(
            time_dimension.CumulVar(prior_index) < time_dimension.CumulVar(next_index)
        )

        # pickup after delivery for internal nodes
        if prior_node in data["dropoff_pickup_alternatives"]:
            pickup_node = data["dropoff_pickup_alternatives"][prior_node]
            pickup_index = manager.NodeToIndex(pickup_node)
            solver.Add(
                time_dimension.CumulVar(prior_index) < time_dimension.CumulVar(pickup_index)
            )
        prior_index = next_index
        prior_node = next_node

As I said at the beginning, this question came up in the forums in the fall of 2023, and I am only writing it up now in the spring of 2024. So as I looked it over again, I didn’t quite understand what was going on. As is often the case, future James (me today) is not as clever as yesterday James, but yesterday James, although smart, was lazy and did not document his work! In my working code I am liberal with my comments for that exact reason. Anyway, it is worth stepping through what is going on with an actual sequence as an example.

Suppose the first sequence is [1, 6, 2, 10]. The first link from 1 to 6 is just a plain pickup and dropoff pair. Let’s look at the second set, from 6 to 2.

The start of the outer loop, with an example sequence

for sequence in data["pickups_deliveries"]:
    # sequence is [1, 6, 2, 10]
    prior_index = None
    prior_node = None
    for next_node in sequence:
        # taking the second pair, from 6 to 2,
        # prior_node is 6,
        # next_node is 2
        next_index = manager.NodeToIndex(next_node)

Note that the trip is not from node 6 to node 2. These are the original locations, and so node 6 is only used for the destination of the trip from node 1. For the trip from 6 to to 2, the loop shown in the preceding section has created a dummy node for that purpose. In this case, the node is actually node 17.

Restating that with the nodes used in practice, the defined sequence [1, 6, 2] is expanded into two trips. The first is from 1 to 6, and the next is from node 17 to node 2. Even though node 6 and node 17 refer to the identical location in the real world, they must be kept separate and unique in the model formulation.

The next line in the code sets the relationship between the dimension at node 6 and node 2. Again, these are not pickup and dropoff pairs. Instead, both of these nodes are dropoff nodes in practice. There is a dummy pickup node in between the two nodes, node 17, that is linked as a pickup and deliver pair with node 2.

The constraint linking the time between two drop off nodes of different segments

        # prior index has to happen before next index in Time, which
        # is common across all vehicles
        solver.Add(
            time_dimension.CumulVar(prior_index) < time_dimension.CumulVar(next_index)
        )

This bit of code brings up another point that is difficult to handle in a one-size-fits-all way. With OR-Tools, it is not the case that dimension values are set to zero when a node is not visited. Taking the above case as an example, because the two dropoff nodes 6 and 2 are not linked by any pickup and deliver API call, it is not the case that they are both active or both inactive. So you cannot just assume that if the prior node is not visited, then its dimension values are zero. I actually think the internal rule in the solver is to let these variables take on any value that is convenient as long as it is greater than or equal to zero (dimension values cannot be negative).

Many times this is just fine. The solver will adjust the value as needed if it is not constrained. Sometimes it is needed, however, for example if a constraint should only be enforced if both nodes are active. In that case, I would add a boolean flag to detect if the node is active or not. Something like:

How to link whether the nodes are active to a constraint. The result is either 0=0 if one or both nodes are inactive, or an actual constraint if both are active.

prior_active = routing.ActiveVar(prior_index)
next_active = routing.ActiveVar(next_index)
both_active = prior_active * next_active
solver.Add(
    both_active * time_dimension.CumulVar(prior_index)
    < both_active * time_dimension.CumulVar(next_index)
)

In this case, I did not bother with adding this. Every time you add these sorts of constraints, it has the potential to slow the solver down. It also can create erratic solutions, because the boolean values can create non-monotonic states that the solver does not anticipate or handle well.

The next part of the loop deals with the dummy pickup nodes. Again, in the running example, the prior node is 6, the next node is 2, and the dummy pickup node is 17.

Adding a constraint to link the time of a dropoff node and the very next pickup node

        # pickup after delivery for internal nodes
        if prior_node in data["dropoff_pickup_alternatives"]:
            # this is true here, as 6 has a paired pickup node of 17

            pickup_node = data["dropoff_pickup_alternatives"][prior_node]
            # the pickup node in this example is 17

            pickup_index = manager.NodeToIndex(pickup_node)

            # link the dimension value at the dropoff node (6) with
            # the pickup node (17)
            solver.Add(
                time_dimension.CumulVar(prior_index) < time_dimension.CumulVar(pickup_index)
            )
            # this is not needed, as by definition pickup and next
            # are visited bythe same vehicle.
            #
            # solver.Add(
            #     time_dimension.CumulVar(pickup_index)
            #     < time_dimension.CumulVar(next_index)
            # )
        prior_index = next_index
        prior_node = next_node

It is important to notice the linking of time here and in the prior constraint. The original question stated that any vehicle can serve any of the links. So it is not necessarily the case that the vehicle that drops off at 6 transitions directly to the pickup at 17. From the perspective of the traveler, this is the case, but dimensions are defined from the perspective of the vehicle, not the items or people the vehicle carries.

Without that link, the solver would be free to schedule the pickup before the dropoff, as was shown in the buggy output. The prior constraint establishes that the dropoff at node 6 must occur before the dropoff at node 2. This constraint is needed to make sure that the pickup at the dummy node 17 (again, really node 6) happens after the dropoff at node 6. Both of these constraints must be set to make sure the right thing happens every time. However, because they are visited by the same vehicle, there is no need to explicitly link node 17 and node 2.

These two constraints, therefore, are exactly the requirement requested by the original question on the forum.

With these constraints in place, the solver does the right thing, as is shown in the final output below.

Test results with correct time constraints

Serving sequence [1, 6, 2, 10]
   1 pickup time 31 -> 6 deliver time 39 on vehicle 0
   6 pickup time 100 -> 2 deliver time 103 on vehicle 3
   2 pickup time 134 -> 10 deliver time 138 on vehicle 3
Serving sequence [1, 6, 1]
   1 pickup time 31 -> 6 deliver time 39 on vehicle 0
   6 pickup time 100 -> 1 deliver time 150 on vehicle 3
Serving sequence [1, 10]
   1 pickup time 31 -> 10 deliver time 43 on vehicle 0
Serving sequence [4, 3, 5, 9]
   4 pickup time 6 -> 3 deliver time 7 on vehicle 0
   3 pickup time 28 -> 5 deliver time 37 on vehicle 0
   5 pickup time 88 -> 9 deliver time 142 on vehicle 3
Serving sequence [7, 8, 15, 11]
   7 pickup time 2 -> 8 deliver time 47 on vehicle 0
   8 pickup time 88 -> 15 deliver time 100 on vehicle 2
   15 pickup time 221 -> 11 deliver time 224 on vehicle 2
Serving sequence [13, 12, 16, 14]
   13 pickup time 95 -> 12 deliver time 225 on vehicle 2
   12 pickup time 286 -> 16 deliver time 294 on vehicle 1
   16 pickup time 385 -> 14 deliver time 387 on vehicle 1
There are 0 dropped nodes.  They are: []

In every case, time increases continuously for all pickup and deliver events in each sequence. The deliver time for a node is less than the pickup time at that same node. Further, all of the sequences (except the single leg sequence) are served by multiple vehicles. So the code accomplishes the original request.

I next expanded the list of sequences and increased the numbers of vehicles, to try to stress the solver and see what happened when it could not serve all of the nodes. I discovered that there are some very interesting issues that occur when I start adding more sequences and more vehicles. The solver has trouble serving all of the nodes with just the usual fixed, very large number disjunction value. The next blog post on this topic will explore this aspect of the problem.