数学建模实战：两阶段随机规划在物流网络应急优化中的应用-平芜编程栈

1. 从赛后总结到可复现的建模实战：一次完整的竞赛复盘

去年带队打完第十三届MathorCup的C题，那份31页的论文和几千行代码躺在硬盘里，总觉得不拿出来聊聊有点可惜。这不仅仅是一份“获奖总结”，更是一次对“电商物流网络应急调运与结构优化”这个经典运筹学问题的深度实践。很多同学在初次接触这类问题时，容易陷入两个极端：要么被复杂的网络流、整数规划模型吓退，觉得无从下手；要么就是照着一些现成的代码和模型硬套，最后论文写得空洞，模型解释力弱，经不起推敲。

我这次想做的，就是拆开这个黑箱。不谈空洞的“获奖感言”，而是聚焦于我们当时是如何一步步把“应急调运”和“结构优化”这两个抽象问题，转化为具体的数学模型、算法和代码的。你会看到我们团队在48小时高压下的真实思考路径、模型构建时的关键取舍、代码实现中遇到的坑，以及最终让论文逻辑自洽的那些“小心思”。无论你是正在备战MathorCup、国赛美赛的新手，还是对物流优化、数学建模感兴趣的朋友，希望这篇超过五千字的复盘，能给你带来比单纯看一篇优秀论文更实在的收获。

2. 赛题核心：如何理解“应急调运”与“结构优化”的耦合关系

拿到C题题目，第一感觉是信息量很大，场景很具体。但核心其实可以提炼为两个环环相扣的子问题，理解它们的耦合关系是建模成败的关键。

2.1 问题一：动态应急调运——在不确定性中寻找可行解

题目通常会给出一个基础的电商物流网络，包括多个仓库（分拨中心）、配送站，以及它们之间的运输线路和成本。然后，抛出一个或多个“应急事件”，比如某个主要仓库因故临时关闭，或某条关键运输线路中断。你的任务是，在满足所有客户点需求的前提下，重新规划包裹的流向，目标是使总运输成本最低，或满足时间限制等。

这里的“应急”二字是精髓。它意味着：

时间紧迫性：调度方案需要在很短时间内（可能是几小时）确定并执行，因此模型的计算效率至关重要，过于复杂的模型可能不实用。
资源约束性：应急状态下，某些线路的运力可能下降，某些仓库的处理能力可能饱和。模型必须严格处理这些容量约束。
需求刚性：客户需求必须被满足，不能像平时一样做延迟或取消，这增加了问题的难度。

我们最初的想法是直接建立一个多商品流网络优化模型。但很快发现，如果网络节点和边很多，商品种类（流向不同目的地的包裹）也多，模型规模会急剧膨胀，求解变得困难。因此，我们做了一个关键的简化：按“目的地区域”对包裹进行聚合。不是追踪每一个包裹，而是将去往同一配送区域的包裹视为一种商品。这大大减少了变量的数量，同时没有丢失问题的本质——我们需要决定的是从哪个仓库发多少货去哪个区域。

2.2 问题二：网络结构优化——为长远稳健性投资

第二个问题往往更宏观。它不再满足于“应急救火”，而是问：如果这类应急事件未来可能再次发生，甚至同时发生多起，我们该如何从长远角度优化物流网络本身？这可能包括：

新增或关闭哪些仓库？
扩建哪些仓库的容量？
加固或新增哪些运输线路？
在何处设置备用库存或缓冲资源？

这是一个典型的设施选址-网络设计混合问题，并且是两阶段随机规划或鲁棒优化的天然场景。第一阶段决策网络结构（哪些仓库开，容量多大），这些决策是“这里和现在”就要投入成本的；第二阶段则是在未来各种可能的应急场景（如仓库A失效、线路B中断等）下，进行运营调度，目标是使“投资成本 + 所有可能场景下的期望运营成本”最小。

我们团队当时最大的争论点在于：如何定义“可能的应急场景”？如果枚举所有仓库和线路的单点故障，场景数会爆炸。我们的解决方案是：基于历史数据或风险评估，选取发生概率最高或影响最大的K个关键节点/边失效作为代表性场景集。这既控制了模型规模，又抓住了主要风险。

