数学/数论——欧几里得算法及扩展，乘法逆元

欧几里得算法

介绍

欧几里得算法是用来求两个数的最大公约数（greatest common divisor(gcd)）。

流程

其计算流程如下：有a,b两个正整数，现在按照以下步骤计算它们的最大公约数，即gcd(a,b)

a/b=q1...r1
b/r1=q2...r2
r1/r2=q3..r3
  ...
rm/rn=qt...0

反复迭代至余数为0，最终得rn为a,b的最大公约数,这就是欧几里得算法。

欧几里得算法那的形式化表示：gcd(a,b)=gcd(b,r1)=gcd(r1,r2)...gcd(rm,rn)

证明

前提公理：a=md,b=nd(m,n,d均为正整数,d>1)，则m,n互质⟺d为，a,b的最大公约数。
假设：gcd(a,b)=d1 即a，b的最大公约数是d1,并且有a = md1,b = nd1(m,n,d1)为正整数。
不妨设a>b,
那么a=bq+r1, b=nd1,则r1=a-bq=md1-(nd1)q=(m-nq)d1

如果n,(m-nq)互质，则gcd(a,b)=gcd(b,r1)=d1(d1为正整数)
采用反证法，假定n,(m-nq)不是互质的，有gcd(n,m-nq)=d2(d2>1)
n=xd2,m-nq=yd2(x,y,d2为正整数)
则 m=nq+yd2=(xq+y)d2
那么 a=md=(xq+y)d1d2,b=nd=xd1d2
a，b的最大公约数为d1d2，不符合假设gcd(a,b)=d1，故n,(m-nq)互质。
则由公理gcd(b,r1)=d1
故有gcd(a,b)=d1=gcd(b,r1)

欧几里得算法扩展

需要证明：若整数a和b的最大公因数为d，即gcd(a,b)=d，那么一定存在整数x，y使得下式成立：
ax+by=gcd(a,b)
等价于证明 gcd(a,b)能够表示为a,b的线性组合。
不妨设a>b,
a=bq1+r1,0<r1<b
b=r1q2+r2,0<r2<r1
r1=r2q3+r3,0<r3<r2
...
rk-2=rk-1qk+rk,0<rk<rk-1
rk-1=rkqk+1

从倒数第二个式子开始,gcd(a,b)=rk=rk-2-rk-1qk,
将倒数第三个式子rk-3=rk-2qk-1+rk-1带入得，gcd(a,b)=rk-2-(rk-3-rk-2qk-1)qk=(1+qkqk-1)rk-2-qkrk-3,即gcd(a,b)为rk-2,rk-3的线性组合
持续带入可以得gcd(a,b)为a,b的线性组合

乘法逆元

介绍

如果两个整数a和b满足(ab-1)mod p =0,即ab≡1(mod p),则称a和b互为乘法逆元。
例如，如果a,b,p满足ab≡1(mod p)，那么b可以取消a和任意整数x在模p时的乘法，即a和b互为模p的乘法逆元。
即xab=x(modp),证明如下；
x mod p =(0+x) mod p = x(ab-1)mod p + x mod p = xa*b mod p

乘法逆元的性质

a与p互质是存在a关于模p的乘法逆元b的充分必要条件
证明：
必要性：
假设a与p互为模p的乘法逆元，且a与p的最大公约数为d，则有
(ab-1)mod p =0 => ab-1=np => xdb -1 = nyd => xb - ny = 1/d
显然xb-ny为整数，当且仅当d为1时，即a与p互质时，上述等式才成立
充分性：
gcd(a,p)=gcd(p,r1)=gcd(r1,r2)=...=gcd(rm-1,rm)=1//因为a与p互质
a=n1p+r1,p=n2r1+r2...rm-3=nm-1rm-2+rm-1,rm-2=nmrm-1,rm-1=0,rm=0
带入或者使用从欧几里得扩展获得的定理
gcd(a,p)= xa+yp=1(x,y为正整数)
即，当a与p互质时，存在x，y互为模p时的乘法逆元

计算方法

根据扩展定理，ax+by=gcd(a,b)=gcd(b,a mod b)=gcd(rk-1,rk)到达边界条件即，rk=0时，rk-1=d=1。即按照如下方法即可最终获得的x即为 a关于模b的乘法逆元

int exgcd(int a,int b,int &x,int &y)
{
    if(b==0) {
        x=1;y=0;
        return a;
    }
    int r = exgcd(b,a%b,x,y);
    int t = y;
    y = x-(a/b)*y;
    x = t;
    return r;
}

费马小定理 & 快速幂

费马小定理

如果p为素数，且gcd（a,p）=1,则满足
a^(p-1)=1(mod p)
a*a^(p-2)=1(mod p)
于是，当a为素数时，a在模p时的乘法逆元为a^(p-2),可用快速幂获得

inline int quick_power(int a, int b, int p) {
  int ans = 1;
  a = (a % p + p) % p;
  for (; b; b >>= 1) {
    if (b & 1) ans = (a * ans) % p;
    a = (a * a) % p;
  }
  return ans;
}
// quick_power(a, p-2, p) 即为 a 的乘法逆元。

线性求解1到n中每个数的逆元

求1到n中每个数i在模p时的乘法逆元
特别的，当i=1时，其乘法逆元i^(-1)为1
当i>1，设k=floor(p/i),j=p mod i,即（p=ki+j）
由 p≡0(mod p)=>(ki+j)≡0(mod p)=>(ki+j)i^(-1)*j(-1)≡0(mod p)=>kj^(-1)+i(-1)≡0(mod p)=> -k*j^(-1)=i(-1)(mod p)
i^(-1)=-floor(p/i)(p mod i)^-1(mod p),i>2;1,i=1

求任意n个整数的乘法逆元

设有n个整数，分别为a1,a2,a3...,对应的逆元为a1^(-1),a2(-1),a3^(-1)...
设Si为a1到ai i个整数的乘积，Si^(-1)为a1(-1)到ai^{(-1)i个数的乘积，那么Si与Si}(-1)互为逆元
因此可以通过扩展欧几里得定理或者快速幂Sn的乘法逆元Sn^(-1)
又有
Si^(-1)*ai+1(-1)ai+1=Si^(-1)(mod p)
Si+1^(-1)*ai+1=Si(-1)(mod p)
由此，可以有ai取消ai^{(-1)乘法运算的性质，由Si}(-1)计算出Si-1^(-1)
这样可以以O(n)计算出所有的Si^(-1)
最后有Si^{(-1)及Si计算出所有的ai}(-1）≡ Si-1*Si^(-1)(mod p)

参考文章

• https://blog.csdn.net/qq_37174526/article/details/90112064
• https://www.cnblogs.com/Plozia/p/16156532.html
• https://blog.csdn.net/Time_Limit/article/details/117390957