混合整数规划建模、求解TSP、VRP问题

宇宙微尘
2025-01-14 11:47:11
运筹优化
技术教程
本帖最后由 宇宙微尘 于 2025-1-23 10:41 编辑


TSP和VRP是在运输领域中常见的两个重要问题。这两个问题在不同的场景中都需要求解最优的路径或路线,以降低运输成本、优化资源利用效率。


MIP(混合整数规模模型)是一种常用的数学建模方法,可以有效地解决这些问题。对于小规模的TSP、VRP问题,可以使用MIP数学建模方法来快速求解。 


本文介绍通过pulp来建立MIP模型,来求解TSP和VRP问题。


安装pulp


pip install pulp

TSP问题


1.导入包


###导入需要的包
import itertools
import numpy as np
import pandas as pd
from scipy.spatial import distance_matrix
import matplotlib
import matplotlib.pylab as plt
import seaborn as sns
import pulp

import warnings
warnings.filterwarnings("ignore")

2.生成模拟数据


# 生成TSP的样例数据
np.random.seed(0)#随机种子,固定数据
n_customer = 9
n_point = n_customer + 1
df = pd.DataFrame({
    'x': np.random.randint(0, 100, n_point),
    'y': np.random.randint(0, 100, n_point),
})

df.iloc[0]['x'] = 0
df.iloc[0]['y'] = 0
df


3.生成距离矩阵


# 生成距离矩阵

distances = pd.DataFrame(distance_matrix(df[['x', 'y']].values, df[['x', 'y']].values), index=df.index, columns=df.index).values

fig, ax = plt.subplots(figsize=(8, 7))
sns.heatmap(distances, ax=ax, cmap='Blues', annot=True, fmt='.0f', cbar=True, cbar_kws={"shrink": .3}, linewidths=.1)
plt.title('distance matrix')
plt.show()


4.查看网点分布图


# 查看网点分布图
plt.figure(figsize=(5, 5))

# draw problem state
for i, row in df.iterrows():
    if i == 0:
        plt.scatter(row['x'], row['y'], c='r')
        plt.text(row['x'] + 1, row['y'] + 1, 'depot')
    else:
        plt.scatter(row['x'], row['y'], c='black')
        plt.text(row['x'] + 1, row['y'] + 1, f'{i}')
        
plt.xlim([-10, 110])
plt.ylim([-10, 110])
plt.title('points: id')
plt.show()


5.建模及求解


%%time

# 建立模型
problem = pulp.LpProblem('tsp_mip', pulp.LpMinimize)

# 设置变量
x = pulp.LpVariable.dicts('x', ((i, j) for i in range(n_point) for j in range(n_point)), lowBound=0, upBound=1, cat='Binary')
# 定义每个网点的名称1-9
u = pulp.LpVariable.dicts('u', (i for i in range(n_point)), lowBound=1, upBound=n_point, cat='Integer')

# 设置目标函数
problem += pulp.lpSum(distances[j] * x[i, j] for i in range(n_point) for j in range(n_point))

# 设置约束
for i in range(n_point):
    problem += x[i, i] == 0

for i in range(n_point):
    problem += pulp.lpSum(x[i, j] for j in range(n_point)) == 1
    problem += pulp.lpSum(x[j, i] for j in range(n_point)) == 1

# 消除子回路
for i in range(n_point):
    for j in range(n_point):
        if i != j and (i != 0 and j != 0):
            problem += u - u[j] <= n_point * (1 - x[i, j]) - 1
            
# 求解
status = problem.solve()

# 输出结果
status, pulp.LpStatus[status], pulp.value(problem.objective)


6.画最优线路图


# 画图

plt.figure(figsize=(5, 5))
for i, row in df.iterrows():
    if i == 0:
        plt.scatter(row['x'], row['y'], c='r')
        plt.text(row['x'] + 1, row['y'] + 1, 'depot')
        
    else:
        plt.scatter(row['x'], row['y'], c='black')
        plt.text(row['x'] + 1, row['y'] + 1, f'{i}')
        
plt.xlim([-10, 110])
plt.ylim([-10, 110])
plt.title('points: id')
# 画线路
routes = [(i, j) for i in range(n_point) for j in range(n_point) if pulp.value(x[i, j]) == 1]
arrowprops = dict(arrowstyle='->', connectionstyle='arc3', edgecolor='blue')
for i, j in routes:
    plt.annotate('', xy=[df.iloc[j]['x'], df.iloc[j]['y']], xytext=[df.iloc['x'], df.iloc['y']], arrowprops=arrowprops)
               