2.3 耦合点：为什么不能分开求解？

这是最容易被忽视的一点。很多人会把问题一和问题二拆成两个独立的模型按顺序求解，这会导致严重的次优甚至不可行。

耦合的核心在于“成本”和“能力”。问题二（结构优化）中新建仓库、扩容线路都需要花钱，这些是固定投资。而问题一（应急调运）的运营成本（运输费）高度依赖于问题二给出的网络结构。一个结构脆弱的网络，即使应急调度算法再高明，运营成本也会居高不下；反之，一个过度投资、非常健壮的网络，虽然应急调度轻松，但总投资成本可能无法接受。

因此，必须建立一个统一的优化框架，将网络结构决策变量和每个场景下的流量决策变量同时纳入模型，目标函数是总投资成本加上所有场景下的期望运营成本。这样才能找到那个在“投资”与“运营”、“平时效率”与“应急韧性”之间的最佳平衡点。我们最终选择了一个基于场景的两阶段混合整数线性规划模型来刻画这种耦合关系。

3. 模型构建：从直觉到严格的数学语言

把思路转化为数学公式，是建模比赛的核心环节。这里分享我们模型的关键部分，以及背后的思考。

3.1 符号定义：清晰是严谨的第一步

我们花了相当时间定义索引、集合、参数和变量，确保无一歧义。例如：

集合：I仓库集合（包括现有和候选），J客户区域集合，S应急场景集合（包含正常场景s0）。
参数：d_js场景s下区域j的需求；c_ij从仓库i到区域j的单位运输成本；f_i开设或扩容仓库i的固定成本；u_i仓库i的最大处理能力；cap_ij线路(i,j)的运输容量。
变量：y_i（0-1变量）是否开设/扩容仓库i；x_ijs（连续变量）场景s下从仓库i运往区域j的货量。

注意：cap_ij这个参数需要仔细处理。在应急场景下，它可能变为0（线路中断）或一个更小的值（运力下降）。这需要在每个场景s的数据中预先定义好，是模型输入的一部分。

3.2 核心模型：两阶段随机规划框架

我们的最终模型骨架如下：

目标函数：Minimize: 总成本 = Σ_i (f_i * y_i) + Σ_s (prob_s * Σ_i Σ_j (c_ij * x_ijs)) （投资成本 + 所有场景的期望运营成本）

约束条件：

需求满足约束（对每个场景s，每个区域j）： Σ_i x_ijs = d_js（每个区域的需求在任何场景下都必须被满足，这是硬约束。）
仓库能力约束（对每个场景s，每个仓库i）： Σ_j x_ijs ≤ u_i * y_i（运出量不能超过该仓库的可用能力。注意y_i在这里的作用：如果仓库i未开设（y_i=0），则其运出量必须为0；如果已开设（y_i=1），则不能超过其最大能力u_i。这个约束将一、二阶段变量耦合在一起。）
线路容量约束（对每个场景s，每条线路(i,j)）： x_ijs ≤ cap_ijs（运量不能超过该线路在当前场景下的容量。cap_ijs是随场景变化的参数。）
非负与整数约束： x_ijs ≥ 0, y_i ∈ {0, 1}

这个模型看起来干净，但暗藏玄机。y_i是第一阶段决策，必须在观察到具体哪个场景发生之前就确定（因为建仓库是长期投资）。x_ijs是第二阶段决策，它依赖于第一阶段决策y_i和实际发生的场景s。目标函数中的prob_s是场景s发生的预估概率，这体现了决策者对风险的态度。

3.3 模型线性化与技巧：让求解器更好工作

我们的模型本质上是混合整数线性规划（MILP），这很棒，因为有很多成熟的求解器（如Gurobi, CPLEX）可以处理。但在实际中，我们仍做了一些处理使其更“友好”：

