多项式牛顿迭代

描述

给定多项式 G(x,y)G\left(x, y\right),已知多项式 f(x)f\left(x\right) 满足:

G(x,f(x))≡0(modxn)G\left(x, f\left(x\right)\right)\equiv 0\pmod{x^{n}}

且存在数值 f1f_1 使 G(x,y)G\left(x, y\right) 满足以下条件:

  • G(0,f1)=0G(0, f_1) = 0;
  • ∂G∂y(0,f1)≠0\dfrac{\partial G}{\partial y}(0, f_1) \neq 0.

求出模 xnx^{n} 意义下的 f(x)f\left(x\right).

Newton's Method

考虑倍增.

首先当 n=1n=1 时,[x0]G(x,f(x))=0\left[x^{0}\right]G\left(x, f\left(x\right)\right)=0 的解需要单独求出,假设中的 f1f_1 即为一个解.

假设现在已经得到了模 x⌈n2⌉x^{\left\lceil\frac{n}{2}\right\rceil} 意义下的解 f⌈n2⌉(x)f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right),要求模 xnx^{n} 意义下的解 f(x)=fn(x)f\left(x\right) = f_n\left(x\right).

将 G(x,f(x))G\left(x, f(x)\right) 对 f(x)f(x) 在 f(x)=f⌈n2⌉(x)f(x)=f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right) 处进行泰勒展开,有:

∑i=0+∞∂iG∂yi(x,f⌈n2⌉(x))i!(f(x)−f⌈n2⌉(x))i≡0(modxn)\sum_{i=0}^{+\infty}\frac{\frac{\partial^i G}{\partial y^i}\left(x, f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right)\right)}{i!}\left(f\left(x\right)-f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right)\right)^{i}\equiv 0\pmod{x^{n}}

因为 f(x)−f⌈n2⌉(x)f\left(x\right)-f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right) 的最低非零项次数最低为 ⌈n2⌉\left\lceil\frac{n}{2}\right\rceil,故有:

∀2⩽i:(f(x)−f⌈n2⌉(x))i≡0(modxn)\forall 2\leqslant i:\left(f\left(x\right)-f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right)\right)^{i}\equiv 0\pmod{x^{n}}

则:

∑i=0+∞∂iG∂yi(x,f⌈n2⌉(x))i!(f(x)−f⌈n2⌉(x))i≡G(x,f⌈n2⌉(x))+∂G∂y(x,f⌈n2⌉(x))[f(x)−f⌈n2⌉(x)]≡0(modxn)\begin{aligned} \sum_{i=0}^{+\infty}\frac{\frac{\partial^i G}{\partial y^i}\left(x, f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right)\right)}{i!}\left(f\left(x\right)-f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right)\right)^{i}&\equiv G\left(x, f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right)\right)+\frac{\partial G}{\partial y}\left(x, f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right)\right)\left[f\left(x\right)-f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right)\right]\\ &\equiv 0\pmod{x^{n}} \end{aligned} fn(x)≡f⌈n2⌉(x)−G(x,f⌈n2⌉(x))∂G∂y(x,f⌈n2⌉(x))(modxn)f_n\left(x\right)\equiv f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right)-\frac{G\left(x, f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right)\right)}{\frac{\partial G}{\partial y}\left(x, f_{\left\lceil\frac{n}{2}\right\rceil}\left(x\right)\right)}\pmod{x^{n}}

或者

f2n(x)≡fn(x)−G(x,fn(x))∂G∂y(x,fn(x))(modx2n)f_{2n}\left(x\right)\equiv f_n\left(x\right)-\frac{G\left(x, f_n\left(x\right)\right)}{\frac{\partial G}{\partial y}\left(x, f_n\left(x\right)\right)}\pmod{x^{2n}}

例题

多项式求逆

设给定函数为 h(x)h\left(x\right),有:

G(x,y)=1y−h(x)G\left(x, y\right)=\frac{1}{y}-h\left(x\right)

应用 Newton's Method 可得:

f2n(x)≡fn(x)−1/fn(x)−h(x)−1/fn2(x)(modx2n)≡2fn(x)−fn2(x)h(x)(modx2n)\begin{aligned} f_{2n}\left(x\right)&\equiv f_{n}\left(x\right)-\frac{1/f_{n}\left(x\right)-h\left(x\right)}{-1/f_{n}^{2}\left(x\right)}&\pmod{x^{2n}}\\ &\equiv 2f_{n}\left(x\right)-f_{n}^{2}\left(x\right)h\left(x\right)&\pmod{x^{2n}} \end{aligned}

