扩展欧几里得算法用于解决二元一次不定方程、一次同余方程、逆元等问题
二元一次不定方程
- 二元一次不定方程的定义:ax+by=c,(a,b,c∈Z且a,b≠0)(1)ax + by = c,(a, b, c \in Z且a, b \neq 0) \tag{1}ax+by=c,(a,b,c∈Z且a,b=0)(1)
- 如果c=0c = 0c=0,则方程的整数解为:x=−bt,y=at,(t∈Z)x = -bt,y = at,(t \in Z)x=−bt,y=at,(t∈Z)
- 如果ccc为负数,则(1)式两边乘以−1-1−1,等号右边就变成正整数;因为所求的xxx,yyy可正可负,所以我们只用讨论aaa,bbb是正整数的情况。
- 设ddd为aaa和bbb的最大公因数,即 d=gcd(a,b),a=a1∗d,b=b1∗d,(a1,b1)=1d = gcd(a, b),a = a_1 * d,b = b_1 * d,(a_1, b_1) = 1d=gcd(a,b),a=a1∗d,b=b1∗d,(a1,b1)=1,则(1)式变为d∗(a1x+b1x)=c(2)d * (a_1x + b_1x) = c \tag{2}d∗(a1x+b1x)=c(2)
结论:当 d∤cd\nmid cd∤c 时,(1)式没有整数解;当 d∣cd\mid cd∣c 时,(1)式一定有整数解,这可以通过贝祖定理(又称裴蜀定理)得到。
贝祖定理:如果aaa和bbb是整数,且 gcd(a,b)=dgcd(a, b) = dgcd(a,b)=d,那么对于任意的整数x,y,ax+byx,y,ax + byx,y,ax+by都一定是ddd的倍数;特别的,一定存在整数x,yx,yx,y,使ax+by=dax + by = dax+by=d成立。
- 贝祖定理的一个重要推论是:aaa和bbb互质的充分必要条件是存在整数x,yx,yx,y,使ax+by=1ax + by = 1ax+by=1
- ax+byax + byax+by能够得到的最小结果就是它们的最小公倍数
我们通过几个例子,解释xxx和yyy的求解步骤,并在此基础上得出扩展欧几里得算法。
例1:17x+5y=117x + 5y = 117x+5y=1
∵gcd(17,5)=1,∴方程有解\because gcd(17, 5) = 1,\therefore 方程有解∵gcd(17,5)=1,∴方程有解
(3∗5+2)x+5y=1(3*5 + 2)x + 5y = 1(3∗5+2)x+5y=1,整理得:
5(3x+y)+2x=15(3x + y) + 2x = 15(3x+y)+2x=1
用换元法,如果我们求出如下方程的解,
5x′+2y′=15x^{'} + 2y^{'} = 15x′+2y′=1,则:x=y′,y=x′−3y′x = y^{'},y = x^{'} - 3y^{'}x=y′,y=x′−3y′
(2∗2+1)x′+2y′=1(2*2 + 1)x^{'} + 2y^{'} = 1(2∗2+1)x′+2y′=1,整理得:
2(2x′+y′)+x′=12(2x^{'} + y^{'}) + x^{'} = 12(2x′+y′)+x′=1
用换元法,如果我们求出如下方程的解,
2x′′+y′′=12x^{''} + y^{''} = 12x′′+y′′=1,则:x′=y′′,y′=x′′−2y′′x^{'} = y^{''},y^{'} = x^{''} - 2y^{''}x′=y′′,y′=x′′−2y′′
(2∗1+0)x′′+y′′=1(2*1 + 0)x^{''} + y^{''} = 1(2∗1+0)x′′+y′′=1,整理得:
1(2x′′+y′′)+0x′′=11(2x^{''} + y^{''}) + 0x^{''} = 11(2x′′+y′′)+0x′′=1
用换元法,如果我们求出如下方程的解,
1x′′′+0y′′′=11x^{'''} + 0y^{'''} = 11x′′′+0y′′′=1,则:x′′=y′′′,y′′=x′′′−2y′′′x^{''} = y^{'''},y^{''} = x^{'''} - 2y^{'''}x′′=y′′′,y′′=x′′′−2y′′′
显而易见,方程 1x′′′+0y′′′=11x^{'''} + 0y^{'''} = 11x′′′+0y′′′=1 有解x′′′=1,y′′′=0x^{'''} = 1,y^{'''} = 0x′′′=1,y′′′=0。
x′′′x^{'''}x′′′的系数111,就是原方程中171717和555的最大公约数
回溯可得,
x′′=0,y′′=1x^{''} = 0,y^{''} = 1x′′=0,y′′=1
x′=1,y′=−2x^{'} = 1,y^{'} = -2x′=1,y′=−2
x=−2,y=7x = -2,y = 7x=−2,y=7
例2:27x+21y=327x + 21y = 327x+21y=3
∵gcd(27,21)=3,∴方程有解\because gcd(27, 21) = 3,\therefore 方程有解∵gcd(27,21)=3,∴方程有解
(1∗21+6)x+21y=3(1*21 + 6)x + 21y = 3(1∗21+6)x+21y=3,整理得:
21(x+y)+6x=321(x + y) + 6x = 321(x+y)+6x=3
用换元法,如果我们求出如下方程的解,
21x′+6y′=321x^{'} + 6y^{'} = 321x′+6y′=3,则:x=y′,y=x′−y′x = y^{'},y = x^{'} - y^{'}x=y′,y=x′−y′
(3∗6+3)x′+6y′=3(3*6 + 3)x^{'} + 6y^{'} = 3(3∗6+3)x′+6y′=3,整理得:
6(3x′+y′)+3x′=36(3x^{'} + y^{'}) + 3x^{'} = 36(3x′+y′)+3x′=3
用换元法,如果我们求出如下方程的解,
6x′′+3y′′=36x^{''} + 3y^{''} = 36x′′+3y′′=3,则:x′=y′′,y′=x′′−3y′′x^{'} = y^{''},y^{'} = x^{''} - 3y^{''}x′=y′′,y′=x′′−3y′′
(2∗3+0)x′′+3y′′=3(2*3 + 0)x^{''} + 3y^{''} = 3(2∗3+0)x′′+3y′′=3,整理得:
3(2x′′+y′′)+0x′′=33(2x^{''} + y^{''}) + 0x^{''} = 33(2x′′+y′′)+0x′′=3
用换元法,如果我们求出如下方程的解,
3x′′′+0y′′′=33x^{'''} + 0y^{'''} = 33x′′′+0y′′′=3,则:x′′=y′′′,y′′=x′′′−2y′′′x^{''} = y^{'''},y^{''} = x^{'''} - 2y^{'''}x′′=y′′′,y′′=x′′′−2y′′′
显而易见,方程 3x′′′+0y′′′=33x^{'''} + 0y^{'''} = 33x′′′+0y′′′=3 有解x′′′=1,y′′′=0x^{'''} = 1,y^{'''} = 0x′′′=1,y′′′=0。
x′′′x^{'''}x′′′的系数333,就是原方程中272727和212121的最大公约数
回溯可得,
x′′=0,y′′=1x^{''} = 0,y^{''} = 1x′′=0,y′′=1
x′=1,y′=−3x^{'} = 1,y^{'} = -3x′=1,y′=−3
x=−3,y=4x = -3,y = 4x=−3,y=4
对于方程:ax+by=gcd(a,b)ax + by = gcd(a, b)ax+by=gcd(a,b)
(⌊a/b⌋∗b+a%b)x+by=gcd(a,b)(\lfloor a/b \rfloor * b + a \% b)x + by = gcd(a, b)(⌊a/b⌋∗b+a%b)x+by=gcd(a,b),整理得:
b(⌊a/b⌋∗x+y)+a%b∗x=gcd(a,b)b(\lfloor a/b \rfloor * x + y) + a \% b * x = gcd(a, b)b(⌊a/b⌋∗x+y)+a%b∗x=gcd(a,b)
用换元法,如果我们求出如下方程的解,
b∗x′+a%b∗y′=gcd(a,b)b * x^{'} + a \% b * y^{'} = gcd(a, b)b∗x′+a%b∗y′=gcd(a,b),则:x=y′,y=x′−⌊a/b⌋∗y′x = y^{'},y = x^{'} - \lfloor a/b \rfloor * y^{'}x=y′,y=x′−⌊a/b⌋∗y′
不断递归,最后我们得到方程:
gcd(a,b)∗x~+0∗y~=gcd(a,b)gcd(a, b) * \tilde x + 0 * \tilde y = gcd(a, b)gcd(a,b)∗x~+0∗y~=gcd(a,b),
可得该方程的一组解:x~=1,y~=0\tilde x = 1,\tilde y = 0x~=1,y~=0。
x~\tilde xx~的系数,就是原方程中aaa和bbb的最大公约数
不断回溯,最后得到原方程的解
扩展欧几里得算法的参考代码:
//二元一次不定方程:ax + by = gcd(a, b)
//函数返回 gcd(a, b);通过引用传递,返回一组解 (x, y)
int exgcd(int a, int b, int &x, int &y) {
if (b == 0) {
x = 1, y = 0;
return a;
}
int gcd = exgcd(b, a % b, x, y);
int t = x;
x = y;
y = t - a / b * y;
return gcd;
}
总结:
对于一般情况,ax+by=cax + by = cax+by=c,(a,b,ca, b, ca,b,c为正整数),设 d=gcd(a,b)d = gcd(a, b)d=gcd(a,b)
- 如果d∤cd \nmid cd∤c,方程无解
- 如果d∣cd \mid cd∣c,
- 先求出方程ax+by=dax + by = dax+by=d的一组解 x0,y0x_0, y_0x0,y0
- 原方程的一组解为 x′=x0∗c/d,y′=y0∗c/dx^{'} = x_0 * c/d,y^{'} = y_0 * c/dx′=x0∗c/d,y′=y0∗c/d
- 原方程的通解为 x=x′+k∗b/d,y=y′−k∗a/dx = x^{'} + k * b / d,y = y^{'} - k * a / dx=x′+k∗b/d,y=y′−k∗a/d,kkk为任意整数
- 设p=b/dp = b / dp=b/d,xxx的最小非负整数解为 (x′% p+p) % p(x^{'} \% \ p + p) \ \% \ p(x′% p+p) % p
- 设q=a/dq = a / dq=a/d,yyy的最小非负整数解为 (y′% q+q) % q(y^{'} \% \ q + q) \ \% \ q(y′% q+q) % q
- 如果aaa或bbb是负数,则函数返回的最大公约数可能为负数,将ppp取绝对值,才能用以上方法得到xxx的最小非负整数解
求解一次同余方程 ax≡b (mod m)ax \equiv b \ (mod \ m)ax≡b (mod m)
- 一次同余方程亦称线性同余方程,指未知数仅出现一次幂的同余方程。
- 该同余方程等价于 ax=m∗(−y)+bax = m * (-y) + bax=m∗(−y)+b,即ax+my=bax + my = bax+my=b
- 如果gcd(a,m)∤bgcd(a, m) \nmid bgcd(a,m)∤b,该方程无解
- 如果gcd(a,m)∣bgcd(a, m) \mid bgcd(a,m)∣b,用上述方法求解
求解乘法逆元
- 对于正整数a,ma, ma,m,如果有ax≡1 (mod m)ax \equiv 1 \ (mod \ m)ax≡1 (mod m),则称xxx的最小正整数解为aaa模mmm的逆元,也叫做数论倒数,记作a−1a^{-1}a−1。aaa和xxx一定都和mmm互质,如果aaa和mmm不互质就不存在逆元。
- 乘法运算中a∗a−1=1a * a^{-1} = 1a∗a−1=1,模运算中,a∗a−1≡1 (mod m)a * a^{-1} \equiv 1 \ (mod \ m)a∗a−1≡1 (mod m)
扩展欧几里得算法求逆元:
- 见求解一次同余方程
费马小定理求逆元:
- 如果模数mmm为质数,可利用费马小定理求逆元:am−1≡1 (mod m)a^{m-1} \equiv 1 \ (mod \ m)am−1≡1 (mod m)
- 因为a∗am−2≡1 (mod m)a * a^{m-2} \equiv 1 \ (mod \ m)a∗am−2≡1 (mod m),所以aaa的逆元是am−2 %ma^{m-2} \ \% mam−2 %m,用快速幂计算
- 例如,mmm为7时,3∗35≡1 (mod 7)3 * 3^5 \equiv 1 \ (mod \ 7)3∗35≡1 (mod 7),35 % 73^5 \ \% \ 735 % 7为333模777的逆元
欧拉函数求逆元:
- 欧拉定理:如果aaa和mmm互质,设ϕ(m)\phi(m)ϕ(m)为小于mmm且与mmm互质的正整数的个数,则有aϕ(m)≡1 (mod m)a^{\phi(m)} \equiv 1 \ (mod \ m)aϕ(m)≡1 (mod m)
- 因为a∗aϕ(m)−1≡1 (mod m)a * a ^ {\phi(m) - 1} \equiv 1 \ (mod \ m)a∗aϕ(m)−1≡1 (mod m),所以aaa的逆元是aϕ(m)−1 %ma^{\phi(m)-1} \ \% maϕ(m)−1 %m,用快速幂计算
- mmm可以不是质数
- 例如,mmm为8时,ϕ(8)=4\phi(8) = 4ϕ(8)=4,3∗33≡1 (mod 8)3 * 3^3 \equiv 1 \ (mod \ 8)3∗33≡1 (mod 8),33 % 83^3 \ \% \ 833 % 8为333模888的逆元
逆元的用处:
- 取余运算有如下性质
- (a+b)%m=(a%m+b%m)%m(a + b) \% m = (a \% m + b \% m) \%m(a+b)%m=(a%m+b%m)%m
- (a−b)%m=(a%m−b%m)%m(a - b) \% m = (a \% m - b \% m) \%m(a−b)%m=(a%m−b%m)%m
- (a∗b)%m=(a%m∗b%m)%m(a * b) \% m = (a \% m * b \% m) \%m(a∗b)%m=(a%m∗b%m)%m
- ab %m=(a%m)b %ma^b \ \% m = (a \% m)^b \ \%mab %m=(a%m)b %m
- 但对除法不成立
- (a / b)%m≠(a%m / b%m)%m(a \ / \ b) \% m \neq (a \% m \ / \ b \% m) \%m(a / b)%m=(a%m / b%m)%m
- 例如,(8/2)%7=4(8 / 2) \% 7 = 4(8/2)%7=4,而 (8%7)/(2%7)%7=0(8 \% 7) / (2 \% 7) \% 7 = 0(8%7)/(2%7)%7=0
- 如果b∣ab \mid ab∣a,则(a / b)%m=((a / b)∗(b∗b−1))%m=(a/b∗b)∗(b−1)%m=(a∗b−1)%m(a \ / \ b) \% m = ((a \ / \ b) * (b * b^{-1})) \% m = (a / b * b) * (b^{-1}) \% m = (a * b^{-1}) \% m(a / b)%m=((a / b)∗(b∗b−1))%m=(a/b∗b)∗(b−1)%m=(a∗b−1)%m。这样我们就把除法取余运算转化成了乘法取余运算。
- 例如,求3的幂的和