Big-M线性化：模型中u_i * y_i项已经是线性的，因为y_i是0-1变量。如果能力u_i本身也是一个决策变量（例如，可以选择不同的扩容等级），那么就会产生y_i * (连续变量)的非线性项。这时就需要引入一个足够大的常数M和额外的约束来进行线性化，这是MILP建模的经典技巧。
对称性破缺：如果问题中有多个完全相同的候选仓库位置，求解器可能会在对称的解之间来回搜索，浪费时间。我们可以添加一些约束，例如强制编号小的候选点优先被考虑，来打破这种对称性，加速求解。
有效不等式：我们根据问题特性，添加了一些显而易见的约束来收紧模型的线性规划松弛，帮助求解器更快找到最优解。例如，所有场景的总需求必须小于等于所有已开设仓库的总能力：Σ_s prob_s * Σ_j d_js ≤ Σ_i (u_i * y_i)。虽然这个约束可能被其他约束隐含，但显式写出能提供更好的下界。

4. 算法实现与代码踩坑实录

模型建好了，把它变成代码并求解才是真正的挑战。我们选择Python + Gurobi的组合，但过程绝非一帆风顺。

4.1 环境搭建与数据预处理

import gurobipy as gp from gurobipy import GRB import pandas as pd import numpy as np # 读取数据 warehouse_df = pd.read_csv('warehouse_data.csv') # 包含坐标、固定成本、能力等 demand_scenario_df = pd.read_csv('demand_scenarios.csv') # 各场景下各区域需求 cost_matrix_df = pd.read_csv('transport_cost.csv') # 运输成本矩阵 capacity_scenario_df = pd.read_csv('capacity_scenarios.csv') # 各场景下线路容量 # 预处理：将数据框转化为字典或列表，方便Gurobi建模 # 例如，运输成本字典 cost[(i,j)] = c_ij cost = {(row['from_wh'], row['to_region']): row['cost'] for _, row in cost_matrix_df.iterrows()}

踩坑点1：数据格式与索引一致性。这是最琐碎也最容易出错的地方。仓库、区域、场景都用唯一的ID（最好是整数或字符串）。确保在所有数据文件中，同一个实体的ID完全一致。我们曾因为一个文件里仓库ID是WH1，另一个文件里是1，导致建模时变量找不到对应的参数，报错信息又很隐晦，调试了半个多小时。

4.2 Gurobi建模核心代码片段

def build_and_solve_model(warehouses, regions, scenarios, prob, demand, cost, capacity, fixed_cost, max_cap): """ 构建并求解两阶段随机规划模型 """ model = gp.Model("Logistics_Network_Design") model.setParam('OutputFlag', 1) # 显示求解过程 model.setParam('TimeLimit', 3600) # 设置时间限制为1小时 # --- 添加变量 --- # 第一阶段变量：是否开设仓库 y = model.addVars(warehouses, vtype=GRB.BINARY, name="y") # 第二阶段变量：各场景下的运输量 x = model.addVars(((i, j, s) for i in warehouses for j in regions for s in scenarios), lb=0.0, name="x") # --- 设置目标函数 --- # 投资成本部分 investment_cost = gp.quicksum(fixed_cost[i] * y[i] for i in warehouses) # 期望运营成本部分 expected_op_cost = gp.quicksum(prob[s] * gp.quicksum(cost[i, j] * x[i, j, s] for i in warehouses for j in regions) for s in scenarios) model.setObjective(investment_cost + expected_op_cost, GRB.MINIMIZE) # --- 添加约束 --- # 1. 需求满足约束 for s in scenarios: for j in regions: model.addConstr(gp.quicksum(x[i, j, s] for i in warehouses) == demand[j, s], f"Demand_{j}_{s}") # 2. 仓库能力约束 (耦合约束！) for s in scenarios: for i in warehouses: model.addConstr(gp.quicksum(x[i, j, s] for j in regions) <= max_cap[i] * y[i], f"Capacity_{i}_{s}") # 3. 线路容量约束 for s in scenarios: for i in warehouses: for j in regions: if (i, j, s) in capacity: # 检查该线路在该场景下是否有容量数据 model.addConstr(x[i, j, s] <= capacity[i, j, s], f"RouteCap_{i}_{j}_{s}") # --- 求解 --- model.optimize() # --- 结果处理与返回 --- if model.status == GRB.OPTIMAL: solution_y = {i: y[i].X for i in warehouses} solution_x = {(i,j,s): x[i,j,s].X for (i,j,s) in x.keys()} return model.ObjVal, solution_y, solution_x else: print(f"求解未达到最优。状态码: {model.status}") return None, None, None

