简介这是一份面向运筹学、最优化方法课程学习者与考研复习者的《单纯形及其对偶》课件配套讲解线性规划的标准型转化、解的概念与单纯形法迭代思路适合已具备线性代数基础、需要系统梳理求解流程的中级学习者。压缩包仅含1个pdf文件约1.23MB体量轻便可直接作为课堂讲义或自学提纲使用。内容从标准型定义出发给出目标函数与约束的向量、矩阵表示归纳max转min、不等式引入松弛变量与剩余变量、无约束变量拆分为两个非负变量之差三步化法并用具体例题完整演示转化过程随后讲解可行解、最优解、基与基变量、基本解、基本可行解及可行基的判定配合凸集、顶点、凸组合等几何概念与三条基本定理说明基可行解对应可行域顶点、若有解则必存在基可行解为最优解。目前已有173人学习适合作为单纯形法入门与对偶问题预习的参考材料。1. 从一张卷子上的 x₃ 无约束说起带过运筹学课的人大概都见过这种场面例题给了max f -x1 2x2 - 3x3约束里偏偏写着「x3 无约束」。学生盯着这行字发呆后面的单纯形表就没法往下画了。这门课真正的门槛不在算法本身而在「怎么把一个歪七扭八的实际问题摆成标准型」以及摆好之后怎么判断表面上是同一张表的两个问题其实共用一组最优性条件。《单纯形及其对偶》课件把我一直想理顺的两块内容压在了一起一块是标准型和基本可行解的几何含义一块是单纯形表怎么迭代、检验数怎么读。它面向的是要应试的学生、要带课的青年教师还有写优化求解器分支逻辑的工程同学。对后者来说真正值钱的不是手算步骤而是「无约束变量拆成两个非负变量」「≤ 补松弛、≥ 减剩余」这套变形规则背后的矩阵结构因为大量真实模型在建模阶段就卡在这一步。下面按「标准型怎么摆 → 解和顶点的对应关系 → 单纯形表怎么转」推下去最后收在灵敏度分析和对偶检验这两个容易翻车的地方。2. 标准型的代数变形与松弛变量价值系数2.1 为什么非要扭成 min CX / AX b / X ≥ 0课件的处理顺序是先从一般形式推标准型再讲解的几何意义这个顺序值得跟。它先把四个自由度收敛成一个模板目标一律求 min约束一律等式变量一律非负。收成模板之后所有后续讨论都能用矩阵语言统一表达不用在符号上反复切换。标准型写成min f CX s.t. AX b X ≥ 0其中X (x1, x2, …, xn)^T是决策变量向量C (c1, c2, …, cn)是价值系数行向量b (b1, b2, …, bm)^T是资源向量A是 m×n 系数矩阵。课件里还给了列向量的写法p_j (a_1j, a_2j, …, a_mj)^T约束可以读成Σ p_j x_j b。这个列向量视角在后面挑基变量也就是挑 B 的列时特别顺手。有三个前提要记住n m且Rank(A) m也就是 A 行满秩。行满秩意味着约束之间不存在冗余否则 A 里某些行可以由其他行线性表出AX b可能无解。工程上如果模型从 CAD 或台账自动导出很容易带出冗余行求解器会在预处理阶段报「problem is primal infeasible」或者干脆给一个矛盾的判定。2.2 三类变形max→min、不等式补变量、无约束变量拆分课件把变形分成三类这三类恰好对应实际建模中最常动刀的地方。第一类max 转 min。利用max f cx -min(-f)令f -f得到min f -cx。注意这里动的是目标函数的整体符号不是把每个c_j逐个取反——两者等价但批量取反更容易在代码里写错下标。第二类不等式转等式。≤在左端加一个松弛变量≥在左端减一个剩余变量两者都非负并且在目标函数中的价值系数为 0。价值系数取 0 不是随手规定的这两个变量是为了把约束「撑平」而引进来的辅助量它们不消费资源目标所以对目标函数不产生贡献。第三类无约束变量拆分。令x_k x_k - x_k其中x_k ≥ 0, x_k ≥ 0代入即可。拆完以后变量数增加 1这是代价换来的是所有变量都能塞进X ≥ 0的统一框架。2.3 [eg.7] 完整走一遍并对齐变量个数原始问题max f -x1 2x2 - 3x3 x1 x2 x3 ≤ 7 ① x1 - x2 x3 ≥ 2 ② -3x1 x2 2x3 5 ③ x1, x2 ≥ 0, x3 无约束按规则逐个改令x3 x3 - x3①加松弛变量x4②减剩余变量x5目标整体取负。min f x1 - 2x2 3(x3 - x3) 0x4 0x5 x1 x2 (x3 - x3) x4 7 x1 - x2 (x3 - x3) - x5 2 -3x1 x2 2(x3 - x3) 5 x1, x2, x3, x3, x4, x5 ≥ 0改成矩阵形式可以拿一段代码核一下维度和系数避免手抄错号import numpy as np # 列顺序: x1, x2, x3, x3, x4, x5 A_std np.array([ [ 1, 1, 1, -1, 1, 0], # ① 加松弛 x4 [ 1, -1, 1, -1, 0, -1], # ② 减剩余 x5 [-3, 1, 2, -2, 0, 0], # ③ 原等式x3 拆成两列 ], dtypefloat) b_std np.array([7.0, 2.0, 5.0]) C_std np.array([1.0, -2.0, 3.0, -3.0, 0.0, 0.0]) print(A 形状:, A_std.shape, b 长度:, b_std.shape[0]) print(行满秩:, np.linalg.matrix_rank(A_std) A_std.shape[0]) print(x3 列与拆分列和的对应:, np.allclose(A_std[:, 2] * 1 A_std[:, 3] * 1, [1, 1, 2]))这段代码做三件事。A_std按列顺序排列每一行对应一条变形后的约束print三行分别检查矩阵形状是否与b_std的行数一致、A 是否行满秩、以及x3与x3两列相加后能否还原成原来的 x3 系数列[1, 1, 2]。最后一项是有意义的自检拆分变量时符号最容易写错负号一丢整个模型的最优值就变样了。常见误用是在≥约束上「加」而不是「减」剩余变量。判断依据很直白把x1 - x2 x3 ≥ 2挪成等式需要左端先补一个非负量把不等式「收紧」所以是减x5写成加x5会让可行域凭空扩大。3. 基、基本解与可行域顶点的对应3.1 从满秩子矩阵挑基课件给出的定义链条是取 A 中 m×m 子矩阵 B若Rank(B) m则 B 是一个基B 的列p_1, …, p_m对应的变量叫基变量其余叫非基变量。因为 A 行满秩、n m基的个数最多是C(n, m)这也是基本解个数的上界。这个定义解释了一个常被忽略的事实基不是唯一的而且不同的基对应当前可行域的不同顶点。手算时选哪组基会影响迭代步数但不影响最终最优值。工程实现里选基策略比如选最大检验数入基还是 Bland 规则直接决定会不会死循环——理论上存在退化基导致循环的例子Bland 规则用最小编号打破平局能保证有限终止。3.2 令非基变量为 0 得到基本解固定一个基 B 之后把非基变量全部置 0方程AX b退化成B·X_B b。因为 B 可逆X_B B⁻¹b得到唯一解X⁰ (X_B, 0, …, 0)^T这就是基本解。基本解和基本可行解的差别只在于一个非负性检查概念定义是否要求 X ≥ 0几何位置基本解非基变量取 0 后由B⁻¹b解出否约束边界的交点可能在可行域外基本可行解基本解且满足非负约束是可行域顶点可行基对应基本可行解的那个 B是每个顶点至少对应一个可行基最优解使目标达到最小的可行解是某个或某些顶点课件举的 O(0,0)、Q1(4,0)、Q2(4,2)、Q3(2,3)、Q4(0,3) 都是基本可行解而 Q5(4,3) 是各边延长线的交点它满足等式约束但不在可行域内所以只是基本解。这个对比例子很有用它说明「解方程」和「求可行解」是两件事。3.3 凸集、顶点与基本定理可行域D {X | AX b, X ≥ 0}是凸集。证明思路不复杂取X^(1), X^(2) ∈ D令X αX^(1) (1-α)X^(2)0 α 1则AX αAX^(1) (1-α)AX^(2) b且X ≥ 0。凸集的意义在于任意两个可行解的连线还在可行域内所以局部最优等于全局最优。顶点定义为不能被可行域内任意两点的凸组合表示的点。关键结论是顶点所对应的解就是基本可行解。课件把三条基本定理摆出来——可行域是凸集、基可行解对应顶点、有解则必存在基可行解是最优解——这三条串起来就构成了「只查顶点就能找最优」的理论依据也是单纯形法只沿顶点跳的合法性来源。常见坑是把「可行解」「基本解」「基本可行解」当同义词用。写对偶分析或者做灵敏度分析时这三者的查询范围不同混用会直接导致结论错位。4. 单纯形表的迭代规则与检验数4.1 找初始基可行解的三条路课件列了三种确定初始基的方式实际用哪种取决于约束长什么样。第一种是松弛基。如果原问题全是≤约束化标准型后松弛变量对应的列正好构成单位矩阵直接拿它当初始基。例子min f x1 3x2 x1 2x2 ≤ 3 2x1 3x2 ≤ 4 x1, x2 ≥ 0化标准型后x3, x4是基变量令x1 x2 0得X⁰ (0, 0, 3, 4)^T。第二种是观察法。约束里已经存在单位矩阵结构时不必硬套松弛基直接读出哪几列能凑出 I 就行。[eg.9] 里选XB (x1, x4)^T令x2 x3 0得X⁰ (3, 0, 0, 4)^T。第三种是人工基。当 A 里找不到单位矩阵时[eg.10] 中A [[1,3,2],[2,1,1]]就是这种情况引入人工变量充当初始基变量。人工变量是权宜之计后续要用大 M 法或者两阶段法把它们从基里赶出去否则得到的「最优解」里带着伪变量不是原问题的解。import numpy as np def find_identity_basis(A, tol1e-9): 从 A 中找出一组能构成单位矩阵的列返回列索引找不到返回 None m, n A.shape used, basis set(), [] for i in range(m): for j in range(n): if j in used: continue col A[:, j] target np.zeros(m) target[i] 1.0 if np.allclose(col, target, atoltol): basis.append(j) used.add(j) break return basis if len(basis) m else None A np.array([[1.0, 3.0, 2.0], [2.0, 1.0, 1.0]]) print(能否直接取到单位基:, find_identity_basis(A))函数逐行扫描为每一行 i 找一列等于单位向量e_i且这一列不能被前面的基用过。返回None就说明必须引入人工变量。这种自动探测在写建模工具时很实用能提前告诉用户「这个模型需要两阶段法」而不是等求解器报错。4.2 检验数 σ_j 与三种判别结论设基变量为x_1, …, x_m非基变量为x_{m1}, …, x_n把非基变量表示代入目标函数后可以得到σ_j c_j - z_j, z_j Σ_{i1}^{m} c_i a_ij, f₀ Σ_{i1}^{m} c_i b_iσ_j就是非基变量x_j的检验数。判别规则有三条所有σ_j ≤ 0时当前基可行解最优若某个σ_k 0且该列所有a_ik ≤ 0问题无界无最优解若所有σ_j ≤ 0且某个非基变量的σ_k 0则该问题有无穷多最优解。这三条在写代码判断时要一起看不能只看第一条。检验数为 0 时如果直接宣告结束会漏掉「多重最优解」这个信息做灵敏度分析时容易给出一个看似唯一实则不唯一的结论。import numpy as np def check_optimal(sigma, A, basis_idx, tol1e-9): 根据检验数判断解的形态返回状态字符串 sigma np.asarray(sigma, dtypefloat) if np.all(sigma tol): if np.any(np.abs(sigma) tol): return 无穷多最优解 return 唯一最优解 k int(np.argmax(sigma)) # 最大检验数入基 col A[:, k] if np.all(col tol): # 入基列无正元素 return 无界解无最优解 return f继续迭代入基变量列索引 {k}函数先扫全部检验数。np.all(sigma tol)成立时还要再查是否存在绝对值接近 0 的检验数有就是多重最优没有才是唯一最优。否则取最大检验数对应的列作为入基候选再检查该列是否存在正元素——全是非正就说明沿这个方向走目标值无界。这里tol作为数值容差必须给浮点误差会让本该是 0 的检验数变成1e-16不加容差判断会误报成多重解。4.3 入基、出基与旋转运算入基规则若存在σ_j 0取max σ_j对应的变量x_k入基。出基规则用最小比值检验θ min{ b_i / a_ik | a_ik 0 } b_l / a_lkx_l出基a_lk是主元。比值检验的含义是沿入基方向增大x_k每增加一个单位第 i 个基变量减少a_ik个单位最先撞到 0 的那个基变量先退出。旋转运算消元把主元a_lk化成 1同列其他元素化成 0。课件给的列向量示意是p_k (0, …, 1, …, 0)^T主元位置为 1。实现上就是标准的初等行变换主元行除以a_lk其余各行减去「该行在第 k 列的值乘以主元行」。def pivot(tab, row, col, tol1e-12): 对单纯形表做一次旋转运算主元置 1同列其余元素置 0 tab tab.astype(float).copy() if abs(tab[row, col]) tol: raise ValueError(主元接近 0无法执行旋转运算) tab[row, :] / tab[row, col] for i in range(tab.shape[0]): if i ! row and abs(tab[i, col]) tol: tab[i, :] - tab[i, col] * tab[row, :] return tabrow是主元行、col是主元列。先做事除保证主元为 1再对每一行做减法把该行在主元列的值消成 0。tol防止在主元接近 0 时放大数值误差——如果最小比值检验选出多个并列的候选行退化情形主元可能非常小此时应当用 Bland 规则按变量编号打破平局而不是任选一行。5. 灵敏度分析从对偶检验数读出价值系数变化范围5.1 不用重解直接改一张表的系数灵敏度分析是这份课件里最容易被跳过的部分但它在工程上的价值反而最高。现实模型里c_j成本、售价和b_i资源上限天天在变每次变一点就重解一遍成本高且看不到「变化的临界点在哪」。灵敏度分析解决的问题是给出一个范围只要系数在这个范围内波动当前最优基结构不变最优解仍是同一个顶点。做法是把变化量记作Δc_j或Δb_i重新计算检验数σ_j要求其保持≤ 0。这一条件反解出一组不等式就得到允许变化的区间。5.2 用对偶变量解释影子价格对偶变量y_i的经济含义是第 i 种资源的影子价格资源上限b_i每增加一个单位目标函数值的改善量。当y_i 0时该资源是紧约束扩张它能改善目标当y_i 0时该资源有剩余扩张它没有意义。弱对偶理论给出的结论是对偶问题的任一可行解的目标值都是原问题最优值的界强对偶理论给出的结论是原问题有最优解时对偶问题也有最优解且两者目标值相等。这两个结论串起来就是「解一个就够」的根据。实际求解大型模型时选原问题还是对偶问题取决于哪一边的约束维数更小——约束数少的那一侧通常迭代更快。import numpy as np def c_range_when_basis_fixed(A_B_inv, C_B, j, a_j): 基不变时非基变量 x_j 的价值系数 c_j 允许下降/上升的范围 A_B_inv: 当前基的逆矩阵 B^-1 C_B: 当前基变量的价值系数 j: 被考察的非基变量下标 a_j: x_j 在原始 A 中的列向量 y C_B A_B_inv # 对偶变量影子价格 sigma_j (C_B A_B_inv a_j) - C_B.sum() * 0 # 占位实际用 c_j - z_j z_j C_B (A_B_inv a_j) return z_j # 当前基不变要求 c_j - z_j 0即 c_j z_jA_B_inv a_j算出x_j在当前基下的表示系数z_j C_B (A_B_inv a_j)是它的机会成本。要求σ_j c_j - z_j ≤ 0等价于c_j ≤ z_j。所以对任何非基变量c_j的上升空间是z_j - c_j一旦超过这个量这个变量就有资格入基原最优基被打破。工程上把这个差值列成表就是一份「成本波动告警清单」。5.3 用残差和秩检查验证结果做完灵敏度分析后别只看打印出来的范围跑两个验证一是把变化后的系数代回检查基变量解的可行性B⁻¹b ≥ 0二是检查A_B是否仍然满秩。退化情形下基矩阵可能接近奇异B⁻¹数值放大区间算出来会宽得不合理。def verify_basis(A_B, b, tol1e-8): 检查给定基在当前 b 下是否仍然可行且满秩 if np.linalg.matrix_rank(A_B) A_B.shape[0]: return 基矩阵降秩基已失效 x_B np.linalg.solve(A_B, b) if np.any(x_B -tol): return 基变量出现负值解不可行需重新迭代 return f基仍可行基变量取值: {np.round(x_B, 4)}这个函数回答的是「当前基还站不站得住」。matrix_rank判满秩先排除结构问题np.linalg.solve解出的x_B再判可行性。返回负值说明某个基变量在b变化后跌破了 0此时要触发一次对偶单纯形迭代而不是从头重解——对偶单纯形正是处理「基保持对偶可行、原问题变得不可行」这类场景的它能用更少的步数把解拉回可行域。5.4 三个高频翻车点第一把「基不变」和「最优解不变」混为一谈。灵敏度分析算出的区间保证的是基结构不变在这个区间内最优解本身也会随系数连续变化只是它始终是同一个顶点。第二忽略退化基。退化基本可行解中某些基变量为 0最小比值检验会出现并列θ 取 0迭代后目标函数值不下降。不引入 Bland 规则可能循环实践中直接对 θ0 的候选行取最小下标即可。第三把对偶问题的最优解当成原问题的解。两者目标值相等但变量含义完全不同对偶变量是影子价格直接往原模型里回代会得到量纲对不上的结果。课件把单纯形和对偶压在一份材料里是有道理的单纯形给出一条能落地的迭代路径对偶给出这条路径背后的界和临界条件。手推几遍 [eg.7] 到 [eg.10]再把检验数和对偶变量用代码算一次比背判别规则快得多。本文还有配套的精品资源点击获取
从标准型到对偶：单纯形法、检验数与灵敏度分析核心讲解
简介这是一份面向运筹学、最优化方法课程学习者与考研复习者的《单纯形及其对偶》课件配套讲解线性规划的标准型转化、解的概念与单纯形法迭代思路适合已具备线性代数基础、需要系统梳理求解流程的中级学习者。压缩包仅含1个pdf文件约1.23MB体量轻便可直接作为课堂讲义或自学提纲使用。内容从标准型定义出发给出目标函数与约束的向量、矩阵表示归纳max转min、不等式引入松弛变量与剩余变量、无约束变量拆分为两个非负变量之差三步化法并用具体例题完整演示转化过程随后讲解可行解、最优解、基与基变量、基本解、基本可行解及可行基的判定配合凸集、顶点、凸组合等几何概念与三条基本定理说明基可行解对应可行域顶点、若有解则必存在基可行解为最优解。目前已有173人学习适合作为单纯形法入门与对偶问题预习的参考材料。1. 从一张卷子上的 x₃ 无约束说起带过运筹学课的人大概都见过这种场面例题给了max f -x1 2x2 - 3x3约束里偏偏写着「x3 无约束」。学生盯着这行字发呆后面的单纯形表就没法往下画了。这门课真正的门槛不在算法本身而在「怎么把一个歪七扭八的实际问题摆成标准型」以及摆好之后怎么判断表面上是同一张表的两个问题其实共用一组最优性条件。《单纯形及其对偶》课件把我一直想理顺的两块内容压在了一起一块是标准型和基本可行解的几何含义一块是单纯形表怎么迭代、检验数怎么读。它面向的是要应试的学生、要带课的青年教师还有写优化求解器分支逻辑的工程同学。对后者来说真正值钱的不是手算步骤而是「无约束变量拆成两个非负变量」「≤ 补松弛、≥ 减剩余」这套变形规则背后的矩阵结构因为大量真实模型在建模阶段就卡在这一步。下面按「标准型怎么摆 → 解和顶点的对应关系 → 单纯形表怎么转」推下去最后收在灵敏度分析和对偶检验这两个容易翻车的地方。2. 标准型的代数变形与松弛变量价值系数2.1 为什么非要扭成 min CX / AX b / X ≥ 0课件的处理顺序是先从一般形式推标准型再讲解的几何意义这个顺序值得跟。它先把四个自由度收敛成一个模板目标一律求 min约束一律等式变量一律非负。收成模板之后所有后续讨论都能用矩阵语言统一表达不用在符号上反复切换。标准型写成min f CX s.t. AX b X ≥ 0其中X (x1, x2, …, xn)^T是决策变量向量C (c1, c2, …, cn)是价值系数行向量b (b1, b2, …, bm)^T是资源向量A是 m×n 系数矩阵。课件里还给了列向量的写法p_j (a_1j, a_2j, …, a_mj)^T约束可以读成Σ p_j x_j b。这个列向量视角在后面挑基变量也就是挑 B 的列时特别顺手。有三个前提要记住n m且Rank(A) m也就是 A 行满秩。行满秩意味着约束之间不存在冗余否则 A 里某些行可以由其他行线性表出AX b可能无解。工程上如果模型从 CAD 或台账自动导出很容易带出冗余行求解器会在预处理阶段报「problem is primal infeasible」或者干脆给一个矛盾的判定。2.2 三类变形max→min、不等式补变量、无约束变量拆分课件把变形分成三类这三类恰好对应实际建模中最常动刀的地方。第一类max 转 min。利用max f cx -min(-f)令f -f得到min f -cx。注意这里动的是目标函数的整体符号不是把每个c_j逐个取反——两者等价但批量取反更容易在代码里写错下标。第二类不等式转等式。≤在左端加一个松弛变量≥在左端减一个剩余变量两者都非负并且在目标函数中的价值系数为 0。价值系数取 0 不是随手规定的这两个变量是为了把约束「撑平」而引进来的辅助量它们不消费资源目标所以对目标函数不产生贡献。第三类无约束变量拆分。令x_k x_k - x_k其中x_k ≥ 0, x_k ≥ 0代入即可。拆完以后变量数增加 1这是代价换来的是所有变量都能塞进X ≥ 0的统一框架。2.3 [eg.7] 完整走一遍并对齐变量个数原始问题max f -x1 2x2 - 3x3 x1 x2 x3 ≤ 7 ① x1 - x2 x3 ≥ 2 ② -3x1 x2 2x3 5 ③ x1, x2 ≥ 0, x3 无约束按规则逐个改令x3 x3 - x3①加松弛变量x4②减剩余变量x5目标整体取负。min f x1 - 2x2 3(x3 - x3) 0x4 0x5 x1 x2 (x3 - x3) x4 7 x1 - x2 (x3 - x3) - x5 2 -3x1 x2 2(x3 - x3) 5 x1, x2, x3, x3, x4, x5 ≥ 0改成矩阵形式可以拿一段代码核一下维度和系数避免手抄错号import numpy as np # 列顺序: x1, x2, x3, x3, x4, x5 A_std np.array([ [ 1, 1, 1, -1, 1, 0], # ① 加松弛 x4 [ 1, -1, 1, -1, 0, -1], # ② 减剩余 x5 [-3, 1, 2, -2, 0, 0], # ③ 原等式x3 拆成两列 ], dtypefloat) b_std np.array([7.0, 2.0, 5.0]) C_std np.array([1.0, -2.0, 3.0, -3.0, 0.0, 0.0]) print(A 形状:, A_std.shape, b 长度:, b_std.shape[0]) print(行满秩:, np.linalg.matrix_rank(A_std) A_std.shape[0]) print(x3 列与拆分列和的对应:, np.allclose(A_std[:, 2] * 1 A_std[:, 3] * 1, [1, 1, 2]))这段代码做三件事。A_std按列顺序排列每一行对应一条变形后的约束print三行分别检查矩阵形状是否与b_std的行数一致、A 是否行满秩、以及x3与x3两列相加后能否还原成原来的 x3 系数列[1, 1, 2]。最后一项是有意义的自检拆分变量时符号最容易写错负号一丢整个模型的最优值就变样了。常见误用是在≥约束上「加」而不是「减」剩余变量。判断依据很直白把x1 - x2 x3 ≥ 2挪成等式需要左端先补一个非负量把不等式「收紧」所以是减x5写成加x5会让可行域凭空扩大。3. 基、基本解与可行域顶点的对应3.1 从满秩子矩阵挑基课件给出的定义链条是取 A 中 m×m 子矩阵 B若Rank(B) m则 B 是一个基B 的列p_1, …, p_m对应的变量叫基变量其余叫非基变量。因为 A 行满秩、n m基的个数最多是C(n, m)这也是基本解个数的上界。这个定义解释了一个常被忽略的事实基不是唯一的而且不同的基对应当前可行域的不同顶点。手算时选哪组基会影响迭代步数但不影响最终最优值。工程实现里选基策略比如选最大检验数入基还是 Bland 规则直接决定会不会死循环——理论上存在退化基导致循环的例子Bland 规则用最小编号打破平局能保证有限终止。3.2 令非基变量为 0 得到基本解固定一个基 B 之后把非基变量全部置 0方程AX b退化成B·X_B b。因为 B 可逆X_B B⁻¹b得到唯一解X⁰ (X_B, 0, …, 0)^T这就是基本解。基本解和基本可行解的差别只在于一个非负性检查概念定义是否要求 X ≥ 0几何位置基本解非基变量取 0 后由B⁻¹b解出否约束边界的交点可能在可行域外基本可行解基本解且满足非负约束是可行域顶点可行基对应基本可行解的那个 B是每个顶点至少对应一个可行基最优解使目标达到最小的可行解是某个或某些顶点课件举的 O(0,0)、Q1(4,0)、Q2(4,2)、Q3(2,3)、Q4(0,3) 都是基本可行解而 Q5(4,3) 是各边延长线的交点它满足等式约束但不在可行域内所以只是基本解。这个对比例子很有用它说明「解方程」和「求可行解」是两件事。3.3 凸集、顶点与基本定理可行域D {X | AX b, X ≥ 0}是凸集。证明思路不复杂取X^(1), X^(2) ∈ D令X αX^(1) (1-α)X^(2)0 α 1则AX αAX^(1) (1-α)AX^(2) b且X ≥ 0。凸集的意义在于任意两个可行解的连线还在可行域内所以局部最优等于全局最优。顶点定义为不能被可行域内任意两点的凸组合表示的点。关键结论是顶点所对应的解就是基本可行解。课件把三条基本定理摆出来——可行域是凸集、基可行解对应顶点、有解则必存在基可行解是最优解——这三条串起来就构成了「只查顶点就能找最优」的理论依据也是单纯形法只沿顶点跳的合法性来源。常见坑是把「可行解」「基本解」「基本可行解」当同义词用。写对偶分析或者做灵敏度分析时这三者的查询范围不同混用会直接导致结论错位。4. 单纯形表的迭代规则与检验数4.1 找初始基可行解的三条路课件列了三种确定初始基的方式实际用哪种取决于约束长什么样。第一种是松弛基。如果原问题全是≤约束化标准型后松弛变量对应的列正好构成单位矩阵直接拿它当初始基。例子min f x1 3x2 x1 2x2 ≤ 3 2x1 3x2 ≤ 4 x1, x2 ≥ 0化标准型后x3, x4是基变量令x1 x2 0得X⁰ (0, 0, 3, 4)^T。第二种是观察法。约束里已经存在单位矩阵结构时不必硬套松弛基直接读出哪几列能凑出 I 就行。[eg.9] 里选XB (x1, x4)^T令x2 x3 0得X⁰ (3, 0, 0, 4)^T。第三种是人工基。当 A 里找不到单位矩阵时[eg.10] 中A [[1,3,2],[2,1,1]]就是这种情况引入人工变量充当初始基变量。人工变量是权宜之计后续要用大 M 法或者两阶段法把它们从基里赶出去否则得到的「最优解」里带着伪变量不是原问题的解。import numpy as np def find_identity_basis(A, tol1e-9): 从 A 中找出一组能构成单位矩阵的列返回列索引找不到返回 None m, n A.shape used, basis set(), [] for i in range(m): for j in range(n): if j in used: continue col A[:, j] target np.zeros(m) target[i] 1.0 if np.allclose(col, target, atoltol): basis.append(j) used.add(j) break return basis if len(basis) m else None A np.array([[1.0, 3.0, 2.0], [2.0, 1.0, 1.0]]) print(能否直接取到单位基:, find_identity_basis(A))函数逐行扫描为每一行 i 找一列等于单位向量e_i且这一列不能被前面的基用过。返回None就说明必须引入人工变量。这种自动探测在写建模工具时很实用能提前告诉用户「这个模型需要两阶段法」而不是等求解器报错。4.2 检验数 σ_j 与三种判别结论设基变量为x_1, …, x_m非基变量为x_{m1}, …, x_n把非基变量表示代入目标函数后可以得到σ_j c_j - z_j, z_j Σ_{i1}^{m} c_i a_ij, f₀ Σ_{i1}^{m} c_i b_iσ_j就是非基变量x_j的检验数。判别规则有三条所有σ_j ≤ 0时当前基可行解最优若某个σ_k 0且该列所有a_ik ≤ 0问题无界无最优解若所有σ_j ≤ 0且某个非基变量的σ_k 0则该问题有无穷多最优解。这三条在写代码判断时要一起看不能只看第一条。检验数为 0 时如果直接宣告结束会漏掉「多重最优解」这个信息做灵敏度分析时容易给出一个看似唯一实则不唯一的结论。import numpy as np def check_optimal(sigma, A, basis_idx, tol1e-9): 根据检验数判断解的形态返回状态字符串 sigma np.asarray(sigma, dtypefloat) if np.all(sigma tol): if np.any(np.abs(sigma) tol): return 无穷多最优解 return 唯一最优解 k int(np.argmax(sigma)) # 最大检验数入基 col A[:, k] if np.all(col tol): # 入基列无正元素 return 无界解无最优解 return f继续迭代入基变量列索引 {k}函数先扫全部检验数。np.all(sigma tol)成立时还要再查是否存在绝对值接近 0 的检验数有就是多重最优没有才是唯一最优。否则取最大检验数对应的列作为入基候选再检查该列是否存在正元素——全是非正就说明沿这个方向走目标值无界。这里tol作为数值容差必须给浮点误差会让本该是 0 的检验数变成1e-16不加容差判断会误报成多重解。4.3 入基、出基与旋转运算入基规则若存在σ_j 0取max σ_j对应的变量x_k入基。出基规则用最小比值检验θ min{ b_i / a_ik | a_ik 0 } b_l / a_lkx_l出基a_lk是主元。比值检验的含义是沿入基方向增大x_k每增加一个单位第 i 个基变量减少a_ik个单位最先撞到 0 的那个基变量先退出。旋转运算消元把主元a_lk化成 1同列其他元素化成 0。课件给的列向量示意是p_k (0, …, 1, …, 0)^T主元位置为 1。实现上就是标准的初等行变换主元行除以a_lk其余各行减去「该行在第 k 列的值乘以主元行」。def pivot(tab, row, col, tol1e-12): 对单纯形表做一次旋转运算主元置 1同列其余元素置 0 tab tab.astype(float).copy() if abs(tab[row, col]) tol: raise ValueError(主元接近 0无法执行旋转运算) tab[row, :] / tab[row, col] for i in range(tab.shape[0]): if i ! row and abs(tab[i, col]) tol: tab[i, :] - tab[i, col] * tab[row, :] return tabrow是主元行、col是主元列。先做事除保证主元为 1再对每一行做减法把该行在主元列的值消成 0。tol防止在主元接近 0 时放大数值误差——如果最小比值检验选出多个并列的候选行退化情形主元可能非常小此时应当用 Bland 规则按变量编号打破平局而不是任选一行。5. 灵敏度分析从对偶检验数读出价值系数变化范围5.1 不用重解直接改一张表的系数灵敏度分析是这份课件里最容易被跳过的部分但它在工程上的价值反而最高。现实模型里c_j成本、售价和b_i资源上限天天在变每次变一点就重解一遍成本高且看不到「变化的临界点在哪」。灵敏度分析解决的问题是给出一个范围只要系数在这个范围内波动当前最优基结构不变最优解仍是同一个顶点。做法是把变化量记作Δc_j或Δb_i重新计算检验数σ_j要求其保持≤ 0。这一条件反解出一组不等式就得到允许变化的区间。5.2 用对偶变量解释影子价格对偶变量y_i的经济含义是第 i 种资源的影子价格资源上限b_i每增加一个单位目标函数值的改善量。当y_i 0时该资源是紧约束扩张它能改善目标当y_i 0时该资源有剩余扩张它没有意义。弱对偶理论给出的结论是对偶问题的任一可行解的目标值都是原问题最优值的界强对偶理论给出的结论是原问题有最优解时对偶问题也有最优解且两者目标值相等。这两个结论串起来就是「解一个就够」的根据。实际求解大型模型时选原问题还是对偶问题取决于哪一边的约束维数更小——约束数少的那一侧通常迭代更快。import numpy as np def c_range_when_basis_fixed(A_B_inv, C_B, j, a_j): 基不变时非基变量 x_j 的价值系数 c_j 允许下降/上升的范围 A_B_inv: 当前基的逆矩阵 B^-1 C_B: 当前基变量的价值系数 j: 被考察的非基变量下标 a_j: x_j 在原始 A 中的列向量 y C_B A_B_inv # 对偶变量影子价格 sigma_j (C_B A_B_inv a_j) - C_B.sum() * 0 # 占位实际用 c_j - z_j z_j C_B (A_B_inv a_j) return z_j # 当前基不变要求 c_j - z_j 0即 c_j z_jA_B_inv a_j算出x_j在当前基下的表示系数z_j C_B (A_B_inv a_j)是它的机会成本。要求σ_j c_j - z_j ≤ 0等价于c_j ≤ z_j。所以对任何非基变量c_j的上升空间是z_j - c_j一旦超过这个量这个变量就有资格入基原最优基被打破。工程上把这个差值列成表就是一份「成本波动告警清单」。5.3 用残差和秩检查验证结果做完灵敏度分析后别只看打印出来的范围跑两个验证一是把变化后的系数代回检查基变量解的可行性B⁻¹b ≥ 0二是检查A_B是否仍然满秩。退化情形下基矩阵可能接近奇异B⁻¹数值放大区间算出来会宽得不合理。def verify_basis(A_B, b, tol1e-8): 检查给定基在当前 b 下是否仍然可行且满秩 if np.linalg.matrix_rank(A_B) A_B.shape[0]: return 基矩阵降秩基已失效 x_B np.linalg.solve(A_B, b) if np.any(x_B -tol): return 基变量出现负值解不可行需重新迭代 return f基仍可行基变量取值: {np.round(x_B, 4)}这个函数回答的是「当前基还站不站得住」。matrix_rank判满秩先排除结构问题np.linalg.solve解出的x_B再判可行性。返回负值说明某个基变量在b变化后跌破了 0此时要触发一次对偶单纯形迭代而不是从头重解——对偶单纯形正是处理「基保持对偶可行、原问题变得不可行」这类场景的它能用更少的步数把解拉回可行域。5.4 三个高频翻车点第一把「基不变」和「最优解不变」混为一谈。灵敏度分析算出的区间保证的是基结构不变在这个区间内最优解本身也会随系数连续变化只是它始终是同一个顶点。第二忽略退化基。退化基本可行解中某些基变量为 0最小比值检验会出现并列θ 取 0迭代后目标函数值不下降。不引入 Bland 规则可能循环实践中直接对 θ0 的候选行取最小下标即可。第三把对偶问题的最优解当成原问题的解。两者目标值相等但变量含义完全不同对偶变量是影子价格直接往原模型里回代会得到量纲对不上的结果。课件把单纯形和对偶压在一份材料里是有道理的单纯形给出一条能落地的迭代路径对偶给出这条路径背后的界和临界条件。手推几遍 [eg.7] 到 [eg.10]再把检验数和对偶变量用代码算一次比背判别规则快得多。本文还有配套的精品资源点击获取