多商品流运输问题——Python实现
多商品流运输问题（Multi-Commodity Flow Problem，MCFP）是网络流理论中的一个重要问题，旨在在给定的网络中同时优化多种商品的流量分配。网络由节点和带容量限制和成本的有向边组成。每种商品都有其特定的源节点和汇节点，以及必须满足的需求量。目标是确定每种商品通过每条边的流量，使得总运输成本最小，同时满足网络的容量限制和各商品的需求。MCFP在许多实际应用中非常重要，如交通网络中的多种货物运输、通信网络中的数据包传输、供应链管理中的多产品配送等。它不仅考虑了如何在网络中传输单一类型的货物，还需要解决多个商品在共享同一网络资源时的协调问题。这一问题通常通过线性规划或整数规划方法来求解，利用约束条件确保在网络中有效地分配资源。
 |
 |
一、多商品运输问题数学模型
多商品流运输问题是一种在网络流理论中的问题，旨在通过一个网络进行多种商品的运输，使得每种商品的需求得到满足，同时最小化运输成本或最大化运输效率。
假设我们有一个图\(G = (V, E)\)，其中\(V\)是节点集，\(E\)是边集。每条边 $(i, j) \in E $ 都有一个容量\(c_{ij}\)，代表其最大运输能力。我们需要运输 \(K\)种商品，每种商品\(k\)有其源节点\(s_k\)和汇节点\(t_k\)，以及需求量\(d_k\)，这些商品可以混装运输。
我们用变量\(f_{ij}^k\)表示通过边\((i, j)\)运输商品\(k\)的流量。
- 目标函数
我们希望最小化总运输成本，假设每条边\((i, j)\)有单位运输成本\(cost_{ij}\)，则目标函数为：
\[\min \sum_{k=1}^K \sum_{(i, j) \in E} cost_{ij} \cdot f_{ij}^k
\]
- 约束条件
容量约束: 对于每条边\((i, j)\)，总流量不能超过其容量：
\[\sum_{k=1}^K f_{ij}^k \leq c_{ij}, \quad \forall (i, j) \in E
\]
流量守恒: 对于每个节点\(i\)和每种商品\(k\)，除源节点和汇节点外，流入等于流出：
\[\sum_{j: (i, j) \in E} f_{ij}^k - \sum_{j: (j, i) \in E} f_{ji}^k = 0, \quad \forall i \in V \setminus \{s_k, t_k\}, \forall k
\]
需求约束: 对于每种商品\(k\)，源节点和汇节点的流量需满足其需求：
\[\sum_{j: (s_k, j) \in E} f_{s_k j}^k - \sum_{j: (j, s_k) \in E} f_{js_k}^k = d_k\\
\sum_{j: (t_k, j) \in E} f_{t_k j}^k - \sum_{j: (j, t_k) \in E} f_{jt_k}^k = -d_k\]
二、Python求解
2.1 案例1
import pulp
import networkx as nx
import matplotlib.pyplot as plt
# 定义网络参数
nodes = ['A', 'B', 'C', 'D']
edges = {('A', 'B'): 10, ('A', 'C'): 10, ('B', 'C'): 5, ('B', 'D'): 15, ('C', 'D'): 10}
costs = {('A', 'B'): 1, ('A', 'C'): 2, ('B', 'C'): 1, ('B', 'D'): 3, ('C', 'D'): 1}
commodities = {
'K1': {'source': 'A', 'sink': 'D', 'demand': 5},
'K2': {'source': 'A', 'sink': 'D', 'demand': 10}
}
# 定义问题
prob = pulp.LpProblem("Multi-Commodity Flow Problem", pulp.LpMinimize)
# 定义变量
flow = {}
for (i, j) in edges:
for k in commodities:
flow[(i, j, k)] = pulp.LpVariable(f'flow_{i}_{j}_{k}', lowBound=0)
# 目标函数
prob += pulp.lpSum([costs[(i, j)] * flow[(i, j, k)] for (i, j) in edges for k in commodities])
# 容量约束
for (i, j) in edges:
prob += pulp.lpSum([flow[(i, j, k)] for k in commodities]) <= edges[(i, j)]
# 流量守恒约束
for k in commodities:
for node in nodes:
if node not in [commodities[k]['source'], commodities[k]['sink']]:
prob += pulp.lpSum([flow[(i, j, k)] for (i, j) in edges if i == node]) - \
pulp.lpSum([flow[(j, i, k)] for (j, i) in edges if i == node]) == 0
# 需求约束
for k in commodities:
source = commodities[k]['source']
sink = commodities[k]['sink']
demand = commodities[k]['demand']
prob += pulp.lpSum([flow[(i, j, k)] for (i, j) in edges if i == source]) - \
pulp.lpSum([flow[(j, i, k)] for (j, i) in edges if i == source]) == demand
prob += pulp.lpSum([flow[(i, j, k)] for (i, j) in edges if i == sink]) - \
pulp.lpSum([flow[(j, i, k)] for (j, i) in edges if i == sink]) == -demand
# 求解
prob.solve()
# 输出结果
for v in prob.variables():
print(v.name, "=", v.varValue)
print("Total Cost = ", pulp.value(prob.objective))
# 创建网络图
G = nx.DiGraph()
for (i, j) in edges:
G.add_edge(i, j, capacity=edges[(i, j)], cost=costs[(i, j)])
# 设置节点层次信息
node_layers = {
'A': 0,
'B': 1,
'C': 1,
'D': 2
}
# 添加节点层次信息
for node, layer in node_layers.items():
G.nodes[node]['subset'] = layer
# 绘制网络图
pos = nx.multipartite_layout(G, subset_key="subset") # 使用多层布局
plt.figure(figsize=(14, 10))
nx.draw(G, pos, with_labels=True, node_size=12000, node_color='lightblue', font_size=40, font_weight='bold', width=6, arrowsize=30)
edge_labels = {(i, j): f"c={edges[(i, j)]}, cost={costs[(i, j)]}" for (i, j) in edges}
nx.draw_networkx_edge_labels(G, pos, edge_labels=edge_labels, font_size=30, label_pos=0.5)
plt.title("Network Graph with Capacities and Costs", fontsize=40)
plt.show()
flow_A_B_K1 = 0.0
flow_A_B_K2 = 10.0
flow_A_C_K1 = 5.0
flow_A_C_K2 = 0.0
flow_B_C_K1 = 0.0
flow_B_C_K2 = 5.0
flow_B_D_K1 = 0.0
flow_B_D_K2 = 5.0
flow_C_D_K1 = 5.0
flow_C_D_K2 = 5.0
Total Cost = 50.0
2.2 案例2
每天早上，三个住宅区的人会去往该区域的两个工业中心和两个商业中心。下表按照出发地和目的地详细列出了各个流的基础信息。例如每天源于住宅节点4的6000次总出行中，有1250次是将工业节点7作为目的地。
住宅区至节点出行次数表
| 起始节点 / 目的地节点 | 1 | 2 | 3 | 4 | 5 | 6 | 7 |
|---|---|---|---|---|---|---|---|
| 1 | - | 900 | 750 | 40 | 10 | 600 | 550 |
| 4 | 100 | 2000 | 1100 | - | 150 | 1400 | 1250 |
| 5 | 110 | 4000 | 2200 | 200 | - | 3300 | 2440 |
目前由于地理因素限制，每种出行都仅有一条单一的路。上图网络的弧数字表示不同点之间的距离(km)。例如从节点1出发、去往节点7的人必须全程绕海港行驶，途中经过节点2~6。三个住宅区每天共产生21100次出行，而人们沿这些路线行驶的总距离则为399250千米/天。
区域规划者考虑采用多种改进办法，通过减少驾车行驶的千米数来减少空气污染。其中一个想法是在上图中用虚线弧(2,6)和(6,2)表示渡船。若引入一艘渡船，它可以在早高峰时段沿每个方向携带2000辆车。想知道这样做能够节省多少千米的行驶距离。
 |
 |
