翻了各课件,发现……
某个英文网站 是这么说的:
There is a way to calculate the GCD and resultants in
O
(
n
log
2
n
)
O(n\log^2n)
O ( n log 2 n ) . To do this you should note that if you consider
a
(
x
)
=
a
0
+
x
k
a
1
a(x)=a_0+x^ka_1
a ( x ) = a 0 + x k a 1 and
b
(
x
)
=
b
0
+
x
k
b
1
b(x)=b_0+x^kb_1
b ( x ) = b 0 + x k b 1 where
k
=
min
(
deg
a
,
deg
b
)
/
2
k=\min(\deg a,\deg b)/2
k = min ( deg a , deg b ) / 2 then basically the first few operations of Euclidean algorithm on
a
(
x
)
a(x)
a ( x ) and
b
(
x
)
b(x)
b ( x ) are defined by the Euclidean algorithm on
a
1
(
x
)
a_1(x)
a 1 ( x ) and
b
1
(
x
)
b_1(x)
b 1 ( x ) for which you may also calculate GCD recursively and then somehow memorize linear transforms you made with them and apply it to
a
(
x
)
a(x)
a ( x ) and
b
(
x
)
b(x)
b ( x ) to lower the degrees of polynomials. Implementation of this algorithm seems pretty tedious and technical thus it’s not considered in this article yet.
你如果脑子一冲,傻乎乎地去直接分治一半算出线性变换,DEBUG 一下会让你冷静过来…… 你发现这个线性变换应用到整体多项式上没有任何效果……
而 Picks 的 Introduction to Polynomials 也写的比较晦涩,因此在这里稍微详细地谈谈,以正视听……
给多项式
a
,
b
∈
F
[
z
]
a, b \in \mathbb F[z]
a , b ∈ F [ z ] ,我们原本的问题是求 GCD 以及扩展欧几里得的解:
x
,
y
∈
F
[
z
]
,
a
x
+
b
y
=
gcd
(
a
,
b
)
x, y\in \mathbb F[z], ax + by = \gcd (a, b)
x , y ∈ F [ z ] , a x + b y = g cd( a , b )
由于多项式取模最坏情况下每次只让一个多项式的度数减小,这种方法的复杂度最好也就是
Θ
(
(
deg
a
)
(
deg
b
)
)
\Theta((\deg a) (\deg b))
Θ ( ( deg a ) ( deg b ) ) 的。
我们考虑求解一个中间情况叫做 HALF-GCD,设
k
=
min
(
deg
a
,
deg
b
)
k = \min(\deg a, \deg b)
k = min ( deg a , deg b ) ,这个东西让我们对于
a
,
b
a, b
a , b 能得到一个矩阵
(
x
1
y
1
x
2
y
2
)
\begin{pmatrix} x_1 & y_1 \\ x_2 & y_2 \end{pmatrix}
( x 1 x 2 y 1 y 2 ) ,其中
deg
x
≤
k
2
+
o
(
k
)
\deg x \le \frac k 2 + o(k)
deg x ≤ 2 k + o ( k ) ,且
deg
(
a
x
1
+
b
y
1
)
,
deg
(
a
x
2
+
b
y
2
)
≤
k
2
+
o
(
k
)
\deg (a x_1 + b y_1), \deg (a x_2 + b y_2) \le \frac k 2 + o(k)
deg ( a x 1 + b y 1 ) , deg ( a x 2 + b y 2 ) ≤ 2 k + o ( k ) (除非它们的
gcd
\gcd
g cd 的
deg
\deg
deg 比这个大)。也就是说这个 HALF-GCD 可以理解成是将欧几里得过程执行到一半的产物。
睿智的你想必已经意识到了这个东西如何得出真正的欧几里得的解。我们只需要调用一次这个算法,然后就规约到了一个
n
2
\frac n2
2 n 的子问题(可能需要额外支付
O
(
M
(
n
)
)
O(M(n))
O ( M ( n ) ) 时间用于调整较大的多项式),那么我们就得到了一个
Θ
(
M
(
n
)
)
+
T
(
n
)
+
T
(
n
/
2
)
+
T
(
n
/
4
)
+
⋯
\Theta(M(n)) + T(n) + T(n/2) + T(n/4) + \cdots
Θ ( M ( n ) ) + T ( n ) + T ( n / 2 ) + T ( n / 4 ) + ⋯ 的方法。
那么接下来考虑通过分治来解决 HALF-GCD:设
k
=
min
(
deg
a
,
deg
b
)
k=\min(\deg a, \deg b)
k = min ( deg a , deg b ) ,我们首先算出
a
/
z
k
/
2
,
b
/
z
k
/
2
a/z^{k/2}, b/z^{k/2}
a / z k / 2 , b / z k / 2 对应的矩阵,这样一个矩阵作用与原多项式有两部分:考虑将原来的表为
a
=
a
0
+
a
1
z
k
/
2
,
b
=
b
0
+
b
1
z
k
/
2
a=a_0 + a_1 z^{k/2}, b=b_0 + b_1 z^{k/2}
a = a 0 + a 1 z k / 2 , b = b 0 + b 1 z k / 2 ,那么
a
x
+
b
y
ax + by
a x + b y 中,
deg
(
(
a
1
x
+
b
1
y
)
z
k
/
2
)
=
3
4
k
+
o
(
k
)
\deg ((a_1 x + b_1 y)z^{k/2}) = \frac 34 k + o(k)
deg ( ( a 1 x + b 1 y ) z k / 2 ) = 4 3 k + o ( k ) ,另外
deg
(
a
0
x
+
b
0
y
)
=
3
4
k
+
o
(
k
)
\deg (a_0 x + b_0 y) = \frac 34 k + o(k)
deg ( a 0 x + b 0 y ) = 4 3 k + o ( k ) 也比较显然。接着我们将得到的两个
3
4
k
\frac 34 k
4 3 k 次多项式除掉前
z
k
/
4
z^{k/4}
z k / 4 项送进去,这样又再次保证了递归返回的又是一个
1
4
k
+
o
(
k
)
\frac 14 k + o(k)
4 1 k + o ( k ) 量级的结果,我们将两个矩阵相乘就得到了变换矩阵。而结果正好满足度数在
k
2
+
o
(
k
)
\frac k 2 + o(k)
2 k + o ( k ) 以内。时间复杂度
T
(
n
)
=
2
T
(
n
/
2
)
+
M
(
n
)
T(n) = 2T(n/2) + M(n)
T ( n ) = 2 T ( n / 2 ) + M ( n ) ,即
M
(
n
)
log
n
M(n)\log n
M ( n ) log n 。消耗了
M
(
n
)
M(n)
M ( n ) 的一方面是进行矩阵乘法,另一方面是可能在一些情况减少两个多项式次数的差,同时也可能发生剪枝。
综上所述,我们在
M
(
n
)
log
n
M(n)\log n
M ( n ) log n 的时间内完成了两个多项式的欧几里得求算,高精度整数的做法也是相似的。
这里 是一份代码,这里处理的时候是令
k
=
max
(
deg
a
,
deg
b
)
k = \max (\deg a, \deg b)
k = max ( deg a , deg b ) ,理论上是差不多的,目前尚且不知道这个代码是否存在一些 edge case 上的问题,欢迎大家指正。