简介围绕拟凸问题转化为凸问题可行性问题的二分算法这份文档面向优化理论、机器学习与运筹学方向的学习者与研究者帮助梳理拟凸函数与凸函数的关系、凸可行性转化思路以及二分求解流程。压缩包内仅含1个docx文档约200KB篇幅集中适合在课程复习、算法梳理或研究入门时快速查阅。目前已有560人学习下载具备一定的参考热度。文档从拟凸函数定义与凸问题特性讲起重点说明如何通过引入新变量和约束将拟凸问题转化为寻找满足阈值条件的凸可行性问题并给出二分算法的初始化范围、中间值计算、区间更新与终止条件等步骤。读者可据此理解转化逻辑、掌握算法执行细节并把它作为优化课程笔记或相关专题的辅助材料对关键判断条件也能形成较完整的认识。1. 一个比值目标让求解器直接无法处理上周有同学拿来一段代码目标是最小化“单位产出成本”写成(c x d) / (a x b)约束全是线性的看着人畜无害。CVXPY 一跑就报 DQCP error提示这个表达式既不是凸的也不是凹的。他以为是写法问题改了三版仍然不行。问题不在写法。这个比值不是凸函数但它满足一个更弱的性质对任意阈值 t集合{x | (cxd)/(axb) t, axb 0}是一条线性不等式刻画出来的凸集。这种“子水平集全是凸集”的函数叫拟凸函数对应的优化叫拟凸优化。它没法直接丢给锥规划求解器但可以换一个问法给定 t问“能不能做到目标不超过 t”这是一个标准的凸可行性问题把 t 从大到小反复问就是二分算法。下面从判据、代码、参数到踩坑一条线走完。2. 子水平集视角拟凸为什么能拆成一串凸可行性2.1 拟凸的定义里藏着什么可用的东西标准定义是对任意 x、y 和 θ ∈ [0,1]有 f(θx (1-θ)y) ≤ max(f(x), f(y))。这个式子读起来别扭等价说法才好用——f 的每个子水平集 S_α {x | f(x) ≤ α} 都是凸集。也就是说把函数图像横着切一刀切口下方那块区域投影回定义域永远是凸的。这条性质直接决定了一件事约束f0(x) ≤ t单独拿出来看是一个凸约束尽管 f0 本身不是凸函数。于是“最小化 f0”这个非凸问题被替换成一族凸可行性问题找 x 使得 f0(x) ≤ t, fi(x) ≤ 0 (i1..m), Ax bt 大可行域大容易可行t 小可行域收缩容易不可行。可行性关于 t 单调这正是二分能成立的唯一前提。有三个地方容易被混淆值得单独点出来。第一拟凸说的是“子水平集凸”不是“函数图像凸”比值函数1/xx0是凸的x/yy0不是凸的但拟凸两者不能混。第二凸函数一定拟凸反过来不成立所以拟凸是严格更宽的类不要看到“能用二分”就以为原问题松弛了。第三拟凸函数的局部极小不一定是全局极小严格拟凸才是所以不能用梯度下降那套直觉去理解二分的行为。还有一条边界要说清楚二分的可行性问题里硬约束fi(x) ≤ 0必须本身就是凸的。如果约束里也混了拟凸或者非凸的东西那这层拆解只是把非凸性从目标挪到了约束问题一点没变简单。2.2 常见拟凸目标的水平集凸化对照实际工程里能碰到的拟凸目标就那么几类先把它们的f0(x) ≤ t长什么样记住比推导定义有用得多。目标 f0(x)定义域要求f0(x) ≤ t 的等价凸约束约束类型(cᵀx d)/(aᵀx b)aᵀx b 0cᵀx d ≤ t(aᵀx b)线性‖Ax b‖₂/(cᵀx d)cᵀx d 0‖Ax b‖₂ ≤ t(cᵀx d)二阶锥maxᵢ (aᵢᵀx bᵢ)/(cᵢᵀx dᵢ)各分母 0对每个 i 列一条线性不等式线性多条‖x − p‖₂/‖x − q‖₂x ≠ q‖x − p‖₂ ≤ t‖x − q‖₂二阶锥xᵀPx/(cᵀx d)P ⪰ 0cᵀx d 0xᵀPx ≤ t(cᵀx d)凸二次 ≤ 仿射exp(aᵀx b)无aᵀx b ≤ log tt 0线性表的读法很关键右两列是固定 t 之后的约束t 只是一个常数所以每条都是凸的。第一行的比值是最常见的一类第三行的 max 形式对应“最坏情况比值”第五行的二次型分母版在信噪比、能量效率这类指标里反复出现。需要提醒的是表里第五行。xᵀPx ≤ t(cᵀx d)等价于xᵀPx − t cᵀx − t d ≤ 0左边是凸二次函数凸函数 ≤ 0 是个凸集所以成立。但如果 P 不定这个集合就不再凸二分立刻失效——这不是算法问题是问题本身不拟凸。2.3 二分能跑的前提一族随 t 单调的凸表示把f0(x) ≤ t写成凸约束时有两种写法要区分开。一种是直接代入型t 只出现在约束右端系数是关于 t 的线性项比如cᵀx d ≤ t(aᵀx b)。这种最省事t 增大时右端增大可行域单调扩张二分严格成立。另一种是引入辅助变量型写成∃z, φ(x, z) ≤ t其中 φ 关于 (x, z) 联合凸。这种写法要额外确认 φ 关于 t 单调——t 变大可行域必须变大。假如某个变换把 t 放到了分母上或者做了1/t之类的替换单调性就断了二分会在两个值之间来回跳最后给出一个既不是可行也不是不可行的结论。我一般会做一次极小的自检取两个阈值 t₁ t₂各解一次可行性确认 t₁ 可行时 t₂ 一定可行。这个检查只要两次求解比事后排查二分震荡便宜太多。再补一个选型上的判断。对线性分式规划这一类其实有精确解法Charnes-Cooper 变换把分式目标整体替换成线性目标一次求解拿到精确最优值不需要二分。二分是通用工具代价是只给到 ε 精度、需要 log 次求解有精确变换的结构就别用二分。选二分的真正理由是目标结构复杂到写不出精确变换比如 max 型、二次型分母、或者分母本身依赖多个变量。3. 最小可复现实现把可行性判定写成求解器能吃的形式3.1 可行性判定的三种写法别用 Minimize(0)判断“约束组是否可行”这件事在数值求解器上不是零成本的有三种写法稳定性差很多。第一种cp.Problem(cp.Minimize(0), cons).solve()然后看prob.status是不是optimal。这是最容易踩坑的写法。可行性问题没有目标求解器只保证在容差内满足约束返回optimal时可能有一条约束违反了 1e-6。二分做几十轮这个误差会累积成一棵歪掉的树最后 u 收敛到一个根本不可行的伪最优值。第二种Minimize(0)但配合optimal_inaccurate状态一起判再手工回代检查约束违反量。能用但等于自己重新实现一遍可行性判定。第三种最小化约束违反量。给目标那条约束加一个非负松弛 s最小化 s。s* ≈ 0 就可行s* 明显为正就不可行。好处是拿到了一个连续量——不只是“行/不行”还有“差多远”后面做插值加速、判断数值地板都靠它。import cvxpy as cp import numpy as np def sublevel_feasible(t, data, delta1e-6, tol1e-7): 判定 f0(x) t 是否可行。f0 (c^T x d) / (a^T x b)。 返回 (是否可行, 最小松弛量 s*, 解 x*)。 a, b, c, d, G, h, n data x cp.Variable(n) s cp.Variable(nonnegTrue) # 松弛量目标约束被违反了多少 cons [ G x h, # 原始硬约束必须是凸的 a x b delta, # 定义域分母必须严格为正 c x d t * (a x b) s, # 子水平集约束的松弛形式 ] prob cp.Problem(cp.Minimize(s), cons) prob.solve() if prob.status not in (optimal, optimal_inaccurate): return False, np.inf, None # 硬约束都不可行直接判否 return float(s.value) tol, float(s.value), x.value代码后三行是判定的核心。prob.status不是 optimal 时说明连硬约束Gx ≤ h都进不去这种情况没必要看松弛量。注意delta和tol是两个不同量纲的东西delta是给分母设的地板防止解跑到定义域外tol是判定可行性的门槛要跟数据量级挂钩。两个都在 1e-6 到 1e-7 量级只是巧合具体怎么定放到第 4 章说。另外实践中我经常把s的系数设成 1 之外的值比如对目标约束的每一行分配不同的权重让数值上更容易违反的那几行被优先修正。这个改动不影响可行性判定但会改善求解器在临界 t 附近的稳定性。3.2 二分主循环区间收缩与解的保存主循环本身只有十几行但有两个细节决定了结果对不对可行时 u 要用回代得到的真实目标值更新而不是用 t不可行时 l 更新到 t 而不是 t 的某个邻域别自作聪明加阻尼。def bisect_quasiconvex(f0, sublevel, lo, hi, eps1e-4, max_iter80): lo: 目标值的下界不可行侧hi: 上界可行侧需事先找到一个可行点。 best_x, best_u None, hi k 0 for k in range(1, max_iter 1): if best_u - lo eps: break t 0.5 * (lo best_u) ok, s, x sublevel(t) if ok: v f0(x) # 回代真实目标值通常 v t if v best_u: best_u, best_x v, x else: lo t # 不可行抬高下界 return best_x, best_u, kv f0(x)这一行不能省。松弛形式下可行只意味着f0(x) ≤ t tol真实值可能比 t 小一截。用 v 而不是 t 更新上界等于每次迭代顺手做了一次免费的目标改进收敛往往比理论迭代数少几次。反过来如果 v 比 best_u 还大数值噪声导致的if v best_u保证不会把上界往回推避免区间不单调。lo的初值必须是严格不可行的下界。常见做法是取 0很多比值目标非负、取解析下界、或者用一次松弛问题的对偶信息推出来。如果 lo 本身是可行的二分会在一个错误的区间上收缩最后返回的“最优值”偏大而且不会报错——这种错最难查。3.3 跑通一次线性分式规划的完整对拍给一份可以立刻跑的数据。目标最小化(cᵀx d)/(aᵀx b)约束是Gx ≤ h维度 n5、约束 8 条随机生成但保证存在内点。rng np.random.default_rng(7) n, m 5, 8 G rng.normal(size(m, n)) x_feas rng.normal(sizen) * 0.3 h G x_feas 1.0 # 保证 x_feas 严格可行 a rng.normal(sizen) b 1.5 - a x_feas # 保证 a^T x_feas b 1.5 0 c rng.normal(sizen) d -c x_feas 0.8 # 保证 f0(x_feas) 为正 data (a, b, c, d, G, h, n) f0 lambda x: (c x d) / (a x b) # 上界拿一个已知可行点算目标值 hi f0(x_feas) # 下界目标非负取 0 x_opt, u_opt, iters bisect_quasiconvex( f0, lambda t: sublevel_feasible(t, data), lo0.0, hihi, eps1e-5 ) print(f二分结果 u{u_opt:.8f}迭代 {iters} 次)同一份数据用 Charnes-Cooper 变换求精确解作为对照def lfp_exact(data): 线性分式规划的精确解令 y x/(a^T x b), s 1/(a^T x b) 0。 a, b, c, d, G, h, n data y cp.Variable(n) s cp.Variable(nonnegTrue) cons [G y h * s, a y b * s 1] # 归一化 缩放后的约束 prob cp.Problem(cp.Minimize(c y d * s), cons) prob.solve() if prob.status ! optimal or s.value 1e-9: return None, None # 无界或数值退化 return float(prob.value), y.value / s.value v_exact, x_exact lfp_exact(data) print(f精确解 {v_exact:.8f}与二分差 {abs(v_exact - u_opt):.2e})两个结果应该在 1e-5 量级内一致。如果不一致先别怀疑二分逻辑去看x_opt是否真的满足Gx ≤ h且aᵀx b 0——九成的问题出在解的可行性上而不是区间收缩上。这套对拍是我每次写新拟凸求解时的固定动作。原因很实际二分的错误通常不表现为崩溃而是表现为“收敛到一个略微偏大的值”没有精确解对照根本发现不了。哪怕只有线性分式这一种结构能对拍也值得写因为它能暴露松弛量容差、定义域地板、区间方向这三类问题。4. 参数怎么定区间、精度、迭代次数与数值地板4.1 上下界从哪来先解一次纯可行性二分的初始区间不是随便给的。hi必须是可达的也就是存在一个 x 既满足硬约束又满足定义域且f0(x) ≤ hi。最省事的构造是先解一次只有硬约束和定义域的可行性问题拿到任意可行点 x₀令hi f0(x₀)。def find_upper_bound(data, delta1e-6): a, b, c, d, G, h, n data x cp.Variable(n) prob cp.Problem(cp.Minimize(0), [G x h, a x b delta]) prob.solve() if prob.status not in (optimal, optimal_inaccurate): return None return (c x.value d) / (a x.value b)返回 None 意味着硬约束本身不可行这时候二分无从谈起应该先修模型而不是调参数。lo的取法有三种按可靠度排序解析下界比如比值目标在可行域上恒非负取 0 松弛问题的对偶下界 拍脑袋取hi - 10。第三种只在快速验证时用正式跑之前必须换掉。如果 lo 与真实最优值差距过大多出来的迭代次数是log2(差距/ε)差 1000 倍也就多 10 次求解代价可以接受但前提是 lo 确实不可行。有一个容易忽略的情况可行域无界导致f0的下确界取不到。这时二分仍然会收敛到一个值但那个值附近不存在最优解回代的 x 会随着 t 变化剧烈漂移。判断方法是观察相邻两次可行解的 x 范数如果 t 变化 1e-3 而 x 范数跳了两个数量级基本可以确定中了这个招。4.2 精度 ε 与迭代次数的换算收缩过程是标准的区间减半从 (lo, hi) 到宽度 ε 需要的迭代次数是ceil(log2((hi − lo)/ε))。这个数只跟区间宽度和精度有关跟问题维度、约束条数完全无关是二分最讨喜的地方。初始区间宽度ε 1e-3ε 1e-5ε 1e-7ε 1e-911017243010142027341e3202734401e630374450表的用法是反查如果你的单次可行性求解要 0.2 秒预算 10 秒那最多 50 次求解ε 就别定到 1e-9。很多场景下 1e-5 的相对精度已经远超下游需求硬压到 1e-9 只是把求解器的数值噪声放大成结果抖动。ε 应该按相对精度定义更稳妥eps_rel eps_abs / max(1, |hi|)。目标值量级是 1e4 时绝对精度 1e-5 相当于相对精度 1e-9这个要求超出了多数一阶/内点求解器在中等规模问题上的能力结果就是最后几轮二分在可行/不可行之间随机翻转。4.3 求解器容差与精度地板所有求解器判断“约束是否满足”都带容差内点法大致在 1e-8 到 1e-9 的相对量级一阶法更松。这意味着两件事。第一sublevel_feasible里的tol不能小于求解器容差乘以数据尺度。经验取值是tol 1e-6 * max(1, ‖c‖, ‖d‖, |t| * (‖a‖ |b|))把目标约束各参与项的量级算进去。第二二分存在一个精度地板当hi − lo小到和tol同量级时t 落在“数值上说不清可行还是不可行”的灰色地带二分的判断结果取决于求解器内部迭代路径同一份数据两次跑可能给出不同答案。这个现象很玄学但只要见过一次就忘不掉——你会看到 u 每次运行差个 1e-6而且怎么调随机种子都复现不了。实用做法是给 ε 定一个下限不要让二分进入灰色地带eps max(eps_target, 100 * tol_rel * scale)。如果业务确实需要 1e-9 的精度那就不是二分能不能跑的问题而是要把问题结构摸清楚、换精确变换或改写成单次凸问题。还有一个跟尺度有关的陷阱。Gx ≤ h里如果某几行系数是 1e5 量级另几行是 1e-2 量级求解器会优先满足大系数行小系数行的违反量可能达到 1e-3。这种数据在二分会表现为“明明可行却判不可行”而且只在特定 t 区间出现。预处理时把每行按行范数归一化能消掉大部分这类问题。def normalize_rows(G, h, eps1e-12): 按行范数归一化约束避免量级悬殊导致的判定偏差。 scale np.linalg.norm(G, axis1, keepdimsTrue) scale np.maximum(scale, eps) return G / scale, h / scale.squeeze()归一化只改数值不改可行域——Gx ≤ h和(G/s)x ≤ h/s是同一个集合s 0。这一步做完前面说的容差取值就可以用统一的常数不必逐行算尺度。5. 避坑与排查二分跑飞、判错可行性、踩到定义域外5.1 现象u 一路降到荒谬的低值回代目标却明显违反原因几乎总是拿Minimize(0)判可行性。可行性问题没有目标求解器返回 optimal 只代表约束在容差内近似满足一条约束违反 1e-7、跑 30 轮二分误差会沿着可行方向一路放大最后 u 收敛到一个根本不成立的阈值。解决办法是换 3.1 里的松弛写法并且每次都回代f0(x) t tol时即使s*报的是小值也按不可行处理。回代是最便宜的一致性检查多写一行能省半天。5.2 现象解出的 x 让分母变成负数u 掉到负无穷附近这是定义域丢失。cᵀx d ≤ t(aᵀx b)与f0(x) ≤ t等价前提是aᵀx b 0。分母为负时两边乘过去不等号方向要反等价性就断了。二分在收缩过程中一旦走到这一侧会看到“目标可以无限小”u 直接崩掉。解决办法是显式加aᵀx b ≥ deltadelta 取数据尺度的 1e-6 到 1e-8。delta 太小起不到保护作用太大会把可行域切掉一块、导致最优值偏大。判断 delta 是否合适跑完看最优解处的分母如果它恰好贴着 delta说明约束是紧的得重新选。另一条路是改用 Charnes-Cooper 归一化aᵀy bs 1从根上避掉这个问题。5.3 现象迭代次数比理论值多一倍最后 u 还是偏大上下界方向搞反了或者初始 hi 对应的点其实不可行。二分要求 hi 侧可达如果 hi 只是“看起来合理”的估计值每次可行判定都会失败l 一路往上爬直到和 hi 撞上。排查办法很简单单独验证 hi 对应的那个 x₀打印G x0 - h的最大值和a x0 b的值。只要有一个为正/为负就去修初始点。我习惯在bisect_quasiconvex开头加一行断言把这个问题挡在循环外面。5.4 现象二分在两个值之间来回震荡不收敛t 不是单调进入约束。某些拟凸表示会做变换比如把 t 放到分母、或者引入t²、1−t之类的项导致可行域随 t 非单调。这时候二分的前提被破坏区间收缩的方向变得不可靠。解决办法是回到原始定义把f0(x) ≤ t重写成 t 只出现在线性位置的形式。写不出来就说明这个目标不是拟凸或者需要换一组变量。震荡的另一个来源是 t 恰好落在精度地板附近表现是多次迭代 t 的变化量小于 1e-9同时判断结果随机——这属于 4.3 说的问题抬高 ε 即可。5.5 现象单独跑每次可行性都对串起来结果就对不上状态污染。CVXPY 的Parameter如果被复用在同一进程里连续求解时上一次的解可能作为热启动进入下一次在临界 t 附近导致判断翻转。用 Python 函数每次重新构造问题和变量是稳妥的代价是编译开销如果确实要复用比如把 t 设成cp.Parameter记得每次求解后清掉x.value或者换用不同的求解器设置。还有一种不易察觉的情况并行跑多组 t 时共用了同一个Variable对象。这个在批量扫描里很常见现象是结果和串行跑不一致且不一致的程度随线程数增加而变大。6. 进阶用松弛量做间隙用最优性条件验解二分只给到 ε 精度但它顺带产出了两样东西很多人扔掉了每轮的松弛量 s*(t)以及可行侧的候选解序列。这两样能把二分从“能用”推到“敢信”。用 s(t) 做割平面加速*。s*(t) 关于 t 单调不增且 s*(t) 0 的临界点就是最优值。把历史 (t, s*) 存下来在 s* 关于 t 近似线性的区段做一次割线插值直接跳到预测的零点附近比纯二分快 3 到 5 次求解。注意插值只能用来“选下一个 t”判定结果仍然靠求解器否则一个不准的插值会把区间推到不可行侧之外。def next_t_by_secant(history, lo, hi, lo_s, hi_s): history: [(t, s), ...]用最近两点做割线预测夹在 (lo, hi) 内。 if len(history) 2 or hi_s 0: return 0.5 * (lo hi) (t1, s1), (t2, s2) history[-2], history[-1] if abs(s2 - s1) 1e-12: return 0.5 * (lo hi) t_pred t2 - s2 * (t2 - t1) / (s2 - s1) # 保守回退不越过区间中点的对称位置避免外推过猛 return min(max(t_pred, lo 0.25 * (hi - lo)), lo 0.75 * (hi - lo))那个min(max(...))不是可有可无的装饰。纯割线外推在 s* 曲线接近水平时会给出远离区间的预测值一次错误的跳跃可能让区间收缩方向反向后面几十轮全废。夹在四分之一和四分之三之间最差退化成二分最好也就快一倍风险收益比是划算的。用最优性条件验解。拟凸问题的最优性条件跟凸问题不一样对可微的拟凸 f0x* 最优当且仅当对所有可行方向 y 都有∇f0(x*)ᵀ(y − x*) 0。落到代码上就是解一个小凸问题min ∇f0(x*)ᵀ(y − x*) s.t. y 满足硬约束与定义域如果最优值 ≥ −tol说明找不到下降方向x* 通过验证如果最优值是明显的负数说明二分给的不是最优解得检查区间和容差。这个验证本身是一次凸求解成本远低于重新跑一遍二分而且它给的是“是/否”的判断比对比两个近似值的差距更可靠。验证项检查方法不通过的典型原因解的可行性max(Gx − h)、aᵀx b − delta松弛容差过松、未归一化目标回代一致f0(x) − u最优性最小化∇f0ᵀ(y − x)区间方向错、ε 过大、非拟凸数值稳定性换求解器重跑一遍精度地板、约束量级悬殊最后一行是我个人最看重的。同一个模型换一个求解器重跑如果最优值的差异超过 ε那不是求解器的问题是你的精度设定已经踩进地板了回去调 ε 比换求解器有用。这个习惯是被坑出来的有一次一批 200 个算例跑完结果看起来都收敛得很好换求解器发现有 11 个的目标值差了 1e-3追下去是分母地板 delta 设得比数据尺度小了三个数量级解被推到了定义域边缘。从那以后我的sublevel_feasible里 delta 都是按数据算的不再写死。希望帮到你。本文还有配套的精品资源点击获取
