Python数学规划案例：路径优化CVRP

容量限制的车辆路径问题(Capacitated Vehicle Routing Problem, CVRP)可描述为从一个配送中心出发，用多辆车对多个需求点进行配送，满足货车容量约束下使配送路程最短。每个需求点都有固定的货物需求，货物不可拆分，只能由一辆车配送；每辆车携带的货物总量不能超过货车的容量。

数学模型

数据集

http://vrp.atd-lab.inf.puc-rio.br/index.php/en/

包含节点坐标、货物需求量和货车容量等数据。

Python代码

0)载入packages

from math import sqrt

from docplex.mp.model import Model

import matplotlib.pyplot as plt

import numpy as np

1)把数据集文件下载下来，本文的例子使用“P-n16-k8.vrp”

n表示节点数；k表示车辆数；配送中心为0号节点。

2)生成class CVRP，指定读写文件目录

class CVRP:

def __init__(self):

self.mDir = "D:\\pycharm\\CVRP\\"

self.mDir_input = "i\\"

3)定义读数据函数read_data

def read_data(self):

ls_fn = self.mDir + self.mDir_input

filename = "P-n16-k8.vrp"

with open(ls_fn + filename, "r") as f:

data_list = []

for i in range(9):

data_list.append(f.readline().strip().split(":"))

K = int(data_list[2][1])

N = int(data_list[5][1])

Q = int(data_list[7][1])

node_x = []

node_y = []

D = []

C = np.zeros([N, N])

for i in range(N):

n, x, y = f.readline().strip().split(" ")

node_x.append(float(x))

node_y.append(float(y))

f.readline()

for i in range(N):

n, d = f.readline().strip().split(" ")

D.append(int(d))

for i in range(N):

for j in range(N):

C[i][j] = sqrt((node_x[i]-node_x[j])**2+(node_y[i]-node_y[j])**2)

C = np.round(C).astype("int32")

return N, K, D, C, Q, node_x, node_y

4)定义求解数学规划模型的函数solvebycolex

其中涉及到cplex的函数或属性，例如 binaryvarcube, binaryvarmatrix, integervarmatrix, minimize, sum, addconstraints, addconstraint, solve, objectivevalue, solutionvalue。

def solve_by_cplex(self, N, K, D, C, Q):

# Model

mdl = Model("CVRP")

# Decision variables

x = mdl.binary_var_cube(N, N, K, name="x")

y = mdl.binary_var_matrix(N, K, name="y")

u = mdl.integer_var_matrix(N, K, name="u")

# Objective

mdl.minimize(mdl.sum(C[i, j] * x[i, j, k] for i in range(N) for j in range(N) for k in range(K))) # (1)

# Constraints

mdl.add_constraints(mdl.sum(y[i, k] for k in range(K)) == 1 for i in range(1, N)) # (2)

mdl.add_constraint(mdl.sum(y[0, k] for k in range(K)) == K) # (3)

mdl.add_constraints(mdl.sum(x[i, j, k] for j in range(N)) == mdl.sum(x[j, i, k] for j in range(N))

for i in range(N) for k in range(K)) # (4)

mdl.add_constraints(mdl.sum(x[i, j, k] for j in range(N)) == y[i, k] for i in range(N) for k in range(K)) # (5)

mdl.add_constraints(mdl.sum(D[i] * y[i, k] for i in range(N)) <= Q for k in range(K)) # (6)

mdl.add_constraints(u[i, k] - u[j, k] + Q * x[i, j, k] <= Q - D[j] for i in range(1, N) for j in range(1, N)

for k in range(K)) # (7)

mdl.add_constraints(D[i] <= u[i, k] for i in range(1, N) for k in range(K)) # (8)

mdl.add_constraints(u[i, k] <= Q for i in range(1, N) for k in range(K)) # (8)

# 约束(9)可加可不加，加了求解速度更快

mdl.add_constraints(x[i, i, k] == 0 for i in range(N) for k in range(K)) # (9)

# solve

mdl.solve(log_output=True)

# output

print("obj:%s" % mdl.objective_value)

tour = {}

for k in range(K):

list_i = []

list_j = []

list_k = [0]

for i in range(N):

for j in range(N):

if x[i, j, k].solution_value > 0:

list_i.append(i)

list_j.append(j)

for l in range(len(list_i)):

ds = list_i.index(list_k[-1])

list_k.append(list_j[ds])

tour[k] = list_k

print("K=", k, ":", list_k)

return tour

5)定义函数plot绘制路径图

def plot(self, tour, node_x, node_y):

for k in range(len(tour)):

list_k = tour[k]

c = (np.random.rand(), np.random.rand(), np.random.rand())

for i in range(0, len(list_k) - 1):

plt.plot([node_x[list_k[i]], node_x[list_k[i + 1]]],

[node_y[list_k[i]], node_y[list_k[i + 1]]], color=c, linewidth=2)

plt.scatter(node_x, node_y)

plt.show()

6)编写main部分代码

if __name__ == "__main__":

a = CVRP()

N, K, D, C, Q, node_x, node_y = a.read_data()

tour = a.solve_by_cplex(N, K, D, C, Q)

a.plot(tour, node_x, node_y)

扩展阅读

Python 数学规划

Python数学规划案例:一维装箱

Python数学规划案例四：资源约束的最短路径

Python数学规划案例三：最短路径

Python数学规划案例二

Python数学规划之Cplex之旅

Python数学规划案例一

Python数学规划案例一模型

Python: 数学规划

Python数学结构

Python数据结构

Python图

Python二叉树

Python排序

Python查找算法

Python递归

Python链表

Python栈

Python基础:没有更简单

几分钟做个Web应用

Python利器

Python: 数据库之SQLite

Python: Pandas