数论代码总结

欧几里得:求a,b的最大公约数。

1
2
3
4
int gcd(int a,int b){
    if(b==0)return a;
    return gcd(b,a%b);
}

扩展欧几里得:(求a,b的最大公约数和不定方程ax+by=d(=gcd(a,b))的一对解)
裴蜀定理
对于任意整数 $a$ , $b$,存在无穷多组整数对 $(x,y)$ 满足不定方程 $ax+by=d$ ,其中 $d=gcd(a,b)$。
在求 $d=gcd(a,b)$的同时,可以求出关于x,y的不定方程$ax+by=d$的一组整数解。
考虑递归计算:假设已经算出了$(b,a \bmod b)$的一对解$(x_0,y_0)$满足$bx_0+(a \bmod b)y_0=d$。
可以得到$bx_0+(a – b \left \lfloor \frac{a}{b} \right \rfloor)y_0=d$。
整理可得$b(x_0-\left \lfloor \frac{a}{b} \right \rfloor y_0) +ay_0=d$。
即$ay_0+b(x_0-\left \lfloor \frac{a}{b} \right \rfloor y_0)=d$。
即$x=y_0,y=x_0-\left \lfloor \frac{a}{b} \right \rfloor y_0$。
假设先行$swap(x_0,y_0)$,则有$x=x_0,y=y_0-\left \lfloor \frac{a}{b} \right \rfloor x_0$。

1
2
3
4
5
6
7
8
9
int exgcd(int a,int b,int &x,int &y){
    if(b==0){
        x=1;y=0;
        return a;
    }
    int d=exgcd(b,a%b,y,x);
    y-=a/b*x;
    return d;
}

乘法逆元
方法1:利用扩展欧几里得求逆元,参见:https://renjikai.com/luogu-p1082/
适用条件:被求数a和模数b互质。
方法2:利用费马小定理(快速幂)求逆元,参见:https://renjikai.com/luogu-p1313/
适用条件:模数为质数。(保证了被求数a和模数b互质。)
方法3:不用记复杂的线性递推公式,就可以用O(N)复杂度推出1到n的所有乘法逆元(模意义下,在后面不再赘述)。除此之外,在此过程中求出$1!$到$n!$的阶乘及其逆元。
前提结论:易知$\prod_{i=1}^{n}\ inv(i)=inv(n!)$。
乘法逆元定义:$inv(i)\cdot i\equiv 1$
令数组$A_i=i!$,容易求出数组A。令数组$B_i=inv(i!)$,那么数组B应该怎么求呢?
用方法1、2容易求出$inv(A_n)$,即$inv(n!)$,即求出了$B_n$。
由前提结论知,$B$数组的递推公式为$B_i=B_{i-1} \cdot inv(i)$,两边同乘$i$,则有$B_i \cdot i =B_{i-1} \cdot inv(i)\cdot i$,由乘法逆元定义有$B_{i-1} = B_i \cdot i $,只需回推,即可求出由$i!$的乘法逆元。
因$B_n=\prod_{i=1}^{n}\ inv(i)=inv(n)\cdot inv(n-1) \cdot inv(n-2) … \cdot inv(1)$,$A_{n-1}=(n-1)!=(n-1)\cdot(n-2)…\cdot(1)$,因此两式相乘则有$inv(n)=B_n\cdot A_{n-1}$,即求出了1到n的所有乘法逆元,是线性复杂度$O(N)$,可能常数大一点,但是好记啊。
C++伪代码如下:(不一定完全正确,请注意识别)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
#define MOD
#define ll long long

ll quickpow(ll a, ll b);
ll fInv(ll x) {
    return quickpow(x, MOD - 2) % MOD;
}
ll n;
ll fact[n + 2];
ll invFact[n + 2];
ll inv[n + 2];
void gen() {
    fact[1] = 1;
    for (ll i = 2; i <= n; i++) {
        fact[i] = (fact[i - 1] * i) % MOD;
    }
    invFact[n] = fInv(fact[n]);
    for (ll i = n - 1; i >= 1; i--) {
        invFact[i] = invFact[i + 1] * (i + 1) % MOD;
    }
    inv[1] = 1;
    for (ll i = 2; i <= n; i++) {
        inv[i] = invFact[i] * fact[i - 1] % MOD;
    }
}

发表评论

电子邮件地址不会被公开。 必填项已用*标注