时间复杂度

T(n)=T(n2)+O(nlog⁡n)=O(nlog⁡n)T\left(n\right)=T\left(\frac{n}{2}\right)+O\left(n\log{n}\right)=O\left(n\log{n}\right)

多项式开方

设给定函数为 h(x)h\left(x\right),有:

G(x,y)=y2−h(x)≡0G\left(x, y\right)=y^{2}-h\left(x\right)\equiv 0

应用 Newton's Method 可得:

f2n(x)≡fn(x)−fn2(x)−h(x)2fn(x)(modx2n)≡fn2(x)+h(x)2fn(x)(modx2n)\begin{aligned} f_{2n}\left(x\right)&\equiv f_{n}\left(x\right)-\frac{f_{n}^{2}\left(x\right)-h\left(x\right)}{2f_{n}\left(x\right)}&\pmod{x^{2n}}\\ &\equiv\frac{f_{n}^{2}\left(x\right)+h\left(x\right)}{2f_{n}\left(x\right)}&\pmod{x^{2n}} \end{aligned}

时间复杂度

T(n)=T(n2)+O(nlog⁡n)=O(nlog⁡n)T\left(n\right)=T\left(\frac{n}{2}\right)+O\left(n\log{n}\right)=O\left(n\log{n}\right)

多项式指数函数

设给定函数为 h(x)h\left(x\right),有:

G(x,y)=ln⁡y−h(x)G\left(x, y\right)=\ln{y}-h\left(x\right)

应用 Newton's Method 可得:

f2n(x)≡fn(x)−ln⁡fn(x)−h(x)1/fn(x)(modx2n)≡fn(x)(1−ln⁡fn(x)+h(x))(modx2n)\begin{aligned} f_{2n}\left(x\right)&\equiv f_{n}\left(x\right)-\frac{\ln{f_{n}\left(x\right)}-h\left(x\right)}{1/f_{n}\left(x\right)}&\pmod{x^{2n}}\\ &\equiv f_{n}\left(x\right)\left(1-\ln{f_{n}\left(x\right)+h\left(x\right)}\right)&\pmod{x^{2n}} \end{aligned}

时间复杂度

T(n)=T(n2)+O(nlog⁡n)=O(nlog⁡n)T\left(n\right)=T\left(\frac{n}{2}\right)+O\left(n\log{n}\right)=O\left(n\log{n}\right)

手算演示

为了方便理解,这里举几个例子演示一下算法流程.

复数多项式模多项式的平方根

假设 hh 是一个不被 xx 整除(有常数项)的复数多项式,求它对模 xnx^n 的平方根.

我们有方程:

G(f(x))=f2(x)−h(x)≡0(modxn)G\left(f(x)\right) = f^2(x)-h(x) \equiv 0\pmod{x^{n}}

Taylor 展开 GG 得到下式.注意这里是对 ff 的展开,所以导数都是对 ff 的偏导数,xx 在这里是当成常数算的.

G(f(x))=∑i=0+∞G(i)(f0(x))i!(f(x)−f0(x))i=G(f0(x))+2f0(x)(f(x)−f0(x))+(f(x)−f0(x))2G(f(x)) = \sum_{i=0}^{+\infty}\frac{G^{\left(i\right)}\left(f_{0}(x)\right)}{i!}\left(f(x)-f_{0}(x)\right)^{i} = G(f_0(x)) + 2f_0(x)(f(x)-f_0(x)) + (f(x)-f_0(x))^2

用倍增计算.假设倍增中的中间结果是 f0(x),f1(x),…,fj(x)f_0(x), f_1(x), \ldots, f_j(x),或者数学严谨地说 fj(x)f_j(x) 是满足 G(fj(x))≡0(modx2j)G(f_j(x))\equiv 0\pmod{x^{2^j}} 的一个复数多项式,且为了唯一性它同时满足以下两个条件:

  • fj(x)f_{j}(x) 次数不超过 x2jx^{2^j};
  • fj+k(x)−fj(x)≡0(modx2j)f_{j+k}(x)-f_j(x)\equiv 0\pmod{x^{2^j}},对所有 kk.

把 fj+1(x)f_{j+1}(x) 和 fj(x)f_j(x) 代入上面的式子就有:

