在网上搜“一元三次方程的解-牛莱公式”翻出来的结果十有八九是跑偏的有人贴一段代码就完事有人把求根公式抄错一个符号还有人把切线迭代和微积分基本定理混成一锅端上来看完依然不知道从哪儿下手。我前段时间接手一个参数拟合的小活中间需要反复求解形如 ax³bx²cxd0 的方程参数每变一次就得重算一遍根正好把这套东西从解析解到数值解完整捋了一遍。这里说的“牛莱公式”在三次方程这个语境下其实指向两块内容一块是以切线逼近为核心的迭代求根公式另一块是微积分基本定理那一套定积分与原函数的关系。前者确实是求三次方程数值根的常用手段后者跟“求根”没有直接关系但在根的隔离、面积计算、结果校验上非常好用。真正能给出闭式解的是古典求根公式它跟这两个名字都不是一回事。这篇就按“解析解打底、数值解兜底、微积分工具做校验”的思路把一元三次方程从系数到根的全流程讲透代码可以直接抄坑也一并标出来。适合正在写数值模块的同学、需要处理拟合与插值的工程岗以及被三次方程卡住过的爱好者。1. 先把这个叫法捋清楚别一上来就写代码1.1 社区里“牛莱公式”的两种指向与各自边界先把名词对齐不然代码写到一半会怀疑人生。在三次方程的讨论区里“牛莱公式”通常被用来指代两样东西。第一样是以切线逼近为核心的迭代公式也就是 x_{n1} x_n - f(x_n)/f(x_n) 这一串它的作用是给定一个初值不断用切线交点去替换当前点快速逼近真实根。这是典型的数值方法给出的是近似值需要收敛判据和迭代上限。第二样是微积分基本定理说的是定积分等于原函数在两端点的差它跟“把根解出来”没有直接因果但在判断根落在哪个区间、计算根与根之间的面积、验证三个根的和是否等于 -b/a 这些事情上非常好用。真正负责给出精确根的是古典求根公式思路是把三次方程通过平移消掉二次项变成缺项形式 t³ptq0然后令 t uv 做一次巧妙的代换把问题降维成一个二次方程解出 u³ 和 v³ 之后再开立方回代。这条路线能一次拿到全部根包括复数根代价是运算量偏大而且在“三个实根”的情形下会撞上实数开立方写出复数中间量的尴尬需要切换到三角解法绕过去。所以这三样东西的关系是古典求根公式负责“一次算全”切线迭代负责“快速精修”微积分工具负责“判断与验证”。工程上最稳的组合拳是先用求根公式拿到近似根再用两三步迭代把精度抛光到机器精度最后用定积分或者根与系数的关系做一次交叉校验。单独用任何一条都容易在边界情形上翻车。1.2 三条求解路线的横向对比下面这张表是我自己踩过一轮之后整理的不是教科书抄来的每条“主要坑”都是实际写代码时真撞上的。路线核心思想输出形式适用场景主要坑古典求根公式平移消二次项令 tuv 降维成二次方程全部根含复数教学推导、需要一次拿全根、参数化批量求解三实根情形需切三角解法浮点开立方负数会炸切线迭代一阶泰勒展开用切线零点逼近单个实根逐次逼近已知初值、需要高精度、实时性要求高初值不当会发散多重根退化成线性收敛因式分解/有理根试探试 ±d/a 的有理数候选有理根精确值系数为整数且根明显绝大多数方程没有有理根纯属碰运气选择逻辑其实很简单如果方程是运行时动态生成的、参数一直在变优先写切线迭代快且代码短如果是一次性求解、要求拿到所有根含复根走古典公式如果系数是整数且大概率有人为设计的根先花几十毫秒试有理根试中了就是白捡的。三者不冲突串起来用才是最省心的。有一点要提前说清楚三次方程在复数域上恒有三个根计重数但在实数域上可能只有一到三个。很多工程场景只关心实根所以后面讨论时我会明确区分“全部根”和“实根个数”。2. 标准化把方程压成最简形式2.1 缺项化处理为什么一定要先做平移拿到 ax³bx²cxd0a≠0第一步永远是两边同除以 a把首项系数变成 1得到 x³ (b/a)x² (c/a)x (d/a) 0。接着做平移 x t - b/(3a)把二次项消掉。这一步的动机跟解二次方程时配方是一个道理二次方程的求根公式之所以简洁是因为配方之后只剩 t² 常数三次方程同理只有把二次项干掉后面的代换才有一致的形式。具体算一下。把 x t - B记 B b/(3a)代入 x³ Ax² Cx D 0展开整理后 t² 项的系数是 3B - A令它为零正好得 B A/3 b/(3a)跟上面的平移量一致。此时一次项系数 p 和常数项 q 分别为p (3ac - b²)/(3a²)q (2b³ - 9abc 27a²d)/(27a³)这两个公式我见过至少三种写错一次的版本建议写代码时直接用上面这两行别自己临时推。标准形式 t³ pt q 0 解出 t 之后再做一次 x t - b/(3a) 就回到原变量。整个流程就这一处容易出错平移量是减号回代时别把符号搞反。举个具体例子x³ - 6x² 11x - 6 0。a1b-6c11d-6。p (33 - 36)/3 -1q (-432 594 - 162)/27 0。所以标准形式是 t³ - t 0即 t(t-1)(t1)0t 取 0、1、-1。平移量 B -b/(3a) 2回代得 x 2、3、1正好是 1、2、3。验算一下三根之和 6 -b/a两两之积和 11 c/a乘积 6 -d/a全对。2.2 一个判别式决定根的类型标准形式拿到之后下一步不是急着套公式而是先算判别式因为它直接决定走哪条分支。常用的判别式有两个等价写法用哪个都行形式一Δ (q/2)² (p/3)³形式二D -4p³ - 27q²这个跟形式一差一个正负号和常数因子我习惯用形式一因为后面的求根公式里正好要算 √Δ衔接顺。Δ 的取值根的情况实数域计算分支Δ 01 个实根 2 个共轭复根直接开立方实数运算Δ 0 且 pq0三重实根t0特殊分支直接返回Δ 0 且 p、q 不全为 02 个不同实根其中一个是二重根开立方后合并Δ 03 个不同实根必须切三角解法这里最容易误判的是 Δ 0 那一行。很多代码只写了一个 Δ0 分支结果遇到三重根时算出 NaN因为此时 √Δ 0、u v 0代换退化公式本身失效。三重根必须单独拎出来处理。另外提醒一句浮点比较的问题。Δ 是浮点算出来的理论上等于零的量在机器里可能是 1e-17。所以不能写if delta 0要写带容差的比较。容差取多少跟系数数量级有关我一般用abs(delta) 1e-12 * max(1.0, abs(p)**3, q*q)这种相对判据比拍一个固定值靠谱。2.3 判别式背后的几何解释判别式不是一个凭空冒出来的数它对应着函数图像跟横轴的交点个数。考虑 f(t) t³ pt q它的导数 f(t) 3t² p 是一个二次函数极值点位置由 p 的符号决定。p ≥ 0 时f(t) ≥ 0 恒成立函数单调递增图像一路向上跟横轴必然只交一次所以只有一个实根。这也解释了为什么 Δ 0 的分支只返回一个根。p 0 时f 有两个零点 t₀ √(-p/3) 和 -t₀函数在这里出现一个局部极大值和一个局部极小值。把 t₀ 代进 ff(t₀) t₀³ p·t₀ q t₀(t₀² p) q t₀·(-p/3 p) q (2p/3)·t₀ q同理 f(-t₀) -(2p/3)·t₀ q。两个极值一正一负异号时函数在极大值和极小值之间必然两次穿过横轴加上左侧那一支一共三个实根。它们的乘积f(t₀)·f(-t₀) q² - (4p²t₀²/9) q² - (4p²/9)·(-p/3) q² 4p³/27三实根的充要条件就是这个乘积小于零即 q² 4p³/27 0两边除以 4 得到 (q/2)² (p/3)³ 0正好是 Δ 0。推导到这儿判别式就不再是需要背的公式了它是“极大值和极小值异号”这个几何事实的代数翻译。这个几何视角还有个直接用途做根的隔离。知道了极值点位置和极值符号就能把三个根分别圈进三个单调区间里每个区间内函数单调迭代法必然收敛不会跑到隔壁区间去。3. 解析解实现从求根公式到三角解法3.1 求根公式的构造思路标准形式 t³ pt q 0 的解法核心是一次代换令 t u v。代进去展开(uv)³ p(uv) q u³ v³ 3uv(uv) p(uv) q 0把含 (uv) 的项合并u³ v³ q (3uv p)(uv) 0这一步的巧妙之处在于我们本来只有一个未知量 t现在人为拆成 u 和 v 两个自由度多了一个可以附加一个约束条件让方程变简单。取 3uv p 0也就是 uv -p/3那么上面整个式子就塌缩成 u³ v³ q 0即 u³ v³ -q。现在有两个关于 u³、v³ 的关系和是 -q积是 (uv)³ -p³/27。由韦达定理u³ 和 v³ 正好是二次方程 z² qz - p³/27 0 的两个根。解这个二次方程z -q/2 ± √((q/2)² (p/3)³) -q/2 ± √Δ所以 u³ -q/2 √Δv³ -q/2 - √Δ各自开立方再相加就得到 t u v。整条路线到这里闭环了三次方程 → 平移消项 → 代换降维 → 二次方程 → 开立方回代。这段推导我建议每个写数值代码的人都亲手推一遍。不是说推完就能背下来而是推过之后你会清楚地知道每一步的假设是什么出问题时能快速定位是哪一步失效了。3.2 不可约情形为什么必须切到三角解法按上面的公式Δ 0 时 √Δ 是虚数u³ 和 v³ 变成一对共轭复数理论上开立方之后虚部会抵消最终得到三个实根。数学上完全成立而且这个历史案例还挺有名——它正是复数概念被认真对待的起点之一。但在浮点计算里这条路非常难走。原因有三个一是你得引入复数类型纯实数代码搞不定二是复数开立方需要先转极坐标再除以 3 次方根中间涉及 atan2 和开方误差会累积三是三个根的虚部理论上是精确的零实际算出来是 1e-16 量级的残渣还得再判断一次该丢还是该留。整套流程又长又脆。三角解法绕开了这一切。当 p 0 且 Δ 0 时令 t 2√(-p/3)·cos θ代入 t³ pt q 0。先算 t³ 那一项t³ 8·(-p/3)^{3/2}·cos³θp·t p·2√(-p/3)·cosθ -2·(-p/3)·3^{1/2}...直接给结论更省事代入化简之后会得到 4cos³θ - 3cosθ cos3θ 这个三倍角公式最终整理成cos 3θ (3q)/(2p) · √(-3/p)右边这个量在 Δ 0 的条件下必然落在 [-1, 1] 区间内这正是判别式的另一个含义所以 arccos 有解θ 有三个分支t_k 2√(-p/3) · cos(θ/3 - 2πk/3)k 0, 1, 2其中 θ arccos((3q)/(2p)·√(-3/p))。三个 k 值对应三个实根一次全出来全程实数运算没有任何复数中间量。这就是三实根情形必须走这条路的根本原因。拿 x³ - 3x 1 0 试一下。p -3q 1Δ 0.25 - 1 -0.75 0走三角分支。2√(-p/3) 2√1 2右边 (3·1)/(2·(-3))·√(3/3) -0.5θ arccos(-0.5) 2π/3。三个根分别是 2cos(2π/9)、2cos(2π/9 - 2π/3)、2cos(2π/9 - 4π/3)算出来约 1.532089、0.347296、-1.879385。代回原方程验证残差都在 1e-16 量级非常干净。3.3 重根处理的细节Δ 0 分两种情况代码里必须分开写。第一种是 p 0 且 q 0此时方程退化成 t³ 0三重根 t 0回到原变量就是 x -b/(3a)。这个分支要放在最前面因为后面的公式会算出 0/0。第二种是 p、q 不全为零。此时由 Δ 0 可以推出 (q/2)² -(p/3)³令 u ∛(-q/2)则三个根为 t 2u单根和 t -u二重根。验证一下把 t 2u 代入8u³ 2pu q。由 u³ -q/2 得 8u³ -4q由 p -3u² 得 2pu -6u³ 3q加起来 -4q 3q q 0成立。再验 t -u-u³ - pu q q/2 3u³ q q/2 - 3q/2 q 0也成立。这里有个实践细节二重根在数值上非常敏感。由于浮点误差Δ 往往不会精确等于 0而是 1e-18 或者 -1e-18。如果落在正的一侧会走 Δ 0 的分支算出来是一个实根如果落在负的一侧会走三角分支算出三个几乎相等的实根。两种结果都不算错但根的形式不一致下游代码如果按“根的数量”做判断就会出问题。我的处理方式是在 Δ 的绝对值小于某个阈值时统一按重根分支处理并在返回值里带上一个标记位让调用方知道“这里有重根”。3.4 浮点开立方的符号陷阱这个坑必须单独占一段因为它太常见了。Python 里对一个负数做分数次幂结果不是负数而是复数 (-8) ** (1/3) (1.00000000000000021.7320508075688772j)原因是**对浮点数走的是pow负底数配分数指数在实数域无定义Python 自动切到复数分支。如果你在求根公式里写(-q/2 s) ** (1/3)遇到负的被开方数就会拿到一个复数后面所有运算全部被污染最终结果的类型会从 float 变成 complex排序、比较、转字符串全炸。正确的写法是用带符号的立方根def _cbrt(x): return math.copysign(abs(x) ** (1.0 / 3.0), x)先取绝对值开立方再把符号贴回去。整个求根公式里凡是要开立方的地方一律走这个函数不要图省事直接写** (1/3)。这个小函数能挡掉后面 80% 的诡异报错。顺便说一句精度问题abs(x) ** (1/3)在 x 很大或很小时精度一般如果要追求极限精度可以用math.exp(math.log(abs(x))/3)做初值再补一步迭代但实际用下来前者对 1e-15 的精度要求已经完全够用没必要过度设计。4. 数值解实现切线迭代的工程细节4.1 迭代公式的推导与几何意义切线迭代的来源很简单把 f(x) 在当前点 x_n 处做一阶泰勒展开得到 f(x) ≈ f(x_n) f(x_n)(x - x_n)令近似式为零解出 x x_n - f(x_n)/f(x_n)把这个 x 作为下一个迭代点。几何上就是在曲线上的点 (x_n, f(x_n)) 处画一条切线切线与横轴的交点就是下一个近似点。对三次函数 f(x) ax³bx²cxd导数是 f(x) 3ax²2bxc两个都能用秦九韶Horner格式快速求值f ((a*x b)*x c)*x d df (3*a*x 2*b)*x c秦九韶格式比逐项乘方快不少而且避免了 x³ 溢出。系数不大时差别不明显但如果在循环里跑几万次差距就出来了。收敛速度方面单根附近切线迭代是二次收敛的也就是每迭代一次误差的有效位数大约翻倍。拿 x³ - 3x 1 0 从 x₀ 1.5 起步实测第一次迭代到 1.5333第二次到 1.53193第三次到 1.5320888第四次就到 1.5320888862跟真值的差已经低于 1e-12。四步非常快。4.2 初值选取为什么你的程序会跑飞切线迭代最大的问题是它不保证全局收敛。f(x_n) 接近零时步长会爆掉直接飞到几万公里外初值离目标根太远时可能收敛到另一个根或者干脆在几个值之间来回震荡。对三次方程来说初值其实很好选因为函数结构是已知的。先算 p再分情况p ≥ 0 时函数单调全场只有一个实根随便找个初值都能收敛。常用的有界估计是 |x| ≤ 1 max(|b/a|, |c/a|, |d/a|)从 0 或者这个界的一半起步都行。p 0 时函数有两个极值点 ±√(-p/3)它们把实数轴切成三段每段内单调。要想把三个根都找出来就从三个区间里各取一个初值左区间取 -2√(-p/3)中间取 0右区间取 2√(-p/3)。这三个点是经典的超界估计能保证分别收敛到三个不同的根如果那个区间里存在根的话。还有一种更稳的做法是做区间二分定位从大范围开始检查端点符号异号就二分直到区间足够小再切到切线迭代做精修。二分保证收敛但只有线性速度切线迭代快但需要好初值两者交替使用是数值分析里的经典组合实测下来几乎没有失败的案例。提示判断某个区间有没有根看 f 在两端点的符号是否相反即可。但要注意如果区间里有两个根端点同号这个方法会漏判。对三次方程来说只要区间是从极值点切开的单调区间就不会出现这种情况。4.3 收敛判据与精度控制只用一个判据容易出问题。我一般用双判据if abs(x_new - x) tol * max(1.0, abs(x_new)) and abs(f_new) ftol: return x_new第一个判据看步长相对量比绝对量稳因为根的量级可能很小也可能很大第二个判据看残差防止步长很小但函数值还很大这通常意味着导数值异常迭代已经失去意义。参数取值上tol 我一般设 1e-14ftol 设 1e-12 乘以系数数量级。迭代上限设 100 足够了正常情况 10 步内必然收敛跑到 100 步还没停说明初值或者其他地方出了问题继续跑也没意义。还有一个必须加的保护导数绝对值过小时直接放弃这一轮。如果 |f(x)| 1e-14说明当前点几乎在驻点上切线几乎平行于横轴下一步会飞出去。这时候应该退回用二分法走一步或者换个初值重来。4.4 多重根收敛退化与修正标准切线迭代在单根附近是二次收敛但在多重根附近会退化成线性收敛。原因不复杂设根为 r重数为 m则 f(x) ≈ C(x-r)^mf(x) ≈ mC(x-r)^{m-1}两者相除 ≈ (x-r)/m所以每步误差只缩小到原来的 (m-1)/m。二重根是 1/2三重根是 2/3都远慢于平方。实测一下f(x) (x-1)³ 从 1.2 起步标准迭代的误差序列大约是 0.2 → 0.133 → 0.089 → 0.059每次缩到 2/3要跑到 1e-12 大约需要六七十步。而如果知道重数 m 3改用 x - m·f/f收敛立刻恢复二阶五六步就到位。问题在于重数通常事先不知道。这时候有两个选择。一是用数值方式估计重数在根附近取两三个点算 f 和 f 的比值看它符合哪个 m。二是换用对多重根更友好的迭代式x_{n1} x_n - 2·f(x_n)·f(x_n) / (2f(x_n)² - f(x_n)·f(x_n))这个式子是二阶修正的产物本质上是 Halley 型迭代单根附近收敛阶是三阶比标准切线迭代还快一档多重根附近仍然是线性但收敛常数更小实测比标准式少花三成左右的步数。对三次函数来说 f(x) 6ax2b成本很低值得换取这点稳定性。还有一个技巧如果已经通过解析解拿到一个近似根 r但又想精修到机器精度不要直接在 r 附近跑标准迭代——如果它是多重根会跑得很慢。正确做法是先判断它是不是重根同时看 |f(r)| 和 |f(r)|两个都接近零就是重根是的话用上面那个修正式不是的话标准式就够。5. 微积分那一套在三次方程上的落脚点5.1 用导数与定积分做根的隔离微积分基本定理说的是 f(b) - f(a) ∫ₐᵇ f(x)dx它的实际意义是函数在区间上的增量等于导数的定积分。这个关系在三次方程的根隔离上有直接用途。三次函数的导数是二次函数二次方程可以精确求解所以 f 的零点也就是 f 的极值点位置是解析可得的x (-b ± √(b²-3ac))/(3a)。拿到这两个点之后实数轴被分成最多三段单调区间。每段里 f 要么单调增要么单调减用端点值判断符号符号翻转就说明这一段里有且只有一个根。这套流程的严密性靠的正是导数信息在单调区间上函数值的增量也就是定积分严格同号不会来回穿越横轴所以最多一个根端点异号则至少一个根两者合起来就是恰好一个。实际代码里我会这样组织先算 b²-3ac小于等于零就直接认定只有一个实根用切线迭代从 0 起步大于零就切两个极值点 x₁ x₂分别检查 f(x₁)、f(x₂) 与两端极限的符号组装出根的个数和位置区间再对每个区间跑迭代。这套逻辑不依赖求根公式也不需要判别式纯靠导数信息写起来直观调试也方便。5.2 面积、体积与加权平均的具体计算定积分在三次方程相关的计算里最常见的用途是求曲线与横轴围成的面积。三次函数的原函数是四次多项式直接写出F(x) a·x⁴/4 b·x³/3 c·x²/2 d·x两根 r₁ r₂ 之间曲线下的面积带符号就是 F(r₂) - F(r₁)取绝对值就是几何面积。拿前面那个例子算一下。f(x) x³ - 6x² 11x - 6三个根 1、2、3。求它在 [1, 2] 上围成的面积。F(x) x⁴/4 - 2x³ 11x²/2 - 6x。F(2) 4 - 16 22 - 12 -2F(1) 0.25 - 2 5.5 - 6 -2.25差的绝对值等于 0.25。手动验证一下令 u x - 2x 从 1 到 2 时 u 从 -1 到 0被积函数变成 (u1)u(u-1) u³ - u积分得 [u⁴/4 - u²/2]在 -1 到 0 之间取值是 0 - (0.25 - 0.5) 0.25完全一致。这个结果有个反直觉的地方面积只有 0.25看着很小但曲线在 1 到 2 之间的峰值高度接近 0.4宽度 1采样式的心算容易高估。这也说明涉及三次函数的面积问题最好老实算积分不要凭感觉估。除了面积还有几个衍生量在工程里常用计算目标公式典型用途曲线 
长河君 · 2026/10/31 23:38:49 · 阅读 1.2k 一元三次方程求解全流程：解析解、切线迭代与微积分校验
在网上搜“一元三次方程的解-牛莱公式”翻出来的结果十有八九是跑偏的有人贴一段代码就完事有人把求根公式抄错一个符号还有人把切线迭代和微积分基本定理混成一锅端上来看完依然不知道从哪儿下手。我前段时间接手一个参数拟合的小活中间需要反复求解形如 ax³bx²cxd0 的方程参数每变一次就得重算一遍根正好把这套东西从解析解到数值解完整捋了一遍。这里说的“牛莱公式”在三次方程这个语境下其实指向两块内容一块是以切线逼近为核心的迭代求根公式另一块是微积分基本定理那一套定积分与原函数的关系。前者确实是求三次方程数值根的常用手段后者跟“求根”没有直接关系但在根的隔离、面积计算、结果校验上非常好用。真正能给出闭式解的是古典求根公式它跟这两个名字都不是一回事。这篇就按“解析解打底、数值解兜底、微积分工具做校验”的思路把一元三次方程从系数到根的全流程讲透代码可以直接抄坑也一并标出来。适合正在写数值模块的同学、需要处理拟合与插值的工程岗以及被三次方程卡住过的爱好者。1. 先把这个叫法捋清楚别一上来就写代码1.1 社区里“牛莱公式”的两种指向与各自边界先把名词对齐不然代码写到一半会怀疑人生。在三次方程的讨论区里“牛莱公式”通常被用来指代两样东西。第一样是以切线逼近为核心的迭代公式也就是 x_{n1} x_n - f(x_n)/f(x_n) 这一串它的作用是给定一个初值不断用切线交点去替换当前点快速逼近真实根。这是典型的数值方法给出的是近似值需要收敛判据和迭代上限。第二样是微积分基本定理说的是定积分等于原函数在两端点的差它跟“把根解出来”没有直接因果但在判断根落在哪个区间、计算根与根之间的面积、验证三个根的和是否等于 -b/a 这些事情上非常好用。真正负责给出精确根的是古典求根公式思路是把三次方程通过平移消掉二次项变成缺项形式 t³ptq0然后令 t uv 做一次巧妙的代换把问题降维成一个二次方程解出 u³ 和 v³ 之后再开立方回代。这条路线能一次拿到全部根包括复数根代价是运算量偏大而且在“三个实根”的情形下会撞上实数开立方写出复数中间量的尴尬需要切换到三角解法绕过去。所以这三样东西的关系是古典求根公式负责“一次算全”切线迭代负责“快速精修”微积分工具负责“判断与验证”。工程上最稳的组合拳是先用求根公式拿到近似根再用两三步迭代把精度抛光到机器精度最后用定积分或者根与系数的关系做一次交叉校验。单独用任何一条都容易在边界情形上翻车。1.2 三条求解路线的横向对比下面这张表是我自己踩过一轮之后整理的不是教科书抄来的每条“主要坑”都是实际写代码时真撞上的。路线核心思想输出形式适用场景主要坑古典求根公式平移消二次项令 tuv 降维成二次方程全部根含复数教学推导、需要一次拿全根、参数化批量求解三实根情形需切三角解法浮点开立方负数会炸切线迭代一阶泰勒展开用切线零点逼近单个实根逐次逼近已知初值、需要高精度、实时性要求高初值不当会发散多重根退化成线性收敛因式分解/有理根试探试 ±d/a 的有理数候选有理根精确值系数为整数且根明显绝大多数方程没有有理根纯属碰运气选择逻辑其实很简单如果方程是运行时动态生成的、参数一直在变优先写切线迭代快且代码短如果是一次性求解、要求拿到所有根含复根走古典公式如果系数是整数且大概率有人为设计的根先花几十毫秒试有理根试中了就是白捡的。三者不冲突串起来用才是最省心的。有一点要提前说清楚三次方程在复数域上恒有三个根计重数但在实数域上可能只有一到三个。很多工程场景只关心实根所以后面讨论时我会明确区分“全部根”和“实根个数”。2. 标准化把方程压成最简形式2.1 缺项化处理为什么一定要先做平移拿到 ax³bx²cxd0a≠0第一步永远是两边同除以 a把首项系数变成 1得到 x³ (b/a)x² (c/a)x (d/a) 0。接着做平移 x t - b/(3a)把二次项消掉。这一步的动机跟解二次方程时配方是一个道理二次方程的求根公式之所以简洁是因为配方之后只剩 t² 常数三次方程同理只有把二次项干掉后面的代换才有一致的形式。具体算一下。把 x t - B记 B b/(3a)代入 x³ Ax² Cx D 0展开整理后 t² 项的系数是 3B - A令它为零正好得 B A/3 b/(3a)跟上面的平移量一致。此时一次项系数 p 和常数项 q 分别为p (3ac - b²)/(3a²)q (2b³ - 9abc 27a²d)/(27a³)这两个公式我见过至少三种写错一次的版本建议写代码时直接用上面这两行别自己临时推。标准形式 t³ pt q 0 解出 t 之后再做一次 x t - b/(3a) 就回到原变量。整个流程就这一处容易出错平移量是减号回代时别把符号搞反。举个具体例子x³ - 6x² 11x - 6 0。a1b-6c11d-6。p (33 - 36)/3 -1q (-432 594 - 162)/27 0。所以标准形式是 t³ - t 0即 t(t-1)(t1)0t 取 0、1、-1。平移量 B -b/(3a) 2回代得 x 2、3、1正好是 1、2、3。验算一下三根之和 6 -b/a两两之积和 11 c/a乘积 6 -d/a全对。2.2 一个判别式决定根的类型标准形式拿到之后下一步不是急着套公式而是先算判别式因为它直接决定走哪条分支。常用的判别式有两个等价写法用哪个都行形式一Δ (q/2)² (p/3)³形式二D -4p³ - 27q²这个跟形式一差一个正负号和常数因子我习惯用形式一因为后面的求根公式里正好要算 √Δ衔接顺。Δ 的取值根的情况实数域计算分支Δ 01 个实根 2 个共轭复根直接开立方实数运算Δ 0 且 pq0三重实根t0特殊分支直接返回Δ 0 且 p、q 不全为 02 个不同实根其中一个是二重根开立方后合并Δ 03 个不同实根必须切三角解法这里最容易误判的是 Δ 0 那一行。很多代码只写了一个 Δ0 分支结果遇到三重根时算出 NaN因为此时 √Δ 0、u v 0代换退化公式本身失效。三重根必须单独拎出来处理。另外提醒一句浮点比较的问题。Δ 是浮点算出来的理论上等于零的量在机器里可能是 1e-17。所以不能写if delta 0要写带容差的比较。容差取多少跟系数数量级有关我一般用abs(delta) 1e-12 * max(1.0, abs(p)**3, q*q)这种相对判据比拍一个固定值靠谱。2.3 判别式背后的几何解释判别式不是一个凭空冒出来的数它对应着函数图像跟横轴的交点个数。考虑 f(t) t³ pt q它的导数 f(t) 3t² p 是一个二次函数极值点位置由 p 的符号决定。p ≥ 0 时f(t) ≥ 0 恒成立函数单调递增图像一路向上跟横轴必然只交一次所以只有一个实根。这也解释了为什么 Δ 0 的分支只返回一个根。p 0 时f 有两个零点 t₀ √(-p/3) 和 -t₀函数在这里出现一个局部极大值和一个局部极小值。把 t₀ 代进 ff(t₀) t₀³ p·t₀ q t₀(t₀² p) q t₀·(-p/3 p) q (2p/3)·t₀ q同理 f(-t₀) -(2p/3)·t₀ q。两个极值一正一负异号时函数在极大值和极小值之间必然两次穿过横轴加上左侧那一支一共三个实根。它们的乘积f(t₀)·f(-t₀) q² - (4p²t₀²/9) q² - (4p²/9)·(-p/3) q² 4p³/27三实根的充要条件就是这个乘积小于零即 q² 4p³/27 0两边除以 4 得到 (q/2)² (p/3)³ 0正好是 Δ 0。推导到这儿判别式就不再是需要背的公式了它是“极大值和极小值异号”这个几何事实的代数翻译。这个几何视角还有个直接用途做根的隔离。知道了极值点位置和极值符号就能把三个根分别圈进三个单调区间里每个区间内函数单调迭代法必然收敛不会跑到隔壁区间去。3. 解析解实现从求根公式到三角解法3.1 求根公式的构造思路标准形式 t³ pt q 0 的解法核心是一次代换令 t u v。代进去展开(uv)³ p(uv) q u³ v³ 3uv(uv) p(uv) q 0把含 (uv) 的项合并u³ v³ q (3uv p)(uv) 0这一步的巧妙之处在于我们本来只有一个未知量 t现在人为拆成 u 和 v 两个自由度多了一个可以附加一个约束条件让方程变简单。取 3uv p 0也就是 uv -p/3那么上面整个式子就塌缩成 u³ v³ q 0即 u³ v³ -q。现在有两个关于 u³、v³ 的关系和是 -q积是 (uv)³ -p³/27。由韦达定理u³ 和 v³ 正好是二次方程 z² qz - p³/27 0 的两个根。解这个二次方程z -q/2 ± √((q/2)² (p/3)³) -q/2 ± √Δ所以 u³ -q/2 √Δv³ -q/2 - √Δ各自开立方再相加就得到 t u v。整条路线到这里闭环了三次方程 → 平移消项 → 代换降维 → 二次方程 → 开立方回代。这段推导我建议每个写数值代码的人都亲手推一遍。不是说推完就能背下来而是推过之后你会清楚地知道每一步的假设是什么出问题时能快速定位是哪一步失效了。3.2 不可约情形为什么必须切到三角解法按上面的公式Δ 0 时 √Δ 是虚数u³ 和 v³ 变成一对共轭复数理论上开立方之后虚部会抵消最终得到三个实根。数学上完全成立而且这个历史案例还挺有名——它正是复数概念被认真对待的起点之一。但在浮点计算里这条路非常难走。原因有三个一是你得引入复数类型纯实数代码搞不定二是复数开立方需要先转极坐标再除以 3 次方根中间涉及 atan2 和开方误差会累积三是三个根的虚部理论上是精确的零实际算出来是 1e-16 量级的残渣还得再判断一次该丢还是该留。整套流程又长又脆。三角解法绕开了这一切。当 p 0 且 Δ 0 时令 t 2√(-p/3)·cos θ代入 t³ pt q 0。先算 t³ 那一项t³ 8·(-p/3)^{3/2}·cos³θp·t p·2√(-p/3)·cosθ -2·(-p/3)·3^{1/2}...直接给结论更省事代入化简之后会得到 4cos³θ - 3cosθ cos3θ 这个三倍角公式最终整理成cos 3θ (3q)/(2p) · √(-3/p)右边这个量在 Δ 0 的条件下必然落在 [-1, 1] 区间内这正是判别式的另一个含义所以 arccos 有解θ 有三个分支t_k 2√(-p/3) · cos(θ/3 - 2πk/3)k 0, 1, 2其中 θ arccos((3q)/(2p)·√(-3/p))。三个 k 值对应三个实根一次全出来全程实数运算没有任何复数中间量。这就是三实根情形必须走这条路的根本原因。拿 x³ - 3x 1 0 试一下。p -3q 1Δ 0.25 - 1 -0.75 0走三角分支。2√(-p/3) 2√1 2右边 (3·1)/(2·(-3))·√(3/3) -0.5θ arccos(-0.5) 2π/3。三个根分别是 2cos(2π/9)、2cos(2π/9 - 2π/3)、2cos(2π/9 - 4π/3)算出来约 1.532089、0.347296、-1.879385。代回原方程验证残差都在 1e-16 量级非常干净。3.3 重根处理的细节Δ 0 分两种情况代码里必须分开写。第一种是 p 0 且 q 0此时方程退化成 t³ 0三重根 t 0回到原变量就是 x -b/(3a)。这个分支要放在最前面因为后面的公式会算出 0/0。第二种是 p、q 不全为零。此时由 Δ 0 可以推出 (q/2)² -(p/3)³令 u ∛(-q/2)则三个根为 t 2u单根和 t -u二重根。验证一下把 t 2u 代入8u³ 2pu q。由 u³ -q/2 得 8u³ -4q由 p -3u² 得 2pu -6u³ 3q加起来 -4q 3q q 0成立。再验 t -u-u³ - pu q q/2 3u³ q q/2 - 3q/2 q 0也成立。这里有个实践细节二重根在数值上非常敏感。由于浮点误差Δ 往往不会精确等于 0而是 1e-18 或者 -1e-18。如果落在正的一侧会走 Δ 0 的分支算出来是一个实根如果落在负的一侧会走三角分支算出三个几乎相等的实根。两种结果都不算错但根的形式不一致下游代码如果按“根的数量”做判断就会出问题。我的处理方式是在 Δ 的绝对值小于某个阈值时统一按重根分支处理并在返回值里带上一个标记位让调用方知道“这里有重根”。3.4 浮点开立方的符号陷阱这个坑必须单独占一段因为它太常见了。Python 里对一个负数做分数次幂结果不是负数而是复数 (-8) ** (1/3) (1.00000000000000021.7320508075688772j)原因是**对浮点数走的是pow负底数配分数指数在实数域无定义Python 自动切到复数分支。如果你在求根公式里写(-q/2 s) ** (1/3)遇到负的被开方数就会拿到一个复数后面所有运算全部被污染最终结果的类型会从 float 变成 complex排序、比较、转字符串全炸。正确的写法是用带符号的立方根def _cbrt(x): return math.copysign(abs(x) ** (1.0 / 3.0), x)先取绝对值开立方再把符号贴回去。整个求根公式里凡是要开立方的地方一律走这个函数不要图省事直接写** (1/3)。这个小函数能挡掉后面 80% 的诡异报错。顺便说一句精度问题abs(x) ** (1/3)在 x 很大或很小时精度一般如果要追求极限精度可以用math.exp(math.log(abs(x))/3)做初值再补一步迭代但实际用下来前者对 1e-15 的精度要求已经完全够用没必要过度设计。4. 数值解实现切线迭代的工程细节4.1 迭代公式的推导与几何意义切线迭代的来源很简单把 f(x) 在当前点 x_n 处做一阶泰勒展开得到 f(x) ≈ f(x_n) f(x_n)(x - x_n)令近似式为零解出 x x_n - f(x_n)/f(x_n)把这个 x 作为下一个迭代点。几何上就是在曲线上的点 (x_n, f(x_n)) 处画一条切线切线与横轴的交点就是下一个近似点。对三次函数 f(x) ax³bx²cxd导数是 f(x) 3ax²2bxc两个都能用秦九韶Horner格式快速求值f ((a*x b)*x c)*x d df (3*a*x 2*b)*x c秦九韶格式比逐项乘方快不少而且避免了 x³ 溢出。系数不大时差别不明显但如果在循环里跑几万次差距就出来了。收敛速度方面单根附近切线迭代是二次收敛的也就是每迭代一次误差的有效位数大约翻倍。拿 x³ - 3x 1 0 从 x₀ 1.5 起步实测第一次迭代到 1.5333第二次到 1.53193第三次到 1.5320888第四次就到 1.5320888862跟真值的差已经低于 1e-12。四步非常快。4.2 初值选取为什么你的程序会跑飞切线迭代最大的问题是它不保证全局收敛。f(x_n) 接近零时步长会爆掉直接飞到几万公里外初值离目标根太远时可能收敛到另一个根或者干脆在几个值之间来回震荡。对三次方程来说初值其实很好选因为函数结构是已知的。先算 p再分情况p ≥ 0 时函数单调全场只有一个实根随便找个初值都能收敛。常用的有界估计是 |x| ≤ 1 max(|b/a|, |c/a|, |d/a|)从 0 或者这个界的一半起步都行。p 0 时函数有两个极值点 ±√(-p/3)它们把实数轴切成三段每段内单调。要想把三个根都找出来就从三个区间里各取一个初值左区间取 -2√(-p/3)中间取 0右区间取 2√(-p/3)。这三个点是经典的超界估计能保证分别收敛到三个不同的根如果那个区间里存在根的话。还有一种更稳的做法是做区间二分定位从大范围开始检查端点符号异号就二分直到区间足够小再切到切线迭代做精修。二分保证收敛但只有线性速度切线迭代快但需要好初值两者交替使用是数值分析里的经典组合实测下来几乎没有失败的案例。提示判断某个区间有没有根看 f 在两端点的符号是否相反即可。但要注意如果区间里有两个根端点同号这个方法会漏判。对三次方程来说只要区间是从极值点切开的单调区间就不会出现这种情况。4.3 收敛判据与精度控制只用一个判据容易出问题。我一般用双判据if abs(x_new - x) tol * max(1.0, abs(x_new)) and abs(f_new) ftol: return x_new第一个判据看步长相对量比绝对量稳因为根的量级可能很小也可能很大第二个判据看残差防止步长很小但函数值还很大这通常意味着导数值异常迭代已经失去意义。参数取值上tol 我一般设 1e-14ftol 设 1e-12 乘以系数数量级。迭代上限设 100 足够了正常情况 10 步内必然收敛跑到 100 步还没停说明初值或者其他地方出了问题继续跑也没意义。还有一个必须加的保护导数绝对值过小时直接放弃这一轮。如果 |f(x)| 1e-14说明当前点几乎在驻点上切线几乎平行于横轴下一步会飞出去。这时候应该退回用二分法走一步或者换个初值重来。4.4 多重根收敛退化与修正标准切线迭代在单根附近是二次收敛但在多重根附近会退化成线性收敛。原因不复杂设根为 r重数为 m则 f(x) ≈ C(x-r)^mf(x) ≈ mC(x-r)^{m-1}两者相除 ≈ (x-r)/m所以每步误差只缩小到原来的 (m-1)/m。二重根是 1/2三重根是 2/3都远慢于平方。实测一下f(x) (x-1)³ 从 1.2 起步标准迭代的误差序列大约是 0.2 → 0.133 → 0.089 → 0.059每次缩到 2/3要跑到 1e-12 大约需要六七十步。而如果知道重数 m 3改用 x - m·f/f收敛立刻恢复二阶五六步就到位。问题在于重数通常事先不知道。这时候有两个选择。一是用数值方式估计重数在根附近取两三个点算 f 和 f 的比值看它符合哪个 m。二是换用对多重根更友好的迭代式x_{n1} x_n - 2·f(x_n)·f(x_n) / (2f(x_n)² - f(x_n)·f(x_n))这个式子是二阶修正的产物本质上是 Halley 型迭代单根附近收敛阶是三阶比标准切线迭代还快一档多重根附近仍然是线性但收敛常数更小实测比标准式少花三成左右的步数。对三次函数来说 f(x) 6ax2b成本很低值得换取这点稳定性。还有一个技巧如果已经通过解析解拿到一个近似根 r但又想精修到机器精度不要直接在 r 附近跑标准迭代——如果它是多重根会跑得很慢。正确做法是先判断它是不是重根同时看 |f(r)| 和 |f(r)|两个都接近零就是重根是的话用上面那个修正式不是的话标准式就够。5. 微积分那一套在三次方程上的落脚点5.1 用导数与定积分做根的隔离微积分基本定理说的是 f(b) - f(a) ∫ₐᵇ f(x)dx它的实际意义是函数在区间上的增量等于导数的定积分。这个关系在三次方程的根隔离上有直接用途。三次函数的导数是二次函数二次方程可以精确求解所以 f 的零点也就是 f 的极值点位置是解析可得的x (-b ± √(b²-3ac))/(3a)。拿到这两个点之后实数轴被分成最多三段单调区间。每段里 f 要么单调增要么单调减用端点值判断符号符号翻转就说明这一段里有且只有一个根。这套流程的严密性靠的正是导数信息在单调区间上函数值的增量也就是定积分严格同号不会来回穿越横轴所以最多一个根端点异号则至少一个根两者合起来就是恰好一个。实际代码里我会这样组织先算 b²-3ac小于等于零就直接认定只有一个实根用切线迭代从 0 起步大于零就切两个极值点 x₁ x₂分别检查 f(x₁)、f(x₂) 与两端极限的符号组装出根的个数和位置区间再对每个区间跑迭代。这套逻辑不依赖求根公式也不需要判别式纯靠导数信息写起来直观调试也方便。5.2 面积、体积与加权平均的具体计算定积分在三次方程相关的计算里最常见的用途是求曲线与横轴围成的面积。三次函数的原函数是四次多项式直接写出F(x) a·x⁴/4 b·x³/3 c·x²/2 d·x两根 r₁ r₂ 之间曲线下的面积带符号就是 F(r₂) - F(r₁)取绝对值就是几何面积。拿前面那个例子算一下。f(x) x³ - 6x² 11x - 6三个根 1、2、3。求它在 [1, 2] 上围成的面积。F(x) x⁴/4 - 2x³ 11x²/2 - 6x。F(2) 4 - 16 22 - 12 -2F(1) 0.25 - 2 5.5 - 6 -2.25差的绝对值等于 0.25。手动验证一下令 u x - 2x 从 1 到 2 时 u 从 -1 到 0被积函数变成 (u1)u(u-1) u³ - u积分得 [u⁴/4 - u²/2]在 -1 到 0 之间取值是 0 - (0.25 - 0.5) 0.25完全一致。这个结果有个反直觉的地方面积只有 0.25看着很小但曲线在 1 到 2 之间的峰值高度接近 0.4宽度 1采样式的心算容易高估。这也说明涉及三次函数的面积问题最好老实算积分不要凭感觉估。除了面积还有几个衍生量在工程里常用计算目标公式典型用途曲线