MVRPDExample
The VRPD (Multi-Period Vehicle Routing Problem with Due Dates) is a problem in which the VMND has achieved excellent results, it has outperformed other heuristics as well, specially in instances with huge networks (100 customers).This problem can be succintly described with the following parameters and characteristics:
- A directed graph G(V, E), in which the set V = {0, 1, ..., n} consists of a depot (the 0 vertex) and the customers or retailers (1, ..., n), denoted as V'. A non-negative cost c_{ij} or a distance is provided for each pair of vertices.
- An amount of periods H.
- A maximum amount of vehicles m that can be used in every period.
- A demand q_{i} for every customer in V'.
- A capacity Q for every vehicle.
- A set of optional customers C, a subset of V'. The demand for the other customers must be met within the periods 1, ..., H.
- Every mandatory customer has asociated to it a release date r_{i} and a due date d_{i}. These are the periods when the products of the i'th customer arrive to the depot and are available to be dispatched and the last period in which his products can be delivered respectively.
- In case the optional customers are not visited, we incur in a penalty of p_{i} in the objective function.
- There is a cost of storing inventory on the depot for every period after the release date of the products of h_{i} for customer i.
As these parameters are considered, MVRPD minimizes the sum of the routing cost, the holding costs in the depot and finally the penalty costs. The optimal solution will be the routes and quantities of each vehicle in each period. Among other type of variables, the formulation of Larrain et al 2019 defines the binary variables x_i_j_t for the decision of traversing the arc (i, j) by some vehicle in the period t and all neighborhoods are built upon them.
In this article two neighborhoods are implemented, namely the period neighborhood and the K-Vicinities Neighborhood. As it can be seen in the illustration down below, the example portrays the period neighborhood for the parameterization t = 2.
Neighborhood definition for the implemented VMND algorithm for the MVRPD.
Also, the K-Vicinities can be seen for k = 3, which freezes the x variables which have at least one vertex inside the three nearest neighbors of i_{p} (V_{3}(i_{p}) with light blue). Let's recall that according to this implementation the parameterizations of the second neighborhood are asociated to each of the vertices, for which the set V_{k}(i_{p}) is calculated.
To start with, we will asume we are located inside the Examples folder (the file MVRPDExample will be contain all this code as well), otherwise, the paths of the imports shall be modified relatively.
import sys
import os
sys.path.append(os.path.pardir)
sys.path.append( os.path.join(os.path.pardir, 'Instances'))
import MVRPD as MVRPD
from sklearn.neighbors import NearestNeighbors
from Neighborhood import Neighborhoods
from VMND import solver
Notice that all these imports are not always necessary, since we are using sklearn only for computing the k-nearest neighbors and MVRPD only for:
- Loading the instance, defining the model, constraints and objective function.
- exporting the guroby model as a mps instance to the MIPLIB folder.
As we can read in the user's guide, the first step to solve an instance is to have the model in a mps format. To do this we have already programmed in the module MVRPD methods for this. As we can see in the code down below, the instance is loaded and written as mps in the MIPLIB folder with the self.exportMPS() method.
inst1 = MVRPD.MVRPD( os.path.join( 'MVRPDInstances' , 'ajs5n25_h_3.dat' ) )
inst1.exportMPS()
path = inst1.pathMPS
Finally we save the path into a string varaible.
Of course, if you wish to solve another instance or model, the code works just the same!
In order to calculate the second neighborhood, it is necessary to compute the k-nearest neighbors of each vertex from the network. To do this, we use the NearestNeighbors method from sklearn. These indices are stored inside the indices variable.
X = inst1.positions
nbrs = NearestNeighbors(n_neighbors=20, algorithm='ball_tree').fit(X)
indices = nbrs.kneighbors(X)[1]
The function to be used is declared (to know more about this function look at the --User's guide--). The important thing about this function is that accepts the parameters associated to the variable's name, the current neighborhood and current parameterization. In this case, it separates the variable's name into its prefix and the indexes and stores these components into the list elements. From this information, the function infers the period and the nodes associated to the variable and will return True if this varaible is to be fixed.
def fNbhs(varName, depth, param):
elements = varName.split('_')
if len(elements) < 4:
return False
else:
tl = int(elements[3])
il = int(elements[1])
jl = int(elements[2])
if depth == 1:
return tl != param
elif depth == 2:
return il not in indices[param] and jl not in indices[param]
else:
print('Error 23 Nbhds Function!! ')
return 0
return True
Another necessary parameter when we are using functional neighborhoods is outerNeighborhoods. In this case, each value of the keys of the dictionary are tuples which contain the names of the parameterizations. In this case they are defined by comprehension.
outer = {
1 : tuple([ tf for tf in range(1, inst1.H + 1) ]),
2 : tuple([ i for i in range(1, inst1.V + 1) ])
}
The list of "keys" or names of varaibles that will be fixed is as well necessary. We should recall here that is it not necessary to put the name of all variables but only those which play a role in Local Search. In the case of MVRPD, these varaibles are the routing variables.
klist = ['x_{}_{}_{}'.format( i, j, t )
for t in range(1, inst1.H + 1) for i in range(inst1.V + 1) for j in range(inst1.V + 1) if i != j]
Finally, the neighborhoods are declared:
functionNeighborhoods = Neighborhoods(
lowest = 1,
highest = 2,
keysList= klist,
randomSet=False,
outerNeighborhoods=outer,
funNeighborhoods= fNbhs,
useFunction=True)
Once we declare successfully the neighborhoods, we just have to invoke the solver method. In this case we don't need lazy constraints, alpha is set to 1 and the time limit are 300 seconds.
modelOutput = solver(
path = path,
addlazy = False,
funlazy= None,
importNeighborhoods= True,
importedNeighborhoods= functionNeighborhoods,
funTest= inst1.genTestFunction(),
alpha = 1,
callback = 'vmnd',
verbose = True,
minBCTime = 0,
timeLimitSeconds= 300
)
The test function inst1.genTestFunction() is used to check subtour feasability of the final solution, in general this is not a necessary parameter (and could be None), but since it is already implemented in the class MVRPPD we put it here.
Of course, all the previous code is already implemented in the class MVRPD with many more features, so you can do the same as before just by invoking the self.run() method.
inst1.run(
outImportedNeighborhoods='function',
writeResult=False,
outVerbose=True,
outCallback = 'vmnd'
)
After running the heuristic with the previous method, the results can be visualized with the following method:
inst1.visualizeRes()