G(fj+1(x))=G(fj(x))+2fj(fj+1(x)−fj(x))+(fj+1(x)−fj(x))2≡0(modx2j+1)G(f_{j+1}(x)) = G(f_j(x)) + 2f_j(f_{j+1}(x)-f_j(x)) + (f_{j+1}(x)-f_{j}(x))^2 \equiv 0 \pmod{x^{2^{j+1}}}

显然 fj+1(x)−fj(x)f_{j+1}(x)-f_j(x) 必然是 x2jx^{2^j} 的倍数.于是得到

fj+1(x)≡fj(x)−fj2(x)−h(x)2fj(x)≡fj(x)2+h(x)2fj(x)(modx2j+1)f_{j+1}(x) \equiv f_j(x) - \frac{f_j^2(x)-h(x)}{2f_j(x)} \equiv \frac{f_j(x)^2 + h(x)}{2f_j(x)} \pmod{x^{2^{j+1}}}

如果 fj(x)f_j(x) 存在,那么 2fj(x)2f_j(x) 不被 xx 整除(有常数项),所以必然有模 x2j+1x^{2^{j+1}} 的逆元.因此数列 f0,f1…,fjf_0,f_1\ldots,f_j 存在当且仅当 f0f_0 存在.不被 xx 整除的复数多项式 h(x)h(x) 模 xx 的平方根是一定存在的,因为 h(x)h(x) 模掉 xx 就是个普通非零复数,一定有两个平方根.所以可以对所有有常数项的 h(x)h(x) 用这个算法.

选 h(x)=x+1h(x)=x+1 举例计算如下:

  • f0(x)=1f_0(x)=1,f1(x)=12+x+12×1mod  x2=12x+1f_1(x)=\dfrac{1^2+x+1}{2\times 1}\mod x^2 = \dfrac{1}{2}x+1,f2(x)=(12x+1)2+x+12×(12x+1)mod  x4=116x3−18x2+12x+1f_2(x)=\dfrac{\left(\dfrac{1}{2}x+1\right)^2+x+1}{2\times \left(\dfrac{1}{2}x+1\right)}\mod x^4 = \dfrac{1}{16}x^3-\dfrac{1}{8}x^2+\dfrac{1}{2}x+1,…\ldots
  • f0(x)=−1f_0(x)=-1,f1(x)=(−1)2+x+12×(−1)mod  x2=−12x−1f_1(x)=\dfrac{(-1)^2+x+1}{2\times (-1)}\mod x^2 = -\dfrac{1}{2}x-1,…\ldots(等于前一个取负)

可以验证两个都是正确的模平方根多项式列.

整数模素数幂的平方根

牛顿迭代算法还可以迁移到整数模素数的幂的情况. 假设 hh 是一个不被 3 整除的「方便的」整数.(「方便」指「必然有解」,具体条件后文再言)假设要算 hh 在模 3n3^n 意义下的平方根 ff.有方程:

G(f)=f2−h≡0(mod3n)G\left(f\right) = f^2-h \equiv 0\pmod{3^{n}}

Taylor 展开 GG 得到:

G(f)=∑i=0+∞G(i)(f0)i!(f−f0)i=G(f0)+2f0(f−f0)+(f−f0)2G(f) = \sum_{i=0}^{+\infty}\frac{G^{\left(i\right)}\left(f_{0}\right)}{i!}\left(f-f_{0}\right)^{i} = G(f_0) + 2f_0(f-f_0) + (f-f_0)^2

用倍增计算.假设倍增中得到的中间结果是 f0,f1,…,fjf_0, f_1, \ldots, f_j,或者严谨地说 fjf_j 是满足 G(fj)≡0(mod32j)G(f_j)\equiv 0\pmod{3^{2^j}} 的一个整数,且为了唯一性它同时满足以下两个条件:

  • 0<fj<32j0 < f_{j} < 3^{2^j};
  • fj+k−fj≡0(mod32j)f_{j+k}-f_j\equiv 0\pmod{3^{2^j}},对所有 kk.

把 fj+1f_{j+1} 和 fjf_j 代入上面的式子就有:

G(fj+1)=G(fj)+2fj(fj+1−fj)+(fj+1−fj)2≡0(mod32j+1)G(f_{j+1}) = G(f_j) + 2f_j(f_{j+1}-f_j) + (f_{j+1}-f_{j})^2 \equiv 0 \pmod{3^{2^{j+1}}}

