Find Mipgap of a Solution Cplex
Using CPLEX and python for finding an exact solution for the CVRP
After seeing the formulation of our main issue which is CVRP ,
( https://medium.com/cmsa-algorithm-for-the-service-of-the-capacitated/state-of-art-of-the-capacitated-vehicle-rooting-problem-30ad4de6c2e9)
Now let's move to the serious work and implement a good program that allow us to generate an exact solution for the CVRP
First of all CPLEX IBM API is an optimizer solves integer programming problems and very large linear programming problems . It has a modeling layer called Concert that provides interfaces to the C++, C, and Java languages. There is a Python language interface based on the C interface. Additionally, connectors to Microsoft Excel and MATLAB are provided. Finally, a stand-alone Interactive Optimizer executable is provided for debugging and other purposes.
The rest of tools that we will be used is described in the table below
After describing our tools that we need to implement an exact solution for the CVRP we will focus on the development by explaining each iteration below:
- Making random numbers with the package Numpy and Creating an object called "rnd" for generating random numbers.
- Initializing our parameters witch are:
n: Number of customers.
Q: Maximum Capacity of the Vehicle.
N: The set of nodes.
V: The set of nodes + the node 0 which is the depot.
d: A collection that contains the demand of each node. - Initializing a random location for each node on a plot with two dimensions (x=200, y=100).
- Plotting the initial view of our random nodes.
- Initializing our set of arcs A and their cost c (c = The set of distances between each two nodes) . The cost in our case is the distance between each two nodes by calculating the differences between their locations in the plot.
- Calling the CPLEX API, initializing our decisions variables, constraints and printing the solution.
- Lastly, plotting our solution.
Coding the exact solution
#Importing the package numpy as np. import numpy as np # "rnd" is an object that generate random numbers. rnd = np.random #"seed(0)" is a method that reset (every time), the same random set of numbers. rnd.seed(0) # Number of customers. n = 10 # Maximum capacity of the vehicle. Q = 15 #The set of nodes without the depot. N = [i for i in range(1,n+1)] # The set of nodes + the depot. V = [0]+ N # A collection that contains the demand of each node. d ={i:rnd.randint(1,10)for i in N} # Generating random numbers between (0 and 15) * 200. loc_x = rnd.rand(len(V))*200 # Generating random numbers between (0 and 15) * 100. loc_y = rnd.rand(len(V))*100 #Importing the package matplotlib.pyplot as plt. import matplotlib.pyplot as plt #Plotting the n nodes without the node 0 (depot) and chose the color blue for each node. plt.scatter(loc_x[1:],loc_y[1:],c='b') # Associating and plotting each demand in the right of each blue node (customer). for i in N: plt.annotate('$d_%d=%d$'%(i,d[i]),(loc_x[i]+2,loc_y[i])) #Ploting the node 0, chosing the red like its color and the square form like a marker. plt.plot(loc_x[0],loc_y[0],c='r' ,marker='s') #Showing the Initial plot. plt.show() #Intializing the set of arcs A. A = [(i,j) for i in V for j in V if i!=j] #Calculating the distance between each node. c= {(i,j):np.hypot(loc_x[i]-loc_x[j],loc_y[i]-loc_y[j]) for i,j in A} #Importing the docplex.mp.model from the CPLEX as Model from docplex.mp.model import Model mdl = Model('CVRP') #Initializing our binary variable x_i,j x=mdl.binary_var_dict (A,name='x') #Initializing our cumulative demand u u=mdl.continuous_var_dict (N,ub=Q ,name = 'u') #Initializing the objectif function mdl.minimize(mdl.sum(c[i,j]*x[i,j]for i,j in A)) #Initialzing the first constraint mdl.add_constraints(mdl.sum(x[i,j]for j in V if j!=i)==1 for i in N) #Initialzing the second constraint mdl.add_constraints(mdl.sum(x[i,j]for i in V if i!=j)==1 for j in N) #Initialzing the third constraint mdl.add_indicator_constraints_(mdl.indicator_constraint(x[i,j],u[i]+d[j]==u[j])for i,j in A if i!=0 and j!=0) #Initialzing the fourth constraint mdl.add_constraints(u[i]>=d[i] for i in N) #Getting the solution solution =mdl.solve(log_output=True) #Printing the solution print(solution) #Identifing the active arcs. active_arcs =[ a for a in A if x[a].solution_value> 0.9] plt.scatter(loc_x[1:],loc_y[1:],c='b') for i in N: plt.annotate('$d_%d=%d$'%(i,d[i]),(loc_x[i]+2,loc_y[i])) for i,j in active_arcs : #Coloring the active arcs plt.plot([loc_x[i],loc_x[j]],[loc_y[i],loc_y[j]],c='g',alpha=0.3) plt.plot(loc_x[0],loc_y[0],c='r' ,marker='s') #Plotting the solution plt.show()
Interpreting and verifying the reliability of the results.
In this part we will verify our constraints if they are satisfied or not form the figure (A plot for the exact solution of the CVRP (10 costumers))that shows that there is no sub-tours so we can conclude that the two first constraints are satisfied with success.
For the rest of the constraints we will use the output below of the program that shows that always the cumulative demand ui + dj = uj and do not exceeds the total capacity Q
Find Mipgap of a Solution Cplex
Source: https://medium.com/cmsa-algorithm-for-the-service-of-the-capacitated/using-cplex-and-python-for-finding-an-exact-solution-for-the-cvrp-ac789ee0d8c4