Binary GCD
代码
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23def binaryGCD(x, y):
x, y = abs(x), abs(y)
if x == 0 or y == 0:
return x + y
if x == y:
return x
cnt = 0
# this cycle is O(N^2)(assume that N = max(lgx, lgy))
while ((x & 1) | (y & 1)) == 0:
cnt += 1
x = x >> 1
y = y >> 1
# the y below is surely odd
# when x-y, x and y are odd, so x will become even
# so the x>>1 will be run every cycles
# so this cycle is O(N^2)
while x != 0:
while (x & 1) == 0:
x = x >> 1
if y > x:
x, y = y, x
x, y = x - y, y
return y * (1 << cnt)复杂度分析：设$n=max(bitlen(x), bitlen(y))$。17行开始的循环中，因为
y必然为奇数，所以x-y为偶数，所以第18行的循环在每次外层循环运行一次时都至少运行一次，所以有$O(n)$次循环，每次循环需要$O(n)$的均摊复杂度——因为第18行的每次都是$O(n)$的复杂度，然后整个算法中这种移位最多有$O(n)$次，所以，整个算法复杂度为$O(n^2)$trivial 版本的GCD如下
1
2
3
4def GCD(x, y):
while y != 0:
x, y = y, x % y
return x需要$O(n^3)$的复杂度
Extend Binary GCD
代码，引用自
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44def extendBinaryGCD(a, b):
"""Extended binary GCD.
Given input a, b the function returns s, t, d
such that gcd(a,b) = d = as + bt."""
if a == 0:
return 0, 1, b
if b == 0:
return 1, 0, a
if a == b:
return 1, 0, a
u, v, s, t, r = 1, 0, 0, 1, 0
while (a % 2 == 0) and (b % 2 == 0):
a, b, r = a // 2, b // 2, r + 1
alpha, beta = a, b
#
# from here on we maintain a = u * alpha + v * beta
# and b = s * alpha + t * beta
#
while a % 2 == 0:
# v is always even
a = a // 2
if (u % 2 == 0) and (v % 2 == 0):
u, v = u // 2, v // 2
else:
u, v = (u + beta) // 2, (v - alpha) // 2
while a != b:
if b % 2 == 0:
b = b // 2
#
# Commentary: note that here, since b is even,
# (i) if s, t are both odd then so are alpha, beta
# (ii) if s is odd and t even then alpha must be even, so beta is odd
# (iii) if t is odd and s even then beta must be even, so alpha is odd
# so for each of (i), (ii) and (iii) s + beta and t - alpha are even
#
if (s % 2 == 0) and (t % 2 == 0):
s, t = s // 2, t // 2
else:
s, t = (s + beta) // 2, (t - alpha) // 2
elif b < a:
a, b, u, v, s, t = b, a, s, t, u, v
else:
b, s, t = b - a, s - u, t - v
return s, t, (2 ** r) * a思路：从19行开始，维护式子b=s*alpha+t*beta 的成立——可以验证，每次s、t更改后，式子还是成立