显然 fj+1−fjf_{j+1}-f_j 必然是 32j3^{2^j} 的倍数.于是得到

fj+1≡fj−fj2−h2fj≡fj2+h2fj(mod32j+1)f_{j+1} \equiv f_j - \frac{f_j^2-h}{2f_j} \equiv \frac{f_j^2 + h}{2f_j} \pmod{3^{2^{j+1}}}

如果 fjf_j 存在,那么 2fj2f_j 不被 33 整除,所以必然有模 32j+13^{2^{j+1}} 的逆元.因此数列 f0,f1…,fjf_0,f_1\ldots,f_j 存在当且仅当 f0f_0 存在.不被 3 整除的整数 hh 模 33 的平方根要么不存在,要么有两个.所以 hh 有模 33 平方根就是整个算法能跑的唯一条件.

这里选 h=46h=46 实际计算示例.

  • f0=1f_0=1,f1=12+462×1mod  9=1f_1=\dfrac{1^2+46}{2\times 1}\mod 9 = 1,f2=12+462×1mod  81=64f_2=\dfrac{1^2+46}{2\times 1}\mod 81 = 64,f3=642+462×64mod  6561=955f_3=\dfrac{64^2+46}{2\times 64}\mod 6561 = 955,…\ldots
  • f0=2f_0=2,f1=22+462×2mod  9=8f_1=\dfrac{2^2+46}{2\times 2}\mod 9 = 8,f2=82+462×8mod  81=17f_2=\dfrac{8^2+46}{2\times 8}\mod 81 = 17,f3=172+462×17mod  6561=5606f_3=\dfrac{17^2+46}{2\times 17}\mod 6561 = 5606,…\ldots(等于前一个取负)

可以验证一下两个都是正确的模平方根数列.

代数证明

这一节对前文进行引申,用抽象代数的语言证明只要 ff 满足初始解条件,牛顿法对所有的 nn 都能给出解,并且可以得到全部的解.

有解的证明

引理 1

设 整环 RR 有多项式或 形式幂级数 f(X)=∑i≥0aiXif(X) = \sum_{i\geq 0}a_iX^i 和 r,p∈Rr,p\in R 使得 f(r)∈Rpf(r)\in Rp(亦即 rr 是 f(X)f(X) 在模 pp 意义下的根)且 f′(r)∈Rf'(r)\in R 在模 pp 意义下是可逆的.这里 f′(X):=∑i≥0(i+1)ai+1Xif'(X) := \sum_{i\geq 0}(i+1)a_{i+1}X^i 是 f(X)f(X) 的 形式导数.那么 f(r−f(r)f′(r))≡0(modp2)f\left(r-\dfrac{f(r)}{f'(r)}\right) \equiv 0\pmod {p^2}.

证明

对所有 s∈Rs\in R,

f(r+sp)=∑i≥0ai(r+sp)i=∑i≥0airi+sp∑i≥1iairi−1+s2p2(…)=f(r)+spf′(r)+s2p2(f′′(r)2!+⋯ ),\begin{aligned} f(r+sp) &= \sum_{i\geq 0}a_i(r+sp)^i \\ &= \sum_{i\geq 0}a_ir^i + sp\sum_{i\geq 1}ia_ir^{i-1} + s^2p^2\left(\ldots\right) \\ &= f(r) + spf'(r) + s^2p^2\left(\frac{f''(r)}{2!} + \cdots\right), \end{aligned}

所以

f(r+sp)∈Rp2  ⟺  f(r)+f′(r)sp∈Rp2f(r+sp) \in Rp^2 \iff f(r)+f'(r)sp \in Rp^2

因为 f(r)∈Rpf(r)\in Rp,且 f′(r)f'(r) 可逆,所以取 sp=−f(r)f′(r)sp = -\dfrac{f(r)}{f'(r)} 即可,这里 1f′(r)\dfrac{1}{f'(r)} 是模 p2p^2 意义下的逆元.因为 f′(r)f'(r) 在模 pp 意义下可逆,所以它在模 p2p^2 意义下也必定存在逆元:设有 a,b,c∈Ra,b,c\in R 使 af′(r)=bp+1af'(r) = bp+1 和 f(r)=cpf(r)=cp,那么 (a2f′(r)−2)f′(r)=b2p2+1\left(a^2f'(r)-2\right)f'(r) = b^2p^2+1,故可以取 s=c(2−a2f′(r))s=c(2-a^2f'(r)).