plt.show()


VRP问题


1.生成模拟数据


与TSP不同的是,生成了每个点的需求量、车辆的最大装载。


# 生成模拟数据
np.random.seed(0)

n_customer = 9
n_point = n_customer + 1
vehicle_capacity = 8

df = pd.DataFrame({
    'x': np.random.randint(0, 100, n_point),
    'y': np.random.randint(0, 100, n_point),
    'demand': np.random.randint(1, 5, n_point),
})

df.iloc[0]['x'] = 0
df.iloc[0]['y'] = 0
df.iloc[0]['demand'] = 0
df


2.生成距离矩阵


与TSP问题生成方式相同。



3.画网点定位图


画图


plt.figure(figsize=(5, 5))
for i, row in df.iterrows():
    if i == 0:
        plt.scatter(row['x'], row['y'], c='r')
        plt.text(row['x'] + 1, row['y'] + 1, 'depot')
    else:
        plt.scatter(row['x'], row['y'], c='black')
        demand = row['demand']
        plt.text(row['x'] + 1, row['y'] + 1, f'{i}({demand})')
        
plt.xlim([-10, 110])
plt.ylim([-10, 110])
plt.title('points: id(demand)')
plt.show()


4.建模和求解


%%time

demands = df['demand'].values

# 建立模型
problem = pulp.LpProblem('cvrp_mip', pulp.LpMinimize)

# 设定变量
x = pulp.LpVariable.dicts('x', ((i, j) for i in range(n_point) for j in range(n_point)), lowBound=0, upBound=1, cat='Binary')
n_vehicle = pulp.LpVariable('n_vehicle', lowBound=0, upBound=100, cat='Integer')

# 设定目标函数
problem += pulp.lpSum([distances[j] * x[i, j] for i in range(n_point) for j in range(n_point)])

# 设置约束
for i in range(n_point):
    problem += x[i, i] == 0
   
for i in range(1, n_point):
    problem += pulp.lpSum(x[j, i] for j in range(n_point)) == 1
    problem += pulp.lpSum(x[i, j] for j in range(n_point)) == 1
        
problem += pulp.lpSum(x[i, 0] for i in range(n_point)) == n_vehicle
problem += pulp.lpSum(x[0, i] for i in range(n_point)) == n_vehicle

# 消除子环路
subtours = []
for length in range(2, n_point):
     subtours += itertools.combinations(range(1, n_point), length)

for st in subtours:
    demand = np.sum([demands for s in st])
    arcs = [x[i, j] for i, j in itertools.permutations(st, 2)]
    problem += pulp.lpSum(arcs) <= np.max([0, len(st) - np.ceil(demand / vehicle_capacity)])

# 模型求解
status = problem.solve()

# 输出结果
status, pulp.LpStatus[status], pulp.value(problem.objective)


5.查看使用了几辆车


pulp.value(n_vehicle)

6.画最终的线路图


# 画图

plt.figure(figsize=(5, 5))

# 画点
for i, row in df.iterrows():
    if i == 0:
        plt.scatter(row['x'], row['y'], c='r')
        plt.text(row['x'] + 1, row['y'] + 1, 'depot')
    else:
        plt.scatter(row['x'], row['y'], c='black')
        demand = row['demand']
        plt.text(row['x'] + 1, row['y'] + 1, f'{i}({demand})')
        
plt.xlim([-10, 110])
plt.ylim([-10, 110])
plt.title('points: id(demand)')

# 画线路
cmap = matplotlib.cm.get_cmap('Dark2')
routes = [(i, j) for i in range(n_point) for j in range(n_point) if pulp.value(x[i, j]) == 1]

for v in range(int(pulp.value(n_vehicle))):
   
    # 定义每一个车辆
    vehicle_route = [routes[v]]
    while vehicle_route[-1][1] != 0:
        for p in routes:
            if p[0] == vehicle_route[-1][1]:
                vehicle_route.append(p)
                break

    # 不同颜色,画每一个车辆的线路
    arrowprops = dict(arrowstyle='->', connectionstyle='arc3', edgecolor=cmap(v))
    for i, j in vehicle_route:
        plt.annotate('', xy=[df.iloc[j]['x'], df.iloc[j]['y']], xytext=[df.iloc['x'], df.iloc['y']], arrowprops=arrowprops)
               
