简介这份数学建模论文围绕公园内部道路优化设计问题面向参加数学建模竞赛的学生及算法爱好者。论文以最小化道路总长度为目标同时满足任意两入口最短路不超过直线距离1.4倍完整呈现Prim算法构建最小生成树、Floyd算法求解最短路径、Fermat点定位交叉点以及斯坦纳树优化等关键方法并按问题一、二、三逐步推进给定交叉点、自由设计交叉点、存在矩形湖障碍三种情形均有详细建模与求解结果。资源为1个doc格式文档压缩包仅495KB下载使用方便目前已有330人学习。文中还包含模型假设、符号说明、数据分析、算法Matlab代码及最终路径与交叉口坐标适合用于赛题复现、算法演练或撰写建模论文时参考能够帮助读者理解图论与几何方法在路径优化中的实际应用。1. 公园内道路设计不是画图题这是一道带约束的离散优化第一次拿到「公园内道路设计」这个数学建模题目时很多人会习惯性把它当成平面几何作图题在公园图纸上拉几条直线、连几个入口看起来顺眼就算完成。真跑起来才发现这题的核心根本不是画图而是一道带约束的路径优化问题公园形状给定、入口位置给定你要修一个道路网络让任意两个入口之间沿着道路走的路程不超过两点空间直线距离的 1.414 倍同时让道路总长尽量短。这里的 1.414 倍是整道题的灵魂。去掉它问题退化为一个普通的加权最短网络加上它道路就可能在三角形内部出现折点、分叉甚至绕行。约束一旦激活答案不是「直连最短」而是「在可行域边界上找最优」。本文不写抽象推导直接给出一套可以照着复现的建模与求解路径坐标建模、Dijkstra 可行性校验、变步长搜索交叉点、常见避坑和交付前的验证方法。适合正在备赛、或者需要处理同类路径优化问题的读者照着跑就能出数换一组入口数据也能复用。2. 把题目翻译成路径约束模型坐标、Dijkstra 可行性与不设交叉口的下界建模的第一步不是画图而是把题面语言翻译成代码能处理的数据结构。公园是一个等边三角形入口分布在三条边上道路只能由直线段组成可以在公园内部任意位置增加交叉点。这样的描述听起来简单但一旦进入求解阶段你会发现坐标系的选取、边的组织方式、约束的判定函数直接决定后续所有代码能不能复用。2.1 建模前先定坐标入口、边界与道路网络的三种数据结构我一般先把问题放到一个统一的直角坐标系里。公园取边长 480 米的等边三角形三个顶点坐标固定为一个左下角、一个右下角、一个顶部顶点。这样三角形内部的判定、距离计算、画图都比较自然。以下是一套模拟题面的入口数据写成了参数化函数方便换题面时只改比例系数import math from itertools import combinations L 480.0 SQ3 math.sqrt(3.0) A (0.0, 0.0) B (L, 0.0) C (L / 2.0, L * SQ3 / 2.0) # 四个入口用“所在边 距起点比例”描述方便替换成你自己的题面数据 def entry_spec(): return { P1: (A, B, 80.0 / L), # AB 边距 A 点 80m P2: (A, B, 400.0 / L), # AB 边距 A 点 400m P3: (B, C, 80.0 / L), # BC 边距 B 点 80m P4: (A, C, 280.0 / L), # AC 边距 A 点 280m } def point_on_seg(p, q, ratio): return (p[0] (q[0] - p[0]) * ratio, p[1] (q[1] - p[1]) * ratio) ENTRANCES {k: point_on_seg(*v) for k, v in entry_spec().items()} def inside_triangle(p): x, y p return (0.0 y SQ3 * x) and (y SQ3 * (L - x))这段代码的核心逻辑是参数化。入口坐标不写死成数值而是用「在哪条边上、离端点多远」来描述换一套入口只需要改entry_spec()里的比例。inside_triangle判断一个点是否落在三角形内部后续搜索交叉点时每走一步都要用它把关。坐标系确定后整道题就变成了三个数据结构的问题入口坐标字典、边的端点对列表、以及每个端点的坐标字典。2.2 用 Dijkstra 判断「任意两入口是否满足 1.414 倍直线距离」道路网络的本质是加权无向图。每条道路是一条直线段权值就是它的欧氏距离交叉点、入口都是图里的节点。约束「任意两入口之间沿道路走的路程不超过直线距离的 1.414 倍」等价于对每一对入口在图里跑一次最短路然后和直线距离乘以 1.414 比较。Dijkstra 在这个场景下足够快因为图的节点数通常只有几个到十几个不需要上复杂的堆优化。下面的函数是整套代码的判定核心def build_adjacency(node_coords, node_pairs): adj {k: {} for k in node_coords} for u, v in node_pairs: w math.dist(node_coords[u], node_coords[v]) adj[u][v] w adj[v][u] w return adj def dijkstra(adj, start, target): dist {start: 0.0} unvisited set(adj.keys()) while unvisited: u min(unvisited, keylambda n: dist.get(n, math.inf)) if u target: return dist[u] unvisited.discard(u) for v, w in adj[u].items(): if dist.get(v, math.inf) dist[u] w: dist[v] dist[u] w return math.inf def feasible(node_coords, node_pairs, ratio1.414): adj build_adjacency(node_coords, node_pairs) for u, v in combinations(ENTRANCES, 2): path_len dijkstra(adj, u, v) direct_len math.dist(node_coords[u], node_coords[v]) if path_len ratio * direct_len: return False return True这个判定函数接受两个输入节点坐标字典和边的端点对列表。注意ratio1.414不是固定值题目给 1.4 或 1.5 时直接把参数传进去即可算法不用改。判定逻辑也很直白有一对入口不满足约束整个方案就判为不可行。现实中很多队伍会在这里犯一个错误只检查相邻入口而不检查所有入口对导致不连通或绕远路的方案混过关。所以combinations(ENTRANCES, 2)这一步是必须的六个入口对一个都不能少。2.3 不设交叉口的下界先跑一个「沿边修正路网」看失败模式不引入交叉点时最简单的方案是把三角形三条边全部修上再补一条入口之间的底边连接。这个方案虽然笨但能用来观察失败模式。把它代入约束判定perimeter_coords dict(ENTRANCES) perimeter_coords.update({A: A, B: B, C: C}) perimeter_edges [ (P1, A), (A, P4), (P4, C), (C, P3), (P3, B), (B, P2), (P1, P2) ] print(feasible(perimeter_coords, perimeter_edges, 1.414))这段网络把道路全部贴着公园边界铺设总长正好是三角形周长 1440 米。跑出来的结果大概率是False原因有两类一是某些入口对沿边界绕路太长超过 1.414 倍直线距离二是就算满足约束总长也远大于合理设计。这说明「不设交叉点、只靠边界道路」并不能解决问题必须在公园内部增加折点或交叉点。更严谨一点的预处理是计算四个入口的最小生成树长度作为任意可行方案的理论下界。如果某条候选方案总长比 MST 还短那一定是算错了——因为任何连通网络的总长都不可能小于最小生成树。这个下界配合约束判定能在早期过滤掉大量不靠谱方案也为后续搜索提供一个基准值。3. 在三角形内部找交叉点单点变步长搜索的 Python 实现第一阶段的判定函数搭好后第二阶段的问题就清晰了在三角形内部选一个或多个交叉点让道路总长最小同时保证所有入口对满足 1.414 倍约束。这个问题的难点在于约束不可导、目标函数又是距离和用解析法很麻烦常见的可靠做法是数值搜索。3.1 为什么交叉点能省路四点连接的斯坦纳直觉四个入口直接两两相连形成完全图时总长非常大只保留最小生成树时总长最小但不满足最短路约束。交叉点的作用是在两者之间找平衡。考虑一个简单情形如果四个入口连到一个公共交叉点 O道路总长就是四个入口到 O 的距离之和。O 选在几何中心附近时这个和通常远小于入口两两直连的总长但可能仍然无法满足任意入口对的最短路约束。于是问题的本质变成了在可行域内找那个能让目标函数最小的点。这个点未必是几何中心因为约束会把解推向特定方向。单交叉点的搜索可以直接在二维平面上做。目标函数是四个距离之和这是一个凸函数但约束区域不是凸的所以不能只做一次梯度下降。我一般先用粗网格扫描全图再在最优位置附近做变步长精修这是最不容易翻车的方案。3.2 单交叉点粗扫全图 局部精修的完整实现先写两个辅助函数一个计算星形拓扑的总长一个按当前交叉点位置构造网络并做可行性校验def total_len(o): return sum(math.dist(o, e) for e in ENTRANCES.values()) def star_check(o, ratio1.414): coords dict(ENTRANCES) coords[O] o pairs [(k, O) for k in ENTRANCES.keys()] return feasible(coords, pairs, ratio)然后是粗扫和精修的主循环step0 20.0 best_o None best_len math.inf # 第一轮以 20m 为步长在整个三角形内做网格扫描 for i in range(0, int(L / step0) 1): for j in range(0, int(L * SQ3 / step0) 1): o (i * step0, j * step0) if not inside_triangle(o): continue if not star_check(o): continue s total_len(o) if s best_len: best_len, best_o s, o # 第二轮以当前最优点为中心步长逐次减半做局部精修 step step0 x0, y0 best_o for _ in range(6): step / 2.0 for dx in (-step, 0.0, step): for dy in (-step, 0.0, step): cand (x0 dx, y0 dy) if not inside_triangle(cand): continue if not star_check(cand): continue s total_len(cand) if s best_len: best_len, best_o s, cand x0, y0 cand print(fbest point: ({best_o[0]:.2f}, {best_o[1]:.2f})) print(froad total: {best_len:.2f} m)粗扫步长取 20 米的逻辑是什么入口之间的距离普遍在百米量级道路总长对交叉点偏移 20 米通常只产生十几米的波动不会漏掉结构性最优位置。网格扫描完成后以最优点为中心做 6 轮精修每轮步长减半最终精度约 0.31 米对建模题完全够用。输出结果时我一般把交叉点坐标和总长同时打印出来方便后续画图验证。3.3 一个点不够时双交叉点的坐标下降扩展单交叉点在某些入口分布下会无解。我跑过几组数据发现当四个入口分成两组、每组内部距离很近但组间距离很远时一个公共点无法同时满足所有入口对的约束。这时候需要引入两个交叉点 O1、O2形成两条分支再加一条连接段。双交叉点的搜索代价比单点高一个数量级直接做二维网格的笛卡尔积会非常慢。常见做法是坐标下降先固定 O1优化 O2再固定 O2优化 O1交替两到三轮收敛速度很快。核心评估函数如下def two_point_cost(o1, o2, ratio1.414): coords dict(ENTRANCES) coords[O1] o1 coords[O2] o2 pairs [ (O1, P1), (O1, P2), (O2, P3), (O2, P4), (O1, O2) ] if not feasible(coords, pairs, ratio): return math.inf return sum(math.dist(coords[u], coords[v]) for u, v in pairs)这里把 P1、P2 分给 O1P3、P4 分给 O2再用 O1-O2 连接两条分支。你可以根据实际入口位置调整分组方式比如把相距最远的两个入口分到不同分支。交替优化的循环就是把上一节里的网格扫描封装成一个函数传入一个固定点和另一个待优化点的初始值返回优化后的点和目标值然后来回调两三次。双点搜索时有一个经验第一轮粗扫不需要全图 20 米步长可以用 40 米甚至 50 米先把 O2 粗定下来再在精修时同步收紧两个点。因为双点目标函数有多个局部极值粗扫太密反而容易陷入第一个局部极小还不自知。先把关键位置都趟一遍再慢慢收步长结果更稳。3.4 步长、收敛精度和性能怎么看这套搜索方案的参数看起来多实际就三组粗扫步长、精修轮数、收敛精度。粗扫步长决定全局探索的完整性精修轮数决定局部极值的收敛精度。我的习惯是粗扫步长取最短入口间距的十分之一到五分之一精修轮数取 5 到 6 轮因为第 6 轮步长已经到亚米级继续细化只会增加计算量对道路总长的影响不到 0.1 米没什么实际意义。性能方面单交叉点的网格扫描约几百个候选点每个点跑六次 Dijkstra总耗时在 1 秒以内。双交叉点的坐标下降会比单点慢十几倍但也在可接受范围。真正拖慢速度的是在feasible里频繁调用math.dist这个问题不大。如果入口数量增加到十几个建议把节点名改成整数索引Dijkstra 里用列表而不是字典性能能再快一个量级。4. 道路设计避坑五个最容易让结果翻车的建模错误这套流程本身不难难的是中间那些看上去不起眼、实则能让结果整体报废的细节。以下五条都是我实际踩过的坑每条按「现象 → 原因 → 解决」写出来比空洞的告诫有用得多。4.1 把「最短路约束」写成了「直线距离约束」现象代码跑出来的方案总长很漂亮画图一看道路网络也不乱但一检查就发现某两个入口之间沿道路走要绕一个大弯路程远大于约束要求方案根本不合法。原因建模时把题目里「任意入口间沿道路的最短路程不超过直线距离的 1.414 倍」误简化成了「交叉点到每个入口的直线距离不超过某值」。这两个条件完全不同前者要求的是图上的最短路径后者只是局部约束完全不保证入口之间的连通性达标。解决所有入口对必须用 Dijkstra 算一遍图上最短路再比较。我在代码里用feasible函数统一收口所有合法性判断任何新增方案都必须过这个函数不允许在目标函数里手写局部距离约束。4.2 交叉点跑出三角形边界检查函数形同虚设现象搜索出的最优交叉点坐标是负的或者 y 坐标远大于三角形高度。把点画出来才发现它落在公园外面但在代码里却被当成合法解输出了。原因只检查了目标函数值没有在每一轮搜索里调用inside_triangle。更隐蔽的情况是检查了边界但检查条件写错比如只判断 x 大于 0漏掉了三角形斜边的约束。解决把inside_triangle放进每一轮候选点筛选的第一步而不是在最后统一过滤。粗扫和精修的每个候选点都必须先过边界检查再进可行性判定。坐标系的符号也要留意三角形顶点在 x 轴上方时y 非负只是一个必要条件两条腰线的约束同样要写全。4.3 图连接关系漏边Dijkstra 算出一个虚假的可行解现象某条入口对明明只能绕很远的路径到达但程序却报「满足约束」打印出的最短路径长度比实际路线短了一截。原因构造边的列表时漏掉了交叉点与某个入口的连接线或者双交叉点方案里忘了加 O1-O2 之间的连接段。Dijkstra 只能使用给定的边漏边之后图可能是不连通的但dijkstra函数在找不到目标时会返回无穷大不会报错更麻烦的是漏边后存在另一条更长的绕路被当成了最短路。解决每次构造完node_pairs后先打印边的总数和端点集合核对入口、交叉点是否全部出现。我习惯在feasible里加一个断言检查所有入口节点都在adj中并且每个入口节点至少有一条边这样漏边问题会在第一时间暴露。4.4 粗扫步长过大最优点从网格缝里溜走现象粗扫结果看起来合理精修后总长只降了一点点但把目标函数画成等值线图发现最优点其实在粗扫网格线的对角线方向上离扫描点差了十几米。原因网格扫描只检查步长的整数倍位置如果最优点恰好落在四个扫描点围成的方格中心而且这个方格的四个角都不可行粗扫会直接漏掉整个区域。精修只在粗扫最优点的附近搜索一旦粗扫选错了区域后续精修怎么修都回不到真正的最优点。解决粗扫步长不要取得太大。20 米对百米尺度的入口布局是安全的但如果入口间距更大建议用最短入口间距的十分之一作为步长上限。另一个补救办法是粗扫时如果某方格内存在不可行点不要急着跳过整个区域可以在方格中心补测一个点很多最优解恰好就落在这种边界区域。4.5 拿到合法解就收手没有区分「可行」和「最优」现象搜索结果总长 780 米满足所有约束于是直接出图交稿。但把约束倍数从 1.414 放宽到 1.5发现总长能降到 690 米说明原来的解离约束边界太远算法根本没有真正逼近最优值。原因搜索只在「可行区域」内找最小目标值如果粗扫阶段没有一个点靠近约束边界精修阶段也会在远离边界的区域里打转。目标函数是距离和通常约束边界附近才有真正的极值点。解决交付前必须做边界敏感性测试。把约束倍数改成 1.3、1.5、1.6 分别跑一遍观察最优交叉点位置和总长的变化。如果 1.5 倍时的总长比 1.414 倍时显著下降说明当前解卡在可行区域内部而不是最优边界上需要缩小粗扫步长或调整搜索策略。5. 用活动约束与等值线圈验证你的解交付前的最后一道检查搜索结束不等于建模完成。最后一步是验证解的质量我一般从三个角度做找活动约束、画等值线图、做敏感性测试。这三件事加起来只需要几十行代码但能避免大部分低级错误。5.1 看比值曲线找出被 1.414 卡住的那对入口最优解通常会有一个或多个入口对恰好让路径长度除以直线距离等于 1.414这就是活动约束。把最优网络代入后打印所有入口对的这个比值最接近 1.414 的那一对就是设计瓶颈。如果没有任何一对接近 1.414说明约束没有激活很可能是约束倍数本身设置过松或者道路网络不是最优结构。def print_ratio(node_coords, node_pairs): adj build_adjacency(node_coords, node_pairs) for u, v in combinations(ENTRANCES, 2): path_len dijkstra(adj, u, v) direct_len math.dist(node_coords[u], node_coords[v]) print(f{u}-{v}: {path_len / direct_len:.3f})这个输出还能帮你判断该增加哪个交叉点。如果瓶颈出现在 P1-P4 这对入口而这两个入口恰好分属两条不同分支下一步优化就优先调整对应分支的结构而不是盲目加第三个点。5.2 把不可行区域画成阴影解的好坏一眼可见数值输出的说服力有限画图是最直观的验证方式。对单交叉点方案可以在三角形范围内遍历候选位置把不可行点标成 NaN可行点按总长着色形成一张带空洞的等值线图。最优解如果落在可行区域和不可行区域的交界带上说明约束确实在起作用如果最优解周围一大片全是可行区域那这个解很可能是「可行但不够优」。import numpy as np import matplotlib.pyplot as plt step 10.0 xs np.arange(0, 481, step) ys np.arange(0, 416, step) X, Y np.meshgrid(xs, ys) Z np.full_like(X, np.nan, dtypefloat) for i in range(X.shape[0]): for j in range(X.shape[1]): o (X[i, j], Y[i, j]) if not inside_triangle(o): continue if not star_check(o): continue Z[i, j] total_len(o) plt.contourf(X, Y, Z, levels24, cmapviridis) plt.colorbar(labelroad total (m)) plt.scatter([e[0] for e in ENTRANCES.values()], [e[1] for e in ENTRANCES.values()], cred, marker*, s120, labelentrance) plt.scatter([best_o[0]], [best_o[1]], cwhite, markero, s80, labelcross point) plt.legend() plt.gca().set_aspect(equal) plt.show()这块代码本身很短但价值很高。等值线图能直接暴露步长选取是否合理如果图上最优区域是细长的狭窄条带说明目标函数在某个方向上很平此时交叉点坐标精度要求就不高如果最优区域是一个尖点坐标差一两米总长就会明显变化交稿时需要把交叉点坐标保留到小数点后两位。5.3 敏感性测试把倍数从 1.3 调到 1.5观察最优交叉点怎么移动最后一道检查是敏感性测试方法很简单把star_check的ratio参数改成 1.3、1.5 分别跑一遍搜索对比最优点的移动距离。如果 1.3 倍时直接无解说明当前拓扑对约束非常敏感可能需要额外增加一个交叉点如果 1.5 倍时最优点大幅向三角形中心移动说明 1.414 倍约束在真实地牵引着解的形状当前结果可信如果 1.5 倍时最优点基本不动说明约束根本不是瓶颈目标函数在当前地形下已经收敛到一个比较稳定的极值点此时可以放心交数。有一次我把 1.5 倍时的结果和 1.414 倍时的结果对比发现总长差了近 60 米而交叉点坐标只偏移了 8 米。后来仔细查才发现是约束边界上一对入口的路程比值恰好卡在 1.415差一点就越界。从那以后我养成一个习惯跑完任何优化问题先打印活动约束再画等值线最后做敏感性测试三步缺一不可。这套流程也许不能帮你证明最优性但能让你交出去的每一个数都经得起追问。希望帮到你。本文还有配套的精品资源点击获取