Recently a question was asked on the OR Tools forum about how to convert a single-day TSP problem into a multi-day TSP problem. The original poster was thorough, and cross posted to Stack Overflow and posted his original code as a gist on GitHub. The problem is interesting in that it takes a away from an idealized lab problem toward integrating real-world constraints. So I responded to the question with my solution, and decided to write it up properly here.
Analyzing the problem and its requirements
The full question is in the links above, but it can be summarized as follows. The code as written finds a TSP path through the nodes (which apparently “represent McDonald’s restaurants in the Minneapolis area”!), with the constraints that the vehicle start at the depot at 6am, and then end up back at the depot node at 6pm. This worked, but the route did not visit all of the nodes. The poster’s question was how to change this to visit all of the nodes over multiple days.
My first thought when reading the problem statement, prior to looking at the code at all, was that the solution would involve creating dummy nodes for the additional days. In the original setup, node 0 is the starting depot, and node 1 is the ending depot. Therefore I would “duplicate” node 1 for the purposes of computing distances to and from all of the other nodes.
Each of these dummy overnight nodes can be used to “reset” the time dimension back to zero. Because there is a hard start and end time and just one vehicle, there is no need to worry about an infinite clock. Time can be treated as “loading up” a single day, and then emptying it out back at the depot.
A similar technique is used in the reload example in the OR Tools source code. The “demand” at the depot nodes are the negative of the vehicle capacity values, as shown here, which remove all possible loads collected at the nodes.
The next question that needed to be answered was how to handle time. Typically in a routing problem, time is modeled using a dimension that “collects” time as the vehicle travels between nodes and handles service events at nodes. Sometimes it is important to make sure that time progresses linearly just like in the real world. In that case time starts with some starting date, and then increases to some ending date, collecting seconds, minutes, or hours along the way. In other cases, time can be handled as an offset. For example, a daily schedule will start time at 00:00 and runs to 23:59, and so time can be treated as integer minutes, with a minimum of 0 and a maximum of 1440. Choosing which way to go is important because it dictates how to handle multiple days. The time as a continuum approach means the dimension just keeps adding, so the second day starts at 1441 seconds. The time as an offset from zero approach means that the value of time in the dimension must be reset to zero at the end of the day.
In most cases, handling time as a continuum is just easier, and it is the approach that I prefer for real-word problems, as it simplifies handling 24hr events, overlapping schedules, events with a period greater than a day, and aperiodic events. However, I had never actually implemented and tested the case of time resetting to zero. Because the problem statement in this case said the day starts at 6am and ends at 6pm, it seemed a good problem to try the reset route.
Reset time using “slack”
In a standard pickup and delivery problem, the load on the vehicle is reduced by visiting locations with a negative demand. Because time in OR Tools is no different than any other dimension, the way to “reset” time to zero is to deliver the collected time to a node with negative demand. The trick is to make sure to reset the dimension all the way to zero. A vehicle traveling for just an hour needs to “deliver” less than a vehicle that has been traveling all day. So if the reset-time node demands a full day of time, it should be able to reset all vehicles to zero. However, unless the vehicle has “enough” quantity stored in a dimension, the solver will not let the vehicle visit a node. The way around this is to use a non-zero slack in the dimension constructor.
In working with others who are new to using OR Tools, the concept of slack is difficult to grasp. I myself cargo-culted the use of slack in my early projects, until I found the definition in the documentation.
The key to understanding slack is understanding how the solver updates dimension quantities as a vehicles moves around the network. According to the documentation, the relationship is as follows.
Given:
i = current node
j = next node
cumul(x) = quantity at node x
Then:
cumul(j) = cumul(i) + transit(i) + slack(i)
There are a few items here that I missed the first time I tried to
understand this. First of all, the quantity of the dimension at node
j, or cumul(j), is only dependent on the values of the preceding
node i. That is, cumul(j) does not reference the value of
transit(j) or slack(j). This is confusing at first, but I like to
think about it as follows:
cumul(j)is the state of the vehicle as it arrives atj.
The second item that I missed when I was starting out is that
transit(i) (which should really be called transit(i,j)) is
associated with activities, loading, services, or other changes to the
dimension performed at i, and then any changes to the dimension
occurring because of the transition from i to j. This is
important when writing code for the function to update a dimension’s
value, especially in a two dimensional callback. In that case, the
callback is passed both a “from” index and a “to” index, and it isn’t
immediately obvious which to use for the node-specific values like
service time, loading an item etc. The “from” index should be used,
as this is node i in the definition. For example:
def time_callback(manager, from_index, to_index):
# Returns the travel plus service time between the two nodes.
# Convert from routing variable Index to time matrix NodeIndex.
from_node = manager.IndexToNode(from_index)
to_node = manager.IndexToNode(to_index)
# get the service time for ... um ...
# which node should it be? Well, the to_node is j while the
# from_node is i, and the documentation says that the
# cumul of j is equal to the transit of i, so it's gotta be from_node Alex
_service_time = node_service_times[from_node]
return travel_time[(from_node,to_node)] + _service_time
In this function the from_index is used to determine the service
time, and both indices are used to determine the travel time between
nodes. Again, these changes produce the state of the vehicle as it
arrives at the destination node.
Which finally brings us to the consideration of the “slack” at node
i. Slack is defined in the dimension constructor, but it can be
adjusted for each node if desired.
An example constructor call looks like this:
# Let the solver know about the callback function
time_callback_index = routing.RegisterTransitCallback(time_callback_fn)
# create time dimension
routing.AddDimension(
time_callback_index,
Slack_Max, # An upper bound for slack (wait times, as well as reset time)
MaxTime, # An upper bound for the total time over each vehicle's route.
False, # Do not fix the starting time at zero.
'Time' # The dimension name
)
So what should be the values for MaxTime and SlackMax? Well, for
this problem, the vehicles need to be back at the depot at 6pm. Using
minutes, that is 1080 minutes (18 hours, 60 minutes to an hour). So
the MaxTime can be set to 1080.
The value for SlackMax is a slightly more interesting
question. Ideally, a vehicle will be working the entire day, and must
leave the depot the next morning at 6am. In that case, if the vehicle
returns to the depot at 1080 minutes, then that value can be “reset”
to 360 minutes (6am) by giving the depot node a “demand” of -720. In
this case, the math works perfectly and no slack is needed to make it
work.
At the other end of the spectrum, if a vehicle does nothing at all,
then it will return to the “next” depot at 6am, meaning that it
arrives with a CumulVar value of 360 minutes. In that case,
subtracting 720 will push the value to zero for the dimension. (But
not negative. I learned the hard way that dimensions have a hard
floor of zero. They are non-negative integers.) This value of zero
needs to be bumped back up to 360, so the maximum possible slack to
balance this equation is SlackMax=360 minutes.
Note that it is required to use a non-zero slack in the dimension constructor if you wish to use slack anywhere in the model. It doesn’t work to initialize with a zero slack and then change one node to have a positive slack. Instead, set the dimension slack to positive, and then reset all other nodes to have zero slack. This sounds less efficient, but it is required to avoid creating a constant, unchanging value of zero for the slack.
Using dummy depot nodes
In order to reset the time dimension using the slack and reset values figured out above, I need to create some dummy “overnight” nodes, one for each day past the default single day case. In an email exchange, it became clear that the original poster wanted to be able to adjust the overnight nodes based on some external logic. For my example, I just chose them to be co-located with node 1 (the ending depot of the original problem). OR Tools also requires the user to track nodes using a numeric index. In this case, I added the dummy nodes on after the original problem nodes, as follows.
import time_details as T
...
num_nodes = T.num_nodes()
# create dummy nodes for returning to the depot every night
night_nodes = range(num_nodes, num_nodes+num_days)
total_nodes = num_nodes + len(night_nodes)
# Create the routing index manager.
manager = pywrapcp.RoutingIndexManager(total_nodes, 1, [0], [1])
These night nodes have non-zero slack, and every other node has zero slack. This is accomplished as follows:
# create time dimension
routing.AddDimension(
time_callback_index,
Slack_Max, # An upper bound for slack (the wait times at the locations).
Capacity, # An upper bound for the total time over each vehicle's route.
False, # Determine whether the cumulative variable is set to zero at the start of the vehicle's route.
'Time')
time_dimension = routing.GetDimensionOrDie('Time')
# get rid of slack for all regular nodes, but keep for depot, night nodes
for node in range(2, num_nodes):
index = manager.NodeToIndex(node)
time_dimension.SlackVar(index).SetValue(0)
I have another file time_details.py that I am using to modularize
everything associated with time, including the distance matrix. The
size of the distance matrix give me the number of nodes. This value
starts the indexing of the night nodes. The range() function
creates one night node for each day, less one. Then the solver
manager is created to be used when translating the solver’s index
values into my node values.
The next issue to handle is to create the illusion that each of these night nodes is located at the same place as node 1. This is done inside the link-to-link transit callback. The original code looked like this:
# Create and register a transit callback.
def transit_callback(manager,
from_index, to_index):
# Returns the travel time between the two nodes.
# Convert from routing variable Index to time matrix NodeIndex.
from_node = manager.IndexToNode(from_index)
to_node = manager.IndexToNode(to_index)
return Matrix[from_node][to_node]
}
This function takes three parameters: the manager, the from index, and
the to index. The idea is to use functools.partial to pass in the
manager parameter, and then pass the curried function to the solver,
which in turn will supply the indices. Now that I’ve created night
nodes, this function will fail as there are no entries for the night
nodes in the distance matrix, but the solver will pass in index
values that correspond to those nodes. Therefore I’ve made the following
changes:
# Create and register a transit callback.
def transit_callback(manager,
night_nodes,
from_index, to_index):
# Returns the travel time between the two nodes.
# Convert from routing variable Index to time matrix NodeIndex.
from_node = manager.IndexToNode(from_index)
to_node = manager.IndexToNode(to_index)
# adjust if either node is a night node
if from_node in night_nodes:
from_node = 1
if to_node in night_nodes:
to_node = 1
return Matrix[from_node][to_node]
}
This expanded version of the transit_callback takes a second
parameter for the functools.partial call: the list of night nodes.
Obviously one can get more creative associating each individual night
node with a either a specific original node value or a completely new
node, but using node 1 works just as well to demonstrate the idea.
Finally, I want to make sure the solver can skip the overnight nodes if it wants to. There is no need to force the solver to visit an overnight node if it can fit all the trips into fewer days. To do this I use disjunctions as follows.
# Allow all overnight nodes to be dropped for free
for node in night_nodes:
routing.AddDisjunction([manager.NodeToIndex(node)], 0)
Initial runs and early mornings
With the dummy overnight nodes, the solver almost gets everything right. For a run with one day, the solver is able to fit in some of the nodes, but certainly not all of them.
solution found. Objective value is 210007252
Dropped
[2, 3, 4, 5, 7, 8, 9, 19, 20, 21, 24, 25, 28, 30, 31, 35, 36, 37, 39, 40, 41]
Scheduled
[node, order, min time, max time]
[0, 0, '06:00', '06:29']
[14, 1, '06:06', '06:35']
[16, 2, '06:40', '07:09']
[26, 3, '07:16', '07:46']
[18, 4, '07:51', '08:20']
[17, 5, '08:28', '08:57']
[6, 6, '09:07', '09:36']
[15, 7, '09:47', '10:16']
[10, 8, '10:24', '10:53']
[12, 9, '10:58', '11:27']
[11, 10, '11:32', '12:01']
[13, 11, '12:07', '12:36']
[23, 12, '12:44', '13:13']
[29, 13, '13:20', '13:49']
[22, 14, '13:55', '14:24']
[27, 15, '14:34', '15:03']
[32, 16, '15:09', '15:38']
[33, 17, '15:47', '16:16']
[34, 18, '16:21', '16:50']
[38, 19, '16:55', '17:24']
[1, 20, '17:30', '18:00']
When given two days to work with, the solver drops just three nodes:
solution found. Objective value is 30019596
Dropped
[5, 9, 37]
Scheduled
[node, order, min time, max time]
[0, 0, '06:00', '06:01']
[22, 1, '06:11', '06:13']
[29, 2, '06:47', '06:49']
[24, 3, '07:21', '07:23']
[8, 4, '08:00', '08:02']
[30, 5, '08:42', '08:44']
[31, 6, '09:24', '09:26']
[27, 7, '10:03', '10:05']
[32, 8, '10:38', '10:40']
[28, 9, '11:21', '11:23']
[41, 10, '12:08', '12:10']
[39, 11, '12:43', '12:45']
[40, 12, '13:22', '13:24']
[36, 13, '14:03', '14:05']
[20, 14, '14:49', '14:51']
[25, 15, '15:28', '15:30']
[19, 16, '16:05', '16:07']
[21, 17, '16:42', '16:44']
[35, 18, '17:17', '17:19']
['Overnight at 42, dummy for 1', 19, '17:58', '18:00']
[16, 20, '06:02', '06:06']
[14, 21, '06:37', '06:40']
[15, 22, '07:15', '07:18']
[10, 23, '07:52', '07:55']
[12, 24, '08:26', '08:29']
[11, 25, '09:00', '09:04']
[13, 26, '09:36', '09:39']
[23, 27, '10:12', '10:16']
[7, 28, '10:53', '10:56']
[3, 29, '11:34', '11:38']
[2, 30, '12:14', '12:17']
[4, 31, '12:54', '12:57']
[6, 32, '13:44', '13:47']
[17, 33, '14:21', '14:25']
[18, 34, '15:00', '15:03']
[26, 35, '15:33', '15:37']
[34, 36, '16:14', '16:17']
[33, 37, '16:47', '16:51']
[38, 38, '17:21', '17:24']
[1, 39, '17:56', '18:00']
And with four days, the solver fits in all the nodes, but does not bother visiting the third overnight node.
solution found. Objective value is 22917
Dropped
[]
Scheduled
[node, order, min time, max time]
[0, 0, '06:00', '07:11']
[14, 1, '06:06', '07:17']
[12, 2, '06:43', '07:54']
[13, 3, '07:18', '08:29']
[23, 4, '07:55', '09:06']
[22, 5, '08:31', '09:42']
[29, 6, '09:07', '10:18']
[24, 7, '09:41', '10:52']
[8, 8, '10:20', '11:31']
[9, 9, '11:06', '12:17']
[30, 10, '11:46', '12:57']
[31, 11, '12:28', '13:39']
[27, 12, '13:07', '14:18']
[32, 13, '13:41', '14:52']
[28, 14, '14:24', '15:35']
[34, 15, '15:06', '16:17']
[33, 16, '15:40', '16:51']
[38, 17, '16:13', '17:24']
['Overnight at 42, dummy for 1', 18, '16:49', '18:00']
[35, 19, '06:00', '09:21']
[21, 20, '06:35', '09:56']
[20, 21, '07:10', '10:32']
[39, 22, '07:51', '11:12']
[41, 23, '08:27', '11:48']
[40, 24, '09:08', '12:29']
[36, 25, '09:49', '13:11']
[37, 26, '10:30', '13:51']
[5, 27, '11:26', '14:47']
[25, 28, '12:11', '15:32']
[19, 29, '12:48', '16:09']
[18, 30, '13:28', '16:49']
[26, 31, '14:01', '17:22']
['Overnight at 43, dummy for 1', 32, '14:38', '18:00']
[16, 33, '06:00', '11:21']
[17, 34, '06:40', '12:02']
[6, 35, '07:18', '12:40']
[4, 36, '08:08', '13:30']
[2, 37, '08:49', '14:11']
[3, 38, '09:28', '14:50']
[7, 39, '10:09', '15:31']
[11, 40, '10:48', '16:10']
[10, 41, '11:23', '16:45']
[15, 42, '12:00', '17:22']
[1, 43, '12:38', '18:00']
These solutions are not optimal, in that I only let the solver run for 30 seconds, but the objective function seemed to stabilize for at least a handful of seconds in the debug output, so they are probably pretty close to optimal. It is also interesting to notice that the sequence of nodes is different for each run depending on how many days are allowed. Because the vehicle has to get back to the overnight node at the end of each day, and because the solver is trying to minimize the total travel time of each vehicle, the routes change so as to minimize the extra travel needed to get back to the overnight nodes, and then every morning get from the overnight node back out to the regular nodes.
However, there is a subtle problem that is not obvious with the two day solution, but is obvious with the 4 day one. The parameters of the problem require that the vehicles leave the depot in the morning at 6:00 am. This is true on the first day, but it is not true on the second and subsequent days. This is due to the way that slack works.
Consider the end of the first day in the 4-day case. The solver says that the vehicle will visit node 38, the overnight node, and then node 35 the next morning.
[node, order, min time, max time]
...
[38, 17, '16:13', '17:24']
['Overnight at 42, dummy for 1', 18, '16:49', '18:00']
[35, 19, '06:00', '09:21']
The calculation for the value of time at the overnight node is:
cumul(42) = cumul(38) + transit(38, 42) + slack(38)
= (16:13) + (310.5s + 30m) + 0s
= 16:13 + 00:35 + 00:00
= 16:48
Note that the slack value at the “regular node” is zero, and that the “transit time” includes both the transit time from 38 to 42, as well as the service time of 30 minutes at node 38. Allowing for rounding error, the value I compute of 16:48 to the solver’s value of 16:49.
The computation for the minimum time at node 35 is similar, but this time the “service time” is the unloading value of -720 minutes (or -43,200 seconds), and the slack is non-zero. The slack can be as large as 360 minutes (or 21,600 seconds), depending on what is needed to fit the dimension constraints at the destination node.
cumul(35) = cumul(42) + transit(42, 35) + slack(42)
= (16:49) + (685.7s - 43200s) + (up to 21,600s)
= 60540s + (-42514.3s) + (up to 21,600s)
= 18025.7s + (up to 21,600s)
= 05:00 + (up to 21,600s)
In words, the earliest arrival time at the night node is 16:49.
Subtracting the 12 hours for the overnight node’s “service time” will
push the CumulVar down to 04:49. The transit time from the night
node (actually node 1) to node 35 is 685.7 seconds, or 11 minutes,
25.7 seconds. So without adding in any slack, the latest possible
“early arrival time” at node 35 is 05:00:25. Remember, the state of
the CumulVar is at the arrival at a node. There is a constraint
that the earliest arrival should be 6:00, so the slack needed to meet
that constraint is about 60 minutes, which is well under the maximum
slack of 360 minutes. So the solver can assign the early arrival time
at node 35 as 06:00.
However, this is wrong. What I really want is for the vehicle to leave the depot at 06:00, and then arrive at node 35 at 06:11. There isn’t any way to do that at the moment. What is needed is another set of dummy nodes, modeling a morning departure from the depot node, with the same time constraints of 06:00 to 18:00. In that case, the solver can assign the correct amount of slack to the overnight node to fix the arrival at the morning node to 06:00. If the morning nodes have a service time of zero, then this means the next node down the line will have the correct arrival time. Something like:
cumul(AM) = cumul(42) + transit(42, AM) + slack(42)
= (16:49) + (0s - 43200s) + (up to 21,600s)
= 60540s + (-42514.3s) + (up to 21,600s)
= 18025.7s + (up to 21,600s)
= 05:00 + (up to 21,600s)
= 06:00
cumul(35) = cumul(AM) + transit(AM, 35) + slack(AM)
= (06:00) + (685.7s + 0s) + 0s
= 06:11
Adding those extra morning nodes is similar to adding the overnight nodes, except with zero service time and zero slack.
# create dummy nodes for returning to the depot every night
night_nodes = range(num_nodes, num_nodes+num_days)
# create dummy nodes linked to night nodes that fix the AM depart time
morning_nodes = range(num_nodes+num_days, num_nodes+num_days+num_days)
total_nodes = num_nodes + len(night_nodes) +len(morning_nodes)
# Create the routing index manager.
manager = pywrapcp.RoutingIndexManager(total_nodes, 1, [0], [1])
As above, the slack at the morning nodes needs to be set to zero.
for node in morning_nodes:
index = manager.NodeToIndex(node)
time_dimension.SlackVar(index).SetValue(0)
The transit callback needs to be adjusted to expect the morning nodes as well.
def transit_callback(manager,
night_nodes,
morning_nodes,
from_index, to_index):
# Returns the travel time between the two nodes.
# Convert from routing variable Index to time matrix NodeIndex.
from_node = manager.IndexToNode(from_index)
to_node = manager.IndexToNode(to_index)
# adjust if either node is a night node
# adjust if either node is a night node
if from_node in night_nodes or from_node in morning_nodes:
from_node = 1
if to_node in night_nodes or to_node in morning_nodes:
to_node = 1
return Matrix[from_node][to_node]
}
And finally, the solver should be allowed to drop the morning nodes without consequences.
# Allow all overnight nodes to be dropped for free
for node in night_nodes:
routing.AddDisjunction([manager.NodeToIndex(node)], 0)
This almost works, but not quite. The solver can sometimes fall into weird solutions where the night nodes transition to night nodes, morning nodes transition to other mornings, and so on. To prevent these transitions, I added “infinite distances” to the transit times for transitions I do not want:
INFINITE_TIME = day_end*1000 # way more than ever possible
# night nodes can *only* transition to morning nodes
if from_node in night_nodes:
if to_node in morning_nodes:
return 0
else:
return INFINITE_TIME
# prevent movement to night nodes from morning nodes
if from_node in morning_nodes:
if to_node in night_nodes:
return INFINITE_TIME
# prevent movement to morning nodes from morning nodes
if from_node in morning_nodes:
if to_node in morning_nodes:
return INFINITE_TIME
# prevent movement to night nodes from night nodes
if from_node in night_nodes:
if to_node in night_nodes:
return INFINITE_TIME
# prevent movement to night nodes, morning nodes from depot nodes
if from_node == 0:
if to_node in night_nodes or to_node in morning_nodes:
return INFINITE_TIME
# prevent movement to start or end nodes from morning nodes
if to_node == 0 or to_node == 1:
if from_node in morning_nodes:
return INFINITE_TIME
In addition, I want to preserve the order of night nodes and morning nodes. To do this, I added a “counting dimension” (which I explained in its own post a while ago.
# make sure the days happen in order. first end day 1, end day 2,
# etc, then node 1
#
# create counting dimension
routing.AddConstantDimension(1, # increment by 1
total_nodes+1, # the max count is visit every node
True, # start count at zero
"Counting")
count_dimension = routing.GetDimensionOrDie('Counting')
In order to use the counting dimension, I need to set some constraints.
# use count dim to enforce ordering of overnight, morning nodes
solver = routing.solver()
for i in range(len(night_nodes)):
inode = night_nodes[i]
iidx = manager.NodeToIndex(inode)
iactive = routing.ActiveVar(iidx)
for j in range(i+1, len(night_nodes)):
# make i come before j using count dimension
jnode = night_nodes[j]
jidx = manager.NodeToIndex(jnode)
jactive = routing.ActiveVar(jidx)
solver.Add(iactive >= jactive)
solver.Add(count_dimension.CumulVar(iidx) * iactive * jactive <=
count_dimension.CumulVar(jidx) * iactive * jactive)
# if night node is active, AND night_node is not the last night,
# must transition to corresponding morning node
if i < len(morning_nodes):
i_morning_idx = manager.NodeToIndex(morning_nodes[i])
i_morning_active = routing.ActiveVar(i_morning_idx)
solver.Add(iactive == i_morning_active)
solver.Add(count_dimension.CumulVar(iidx) + 1 ==
count_dimension.CumulVar(i_morning_idx))
for i in range(len(morning_nodes)):
inode = morning_nodes[i]
iidx = manager.NodeToIndex(inode)
iactive = routing.ActiveVar(iidx)
for j in range(i+1, len(morning_nodes)):
# make i come before j using count dimension
jnode = morning_nodes[j]
jidx = manager.NodeToIndex(jnode)
jactive = routing.ActiveVar(jidx)
solver.Add(iactive >= jactive)
solver.Add(count_dimension.CumulVar(iidx) * iactive * jactive <=
count_dimension.CumulVar(jidx) * iactive * jactive)
These are a bit tricky, but this post is already much longer than I
want. The key is the use of ActiveVar() calls to make sure that a
node is active before constraining it. The ActiveVar() calls
essentially return a zero-one variable. This can be multiplied by the
constraint values to zero them out if the nodes are not active. If
either or both are nodes inactive, then the
conditions
solver.Add(count_dimension.CumulVar(iidx) * iactive * jactive <=
count_dimension.CumulVar(jidx) * iactive * jactive)
will degenerate to 0 <= 0.
With these additions, the multi-day solver run results finally look like what I want.
solution found. Objective value is 22917
Dropped
[]
Scheduled
[node, order, min time, max time]
[0, 0, '06:00', '07:11']
[14, 1, '06:06', '07:17']
[12, 2, '06:43', '07:54']
[13, 3, '07:18', '08:29']
[23, 4, '07:55', '09:06']
[22, 5, '08:31', '09:42']
[29, 6, '09:07', '10:18']
[24, 7, '09:41', '10:52']
[8, 8, '10:20', '11:31']
[9, 9, '11:06', '12:17']
[30, 10, '11:46', '12:57']
[31, 11, '12:28', '13:39']
[27, 12, '13:07', '14:18']
[32, 13, '13:41', '14:52']
[28, 14, '14:24', '15:35']
[34, 15, '15:06', '16:17']
[33, 16, '15:40', '16:51']
[38, 17, '16:13', '17:24']
['Overnight at 42, dummy for 1', 18, '16:49', '18:00']
[35, 19, '06:00', '09:21']
[21, 20, '06:35', '09:56']
[20, 21, '07:10', '10:32']
[39, 22, '07:51', '11:12']
[41, 23, '08:27', '11:48']
[40, 24, '09:08', '12:29']
[36, 25, '09:49', '13:11']
[37, 26, '10:30', '13:51']
[5, 27, '11:26', '14:47']
[25, 28, '12:11', '15:32']
[19, 29, '12:48', '16:09']
[18, 30, '13:28', '16:49']
[26, 31, '14:01', '17:22']
['Overnight at 43, dummy for 1', 32, '14:38', '18:00']
[16, 33, '06:00', '11:21']
[17, 34, '06:40', '12:02']
[6, 35, '07:18', '12:40']
[4, 36, '08:08', '13:30']
[2, 37, '08:49', '14:11']
[3, 38, '09:28', '14:50']
[7, 39, '10:09', '15:31']
[11, 40, '10:48', '16:10']
[10, 41, '11:23', '16:45']
[15, 42, '12:00', '17:22']
[1, 43, '12:38', '18:00']
One last item to notice. The final solution of the 3+ day run is the same as before: 22917. However, the solution of the 2 day run is different! The addition of the morning nodes fixing the departure time at 06:00 means that the best solution found drops 4 nodes, not 3.
solution found. Objective value is 40019479
Dropped
[5, 9, 28, 37]
Scheduled
[node, order, min time, max time]
...
Link to code
This post is a bit more complicated than can be digested in one go, so I’ve posted the code to github at https://github.com/jmarca/tsp_multiple_days. The original state is copied from the mailing list post, and all of my changes are detailed in the commit log (including the edits I made as I was writing this post). I added a command line option, so that I could test the differences between using morning nodes, and the effect of using the guided local heuristic (which I did not talk about here, but which is required to break out of local minima).
The list of options is:
python tsp_multiple_days.py -h
usage: tsp_multiple_days.py [-h] [--days DAYS] [--start START] [--end END] [--waittime SERVICE] [-t, --timelimit TIMELIMIT] [--debug] [--no_guided_local]
[--skip_mornings]
Solve routing problem to a fixed set of destinations
optional arguments:
-h, --help show this help message and exit
--days DAYS Number of days to schedule. Default is 2 days
--start START The earliest any trip can start on any day, in hours. Default 6
--end END The earliest any trip can end on any day, in hours. Default 18 (which is 6pm)
--waittime SERVICE The time required to wait at each visited node, in minutes. Default is 30
-t, --timelimit TIMELIMIT
Maximum run time for solver, in seconds. Default is 10 seconds.
--debug Turn on solver logging.
--no_guided_local Whether or not to use the guided local search metaheuristic
--skip_mornings Whether or not to use dummy morning nodes. Default is true