plt.show()


本文详细介绍了使用pulp通过建立MIP模型来解决TSP、VRP问题的步骤,代码均直接提供,感兴趣的朋友可直接复制本文代码运行。




本文转载自微信公众号:Python学习杂记

1103
0
0
0
关于作者
相关文章
  • 量子计算赋能能源材料革命:从分子模拟到材料设计,开启高效研发 ...
    国际团队在《ACS Energy Letters》发表《Harnessing Quantum Computing for Energy Materials: O ...
    了解详情 
  • 量子力学百年之际的 MLIP:化学精度、效率与通用泛化的演进与展 ...
    随着量子力学在2025年迎来其百年纪念,机器学习原子间势(MLIP)已逐渐发展为分子建模领域的变革性 ...
    了解详情 
  • 具身智能再迎突破!俄罗斯团队用量子退火技术破解机器人运动控制 ...
    内容提要一项由俄罗斯多家科研机构联合完成的最新研究表明,量子退火(Quantum Annealing)技术 ...
    了解详情 
  • Transformer 的尽头是 Ising 机
    华人学者闪耀2026元旦,前有DeepSeek mHC:一次将 Transformer 残差流拉回重整化轨道的重大升级 ...
    了解详情 
  • 机器学习赋能靶向蛋白降解药物设计:PROTACs与分子胶的技术综述 ...
    靶向蛋白降解(TPD)通过利用泛素–蛋白酶体系统(UPS)实现对致病相关蛋白的催化性去除。蛋白 ...
    了解详情 
领取成功
本月5个550bit真机配额已发放给您,配额将在2个月后到期,请及时使用哦~
活动中心
联系我们
二维码
返回顶部
返回
活动中心

完成任务,轻松获取真机配额

×
每日必做
新手任务
长期任务
其他任务
快速回复 返回顶部 返回列表
玻色有奖小调研
填写问卷,将免费赠送您1个1000bit真机配额
(单选) 您是从哪个渠道得知我们的?*
您是从哪个社交媒体得知我们的?*
您是通过哪个学校的校园宣讲得知我们的呢?
取消

提交成功

真机配额已发放到您的账户,可前往【云平台】查看

量子AI开发者认证

考核目标

开发者能够成功搭建Kaiwu-PyTorch-Plugin项目基础环境,并成功运行示例代码,根据示例提示,输出指定的值并填写至相应的输入框中。

通过奖励

5个一年效期的1000量子比特真机配额

专属「量子AI开发者」社区认证标识

开发者权益

每月固定权益:5个550量子比特真机配额
前往考核

第一步

按照README提示成功安装Kaiwu-PyTorch-Plugin库环境依赖
前往GitHub

第二步

运行 community-assessment 分支下的 run_rbm.py 代码示例

第三步

理解示例代码,手动打印并填写如下数值:

正相采样的状态

负相采样的状态

正相的能量值

负相的能量值

*

提交答案

开发者权益

每月固定权益:5个550量子比特的真机配额

恭喜您完成考核

您将获得量子AI开发者认证标识及考核奖励

1000 bit*5

配额

Quantum AI Developer Certification

Assessment Objectives

Developers should successfully set up the basic environment for the Kaiwu-PyTorch-Plugin project, run the QBM-VAE sample code, and calculate the correct FID value based on the random seed value provided by the system.

Pass Rewards

10 quotas for 550-qubit real quantum machines with a one-year validity period

Exclusive "Quantum AI Developer" Community Certification Badge

Developer Benefits

Fixed Monthly Benefits: 5 quotas for 550-qubit real quantum machines
Proceed to Assessment

Step 1

Install the environment dependencies for the Kaiwu-PyTorch-Plugin library according to the README instructions
Go to GitHub

Step 2

Replace the Seed Value

Your seed value is

Step 3

Enter the FID Value You Calculated

*

Submit Answer

Developer Benefits

Fixed Monthly Benefits: 5 quotas of 550-qubit real machines

Congratulations on Completing the Assessment

You will receive the Quantum AI Developer Certification Badge and Assessment Rewards

550bit*10

Quotas