【数论相关】拓展欧几里得解线性同余方程和不定方程
前置知识：欧几里得算法（辗转相除法）
裴蜀定理
裴蜀定理（或贝祖定理，Bézout's identity）得名于法国数学家艾蒂安·裴蜀，说明了对任何整数a、b和它们的最大公约数d，关于未知数x和y的线性不定方程（称为裴蜀等式）：若a,b是整数,且gcd(a,b)=d，那么对于任意的整数x,y,ax+by都一定是d的倍数，特别地，一定存在整数x,y，使ax+by=d成立。
它的一个重要推论是：a,b互质的充要条件是存在整数x,y使ax+by=1.
（来自百度百科）
相关证明可以查看百度百科的内容
百科链接
拓展欧几里得
顾名思义，就是欧几里得算法的拓展。可以在高速求出\(gcd(a,b)\)的同时，求出裴蜀定理中不定方程\(ax+by=gcd(a,b)\)的中\(x\)和\(y\)的值。
下面是推导过程：
根据欧几里得算法可以得出\(gcd(a,b)=gcd(a,a \mod b)\)
由此可以列出下面两个方程
\(ax_0+by_0=gcd(a,b)\)
\(bx_1+(a \mod b) y_1=gcd(a,b)\)
因为\(a \mod b=a-(a/b)*b\)
所以原式可以写成：
\(bx_1+(a-(a/b)*b)y_1=gcd(a,b)\)
\(bx_1+ay_1-(a/b)*by_1=gcd(a,b)\)
\(ay_1+b(x_1-a/b*y_1)=gcd(a,b)\)
结合第一条方程，得出\(x_0=y_1,y_0=(x_1-a/b*y_1)\)
这样，拓展欧几里得的递归模式就显示出来了。
那么边界是什么呢？
显然是\(b=0\)时，此时\(a=gcd(a,b)\)，显然有解\(x=1,y=0\)
练习
Luogu P1082
\(ax\equiv1\pmod {b}\)
对于这样一个同余方程，我们可以通过化为不定方程利用拓展欧几里得来求解
\(ax\equiv1\pmod {b}\)等价于\(ax\mod b=1 mod b\)
即\(ax\mod b=1 mod b\)
可以得出\(b|ax-1\)（即\(ax-1\)是\(b\)的倍数）,不妨假设为-y倍，则原方程可以写成
\((ax-1)/b=-y\)
化简并移项整理得到\(ax+by=1\)
根据裴蜀定理，此时\(gcd(a,b)=1\)
所以我们就可以利用拓展欧几里得来求解了。
#include<cstdio>
using namespace std;
long long a,b,x,y;
void exgcd(int a,int b)
{
if (b==0)
{
x=1;
y=0;
return ;
}
exgcd(b,a%b);
int tmp=x;
x=y;
y=tmp-(a/b)*y;
}
int main()
{
scanf("%lld%lld",&a,&b);
exgcd(a,b);
printf("%lld",(x%b+b)%b);//防止x出现负数的情况
return 0;
}
乘法逆元
乘法逆元，是指数学领域群G中任意一个元素a，都在G中有唯一的逆元a‘，具有性质a×a'=a'×a=e，其中e为该群的单位元。
（来自百度百科）
百科的定义十分难懂，我们用比较容易理解的方式来讲述：
乘法逆元就是一个数在模\(p\)意义下的倒数。
即\(a/b \mod p=k\)且\(a*c\mod p=k\) 此时称\(c\)为\(b\)在模\(p\)意义下的乘法逆元。
举个例子 \(3/3 \mod 11=1\)
\(3*4 \mod 11=1\)
4就是3在模11意义下的乘法逆元。
那么就可以根据倒数的性质，推出原数和乘法逆元的关系\(x*inv(x)\mod p=1\)
从而得到同余方程 \(x*inv(x)\equiv1\pmod {b}\)
于是就能够利用拓展欧几里得算法求解了。
#include<cstdio>
using namespace std;
long long p,n,x,y;
void exgcd(int a,int b)
{
if (b==0)
{
x=1;
y=0;
return ;
}
exgcd(b,a%b);
int tmp=x;
x=y;
y=tmp-(a/b)*y;
}
int main()
{
scanf("%lld%lld",&n,&p);
for (int i=1;i<=n;i++)
{
exgcd(i,p);
printf("%lld\n",(x%p+p)%p);
}
return 0;
}
但是这种做法对于这道题只有80分。
因为数据量太大了……
于是又有神犇发明了线性求逆元的方式。
证明方式：不会略
#include<iostream>
using namespace std;
long long p,n,inv[20000530];
int main()
{
scanf("%lld%lld",&n,&p);
inv[1]=1;
for (int i=2;i<=n;i++)
inv[i]=-(p/i)*inv[p%i]%p;
for (int i=1;i<=n;i++) printf("%lld\n",((inv[i]%p+p))%p);
return 0;
}
