1. 从赛题到方案:一次完整的数学建模实战复盘
又到了一年一度的MathorCup,看着今年A题《新能源城市配送优化》的题目,是不是感觉既熟悉又头疼?熟悉的是,这又是一个典型的运筹优化问题,核心离不开车辆路径规划、资源调度和成本控制;头疼的是,题目里那些“动态需求”、“充电策略”、“时间窗”的约束,还有“新能源”这个背景带来的额外复杂度,怎么把它们揉成一个逻辑自洽、可求解的模型,再变成能跑出结果的代码,每一步都够喝一壶的。我指导过不少队伍,也自己上手做过很多次,发现大家最大的瓶颈往往不是某个具体的算法,而是从拿到题目到交出论文这个完整链条的断裂。很多人模型建得天花乱坠,但代码实现不了;或者代码跑出来了,但论文里讲不清楚逻辑。今天,我就以2024年MathorCup A题(虽然具体题目细节未公开,但结合“新能源城市配送优化”这一高度相关的热词,我们可以构建一个极具代表性的分析框架)为引子,抛开那些华而不实的理论,从头到尾拆解一遍数学建模的完整实战过程。我会重点分享我们是怎么想的、为什么这么选、以及实现过程中那些教科书里不会写的“坑”。目标很简单:让你看完之后,不仅能复现出一个完整的项目,更能掌握一套应对这类优化问题的通用思考框架和工程化方法。
2. 破题与抽象:把现实问题翻译成数学语言
拿到一个像“新能源城市配送优化”这样的题目,第一步绝对不是打开MATLAB或者Python就开始写代码。那等于蒙着眼睛跑步,大概率会撞墙。我们首先要做的,是当好一个“翻译官”,把充满现实细节的题目描述,精准地翻译成严谨的数学语言。
2.1 核心要素提取与定义
题目通常会给出一个故事背景,比如“某物流公司拥有若干辆新能源货车,需要服务分布在城市的多个客户点,客户有确定的需求量和服务时间窗,车辆从配送中心出发,电量有限,途中可能需要充电...”。我们的任务就是从这段话里,把所有的“零件”都拆解出来,并给它们起好数学名字。
实体定义:
- 配送中心 (Depot):通常记为节点
0,是所有车辆的起点和终点。 - 客户点 (Customers):记为节点集合
C = {1, 2, ..., n}。每个客户i都有一系列属性,这是后续建模的基石。 - 车辆 (Vehicles):记为集合
K = {1, 2, ..., m}。这里的关键是“新能源”属性,意味着每辆车k有一个初始电量Battery_k,一个最大容量Q_k(载重),以及一个耗电率(通常与行驶距离和载重相关)。
- 配送中心 (Depot):通常记为节点
参数与变量梳理:
- 客户参数:对于每个客户
i,我们需要明确其需求量d_i(重量或体积),服务时间s_i(卸货耗时),以及最重要的——时间窗[e_i, l_i]。这意味着车辆到达客户i的时间必须在e_i之后,且在l_i之前开始服务。这是硬约束,违反则方案不可行。 - 网络参数:任意两点
i和j之间的距离d_ij或行驶时间t_ij。这里要注意,城市配送中,时间往往比单纯的欧氏距离更重要,可能需要考虑路网、拥堵等因素,题目数据有时会直接给出时间矩阵。 - 新能源特性参数:这是区别于传统车辆路径问题 (VRP) 的关键。包括:
- 车辆百公里耗电量(或单位距离耗电函数)。
- 电池容量。
- 充电站信息(如果有的话):位置、充电功率(快充/慢充)、充电成本。
- 充电策略:是充满电再走,还是可以充到任意电量?充电时间如何计算?(通常是充电量除以充电功率)。
- 客户参数:对于每个客户
决策变量定义:这是模型的“心脏”,它描述了我们的方案到底是什么。
- 最核心的是一组二进制变量
x_{ijk}:如果车辆k从节点i行驶到节点j,则为1,否则为0。 - 为了处理时间窗和电量,我们还需要引入一系列辅助变量,如车辆
k到达节点i的时间A_{ik},离开节点i的时间D_{ik},以及在节点i服务后的剩余电量E_{ik}。
- 最核心的是一组二进制变量
为什么这么定义?因为x_{ijk}这种流变量的定义方式,是描述路径问题最经典和强大的工具,它能自然地与图论、网络流理论结合,方便我们写出流平衡约束(每个客户点只能被访问一次,车辆从仓库出发最后返回仓库)。而时间A_{ik}和电量E_{ik}作为状态变量,是串联起整个顺序决策过程的关键,它们的前后依赖关系构成了模型的核心动态。
2.2 目标函数确立:我们到底要优化什么?
题目问“优化”,优化哪个指标?成本最低?时间最短?车辆最少?还是综合考量?MathorCup的题目通常会有明确导向,也可能需要你自行定义一个合理的目标。对于城市配送,常见的目标有:
总成本最小化:这是最务实的。成本可能包括:固定成本(使用的车辆数 * 每车固定成本)、运输成本(与总行驶距离成正比)、时间成本(司机工时,或早到/晚到的惩罚)、充电成本。
Minimize: α * (车辆数) + β * (总行驶距离) + γ * (总行驶时间) + δ * (总充电成本) + ζ * (时间窗违反惩罚)- 这里的系数 α, β, γ, δ, ζ 可能需要根据题目赋予权重,或者通过层次分析法等确定。
总行驶距离/时间最小化:在车辆数固定的情况下,这是一个很直接的目标。
使用的车辆数最小化:这有助于降低固定资产和司机成本。
在我们的模拟分析中,我们假设一个以总成本最小化为核心的复合目标,因为它最贴近商业实际。我们需要在模型中为每一种成本找到对应的数学表达式。例如,时间窗违反惩罚,可以设计为软约束:如果早于e_i到达,产生等待成本;如果晚于l_i到达,产生延迟惩罚,惩罚函数可以是线性的,也可以是阶梯式的。
2.3 约束条件形式化:给解决方案戴上“紧箍咒”
约束是把我们的数学方案拉回现实的缰绳。必须一条条罗列清楚:
流平衡约束:确保每个客户点被且仅被一辆车访问一次,车辆从仓库出发并返回仓库。
∑_k ∑_j x_{ijk} = 1, 对于所有客户点i。(每个客户只被服务一次)∑_j x_{0jk} = 1, 对于所有车辆k。(每辆车从仓库出发)∑_i x_{ihk} - ∑_j x_{hjk} = 0, 对于所有节点h和车辆k。(路径连续性,即到达某点后必须离开)∑_i x_{i0k} = 1, 对于所有车辆k。(每辆车返回仓库)
载重量约束:车辆在任何路段上的累计载重不能超过其最大容量。
- 这需要引入新的辅助变量
Load_{ik}表示车辆k离开节点i时的载重。约束为:Load_{jk} >= Load_{ik} + d_j - M*(1 - x_{ijk}),其中M是一个很大的数。这确保了如果车辆k从i走到j,那么在j点的载重必须等于i点载重加上j点的需求量。
- 这需要引入新的辅助变量
时间窗约束:这是难点之一。
- 首先有时间传递关系:
A_{jk} >= D_{ik} + t_{ij} - M*(1 - x_{ijk})。即,如果从i到j,那么到达j的时间至少是离开i的时间加上路程时间。 - 服务时间:
D_{ik} = A_{ik} + s_i(假设服务立即开始)。 - 硬时间窗要求:
e_i <= A_{ik} <= l_i。如果是软时间窗,则将其放入目标函数作为惩罚项。
- 首先有时间传递关系:
电量约束(核心难点):这是新能源VRP (E-VRP) 的灵魂。
- 电量消耗模型:最简单的假设是耗电量与距离成正比
e_ij = c * d_ij。更复杂的模型可能考虑载重影响e_ij = (c0 + c1 * Load) * d_ij。 - 电量传递关系:
E_{jk} <= E_{ik} - e_ij + M*(1 - x_{ijk})。即,到达j点的电量小于等于离开i点的电量减去路途消耗。 - 电量边界约束:
0 <= E_{ik} <= Battery_max。电量不能为负,也不能超过电池容量。 - 充电决策整合:如果模型包含充电站(记为节点集合
F),复杂度激增。我们需要决定是否去充电站、去哪个、充多少。这需要引入新的二进制变量y_{ifk}表示是否在节点i(可能是客户点或仓库)后去充电站f,以及充电量变量Ch_{ik}。电量约束将变为分段函数,在普通节点间消耗,在充电站节点补充。
- 电量消耗模型:最简单的假设是耗电量与距离成正比
实操心得:在纸上或白板上画出所有这些约束的依赖图,理解变量之间如何相互影响,是构建正确模型的关键。很多同学模型报错“不可行”,根源就是约束之间存在隐藏的矛盾,比如时间窗太紧,而车辆速度太慢,导致无论如何都无法满足。这时可能需要回头检查问题假设,或者引入软约束。
3. 模型求解策略:精确解与启发式的权衡
当我们把目标函数和所有约束都用数学公式写出来后,一个完整的混合整数线性规划 (MILP) 模型就诞生了。你可以用CPLEX,Gurobi,OR-Tools等求解器直接去求解。对于小规模问题(比如客户点<50),这可能行得通。但对于稍具规模的实际问题(客户点>100),MILP模型可能会因为计算复杂度过高而无法在比赛时间内得到最优解,甚至得不到可行解。
这时,就必须诉诸启发式算法。我们的策略通常是“分而治之,逐步优化”。
3.1 初始解构造:先有一个“能用”的方案
“从零到一”比“从一到优”更重要。一个快速的初始解构造算法能为后续优化提供起点。常用方法:
- 最近邻法:从仓库出发,总是选择距离当前点最近且满足所有约束(时间窗、电量、载重)的未服务客户点加入路径。如果无法加入,则派出一辆新车。这种方法速度快,但解的质量一般。
- 节约算法:这是解决VRP的经典启发式。其核心思想是,比较将两个客户点分别用两辆车服务与合并到同一辆车服务所带来的距离“节约值”,优先合并节约值大的点对。我们需要在合并时实时检查约束是否满足。
- 针对新能源的适配:在构造初始解时,就必须考虑电量。一个简单的策略是,在最近邻选择时,不仅看距离,还看前往该客户后,剩余电量是否足够返回最近的可能充电点(或仓库)。这相当于增加了一个“电量安全边际”的检查。
代码片段示意(最近邻法思路):
def construct_initial_solution(customers, vehicle_capacity, battery_capacity, distance_matrix): routes = [] unserved = customers.copy() while unserved: current_route = [0] # 从仓库开始 current_load = 0 current_battery = battery_capacity current_time = 0 while True: # 找出所有未服务且满足约束的候选客户 candidates = [] for cust in unserved: if (current_load + cust.demand <= vehicle_capacity and current_battery >= distance_matrix[current_route[-1]][cust.id] * energy_rate and current_time + travel_time <= cust.time_window_end): # 简化检查 # 计算一个综合得分,如:距离倒数 + 紧急程度 score = 1.0 / distance_matrix[current_route[-1]][cust.id] + (1.0 / (cust.time_window_end - current_time)) candidates.append((score, cust)) if not candidates: break # 当前车无法再添加客户,结束路径 # 选择得分最高的客户 _, next_cust = max(candidates, key=lambda x: x[0]) current_route.append(next_cust.id) current_load += next_cust.demand current_battery -= distance_matrix[current_route[-2]][next_cust.id] * energy_rate current_time += travel_time + next_cust.service_time unserved.remove(next_cust) # 简单电量检查:如果电量不足以返回仓库,则提前结束路径(实际中应考虑充电) if current_battery < distance_matrix[current_route[-1]][0] * energy_rate: break current_route.append(0) # 返回仓库 routes.append(current_route) return routes注意:这是一个极度简化的示意,忽略了时间窗的精确计算、充电决策等复杂情况,但它展示了如何将约束检查融入解构造过程。
3.2 局部搜索与元启发式优化:让方案变得“更好”
有了初始解,我们就可以开始“折腾”它,寻找更好的方案。核心是定义一系列邻域动作,然后搜索这些动作。
常用邻域动作:
- Relocate:将一个客户点从一条路径移动到另一条路径的某个位置。
- Swap:交换两条路径中的两个客户点。
- 2-opt:在一条路径内部,反转一段子路径的顺序。这对优化单条路径的行驶距离非常有效。
- Cross-exchange:交换两条路径中的两段子路径。
优化框架:
- 局部搜索:遍历当前解的所有可能邻域动作,接受第一个能使目标函数改进的动作,然后在新解上重复此过程,直到找不到改进动作为止。容易陷入局部最优。
- 模拟退火:允许以一定的概率接受恶化解,从而有机会跳出局部最优。核心参数是初始温度、降温速率和终止温度。
- 变邻域搜索:系统性地切换不同的邻域结构进行搜索,当在一个邻域中找不到更好解时,就切换到另一个更大的邻域。
- 遗传算法:将路径编码为染色体,通过选择、交叉、变异等操作模拟进化过程。编码和交叉算子的设计需要技巧,要保证生成的后代仍是可行解。
为什么选择VNS或模拟退火?对于数学建模竞赛,变邻域搜索是一个非常好的选择。它结构清晰,易于实现和调试,而且通常能取得不错的效果。模拟退火需要对温度参数有较好的把握。遗传算法则更复杂,调整参数多,在有限时间内不一定能调出最佳性能。
3.3 针对新能源特性的特殊优化算子
除了通用的路径优化算子,我们必须设计专门针对电量约束的算子:
- 充电站插入/删除/移动:在路径中智能地插入充电站访问。策略可以是:当检测到某段行程后电量可能不足时,在沿途寻找最近的充电站插入;或者评估移除一个充电站是否仍能满足全程电量要求。
- 电量可行的路径片段交换:在进行
Swap或Cross-exchange时,交换后的新路径片段必须重新进行电量模拟,确保从片段起点到终点的电量消耗是可行的。如果不可行,则需要尝试在片段内部或附近插入充电站。 - 基于电量的路径分割与合并:对于一条过长的路径,可以基于电量约束将其在合适的位置(如充电站附近)分割成两条;反之,也可以尝试合并两条电量允许的短路径。
实操心得:在实现这些优化算子时,可行性检查函数的效率至关重要。你需要一个快速函数,给定一条路径(包含客户点和可能的充电站序列),能迅速判断其是否满足载重、时间窗和电量约束。这个函数会被调用成千上万次,它的效率直接决定了整个优化过程的快慢。建议将路径的累积载重、累积时间、累积耗电量等状态预先计算并缓存,避免每次检查都从头计算。
4. 代码实现与工程化要点
模型和算法设计得再好,最终都要落地为代码。这里我用Python为例,分享几个关键的实现模块和踩坑点。
4.1 数据结构设计:一切效率的基础
不要用简单的列表套列表来管理一切。设计清晰的数据结构能让你的代码更健壮、更高效。
class Customer: def __init__(self, id, x, y, demand, service_time, time_window_start, time_window_end): self.id = id self.x = x self.y = y self.demand = demand self.service_time = service_time self.tw_start = time_window_start self.tw_end = time_window_end class Vehicle: def __init__(self, id, capacity, battery_capacity, energy_rate): self.id = id self.capacity = capacity self.battery_capacity = battery_capacity self.energy_rate = energy_rate # 单位距离能耗 class ProblemInstance: def __init__(self): self.customers = [] # 索引0通常是仓库 self.vehicles = [] self.distance_matrix = None self.time_matrix = None self.charging_stations = [] def calculate_matrices(self): # 预计算所有点之间的距离和时间矩阵 # 这是一个O(n^2)的操作,但只需做一次,可以节省大量后续计算 pass4.2 核心算法模块实现
- 初始解生成器:实现前面提到的节约算法或最近邻法,并确保返回的是一个包含多条路径的可行解。
- 邻域动作生成器:实现
relocate,swap,2-opt等操作。关键是要生成增量变化,而不是每次都对整个解进行深拷贝和全量评估。def generate_relocate_moves(solution): moves = [] for route_from_idx, route_from in enumerate(solution.routes): for i, node_i in enumerate(route_from[1:-1]): # 跳过首尾的仓库 for route_to_idx, route_to in enumerate(solution.routes): if route_from_idx == route_to_idx: continue for j in range(1, len(route_to)): # 插入到路径中的位置 # 创建一个“移动”对象,记录从哪条路、哪个位置、移到哪条路、哪个位置 move = RelocateMove(route_from_idx, i, route_to_idx, j) # 快速评估这个移动的成本变化 delta_cost delta_cost = evaluate_relocate_delta(solution, move) moves.append((delta_cost, move)) return moves - 可行性检查与成本评估:这是最核心的模块。对于一条给定的路径序列,需要模拟车辆行驶过程,计算到达每个点的时间、剩余电量,并检查约束。
def evaluate_route(route, problem, vehicle): current_time = 0 current_load = 0 current_battery = vehicle.battery_capacity total_cost = 0 feasible = True prev_node = problem.depot for node in route[1:]: # route[0]是仓库 # 计算距离和时间 dist = problem.distance_matrix[prev_node.id][node.id] travel_time = problem.time_matrix[prev_node.id][node.id] # 检查电量 energy_consumed = dist * vehicle.energy_rate if current_battery < energy_consumed: feasible = False break current_battery -= energy_consumed # 检查时间窗(假设为硬约束) arrival_time = current_time + travel_time if arrival_time < node.tw_start: current_time = node.tw_start # 等待 elif arrival_time > node.tw_end: feasible = False break else: current_time = arrival_time # 服务 current_time += node.service_time current_load += node.demand if current_load > vehicle.capacity: feasible = False break # 累加成本(例如距离成本) total_cost += dist * cost_per_km prev_node = node # 最后返回仓库 dist_back = problem.distance_matrix[prev_node.id][problem.depot.id] if current_battery < dist_back * vehicle.energy_rate: feasible = False total_cost += dist_back * cost_per_km return feasible, total_cost, current_time注意:这个函数在优化循环中会被调用无数次。任何微小的优化都能带来显著的加速,比如使用
numpy数组进行矩阵运算,缓存部分计算结果。
4.3 调试与可视化:相信你的眼睛
在复杂的优化算法中,Bug难以避免。除了看日志,可视化是终极武器。
- 路径可视化:使用
matplotlib绘制所有车辆的行驶路径,用不同颜色区分。一眼就能看出路径是否交叉、是否合理、是否有车辆空跑。import matplotlib.pyplot as plt def plot_solution(solution, problem): plt.figure(figsize=(10, 8)) colors = plt.cm.tab10(np.linspace(0, 1, len(solution.routes))) for route, color in zip(solution.routes, colors): x_coords = [problem.customers[node_id].x for node_id in route] y_coords = [problem.customers[node_id].y for node_id in route] plt.plot(x_coords, y_coords, 'o-', color=color, linewidth=2, markersize=8) # 标注路径顺序 for i, node_id in enumerate(route): plt.annotate(str(node_id), (x_coords[i], y_coords[i]), fontsize=9) # 标出仓库 depot = problem.depot plt.plot(depot.x, depot.y, 'ks', markersize=12, label='Depot') plt.legend() plt.grid(True, alpha=0.3) plt.title('Vehicle Routes') plt.show() - 收敛曲线:绘制每次迭代后最优解的成本变化曲线。这能帮你判断算法是否在有效搜索、是否已收敛、模拟退火的温度设置是否合适。
- 关键变量跟踪:在日志中输出每次重要邻域动作后的成本变化、可行性状态,帮助你定位是哪个算子或哪个检查函数出了问题。
踩坑实录:我曾遇到一个Bug,2-opt算子优化后,路径的总距离确实缩短了,但总成本却增加了。排查了很久才发现,我的成本计算函数只考虑了距离成本,但2-opt改变了客户访问顺序,导致某些客户的时间窗被违反,产生了高额惩罚,而我在评估动作时,只快速计算了距离变化,没有重新计算时间窗惩罚。教训是:任何邻域动作的增量评估,必须覆盖目标函数的所有组成部分。
5. 从结果到论文:如何讲好你的建模故事
代码跑出了结果,只算成功了一半。另一半是把你的工作清晰、有说服力地呈现在论文里。论文不是代码的流水账,而是一个有逻辑的故事。
5.1 模型阐述:清晰性与严谨性并重
- 符号说明表:这是论文的门面,务必清晰、完整。按照实体、集合、参数、决策变量的顺序列出,并给出单位(如距离:km,时间:minute)。
- 模型公式:将你在“破题”阶段写出的目标函数和约束,用专业的数学公式排版。建议使用
LaTeX环境(如Overleaf)撰写论文,这是数学建模的标配。公式要编号,并在文中引用。 - 模型解释:不要假设评委能一眼看懂你的公式。在公式下方,用一两句话解释这个约束的物理或商业含义。例如,在写出电量约束后,可以说明“该约束确保了车辆在任意路段的电量消耗不会导致其剩余电量为负,模拟了电池的物理限制”。
5.2 算法描述:突出你的创新与适配
- 流程图是必备的:用清晰的流程图展示你算法的整体框架,比如“初始解构造 -> 变邻域搜索主循环(包含几个邻域算子) -> 接受准则 -> 终止条件”。这比大段文字描述直观得多。
- 伪代码:对于核心算法(如你的变邻域搜索主循环、关键的邻域算子),给出伪代码。伪代码应介于编程语言和自然语言之间,突出逻辑而非语法细节。
- 强调针对性的设计:专门用一小节说明你是如何适配“新能源”和“时间窗”特性的。例如,“在初始解构造中,我们加入了电量安全边际检查”;“在邻域动作的可行性验证中,我们设计了一个快速的电量模拟函数,其时间复杂度为O(L),其中L为路径长度”。
5.3 实验设计与结果分析:用数据说话
- 基准数据集:如果题目没有给数据,你需要自己生成或寻找公开的基准数据集(如Solomon的VRPTW数据集)。说明你生成数据的规则(客户点分布、时间窗宽度、需求量范围等),以保证实验的可复现性。
- 对比实验:这是体现你算法价值的关键。至少对比:
- 不同初始解策略的效果:对比最近邻法、节约算法等对最终结果的影响。
- 不同优化算法的效果:对比纯局部搜索、模拟退火、你的VNS算法在求解质量和时间上的差异。
- 关键参数敏感性分析:例如,分析电池容量大小对总成本、车辆使用数的影响;或者分析时间窗严格程度对路径规划可行性的影响。这能体现你对问题深度的理解。
- 结果可视化:
- 表格:清晰列出不同算法/参数下的目标函数值、计算时间、车辆使用数等关键指标。
- 图表:除了前面提到的路径图、收敛曲线,还可以绘制柱状图对比不同方案的成本构成(固定成本、运输成本、惩罚成本各占多少),绘制折线图展示参数敏感性。
5.4 模型检验与鲁棒性讨论
这是很多论文的薄弱环节,但却是加分项。你需要证明你的模型和方案是可靠的。
- 极端情况测试:设计一些极端场景,比如某个客户的需求量突然极大,或者某个路段临时封闭(距离设为无穷大),看你的算法能否给出合理的应对方案(如启用备用车辆、路径重规划)。
- 数据扰动分析:对客户需求量、服务时间等参数加入微小随机扰动(例如±10%),重新运行算法,观察结果的变化是否在可接受范围内。这体现了模型对数据误差的鲁棒性。
- 灵敏度分析报告:系统地汇报某个关键参数(如单位运输成本、时间窗惩罚系数)在一定范围内变动时,最优解如何变化。这能为决策者提供有价值的参考。
个人体会:写论文时,一定要站在评委的角度思考。评委时间有限,他们希望快速抓住你的核心工作。因此,摘要要精炼,突出问题、方法、主要结果和结论;图表要自明,不看正文也能懂个七八分;模型和算法部分要逻辑连贯,像讲故事一样层层递进。最后,检查全文的格式、编号、参考文献引用,这些细节决定了论文的专业第一印象。数学建模竞赛,既是智力的比拼,也是工程实现和学术表达的全面考验。从理解问题到编码实现,再到撰写报告,每一个环节的扎实程度,最终都会体现在你的论文里。