4.3 求解过程中的“惊险时刻”与调优

即使模型正确，求解大规模MILP也可能遇到麻烦。

问题规模与求解时间：当仓库、区域、场景数量达到几十、上百时，变量和约束的数量是乘积级增长。我们的第一个完整模型有近十万个变量，直接求解在时限内无法找到可行解。
- 我们的策略：
  1. 场景削减：采用K-means等方法对大量随机生成的场景进行聚类，用几个典型场景（聚类中心）及其权重来近似原分布，将场景数从上百个减少到10个以内。
  2. 启发式初始解：我们先快速求解一个简化模型（例如，只考虑正常场景，或放松整数约束），用得到的y_i解作为完整模型的初始解（model.setStart()），这能显著提升求解速度。
  3. Gurobi参数调优：我们调整了MIPGap（允许的间隙）、MIPFocus（侧重寻找可行解还是提升下界）等参数。对于我们的问题，将MIPFocus设为1（快速找到可行解）初期效果更好。
内存溢出：在添加约束的循环中，如果直接使用model.addConstrs()并传入一个巨大的生成器表达式，Gurobi在构建模型时可能会消耗大量内存。
- 我们的策略：改为在循环内逐条添加约束，虽然代码稍长，但内存控制更精细。对于特别大的约束集（如线路容量约束），我们评估后发现很多线路容量无限大（cap_ijs为一个很大的数），这些约束是无效的，在添加前就通过if语句过滤掉，减少了约30%的约束数量。
结果分析与验证：求解器输出最优解后，必须进行“常识验证”。
- 检查是否所有需求都被满足（约束1）。
- 检查开设的仓库是否真的被用到（y_i=1的仓库，其总运量sum(x_ijs)应大于0）。
- 检查在应急场景下，方案是否真的绕开了失效的节点或边。我们编写了可视化脚本，用网络图的形式展示正常和应急场景下的货流，一目了然。

5. 论文写作：如何将代码和结果转化为31页的故事

数学建模竞赛，论文是最终的交付物。模型再精巧，求解再成功，如果论文讲不清楚，一切白费。我们的31页论文（含附录）结构大致如下，重点分享几个核心部分的写法。

5.1 摘要：浓缩的精华，决胜的关键

摘要必须在有限字数内讲清所有关键点。我们采用“问题概述-模型思路-方法特色-主要结果”的四段式。

首句破题：直接点明研究的是“电商物流网络在突发事件下的应急调运与长期结构优化协同决策问题”。
模型亮点：明确指出我们建立了“基于场景的两阶段随机规划模型”，以“最小化固定投资与期望应急运营总成本”为目标，并强调了“按目的区域聚合需求”和“基于关键节点失效构建代表性场景集”两个简化技巧。
方法简述：说明采用混合整数线性规划（MILP）框架，并使用Gurobi优化器求解。
核心结论：用数据说话。例如：“结果表明，在给定投资预算下，优化后的网络能将重大应急事件下的平均运输成本提升控制在15%以内，相较于无优化网络降低约40%的期望总成本。” 并简要提及灵敏度分析的主要发现。

5.2 模型假设与符号说明：体现严谨性

这部分容易被轻视，实则重要。合理的假设能简化问题，聚焦核心。

我们明确假设“运输成本与货量成线性关系”、“单个场景内需求确定已知”、“应急事件发生时网络拓扑变化已知”。
符号说明采用三线表，清晰列出所有集合、参数、变量，并注明单位（如成本：元/件，能力：件/日）。

5.3 模型建立与求解：图文并茂，层层递进

这是论文的主体。我们避免了大段堆砌公式。

问题分析图：首先画了一张图，展示正常网络和应急场景下网络的对比，直观引出“调运”和“优化”两个问题。
模型推导：从最简单的单场景确定性模型讲起，逐步引入多场景、概率，最终导出两阶段随机规划模型。每一步都说明为什么这样扩展，让评审老师跟上思路。
算法流程图：画出了从数据输入、预处理、模型构建、求解到结果输出的完整流程图，并在旁边附上对应的代码文件名（如data_preprocess.py,main_model.py），体现工作的系统性。
关键代码展示：在正文中只贴最核心的代码片段（如目标函数和核心约束的添加代码），完整代码放附录。在代码旁用注释解释关键行。

