文章参考了ZhaoDongqiang（C++调用cplex求解VRPTW模型），补充了在Jyputer编辑器中用Python调用Cplex求解VRPTW问题，修正了模型中的小错误，具体如下：
目录
1 VRPTW数学模型
2 Python调用Cplex求解
2.1 调用所需的库
import pandas as pd
import numpy as np
from docplex.mp.model import Model
import matplotlib.pyplot as plt2.2 初始化参数
data.txt 内容如下：
40 50 0 0 1236 0
45 68 10 912 967 90
45 70 30 825 870 90
42 66 10 65 146 90
42 68 10 727 782 90
42 65 10 15 67 90
40 69 20 621 702 90
40 66 20 170 225 90
38 68 20 255 324 90
38 70 10 534 605 90
35 66 10 357 410 90
35 69 10 448 505 90
25 85 20 652 721 90
22 75 30 30 92 90
22 85 10 567 620 90
20 80 40 384 429 90
20 85 40 475 528 90
18 75 20 99 148 90
15 75 20 179 254 90
15 80 10 278 345 90
30 50 10 10 73 90
30 52 20 914 965 90
28 52 20 812 883 90
28 55 10 732 777 90
25 50 10 65 144 90
25 52 40 169 224 90
40 50 0 0 1236 0df = pd.read_csv('data.txt',sep='\t',header=None)
M=10000
c_num = 25 #顾客数
v_num = 3 #卡车数
q = 500 #卡车的容量
CM = [i for i in range (1,c_num+1)] #顾客列表
V = [i for i in range(1,v_num+1)] #卡车列表
N = [0]+CM+[26] #结点列表
#读取txt文件
x_loc = list(df.iloc[:,0])
y_loc = list(df.iloc[:,1])
d = list(df.iloc[:,2])
lbtw = list(df.iloc[:,3])
ubtw = list(df.iloc[:,4])
t_service = list(df.iloc[:,5])
t_service[0]
#各个节点之间的弧
A = [(i,j) for i in N for j in N if j!=i]
#各节点之间的距离
C={(i,j):np.floor(np.hypot((x_loc[i]-x_loc[j]),(y_loc[i]-y_loc[j]))*10)/10 for i,j in A}
#将两点距离作为时间
t={(i,j):0 for i,j in A}
for i,j in C:
t[i,j]=C[i,j]+t_service[i]
#画出坐标散点图
import matplotlib.pyplot as plt
plt.scatter(x_loc[1:25],y_loc[1:25],c='b')
plt.plot(x_loc[0],y_loc[0],c='r',marker='s')
for i in CM:
plt.annotate('%d'%i,(x_loc[i]+0.5,y_loc[i]))
plt.annotate("%d"%0,(x_loc[0]+0.5,y_loc[0]))
#初始化A_car 和 S_time，作为决策变量s和x字典的键值
A_car=[(k,i,j) for k in V for i,j in A ]
S_time=[(k,i) for k in V for i in N]
图1 各点的坐标散点图
2.3 调用docplex建模
mdl=Model("VRPTW")
#决策变量
x=mdl.binary_var_dict(A_car,name='x')
s=mdl.continuous_var_dict(S_time,name='s')
#目标函数
mdl.minimize(mdl.sum(C[i,j]*x[k,i,j] for k,i,j in A_car))
#约束1-8：
#约束1
mdl.add_constraints(mdl.sum(x[k,i,j] for k in V for j in N if j!=i\
)==1 for i in CM)
#约束2
mdl.add_constraints(mdl.sum(d[i]*mdl.sum(x[k,i,j] for j in N if j!=i)\
for i in CM)<=q for k in V)
#约束3
mdl.add_constraints(mdl.sum(x[k,0,j] for j in N if j!=0)==1 for k in V)
#约束4
mdl.add_constraints(mdl.sum(x[k,i,h] for i in N if i!=h)==mdl.sum(x[k,h,j] for j in N if j!=h) for h in CM for k in V)
#约束5
mdl.add_constraints(mdl.sum(x[k,i,c_num+1] for i in N if i!=(c_num+1))==1 for k in V)
#约束6
mdl.add_constraint(mdl.sum(x[k,0,j] for k in V for j in N if j!=0)<=v_num)
#约束7
mdl.add_constraints(s[k,i]+t[i,j]-M*(1-x[k,i,j])<=s[k,j] for k in V for i,j in A)
#约束8
mdl.add_constraints(s[k,i]<=ubtw[i] for k in V for i in N)
mdl.add_constraints(s[k,i]>=lbtw[i] for k in V for i in N);
#求解
solution=mdl.solve()
print(solution)
solution for: VRPTW objective: 191.3 x_1_0_20=1 x_1_20_24=1 x_1_21_26=1 x_1_22_21=1 x_1_23_22=1 x_1_24_25=1 x_1_25_23=1 x_2_0_13=1 x_2_12_26=1 x_2_13_17=1 x_2_14_12=1 x_2_15_16=1 x_2_16_14=1 x_2_17_18=1 x_2_18_19=1 x_2_19_15=1 x_3_0_5=1 x_3_1_26=1 x_3_2_1=1 x_3_3_7=1 x_3_4_2=1 x_3_5_3=1 x_3_6_4=1 x_3_7_8=1 x_3_8_10=1 x_3_9_6=1 x_3_10_11=1 x_3_11_9=1 s_1_1=967.000 s_1_2=825.000 s_1_3=65.000 s_1_4=727.000 s_1_5=15.000 s_1_6=621.000 s_1_7=170.000 s_1_8=255.000 s_1_9=534.000 s_1_10=357.000 s_1_11=448.000 s_1_12=652.000 s_1_13=30.000 s_1_14=567.000 s_1_15=384.000 s_1_16=475.000 s_1_17=99.000 s_1_18=179.000 s_1_19=278.000 s_1_20=10.000 s_1_21=917.000 s_1_22=825.000 s_1_23=732.000 s_1_24=132.000 s_1_25=224.000 s_1_26=1075.600 s_2_1=912.000 s_2_2=825.000 s_2_3=65.000 s_2_4=727.000 s_2_5=15.000 s_2_6=621.000 s_2_7=170.000 s_2_8=255.000 s_2_9=534.000 s_2_10=357.000 s_2_11=448.000 s_2_12=687.800 s_2_13=30.800 s_2_14=594.800 s_2_15=407.800 s_2_16=502.800 s_2_17=124.800 s_2_18=217.800 s_2_19=312.800 s_2_20=10.000 s_2_21=914.000 s_2_22=812.000 s_2_23=732.000 s_2_24=65.000 s_2_25=169.000 s_2_26=815.800 s_3_1=940.600 s_3_2=848.600 s_3_3=106.100 s_3_4=755.000 s_3_5=15.100 s_3_6=662.800 s_3_7=198.100 s_3_8=290.900 s_3_9=570.600 s_3_10=384.500 s_3_11=477.500 s_3_12=652.000 s_3_13=30.000 s_3_14=567.000 s_3_15=384.000 s_3_16=475.000 s_3_17=99.000 s_3_18=179.000 s_3_19=278.000 s_3_20=10.000 s_3_21=914.000 s_3_22=812.000 s_3_23=732.000 s_3_24=65.000 s_3_25=169.000 s_3_26=1049.200
2.4 求解结果
#选出那些决策变量为1时所经过的弧
active_arcs=[ a for a in A_car if x[a].solution_value > 0.9]
#画出VRPTW的结果
import matplotlib.pyplot as plt
plt.scatter(x_loc[1:25],y_loc[1:25])
plt.plot(x_loc[0],y_loc[0],c='r',marker='s')
for i in CM:
plt.annotate('%d'%i,(x_loc[i]+0.5,y_loc[i]))
plt.annotate("%d"%0,(x_loc[0]+0.5,y_loc[0]))
for k,i,j in active_arcs:
if k==1:
plt.plot((x_loc[i],x_loc[j]),(y_loc[i],y_loc[j]),c='r')
if k==2:
plt.plot((x_loc[i],x_loc[j]),(y_loc[i],y_loc[j]),c='g')
if k==3:
plt.plot((x_loc[i],x_loc[j]),(y_loc[i],y_loc[j]),c='b')
3 结语
目前Python调用Cplex的相关文章较少，本文在ZhaoDongqiang（C++调用cplex求解VRPTW模型）的基础之上，利用Python调用Cplex解决了VRPTW问题。