import networkx as nx
import matplotlib.pyplot as plt
# 出行数表格
flows = {
(1, 2): 900, (1, 3): 750, (1, 4): 40, (1, 5): 10, (1, 6): 600, (1, 7): 550,
(4, 1): 100, (4, 2): 2000, (4, 3): 1100, (4, 5): 150, (4, 6): 1400, (4, 7): 1250,
(5, 1): 110, (5, 2): 4000, (5, 3): 2200, (5, 4): 200, (5, 6): 3300, (5, 7): 2440,
}
# 各路径距离
distances = {
(1, 2): 3.5, (1, 3): 6.5, (1, 4): 11.5, (1, 5): 11.5, (1, 6): 8.5, (1, 7): 8.5,
(4, 1): 11.5, (4, 2): 8.0, (4, 3): 5.0, (4, 5): 14.5, (4, 6): 12.0, (4, 7): 14.5,
(5, 1): 11.5, (5, 2): 8.0, (5, 3): 11.0, (5, 4): 25.0, (5, 6): 4.0, (5, 7): 6.5,
(3, 2): 3.0, # 添加缺失的路径距离
}
# 计算引入渡船前的总距离
total_distance = sum(flows[route] * distances[route] for route in flows)
# 引入渡船后的距离
ferry_distance = 0.0
distances_with_ferry = distances.copy()
distances_with_ferry.update({
(2, 6): ferry_distance,
(6, 2): ferry_distance
})
# 更新路径距离
distances_with_ferry[(1, 6)] = distances[(1, 2)] + ferry_distance
distances_with_ferry[(4, 6)] = distances[(4, 3)] + distances[(3, 2)] + ferry_distance
distances_with_ferry[(5, 2)] = distances[(5, 6)] + ferry_distance
# 计算使用渡船的流量限制
ferry_capacity = 2000
# 调整流量表以考虑渡船容量
adjusted_flows = flows.copy()
adjusted_flows[(2, 6)] = min(ferry_capacity, flows.get((1, 6), 0) + flows.get((4, 6), 0))
adjusted_flows[(6, 2)] = min(ferry_capacity, flows.get((5, 2), 0))
# 分配流量：优先考虑使用摆渡，然后考虑原路线
remaining_flows = flows.copy()
if (2, 6) in adjusted_flows:
total_ferry_2_6 = min(adjusted_flows[(2, 6)], ferry_capacity)
remaining_flows[(1, 6)] -= min(remaining_flows[(1, 6)], total_ferry_2_6)
total_ferry_2_6 -= min(remaining_flows[(1, 6)], total_ferry_2_6)
remaining_flows[(4, 6)] -= total_ferry_2_6
if (6, 2) in adjusted_flows and (5, 2) in remaining_flows:
remaining_flows[(5, 2)] -= adjusted_flows[(6, 2)]
# 计算引入渡船后的总距离
total_distance_with_ferry = sum(remaining_flows[route] * distances[route] for route in remaining_flows)
total_distance_with_ferry += adjusted_flows[(2, 6)] * ferry_distance + adjusted_flows[(6, 2)] * ferry_distance
# 计算节省的距离
distance_saved = total_distance - total_distance_with_ferry
# 输出结果
print(f"引入渡船前的总行驶距离: {total_distance:.1f} 千米/天")
print(f"引入渡船后的总行驶距离: {total_distance_with_ferry:.1f} 千米/天")
print(f"节省的总行驶距离: {distance_saved:.1f} 千米/天")
# 绘制网络图
def draw_graph(edges, title, node_layers, ferry_edges=None):
G = nx.DiGraph()
for (u, v), d in edges.items():
G.add_edge(u, v, weight=d)
# 添加节点层次信息
for node, layer in node_layers.items():
G.nodes[node]['subset'] = layer
pos = nx.multipartite_layout(G, subset_key="subset")
edge_labels = {(u, v): f'{d}' for (u, v), d in edges.items()}
plt.figure(figsize=(14, 10))
nx.draw(G, pos, with_labels=True, node_size=12000, node_color='lightblue', font_size=60, font_weight='bold', width=6, arrowsize=30)
nx.draw_networkx_edge_labels(G, pos, edge_labels=edge_labels, font_size=40)
if ferry_edges:
nx.draw_networkx_edges(G, pos, edgelist=ferry_edges, edge_color='red', width=7.5, style='dashed')
plt.title(title, fontsize=60)
plt.show()
# 节点分层信息
node_layers = {
1: 0, 2: 1, 3: 2, 4: 3, 5: 2, 6: 1, 7: 0
}
# 设置中文字体
plt.rcParams['font.sans-serif'] = ['SimHei'] # 设置默认字体为黑体
plt.rcParams['axes.unicode_minus'] = False # 解决保存图像时负号'-'显示为方块的问题
# 绘制不加摆渡的网络图
draw_graph(distances, "不加摆渡的网络图", node_layers)
# 绘制加摆渡的网络图
ferry_edges = [(2, 6), (6, 2)]
draw_graph(distances_with_ferry, "加摆渡的网络图", node_layers, ferry_edges=ferry_edges)
引入渡船前的总行驶距离: 169650.0 千米/天
引入渡船后的总行驶距离: 124550.0 千米/天
节省的总行驶距离: 45100.0 千米/天
| 没摆渡 | 有摆渡 |
|---|---|
 |
 |