5.4 结果分析与可视化：用图说话，深入讨论

这是展示工作深度的部分。

方案对比：我们设计了多个对比方案：① 仅优化当前网络（不投资）；② 均匀投资方案；③ 我们的优化方案。用表格对比它们的投资成本、正常运营成本、应急期望成本、总成本。
丰富的可视化：
- 网络拓扑图：用不同颜色和大小表示仓库的开设状态和能力。
- 货流桑基图：展示正常和主要应急场景下货物流向的变化，非常直观。
- 成本构成饼图/柱状图：展示总成本中投资和运营的占比。
- 灵敏度分析图：分析关键参数（如应急事件发生概率、投资预算）变化时，最优方案和总成本如何变化。这体现了模型的稳健性和洞察力。
管理启示：根据结果，提炼出几条给物流管理者的建议。例如：“投资应优先集中于网络中的关键枢纽节点，而非均匀分配”、“对于低概率高影响的‘黑天鹅’事件，建立少量战略性备用链路比普遍加固更经济”。

5.5 附录与代码规范

附录不是垃圾堆。我们将完整的代码、大规模的数据表格、额外的灵敏度分析结果放在这里。代码文件结构清晰，有详细的README.md说明运行环境（Python 3.8+, Gurobi 9.5+）和依赖库。重要函数都有文档字符串。这体现了良好的科研习惯。

6. 参赛心得：那些比获奖更重要的收获

回顾整个参赛过程，从最初的茫然到最后的豁然开朗，有几个体会特别深刻，可能比奖项本身更有价值。

第一，团队协作高于个人英雄主义。我们三人分工明确：一人主攻模型推导与论文写作（理论担当），一人负责算法实现与代码调试（编程担当），一人专注于数据整理、可视化与结果分析（数据担当）。但分工不意味着割裂。每天至少两次集中讨论，随时在白板上同步思路。经常出现的情况是，编程的同学在实现时发现了模型定义的一个模糊点，立刻提出来，大家共同修正。这种紧密的协作是能在高压下完成高质量作品的基础。

第二，“先求可行，再求优秀”。比赛时间有限，不要一开始就追求最复杂、最前沿的模型。我们首先用了一个下午，建立了一个最简化的单场景确定性模型并跑通，得到了一个基准解和完整的求解流程。这给了我们巨大的信心。在此基础上，我们再逐步增加随机场景、整数变量等复杂性。如果一开始就搞复杂的随机规划，很可能在调试中耗尽时间。

第三，重视“可解释性”和“讲故事”。数学模型不是炫技。每一个约束、每一个变量都应对应实际问题中的一个逻辑。在论文中，要能用通俗的语言把这种对应关系讲清楚。同样，结果分析不能只是罗列数字，要解读数字背后的业务含义。例如，我们发现优化方案建议关闭某个边缘仓库，不是因为它的运营成本高，而是因为它只在极少数特殊应急场景下才有用，为其支付固定投资不划算。这样的洞察才是论文的亮点。

第四，工具要熟，但思维更重要。熟悉Gurobi/Python/LaTeX当然重要，但最重要的是你如何定义问题、抽象模型、分析结果。比赛后期，我们大部分时间不是在写代码，而是在争论：“用这个代表性场景集是否足以反映风险？”“我们模型中的概率是客观数据还是主观估计？如果是主观的，如何进行鲁棒性处理？”这些思维层面的碰撞，才是数学建模训练的核心。

最后，那份31页的论文和代码，对我们而言，已不仅仅是一个竞赛作品。它更像一个完整的项目原型，清晰地记录了我们如何将一个复杂的现实问题，通过数学建模、算法设计和软件实现，最终转化为一个可量化、可分析、可支持决策的方案。这个过程本身，就是最大的收获。希望这篇复盘，能帮你少走一些我们曾经走过的弯路。