一、欧几里得算法(辗转相除法)
首先我们要知道,初等数论中:
若**(p,q)=1**,存在整数s,t 使得sp+ tq=1
(x,y)表示x,y的最大公约数
[x,y]表示x,y的最小公倍数
先举个例子:
gcd(6,16) = gcd(16,6%16) = gcd(16, 6) = gcd(6,16%6) = gcd(6,4) = gcd(4,6%4) = gcd(4,2) = gcd(2,4%2) = gcd(2,0) = 2
由此可得,欧几里得算法就是计算最大公约数的一个算法,递推式为:
gcd(a , b)= gcd(b , a % b)
直到a,b中有一个数为0,可得到他们的最大公约数。
怎么证明上述式子呢?`
令(A , B) = R,证明(A,B) = (B,R)
因 (A , B) = R
则 A = B * q + R
令 A = a * u, B = b * u
则 R = A - B*q = a * u - b * u * q = u * (a - b * q) = u * n
即 u 是A,B的最大公因数则 u 是B,R的最大公因数.
同理:
令 B = b’ * v
令 R = r’ * v
得 A = a’ * v
即 v 是B,R的最大公因数则 v 是A,B的最大公因数.
证明成立
二、裴蜀定理
内容如下:
若 a , b 是 任 意 整 数 , 且 g c d ( a , b ) = d , 那 么 对 于 任 意 的 整 数 x , y , a x + b y 都 一 定 是 d 的 倍 数 , 特 别 地 , 一 定 存 在 整 数 x , y , 使 a x + b y = d 成 立 。 若 a , b 是任意整数,且 gcd ( a , b ) = d ,那么对于任意的整数 x , y , a x + b y 都一定是 d 的倍数,特别地,一定存在整数 x , y,使 a x + b y = d 成立。 若a,b是任意整数,且gcd(a,b)=d,那么对于任意的整数x,y,ax+by都一定是d的倍数,特别地,一定存在整数x,y,使ax+by=d成立。
推论:
a,b互质的充分必要条件是存在整数x,y使ax+by=1
假设我们有一个关于x和y的线性方程 ax + by = d,现在要求判断这个方程是否存在整数解。
裴蜀定理告诉我们,ax + by = d存在整数解当前仅当 gcd ( a , b ) | d 。例如说3x + 6y = 2就不存在整数解,3x + 6y = 3就存在整数解 x =1 ,y = 0 。
如何证明?显然gcd( a , b) | (ax + by),如果存在整数的解的话必然有 gcd ( a , b ) | d 。
那么在gcd( a, b) | d的时候一定都存在整数解吗?我们来通过扩展欧几里得算法解释这个问题,实际上,扩展欧几里得的算法不仅回答了这个问题,而且还把解构造了出来。
n个整数的裴蜀定理:
设a1,a2,a3…an为n个整数,d是它们的最大公约数,那么存在整数x1…xn使得x1a1+x2a2+…xn*an=d。
特别来说,如果a1…an存在任意两个数是互质的(不必满足两两互质),那么存在整数x1…xn使得x1 * a1 + x2 * a2+…xn * an=1。证法类似两个数的情况。
三、扩展欧几里得算法
定义:在得到整数a,b的最大公约数后,还希望得到整数x,y:使得ax+by=gcd(a,b)
1、
对
于
整
数
a
>
b
显
然
b
=
0
时
,
g
c
d
(
a
,
b
)
=
a
;
此
时
x
=
1
,
y
=
0
对于整数a>b显然b=0时,gcd(a,b) = a;此时x=1,y=0
对于整数a>b显然b=0时,gcd(a,b)=a;此时x=1,y=0
2、
设
a
x
1
+
b
y
1
=
g
c
d
(
a
,
b
)
设ax_1+by_1=gcd(a,b)
设ax1+by1=gcd(a,b)
3、
有
b
x
2
+
(
a
%
b
)
y
2
=
g
c
d
(
b
,
a
%
b
)
有bx_2 + (a\%b)y_2=gcd(b,a\%b)
有bx2+(a%b)y2=gcd(b,a%b)
4、
由
于
g
c
d
(
a
,
b
)
=
g
c
d
(
b
,
a
%
b
)
,
那
么
a
x
1
+
b
y
1
=
b
x
2
+
(
a
%
b
)
y
2
由于gcd(a,b)=gcd(b,a\%b),那么ax_1+by_1=bx_2+(a\%b)y_2
由于gcd(a,b)=gcd(b,a%b),那么ax1+by1=bx2+(a%b)y2
5、
即
a
x
1
+
b
y
1
=
b
x
2
+
(
a
−
[
a
/
b
]
∗
b
)
y
2
=
a
y
2
+
b
x
2
−
[
a
/
b
]
∗
b
y
2
即ax_1+by_1=bx_2+(a-[a/b]*b)y_2=ay_2+bx_2-[a/b]*by_2
即ax1+by1=bx2+(a−[a/b]∗b)y2=ay2+bx2−[a/b]∗by2
6、
也
就
是
a
x
1
+
b
y
1
=
=
a
y
2
+
b
(
x
2
−
[
a
/
b
]
∗
y
2
)
也就是ax_1+by_1==ay_2+b(x_2-[a/b]*y_2)
也就是ax1+by1==ay2+b(x2−[a/b]∗y2)
7、
根
据
恒
等
式
定
理
得
:
x
1
=
y
2
;
y
1
=
x
2
−
[
a
/
b
]
∗
y
2
根据恒等式定理得:x_1=y_2;y_1=x_2-[a/b]*y_2
根据恒等式定理得:x1=y2;y1=x2−[a/b]∗y2
8、
这
样
我
们
就
得
到
了
求
解
x
1
,
y
1
的
方
法
:
x
1
,
y
1
的
值
基
于
x
2
,
y
2
这样我们就得到了求解x_1,y_1的方法:x_1,y_1的值基于x_2,y_2
这样我们就得到了求解x1,y1的方法:x1,y1的值基于x2,y2
9、上面的思想是以递归定义的,gcd不断地递归求解一定会有b=0的时候,所以递归可以结束.
四、扩展欧几里得算法求逆元
乘法逆元的定义:
若存在正整数a,b,p, 满足ab = 1(mod p), 则称a 是b 的乘法逆元, 或称b 是a 的乘法逆元。b ≡ a-1 (mod p),a ≡ b-1 (mod p)
比如说, 在模7 意义下,3 的乘法逆元是5, 也可以说模7 意义下5的乘法逆元是3。模13意义下5的逆元是8……
存在性
看起来和同余方程很相似(其实下面真的可以用exgcd求的!),在同余方程中, ab ≡ 1(mod p)
若a 与p 互质, 则一定存在一个正整数解b, 满足b < p,若a 与p 不互质, 则一定不存在正整数解b. 所以逆元要求a与p互质
定义函数
def exgcd(a,b):
global x,y
if b == 0:
x,y = 1,0
return a
d = exgcd(b,a%b)
x,y = y,x
y -= (a//b)*x
return d
用gmpy2库
gmpy2.is_prime(n) # 判断n是不是素数
gmpy2.gcd(a,b) # 欧几里得算法
gmpy2.gcdext(a,b) # 扩展欧几里得算法
说明:gmpy2.gcdext(a,b)返回一个元组(g,s,t) , g =gcd(a,b) , g=as+bt
import gmpy2
print(*map(int,gmpy2.gcdext(15,10)))
输出:
5 1 -1 #即gcd(15,10) = 5 = 1*15+(-1)*10
五、中国剩余定理
例1、求满足下同余方程组的x
{
x
≡
1
(
m
o
d
3
)
x
≡
3
(
m
o
d
5
)
x
≡
5
(
m
o
d
7
)
\begin{cases} x \equiv 1(mod\quad3)\\ x \equiv 3(mod\quad5)\\ x \equiv 5(mod\quad7)\\ \end{cases}
⎩⎪⎨⎪⎧x≡1(mod3)x≡3(mod5)x≡5(mod7)
⇕
\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\Updownarrow
⇕等价
{
x
≡
−
2
(
m
o
d
3
)
x
≡
−
2
(
m
o
d
5
)
x
≡
−
2
(
m
o
d
7
)
\begin{cases} x \equiv -2(mod\quad3)\\ x \equiv -2(mod\quad5)\\ x \equiv -2(mod\quad7)\\ \end{cases}
⎩⎪⎨⎪⎧x≡−2(mod3)x≡−2(mod5)x≡−2(mod7)
则可以知道 3,5,7|(x+2)
⇒
\Rightarrow
⇒ 105|(x+2)
⇒
\Rightarrow
⇒ x=103、103+105
⋯
\cdots
⋯ 即x最小值为103。
例2、
{
k
≡
a
(
m
o
d
3
)
k
≡
b
(
m
o
d
5
)
k
≡
c
(
m
o
d
7
)
(2.1)
\begin{cases} k \equiv a(mod\quad3)\\ k \equiv b(mod\quad5) \tag{2.1} \\ k \equiv c(mod\quad7)\\ \end{cases}
⎩⎪⎨⎪⎧k≡a(mod3)k≡b(mod5)k≡c(mod7)(2.1)
可以先寻找满足以下条件的方程组 (可看作基础解系)
{
x
≡
1
(
m
o
d
3
)
x
≡
0
(
m
o
d
5
)
x
≡
0
(
m
o
d
7
)
\begin{cases} x \equiv 1(mod\quad3)\\ x \equiv 0(mod\quad5)\\ x \equiv 0(mod\quad7) \end{cases}
⎩⎪⎨⎪⎧x≡1(mod3)x≡0(mod5)x≡0(mod7)
{
y
≡
0
(
m
o
d
3
)
y
≡
1
(
m
o
d
5
)
y
≡
0
(
m
o
d
7
)
\begin{cases} y \equiv 0(mod\quad3)\\ y \equiv 1(mod\quad5)\\ y \equiv 0(mod\quad7) \end{cases}
⎩⎪⎨⎪⎧y≡0(mod3)y≡1(mod5)y≡0(mod7)
{
z
≡
0
(
m
o
d
3
)
z
≡
0
(
m
o
d
5
)
z
≡
1
(
m
o
d
7
)
\begin{cases} z \equiv 0(mod\quad3)\\ z \equiv 0(mod\quad5)\\ z \equiv 1(mod\quad7) \end{cases}
⎩⎪⎨⎪⎧z≡0(mod3)z≡0(mod5)z≡1(mod7)
则式2.1记为
k
=
a
x
+
b
y
+
c
z
(
m
o
d
105
)
k=ax+by+cz \quad (mod \quad 105)
k=ax+by+cz(mod105)
由基础解系知道:35|x,21|y,15|z
x = 35t-1=3s,即存在整数s,t使得35t-3s=1
可以得出一组解:
{
x
=
70
y
=
21
z
=
15
\begin{cases} x=70\\ y=21\\ z=15\\ \end{cases}
⎩⎪⎨⎪⎧x=70y=21z=15
即:
k
=
70
a
+
21
b
+
15
c
(
m
o
d
105
)
k=70a+21b+15c(mod \quad105)
k=70a+21b+15c(mod105),余数就是k的值
根据上述再求例1
x
=
70
∗
1
+
21
∗
3
+
15
∗
5
(
m
o
d
105
)
=
208
m
o
d
105
=
103
。
x = 70*1+21*3+15*5 (mod\quad105)=208\quad mod \quad105=103。
x=70∗1+21∗3+15∗5(mod105)=208mod105=103。
用现代数学的语言说明,中国剩余定理给出了以下的一元线性同余方程组:
(
S
)
:
{
x
≡
a
1
(
m
o
d
m
1
)
x
≡
a
2
(
m
o
d
m
2
)
⋮
x
≡
a
n
(
m
o
d
m
n
)
(S):\begin{cases} x \equiv a_1 \quad (mod \quad m_1)\\ x \equiv a_2 \quad (mod \quad m_2)\\ \quad\vdots \\ x \equiv a_n \quad (mod \quad m_n)\\ \end{cases}
(S):⎩⎪⎪⎪⎪⎨⎪⎪⎪⎪⎧x≡a1(modm1)x≡a2(modm2)⋮x≡an(modmn)
整数
m
1
m_1
m1、
m
2
m_2
m2
⋯
\cdots
⋯
m
n
m_n
mn两两互质,则对任意的整数
a
1
a_1
a1、
a
2
a_2
a2
⋯
\cdots
⋯
a
n
a_n
an方程组(S)有解;
x
=
a
1
y
1
+
a
2
y
2
+
⋯
+
a
n
y
n
(
m
o
d
m
1
m
2
⋯
m
n
)
x=a_1y_1+a_2y_2+\cdots+a_ny_n(mod \quad m_1m_2\cdots m_n)
x=a1y1+a2y2+⋯+anyn(modm1m2⋯mn)
其中
y
1
=
(
m
2
⋯
m
n
)
t
−
1
=
m
1
s
y_1=(m_2\cdots m_n)t-1=m_1s
y1=(m2⋯mn)t−1=m1s,因为
(
m
1
,
m
2
⋯
m
n
)
=
1
(m1,m_2\cdots m_n)=1
(m1,m2⋯mn)=1,满足条件的整数s,t一定存在。
六、中国剩余定理拓展
中国剩余定理要求模数两两互质,而一般情况下我们要解决的都是 m i m_i mi不互质问题,这就推广到中国剩余定理的拓展!
从只含有两个的方程的同余方程组开始
{
x
≡
a
1
(
m
o
d
m
1
)
x
≡
a
2
(
m
o
d
m
2
)
\begin{cases} x \equiv a_1(mod \quad m_1)\\ x \equiv a_2(mod \quad m_2) \end{cases}
{x≡a1(modm1)x≡a2(modm2)
等价于
{
x
=
a
1
+
k
1
m
1
x
=
a
2
+
k
2
m
2
(1)
\begin{cases} x = a_1+k_1m_1\\ x =a_2+k_2m_2 \tag{1} \end{cases}
{x=a1+k1m1x=a2+k2m2(1)
消去x后得
a
1
+
k
1
m
1
=
a
2
+
k
2
m
2
(2)
a_1+k_1m_1=a_2+k_2m_2 \tag{2}
a1+k1m1=a2+k2m2(2)
因为模数不是两两互质,这个方程是否存在整数解就需要考虑一下。根据裴蜀定理,我们可以知道如果
g
c
d
(
m
1
,
m
2
)
∣
(
a
1
,
a
2
)
gcd(m_1,m_2)|(a_1,a_2)
gcd(m1,m2)∣(a1,a2)那么这个方程一定有整数解,否则不存在整数解。
当它存在整数解的时,假设我们已经根据扩展欧几里得算法求得一组特殊解
(
k
1
′
,
k
2
′
)
(k_1',k_2')
(k1′,k2′),将任意一个解回代到(1)中,就会得到满足同余方程组的解,记
x
0
x_0
x0。
由此推出所有解:
我们设
g
=
g
c
d
(
m
1
,
m
2
)
g = gcd(m_1,m_2)
g=gcd(m1,m2),又由于
a
2
−
a
1
a_2-a_1
a2−a1是
g
g
g的倍数,那么可以把方程(1)改写为
k
1
m
1
g
−
k
2
m
2
g
=
a
2
−
a
1
g
(3)
k_1\frac{m_1}{g}-k_2\frac{m_2}{g}=\frac{a_2-a_1}{g}\tag3
k1gm1−k2gm2=ga2−a1(3)
这样的话,
m
1
g
和
m
2
g
就
互
质
了
\frac{m_1}{g}和\frac{m_2}{g}就互质了
gm1和gm2就互质了,那么这个方程所有整数解就是
{
k
1
=
m
2
g
t
+
k
1
′
k
2
=
m
1
g
t
+
k
2
′
\begin{cases} k_1= \frac{m_2}{g}t+k_1'\\ k_2= \frac{m_1}{g}t+k_2' \end{cases}
{k1=gm2t+k1′k2=gm1t+k2′
这里的t是任意整数。
往回代入可以得到
x
=
a
1
+
k
1
m
1
=
x
0
+
m
1
m
2
g
t
=
x
0
+
l
c
m
(
m
1
,
m
2
)
t
x=a_1+k_1m_1 \\ =x_0+ \frac{m_1m_2}{g}t \\ =x_0+lcm(m_1,m_2)t
x=a1+k1m1=x0+gm1m2t=x0+lcm(m1,m2)t
这个解等价于同余方程
z
≡
x
0
(
m
o
d
l
c
m
(
m
1
,
m
2
)
)
z\equiv x_0(mod \quad lcm(m_1,m_2))
z≡x0(modlcm(m1,m2))
由此得出,更完整的中国剩余定理
对于一个同余方程组
{
x
≡
a
1
(
m
o
d
m
1
)
x
≡
a
2
(
m
o
d
m
2
)
⋮
x
≡
a
n
(
m
o
d
m
n
)
\begin{cases} x\equiv a_1(mod \quad m_1)\\ x\equiv a_2(mod \quad m_2)\\ \quad\vdots \\ x\equiv a_n(mod \quad m_n) \end{cases}
⎩⎪⎪⎪⎪⎨⎪⎪⎪⎪⎧x≡a1(modm1)x≡a2(modm2)⋮x≡an(modmn)
它等价于
x
≡
x
0
(
m
o
d
l
c
m
(
m
1
,
m
2
,
⋯
,
m
n
)
)
x\equiv x_0(mod\quad lcm(m_1,m_2,\cdots,m_n))
x≡x0(modlcm(m1,m2,⋯,mn))
这里
x
0
x_0
x0是这个方程的任意解,可以通过刚才的合并计算得到。
十三届蓝桥杯习题:
方法1(暴力): 用后面的数找规律,得等差数列,再遍历满足前面条件参考
方法2:中国剩余定理
math.lcm求最小公倍数,math.gcd求最大公约数
print(math.lcm(8,41,16))
最大公约数记为(a,b),最小公倍数记为[a,b]
import math
from functools import reduce
from sympy import rem
import gmpy2
a = [1,2,1,4,5,4,1,2,9,0,5,10,11,14,9,0,11,18,9,11,11,15,17,9,23,20,25,16,29,27,25,11,17,4,29,22,37,23,9,1,11,11,33,29,15,5,41,46]
l = []
for i in range(len(a)):
l.append(i+2)
def zhi(x):
for i in range(2,int(x**0.5+1)):
if x%i==0:
return False
else:
return True
mod_num = []
rem_num = []
for i in range(len(a)):
if zhi(l[i]):
mod_num.append(l[i])
rem_num.append(a[i])
else:
continue
print(mod_num)
print(rem_num)
M = reduce(lambda x,y:x*y,mod_num)
Mm = []
y = 0
for i in range(len(mod_num)):
t = M//mod_num[i]
_,a,_ = gmpy2.gcdext(t,mod_num[i])
Mm.append(int(a%mod_num[i]))
y += (t*Mm[i]*rem_num[i])
print(y%M)
输出:
[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47]
[1, 2, 4, 4, 0, 10, 0, 18, 15, 16, 27, 22, 1, 11, 5]
2022040920220409
方法3:中国剩余定理拓展
a = [1,2,1,4,5,4,1,2,9,0,5,10,11,14,9,0,11,18,9,11,11,15,17,9,23,20,25,16,29,27,25,11,17,4,29,22,37,23,9,1,11,11,33,29,15,5,41,46]
l = []
for i in range(len(a)):
l.append((i+2,a[i]))
def exgcd(a, b):
global k1, k2
if b == 0:
k1, k2 = 1, 0
return a
d = exgcd(b, a % b)
k1, k2 = k2, k1
k2 -= (a // b) * k1
return d
m1 = l[0][0]
b1 = l[0][1]
flag = 0
for i in range(1,len(l)):
m2 = l[i][0]
b2 = l[i][1]
k1,k2 = 0,0
d = exgcd(m1,m2)
if (b2-b1)%d:
flag = 1
break
k1 *= (b2-b1)//d
# k1' = k1 + k * (m2 // d) , k取任意整数
t = m2//d
k1 = k1%t #取最小的k1
# x = b + km
b1 = k1*m1+b1
m1 = m1//d*m2
if flag:
print(-1)
else:
print(b1%m1) #x的最小整数解
输出:
2022040920220409