2.3 案例3
将城市 v0 生产的商品 A、B 通过铁路网分别运送到 城市 v6、v5，商品 A、B 的需求量分别为 6.0、5.0 单位。铁路网的结构如图所示，每条线路（边）上的数值表示该线路的容量和单位运费。求完成商品运输任务的最小费用的运输方案。这是一个多商品流问题，具有 1 个供应点 V0 和 2个需求点 V5、V6。对于多源多汇的多商品流最小费用流问题，构造一个费用网络 G，网络 G 的各条边具有容量限制，单位流量的运输成本为边的权值 w。源点（供应点）为流量净输出，汇点（需求点）为流量净输入；中间节点不带存储，每种商品的净流量都为 0。所有线路的商品总流量不超过边的容量限制。
import numpy as np
import matplotlib.pyplot as plt # 导入 Matplotlib 工具包
import networkx as nx # 导入 NetworkX 工具包
import pulp # 导入 pulp库
#5.1 创建有向图 (Multi-commodity Flow Problem)
G2 = nx.DiGraph() # 创建一个有向图 DiGraph
G2.add_edges_from([('v0','v1',{'capacity': 7, 'weight': 4}),
('v0','v2',{'capacity': 8, 'weight': 4}),
('v1','v3',{'capacity': 9, 'weight': 1}),
('v2','v1',{'capacity': 5, 'weight': 5}),
('v2','v4',{'capacity': 9, 'weight': 4}),
('v3','v4',{'capacity': 6, 'weight': 2}),
('v3','v5',{'capacity': 10, 'weight': 6}),
('v3','v6',{'capacity': 10, 'weight': 6}),
('v4','v1',{'capacity': 2, 'weight': 1}),
('v4','v5',{'capacity': 5, 'weight': 2}),
('v4','v6',{'capacity': 5, 'weight': 2})]) # 添加边的属性 'capacity', 'weight'
pos={'v0':(0,5),'v1':(4,2),'v2':(4,8),'v3':(10,2),'v4':(10,8),'v5':(15,3),'v6':(15,7)} # 指定顶点绘图位置
supplyA = 6.0 # 商品 A 供应量
supplyB = 5.0 # 商品 B 供应量
print("Supply A:{}\tSupply B:{}".format(supplyA,supplyB))
# 整理边的标签
edgeLabel1 = nx.get_edge_attributes(G2, 'capacity')
edgeLabel2 = nx.get_edge_attributes(G2, 'weight')
edgeLabel = {}
for i in edgeLabel1.keys():
edgeLabel[i] = f'({edgeLabel1[i]:},{edgeLabel2[i]:})' # 边的(容量，成本)
# 5.2 绘制有向网络图
fig, ax = plt.subplots(figsize=(8,6))
nx.draw(G2,pos,with_labels=True,node_color='skyblue',node_size=400,font_size=10) # 绘制有向图，显示顶点标签
nx.draw_networkx_nodes(G2, pos, nodelist=['v0'], node_color='orange',node_size=400) # 设置指定顶点的颜色、宽度
nx.draw_networkx_nodes(G2, pos, nodelist=['v5','v6'], node_color='c',node_size=400) # 设置指定顶点的颜色、宽度
nx.draw_networkx_edge_labels(G2,pos,edgeLabel,font_size=10) # 显示边的标签：'capacity','weight'
ax.set_title("Multi-commodity Flow Problem")
ax.text(-1.8,5.2,"A:{}".format(supplyA),color='m')
ax.text(-1.8,4.6,"B:{}".format(supplyB),color='navy')
ax.text(15.8,7.0,"A:{}".format(supplyA),color='m')
ax.text(15.8,2.8,"B:{}".format(supplyB),color='navy')
plt.xlim(-3, 18)
plt.ylim(1, 9)
plt.axis('on')
plt.show()
# 5.3 用 PuLP 求解多商品流最小费用问题 (Multi-commodity Flow Problem YouCans-XUPT)
edgeWeight = nx.get_edge_attributes(G2, 'weight') # 'weight', 单位流量的成本
edgeCapacity = nx.get_edge_attributes(G2, 'capacity') # 'capacity', 边的容量
maxCapacity = max(edgeCapacity.values()) # 边的容量的最大值
print(edgeWeight)
print("max(Weight)",max(edgeWeight.values()))
print("max(Capacity)",max(edgeCapacity.values()))
# (1) 建立优化问题 MCFproblem: 求最小值(LpMinimize)
MCFproblem = pulp.LpProblem("MultiCommodityFlowProb", sense=pulp.LpMinimize) # 定义问题，求最小值
# (2) 定义决策变量 fA(edges), fB(edges)
# itemsG2 = ["A","B"] # 商品种类
fA = pulp.LpVariable.dicts("FlowA", G2.edges(), lowBound=0.0, upBound=maxCapacity, cat='Continuous')
fB = pulp.LpVariable.dicts("FlowB", G2.edges(), lowBound=0.0, upBound=maxCapacity, cat='Continuous')
print(fA)
print(fB)
# (3). 设置目标函数
MCFproblem += pulp.lpSum([edgeWeight[edge] * (fA[edge]+fB[edge]) for edge in G2.edges]) # 总运输费用
# (4) 设置约束条件
for edge in G2.edges: # 边的最大流量约束
MCFproblem += (fA[edge] + fB[edge] - edgeCapacity[edge] <= 0) # edgeCapacity[edge], 边的容量
for node in G2.nodes: # 顶点的净流量约束
if G2.in_degree(node) == 0: # 入度为 0，判断是否为源点
print("源点：{}, 出边：{}".format(node,G2.out_edges(nbunch=node)))
MCFproblem += (sum(fA[edge] for edge in G2.out_edges(nbunch=node)) <= supplyA) # A 供应量约束
MCFproblem += (sum(fB[edge] for edge in G2.out_edges(nbunch=node)) <= supplyB) # B 供应量约束
elif G2.out_degree(node) == 0: # 出度为 0，判断是否为汇点
print("汇点：{}, 入边：{}".format(node,G2.in_edges(nbunch=node)))
if node=='v6': # 题目条件, v6 需求为 B
MCFproblem += (sum(fA[edge] for edge in G2.in_edges(nbunch=node)) == supplyA) # A 需求量约束
if node=='v5': # 题目条件, v5 需求为 A
MCFproblem += (sum(fB[edge] for edge in G2.in_edges(nbunch=node)) == supplyB) # B 需求量约束
else: # 中间节点，每种商品都是流量平衡
print("中间点：{}, 入边：{}, 出边：{}".format(node,G2.in_edges(nbunch=node),G2.out_edges(nbunch=node)))
MCFproblem += (sum(fA[edge1] for edge1 in G2.out_edges(nbunch=node))
- sum(fA[edge2] for edge2 in G2.in_edges(nbunch=node)) == 0) # 总流出 = 总流入
MCFproblem += (sum(fB[edge1] for edge1 in G2.out_edges(nbunch=node))
- sum(fB[edge2] for edge2 in G2.in_edges(nbunch=node)) == 0) # 总流出 = 总流入
# (5) 求解线性规划问题
MCFproblem.solve()
# (6) 输出优化结果
print("2. PuLP 线性规划（最小费用的优化结果）：") # PuLP 工具包
print(MCFproblem) # 输出问题设定参数和条件
print("Optimization result") # PuLP 工具包
for v in MCFproblem.variables():
print(v.name, " = ", v.varValue) # 输出每个变量的最优值
print("\n最小运输费用 = ", pulp.value(MCFproblem.objective)) # 输出最优解的目标函数值
plt.show()
#最小运输费用 = 105
多商品流运输问题（MCFP）涉及在一个带有容量限制和成本的有向图中，优化多种商品的流动。其核心在于确定如何在网络中分配不同商品的流量，以满足每种商品的需求并使总运输成本最小。该问题需要解决的关键挑战包括：网络容量限制，即每条边的流量不能超过其容量；流量守恒，即除源节点和汇节点外，每个节点的流入流量必须等于流出流量；需求约束，即每种商品从其源节点到汇节点的需求量必须得到满足。
在实际应用中，MCFP可以应用于交通系统、通信网络、物流配送等领域。例如，在一个城市的交通网络中，不同的商品（如食品、燃料、电子产品）需要通过道路网络从多个仓库运送到各个商店。MCFP的解决方案能够帮助制定最优的运输计划，确保各条道路不超载，并且总运输成本最低。MCFP通常通过数学优化技术求解，包括线性规划（LP）和整数规划（IP）。这些方法能够处理复杂的约束条件和多个目标函数，提供有效的流量分配方案。然而，由于问题的复杂性，求解大规模MCFP可能需要大量计算资源，因此常常借助启发式算法和分解技术来提高求解效率。
