模乘逆元简介
所谓模乘逆元(modular multiplicative inverse)就是对于数字a,和与它互质的一个正整数m，找出一个x，满足以下关系:
a
x
≡
1
(
m
o
d
m
)
ax \equiv 1\;(mod\;m)
a
x
≡
1
(
m
o
d
m
)
上面这个公式的意思就是ax除于m的余数是1。公式里的写法叫做同余，就是≡两边的数除于m的余数相同，而1除于其他正整数的余数一定是1。所以x是a对于m的模乘逆元，记为
a
−
1
a^{-1}
a
−
1
。这个写法其实有bug，因为里面没有m，哈哈。其实这个问题，等同于这个问题：
a
x
+
m
y
=
1
ax+my=1
a
x
+
m
y
=
1
解这个问题，需要用到扩展欧几里得法。
扩展欧几里得法
欧几里得法仅仅是求出两个数之间的最大公约数。但是有时候（比如在线性丢番图方程中）需要用这两个数表示最大公约数。也就是要求以下方程的解，公式里的g是a和b的最大公约数。
a
x
+
b
y
=
g
ax+by=g
a
x
+
b
y
=
g
所以模乘逆元不过是扩展欧几里得法中gcd为1的特殊情况而已。实现了扩展欧几里得法便可以直接拿来计算模乘逆元。
x和y，必须为整数，可以为负数。但是这是一个二元一次方程，只有一个方程的话会有N个解。如何保证只有一个解呢？先看看递归过程和例子，我再讲如何保证只有一个解。这是一个利用欧几里得法循环解决问题的过程。利用欧几里得法，把方程简化为以下形式：
b
x
1
+
(
a
m
o
d
b
)
y
1
=
g
bx_1+(a\;mod\;b)y_1=g
b
x
1
+
(
a
m
o
d
b
)
y
1
=
g
因为递归到最后，a%b会变成0，结束递归。既然是递归就要找关系，我们要找的是x,y与x1,y1之间的关系。根据取余的定义，我们可以得到一个公式：
(
a
m
o
d
b
)
=
a
−
⌊
a
b
⌋
b
(a\,mod\,b)=a-\lfloor\frac{a}{b}\rfloor b
(
a
m
o
d
b
)
=
a
−
⌊
b
a
⌋
b
需要注意的是
⌊
a
b
⌋
\lfloor\frac{a}{b}\rfloor
⌊
b
a
⌋
代表的是整数除法
把
(
a
m
o
d
b
)
=
a
−
⌊
a
b
⌋
b
(a\,mod\,b)=a-\lfloor\frac{a}{b}\rfloor b
(
a
m
o
d
b
)
=
a
−
⌊
b
a
⌋
b
代入
b
x
1
+
(
a
m
o
d
b
)
y
1
=
g
bx_1+(a\;mod\;b)y_1=g
b
x
1
+
(
a
m
o
d
b
)
y
1
=
g
得到
b
x
1
+
(
a
−
⌊
a
b
⌋
b
)
y
1
=
g
bx_1+(a-\lfloor\frac{a}{b}\rfloor b)y_1=g
b
x
1
+
(
a
−
⌊
b
a
⌋
b
)
y
1
=
g
化简，得到
a
y
1
+
b
(
x
1
−
y
1
⌊
a
b
⌋
)
=
g
ay_1+b(x_1-y_1\lfloor\frac{a}{b}\rfloor)=g
a
y
1
+
b
(
x
1
−
y
1
⌊
b
a
⌋
)
=
g
与下面公式对比
a
x
+
b
y
=
g
ax+by=g
a
x
+
b
y
=
g
我们就找到了关系
x
=
y
1
y
=
x
1
−
y
1
⌊
a
b
⌋
x=y_1\\ y= x_1-y_1\lfloor\frac{a}{b}\rfloor
x
=
y
1
y
=
x
1
−
y
1
⌊
b
a
⌋
得到递归公式之后，代码就非常容易写了啊。
实例
以a=10,b=6为例子，详细分解这一过程：
第一次化简
a
=
10
b
=
6
6
x
1
+
(
10
m
o
d
6
)
y
1
=
6
x
1
+
4
y
1
=
2
x
=
y
1
y
=
x
1
−
y
1
⌊
10
6
⌋
=
x
1
−
y
1
⋅
1
=
x
1
−
y
1
a=10\\ b=6\\ 6x_1+ (10\,mod\,6)y_1=6x_1+ 4y_1=2\\ x=y_1\\ y=x_1-y_1\lfloor\frac{10}{6}\rfloor=x_1-y_1\cdot1=x_1-y_1\\
a
=
1
0
b
=
6
6
x
1
+
(
1
0
m
o
d
6
)
y
1
=
6
x
1
+
4
y
1
=
2
x
=
y
1
y
=
x
1
−
y
1
⌊
6
1
0
⌋
=
x
1
−
y
1
⋅
1
=
x
1
−
y
1
得出问题为6*x
1
+ 2*y
1
=2,再进行第二次化简
a
1
=
6
b
1
=
4
4
x
2
+
(
6
m
o
d
4
)
y
2
=
=
4
x
2
+
2
y
2
=
2
x
1
=
y
2
y
1
=
x
2
−
y
2
⌊
6
4
⌋
=
x
2
−
y
2
⋅
1
=
x
2
−
y
2
a_1=6\\ b_1=4\\ 4x_2+(6\,mod\,4)y_2==4x_2+2y_2=2\\ x_1=y_2\\ y_1=x_2-y_2\lfloor\frac{6}{4}\rfloor=x_2-y_2\cdot1=x_2-y_2
a
1
=
6
b
1
=
4
4
x
2
+
(
6
m
o
d
4
)
y
2
=
=
4
x
2
+
2
y
2
=
2
x
1
=
y
2
y
1
=
x
2
−
y
2
⌊
4
6
⌋
=
x
2
−
y
2
⋅
1
=
x
2
−
y
2
得出问题为4x
2
+2y
2
=2,再进行第三次化简
a
2
=
4
b
2
=
2
2
x
3
+
(
4
m
o
d
2
)
y
3
=
2
x
3
+
0
y
3
=
2
x
2
=
y
3
y
2
=
x
3
−
y
3
⌊
4
2
⌋
=
x
3
−
y
3
⋅
2
a_2=4\\ b_2=2\\ 2x_3+(4\,mod\,2)y_3=2x_3+0y_3=2\\ x_2=y_3\\ y_2=x_3-y_3\lfloor\frac{4}{2}\rfloor=x_3-y_3\cdot2
a
2
=
4
b
2
=
2
2
x
3
+
(
4
m
o
d
2
)
y
3
=
2
x
3
+
0
y
3
=
2
x
2
=
y
3
y
2
=
x
3
−
y
3
⌊
2
4
⌋
=
x
3
−
y
3
⋅
2
到了这一步，出现了系数0，y
3
有无数种取值，只取y
3
=0。我们前面学过欧几里得公约数算法原理，肯定知道到了这一步，x
3
一定等于1.
∴
y
3
=
0
,
x
3
=
1
\therefore y_3 = 0,x_3=1\\
∴
y
3
=
0
,
x
3
=
1
然后按顺序代回去：
∵
x
2
=
y
3
y
2
=
x
3
−
y
3
⋅
2
∴
x
2
=
0
y
2
=
1
∵
x
1
=
y
2
y
1
=
x
2
−
y
2
∴
x
1
=
1
y
1
=
−
1
∵
x
=
y
1
y
=
x
1
−
y
1
∴
x
=
−
1
y
=
2
\because x_2=y_3\\ y_2=x_3-y_3\cdot2\\ \therefore x_2=0\\ y_2= 1\\ \because x_1=y_2\\ y_1=x_2-y_2\\ \therefore x_1 = 1\\ y_1=-1\\ \because x=y_1\\ y=x_1-y_1\\ \therefore x=-1\\ y= 2\\
∵
x
2
=
y
3
y
2
=
x
3
−
y
3
⋅
2
∴
x
2
=
0
y
2
=
1
∵
x
1
=
y
2
y
1
=
x
2
−
y
2
∴
x
1
=
1
y
1
=
−
1
∵
x
=
y
1
y
=
x
1
−
y
1
∴
x
=
−
1
y
=
2
所以结果为
−
1
⋅
10
+
2
⋅
6
=
2
-1\cdot10+2\cdot6=2
−
1
⋅
1
0
+
2
⋅
6
=
2
。可以看到结果会可能产生负数，所以需要加上N个m，使得它变成正数。
python实现
整个过程可以用python实现，算法首先是要保存(x
0
,y
0
)->(x
1
,y
1
)的递归关系。
而这个关系，只有一个未知数，那就是
⌊
a
b
⌋
\lfloor\frac{a}{b}\rfloor
⌊
b
a
⌋
所以只需要保存这个数就好了。代码如下：
# _*_ coding:utf-8 _*_
def extend_euclid(a: int, b: int):
if a < b:
temp = extend_euclid(b, a)
return reversed(temp)
# 保存a/b的系数
coefficients = []
while b != 0:
coefficients.append(a // b)
a, b = b, a % b
# 循环代入
x, y = 1, 0
coefficients.reverse()
for i in coefficients:
x, y = y, x - y * i
return x, y
def modular_inverse(a: int, m: int):
x, y = extend_euclid(a, m)
while x < 0:
x += m
return x
if __name__ == '__main__':
a, m = 3, 7
x = modular_inverse(3, 7)
print(f'{a} * {x} % {m} = {a * x % m}')
测试结果：
3 * 5 % 7 = 1