数论 —— 逆元（费马小定理/扩展欧几里得）

拓展欧几里得：

(证明过程参考自百度百科)

原式: ax+by=gcd(a,b)(假设a≥b)


当b=0时有gcd(a,b)==a,此时x=1,y=0


当b不为0时,根据欧几里得定理gcd(a,b)=gcd(b,amodb)可得ax+by=gcd(a,b)=gcd(b,amodb)=bx′+(amodb)y′,即 


ax+by=bx′+(amodb)y′=bx′+(a−b∗⌊a/b⌋)y′


移项得 
ax+by=bx′+(amodb)y′=ay′+b(x′−⌊a/b⌋y′)


根据恒等定理,有 
{x=y′y=x′−⌊a/b⌋y′


这有什么用呢?x′和y′还是不知道呀. 
重新来看看我们得到的两个等式.x和y是gcd(a,b)=ax+by的解,而x’和y’是在对gcd(a,b)按欧几里德算法进行一步后的结果对应的贝祖等式gcd(b,amodb)=bx′+(amodb)y′的解.也就是说,gcd(a,b)对应的贝祖等式的解x,y可以由gcd(b,amodb)对应等式的解x’,y’计算得出 
由于欧几里德算法最后一步为gcd(d,0)=d,此时对应的等式的解为x=1,y=0,因此只要如上述代码,从gcd(d,0)往前处理,在进行欧几里德算法的递归的时候根据相邻两次调用间x,y和x’,y’的关系计算即可求出ax+by=gcd(a,b)的解.
更进一步,对于任意不定式ax′+by′=c,只需要在等式ax+by=gcd(a,b)=d两边乘上c/d即可得到解为x′=x∗c/d,y′=y∗c/d


如何得到所有解? 
实际上在之前的计算和证明中我们得到的只是不定方程的一组解,那么怎样得到所有解呢?对于一般形式ax+by=c有通解x=p+kb,y=q−ka(k为任意整数).(证明略,只要代入一下就知道为什么通解是这个了)

int extend_GCD(int a,int b,int &x,int &y){
    int tmp,d;
    if(b==0){
        x=1;
        y=0;
        return a;
    }
    int res=extend_GCD(b,a%b,x,y);
    tmp = x;
    x = y;
    y = tmp - (a/b)*y;
    return res;
}

void Extend_Gcd(int a,int b,int &n,int &x,int &y){
    if(!b){
        n=a,x=1,y=0;
        printf("1:%d %d\n",x,y);
    }else{
        Extend_Gcd(b,a%b,n,y,x);
        y-=x*(a/b);
        printf("---%d %d\n",x,y);
    }
}

费马小定理：

x与p互质:

若gcd(x, p) == 1

x^(p-1) == 1(mod p)

快速幂求解即可

逆元：

a*x  == 1（mod m）

x是相对于m的a的逆元

a是相对于m的x的逆元

1.费马小定理：

若gcd(x,p) == 1

x^(p-1) == 1 (mod p)

x*x^(p-2) == 1(mod p)

所以x是x^(p-2)关于p的逆元

LL quick(LL x,LL n){
    LL res=1;
    while(n){
        if(n%2)res=x*x%Mod;
        x*=x%Mod;
        n>>=1;
    }
    return res;
}
LL inv(LL x){
    return quick(x,Mod-2);
}

2.扩展欧几里得求逆元：

    a*x == 1(mod p)

    a*x + kp == 1

扩展欧几里得解得extend_gcd(a,p)

求解得 a*x0+p*y0 = 1

x0是a的一个特解 

通解 x == x0+pk 

3.

当p是个质数的时候有

inv(a) = (p - p / a) * inv(p % a) % p

证明：

设x = p % a,y = p / a
于是有 x + y * a = p
(x + y * a) % p = 0
移项得 x % p = (-y) * a % p
x * inv(a) % p = (-y) % p
inv(a) = (p - y) * inv(x) % p
于是 inv(a) = (p - p / a) * inv(p % a) % p

然后一直递归到1为止，因为1的逆元就是1

 

LL inv(LL t, LL p) {//求t关于p的逆元，注意:t要小于p，最好传参前先把t%p一下 
    return t == 1 ? 1 : (p - p / t) * inv(p % t, p) % p;
}

可以预处理出前N项逆元：

#include<cstdio>
const int N = 200000 + 5;
const int MOD = (int)1e9 + 7;
int inv[N];
int init(){
    inv[1] = 1;
    for(int i = 2; i < N; i ++){
        inv[i] = (MOD - MOD / i) * 1ll * inv[MOD % i] % MOD;
    }
}
int main(){
    init();
}