对于域 kk 上的多项式环 k[X]k[X],设有 G(X,Y)∈k[X,Y]G(X, Y)\in k[X, Y] 和 fn∈k[X]f_n\in k[X] 使 G(X,fn(X))∈k[X]XnG(X, f_n(X))\in k[X]X^n,那么应用引理 1 就可得到

G(X,fn(X)−G(X,fn(X))∂G∂Y(X,fn(X)))≡0(modX2n)G\left(X, f_n(X) - \frac{G(X, f_n(X))}{\frac{\partial G}{\partial Y}(X, f_n(X))} \right)\equiv 0 \pmod {X^{2n}}

而倍增的初始条件只要有 f1∈kf_1\in k 使得 G(X,f1)≡0(modX)G(X, f_1)\equiv 0\pmod X 和 ∂G∂Y(X,f1)≢0(modX)\dfrac{\partial G}{\partial Y}(X, f_1)\not\equiv 0\pmod X.后一个条件保证了 ∂G∂Y\dfrac{\partial G}{\partial Y} 有非零常数项,同时因为 X∣G(X,fn(X))∂G∂Y(X,fn(X))X\left| \dfrac{G(X, f_n(X))}{\frac{\partial G}{\partial Y}(X, f_n(X))} \right.,故而对所有的 nn,∂G∂Y(X,fn)\dfrac{\partial G}{\partial Y}(X, f_n) 总是模 XnX^n 意义下可逆的,也就满足了下一次迭代的条件.

得到全部解的证明

引理 2

若 RR 为 UFD,f,r,pf,r,p 定义同引理 1.则引理 1 给出的 r−f(r)f′(r)r-\dfrac{f(r)}{f'(r)} 是模 p2p^{2} 意义下唯一满足以下两条件的 xx 的值:

  • f(x)∈Rp2f(x)\in Rp^{2}
  • x−r∈Rpx-r\in Rp

亦即

∀x∈R,p2∣f(x)∧p∣(x−r)  ⟹  x≡r−f(r)f′(r)(modp2)\forall x\in R,\qquad p^2\mid f(x)\wedge p\mid (x-r) \implies x\equiv r-\dfrac{f(r)}{f'(r)} \pmod {p^2}
证明

令 s=−f(r)f′(r)ps = -\dfrac{f(r)}{f'(r)p} 和 u=r+spu = r+sp,引理 1 保证 uu 满足两个条件,且 f(r)+f′(r)sp∈Rp2f(r) + f'(r)sp \in Rp^{2}. 设 vv 是满足上述条件的值,则有 v=r+tpv = r+tp 和 f(r)+f′(r)tp∈Rp2f(r) + f'(r)tp \in Rp^{2}. 于是有 f′(r)(t−s)p∈Rp2f'(r)(t-s)p\in Rp^{2} 和 v−u∈Rp2v-u\in Rp^{2}.

牛顿法可以保证得到模 X2nX^{2^n} 的全部解.假设 G(X,h)≡0(modX2n)G(X, h)\equiv 0\pmod {X^{2^n}},那么令 h2i:=h(modX2i)h_{2^i} := h\pmod {X^{2^i}},然后取 f1=h1f_1 = h_1 并用牛顿法,根据引理 2 可得 f2i≡h2i(modX2i)f_{2^i} \equiv h_{2^i}\pmod {X^{2^i}},所以一定有 f2n=hf_{2^n} = h.

上面的论证也说明了,在 ∂G∂y(0,y)\dfrac{\partial G}{\partial y}(0, y) 永远可逆时,G(X,f)≡0(modXn)G(X, f)\equiv 0\pmod {X^n} 的解的个数等于 G(0,f)≡0(modX)G(0, f)\equiv 0\pmod X 的解的个数.这个结论并非平凡.请看下面的例子.

牛顿法无效时解的个数随次数而变多的例子

模 XX 意义下 X2X^2 的平方根只有 00,但是模 X4X^4 意义下 X2X^2 的平方根有 X,−X,X3+X,…X, -X, X^3+X, \ldots.

本文内容来自 OI Wiki《多项式牛顿迭代》;作者与修订历史见原文。AIOJ 对排版、站内链接和资源路径做了格式适配。

内容依 CC BY-SA 4.0 与 SATA 许可共享;附加条款以原文声明为准。感谢 OI Wiki 作者与贡献者;请按 SATA 约定前往 OI Wiki 仓库 Star 支持。