數論基礎

數論是純數學的一個研究分支，主要研究整數的性質，初等數論包括整除理論、同余理論、連分數理論，這一篇主要記錄的是同余相關的基礎知識，

取模

取模是一種運算，本質就是帶余除法，運算結果就是余數，取模運算結果的符號由被模數（被除數）決定，

\[7\%4=3;\space7\%(-4)=3;\\ (-7)\%4=-3;\space(-7)\%(-4)=-3 \]

取模運算的性質

\[設a>b>0，有：\\ (a+b)\%c=(a\%c+b\%c)\%c\\ (a-b)\%c=(a\%c-b\%c+c)\%c \\ (a\times b)\%c=(a\%c\times b\%c)\%c \]

例題1

\[給定兩個非負整數a，b和正整數n，計算f(a^b)\%n，\\ 其中，f(0)=f(1)=1,且對于\forall i>0，f(i+2)=f(i+1)+f(i)\\ (0\le a,b\le 2^{64},1\le n\le 1000) \]

找規律，對于mod n來說，最多n2項就會出現重復，找出重復項即可得到回圈周期，然后找對應的項即可，

例題2：hdu 1212

題意是給一個位數不超過1,000的正整數A和一個大小不超過100,000的正整數B，求A%B，

這里只需要用到上面的兩條性質，把A截斷，從高到低模擬做除法即可，

#include <iostream>
string s; int mod, ans;
int main() {
    while(std::cin >> str >> mod) {
        ans = 0;
        for(int i = 0; s[i]; i++) ans = (ans * 10 + (s[i] - '0')) % mod;
        std::cout << ans << std::endl;
    }
    return 0;
}

GCD

GCD即Greatest Common Divisor，最大公約數，最大公約數即能同時整除給定兩個整數的最大正整數，特殊地：gcd（a，0）=a，求解最大公約數通常使用輾轉相除法，下面是輾轉相除法的代碼實作：

typedef long long LL;
LL gcd(LL a, LL b) {
    return b ? gcd(b, a % b) : a;
}

輾轉相除法得到的余數序列增長速度比斐波那契數列更快，已知斐波那契的增長是指數級別，則輾轉相除法的復雜度是對數級別，

輾轉相除法的證明：

1.設兩個數a、b（a＞b），他們的最大公約數為gcd（a，b），r = a % b，k = （a - a % b）/ b，那么證明輾轉相除法，即是要證明：gcd（a，b）=gcd（b，r），首先我們令c = gcd（a，b），那么肯定存在互質的整數m，n，使得a = mc，b = nc，那么有r = a - kb = mc - knc = (m - kn)c，根據這條式子，c也是r的因數，回過頭看，如果能證明m-kn和n互質，那么就可以說明c=gcd（（m-kn）c，nc）=gcd（b，r），所以把問題再次轉化為：求證m-kn和n互質，

2.反證法證明m-kn和n互質：假設m-kn和n不互質，用數學語言描述為：假設存在整數x，y，d，其中d＞1，使得m-kn=xd，n=yd，那么有m=kn+xd=kyd+xd=(ky+x)d，從而a=mc=(ky+x)cd，b=nc=ycd，則gcd（a，b）=cd≠c，與前設矛盾，故m-kn和n互質得證，也就證明了gcd（a，b）=gcd（b，r），

唯一分解定理

對于任意大于2的正整數X，它總能寫成如下形式：

\[X=p_1^{a_1}\times p_2^{a_2} \times ......\times p_k^{a_k}，其中p是各不相同的質數 \]

即質因數分解，且這個分解唯一，

LCM

LCM即Least Common Multiple，最小公倍數，

\[\begin{align} &根據唯一分解定理，對于給定的a，b：\\ &a=p_1^{a_1}\times p_2^{a_2} \times......\times p_k^{a_k}\\ &b=p_1^{b_1}\times p_2^{b_2} \times......\times p_k^{b_k}\\ &那么：\\ &gcd(a,b)=p_1^{min(a_1,b_1)}\times p_2^{min(a_2,b_2)}\times ......p_1^{min(a_k,b_k)}\\ &lcm(a,b)=p_1^{max(a_1,b_1)}\times p_2^{max(a_2,b_2)}\times ......p_1^{max(a_k,b_k)}\\ &因此，gcd(a,b)\times lcm(a,b)=a\times b \end{align} \]

LL lcm(LL a, LL b) {
    return a / gcd(a, b) * b;
}

例題：hdu1788

題意是求最小的正整數N，滿足除以每一個給定的數Mi，其結果均為Mi-a，作如下轉化：

\[N\%M_i\equiv(M_i-a)\%M_i\Rightarrow (N+a)\%M_i=0（a<M_i<100） \]

即N是Mi的倍數，那么題意也就是求所有M的最小公倍數，注意數字范圍即可，

同余

對于兩個不同的數a，b，如果有a % p = b % p（p＞1），那我們就稱a和b模p同余，記作：

\[a\equiv b(mod\space p) \]

同余的性質

\[\begin{aligned} &1.自反性：a\equiv a(mod\space m)\\ &2.對稱性：a\equiv b(mod\space m)，則b\equiv a(mod\space m)\\ &3.傳遞性：若a\equiv b(mod\space m)，b\equiv c(mod\space m)，則a\equiv c(mod\space m)\\ &4.線性運算：若a\equiv b(mod\space m)，c\equiv d(mod\space m)\\ &\space\space\space則：a\pm c\equiv b\pm d(mod\space m);\space a\times c\equiv b\times d(mod\space m)\\ &5.冪運算：若a\equiv b(mod\space m)，那么a^n\equiv b^n(mod\space m)\\ &6.若ac\equiv bc(mod\space m)，且c≠0，那么a\equiv b(mod\space {m\over gcd(c,m)})\\ &\space\space\space特殊地，若gcd(c,m)=1，則a\equiv b(mod\space m) \end{aligned} \]

EXGCD

貝祖定理（Bezouts Identity）：若設a,b是不全為0的整數，則存在整數x,y，使得ax+by=gcd（a,b），(a,b)代表最大公因數，則設a,b是不全為零的整數，則肯定存在整數x,y，使得ax+by=(a，b)，這個式子稱為貝祖等式，

EXtend GCD 即擴展歐幾里得演算法，它可以用于解關于x，y的貝祖等式，

EXGCD的可行性

當a=0時，ax+by=gcd（a，b）=gcd（0，b）=b，這時可以解出y=1而x為任意，同理可得b=0時，解得x=1而y為任意，當a，b都不為0時：

\[\begin{aligned} &設存在一組解為x_1，y_1，即ax_1+by_1=gcd(a,b)\\ &由歐幾里得定理：gcd(a,b)=gcd(b,a\%b)得：\\ &ax_1+by_1=gcd(a,b)=gcd(b,a\%b)\\ &那么,bx_2+(a\%b)y_2=gcd(b,a\%b)=ax_1+by_1\\ &又a\%b=a-a/b\times b\\ &\begin {align} 那么，&ax_1+by_1\\ =&bx_2+(a-a/b\times b)y_2\\ =&bx_2+ay_2-a/b*b*y_2\\ =&ay_2+b(x_2-a/b*y_2) \end{align}\\ &故可以得到：x_1=y_2,\space y_1=x_2-a/b*y_2 \end{aligned} \]

a或b為0即迭代求解的出口，迭代程序用代碼表示如下：

typedef long long LL;
// ver 1
// 回傳值為gcd（a，b），x和y即求出的特解
LL exgcd(LL a, LL b, LL &x, LL &y) {
    if (b == 0) {
        x = 1; y = 0;
        return a;
    }
    LL ans = exgcd(b, a % b, x, y);
    //這里通過迭代求出了一組解x，y，這組解對應上面推導程序中的x2和y2，
    LL temp = x;
    x = y;
    y = temp - a / b * y;
    return ans;
}

// ver 2
// 在熟悉了推導程序之后，我們可以利用參考的性質得到更為簡潔的ver2寫法
LL exgcd(LL a, LL b, LL &x, LL &y) {
    if (b == 0) {
        x = 1; y = 0;
        return a;
    }
    LL ans = exgcd(b, a % b, y, x);
    y -= a / b * x;
    return ans;
}

EXGCD的應用

求關于x，y的方程ax+by=c的一組解

這里只需要對貝祖等式稍作變換即可，假設EXGCD求出的方程ax+by=gcd（a，b）的一組解為x1和y1，在方程兩邊同時乘上一個數m，使得c=m*gcd（a，b），這里也就要求c是gcd（a，b）的倍數，即c%gcd（a，b）=0，這也是方程有解的條件，此時方程為：

\[\begin{align} &ax_1\times m+by_1\times m=gcd(a,b)\times m=c \end{align} \]

對于ax+by=gcd（a，b），其引數解為：

\[\begin{align} &x=x_1+{b\over gcd(a,b)}\times t,y=y_1-{a\over gcd(a,b)}\times t,t為引數\\ &這里選用{b\over gcd(a,b)}和{a\over gcd(a,b)}是為了保證結果均為整數 \end{align} \]

那么對比方程ax+by=c，可以得到其特解為：

\[\begin{align} &x_0=x_1\times m=x_1\times {c\over gcd(a,b)}\\ &y_0=y_1\times m=y_1\times {c\over gcd(a,b)} \end{align} \]

那么其通解為：

\[\begin{align} &X=x_0+{b\over gcd(a,b)}\times t,Y=y_0-{a\over gcd(a,b)}\times t,t為引數\\ &這里選用{b\over gcd(a,b)}和{a\over gcd(a,b)}作為系數是為了確保得到均為整數解， \end{align} \]

對于任意一個確定的t，都有一組確定的解與之對應，只需要根據需要找出對應的解即可，

例如求滿足ax+by=c的X的最小非負整數解，那么在得到X的通解之后只需調整t，調整程序用偽代碼表示為：

S = b / gcd(a, b);
X = (x0 % S + S) % S;
// 由于x0可能是負數，故在模之后還要加上S在取一次模，
// 同理可得Y的最小非負整數解
T = a / gcd(a, b);
Y = (y0 % T + T) % T;

完整的exgcd求解方程ax+by=c的X和Y的最小非負整數解得代碼如下：

LL exgcd(LL a, LL b, LL c, LL& x, LL& y) {
    if (b == 0) {
        x = 1; y = 0;
        return a;
    }
    LL ans = exgcd(a, b, c, y, x);
    y -= a / b * x;
    LL S = b / ans, T = a / ans;
    x = (x % S + S) % S;
    y = (y % T + T) % T;
    return ans;
}

求關于x的方程ax≡b（mod m）的最小非負整數解

ax≡b（mod m）即ax%m=b%m，即求ax+my=b%m（取模其實就相當于減去（加上）了若干個模數），但只有a和m互質時有唯一解，否則無解，轉化為ax+my=b（假設b<m）后，套用exgcd即可，

求a關于模數p的乘法逆元

設x是a關于p的乘法逆元，那么有ax≡1（mod p），這是第二個應用的特殊情況，顯然也可以通過類似方法求得，

逆元

Inverse Element，逆元，推廣了加法中的加法逆元和乘法中的倒數，直觀地說，它是一個可以取消另一給定元素運算的元素，a關于模p的逆元存在的條件是gcd（a，p）=1，

\[\begin{align} &在模p意義下，設A的逆元為A^{-1}，那么有A\times A^{-1}\equiv 1(mod\space p) \end{align} \]

為什么需要逆元呢？

\[\begin{align} &在模意義下，{A\over B}\%p = {A\%p\over B\%p}\%p并不成立，\\ &那么設B在模p意義下的逆元表示為B^{-1}，根據逆元的定義，有{A\over B}=A\times B^{-1}(mod \space p)， \end{align} \]

這里，我們把除法轉化為乘法，就可以運用取模運算的性質：(a * b) % c = (a % c * b % c) % c，優化演算法，

逆元的求法

EXGCD求法

給定模數m，求a的逆元相當于求解關于x的方程ax≡1(mod m)，這個方程可以轉化為ax-my=1 ，用EXGCD可以求得一組x，y和gcd，檢查gcd是否為1（因為EXGCD是把問題轉化為求解方程ax-my=gcd（a，m），顯然只有gcd（a，m）=1時才能轉化），gcd不為1則說明逆元不存在，在能夠成功求解的情況下，把解調整x到0~m-1的范圍中即可，

費馬小定理

在模p為質數的情況下，有費馬小定理：

\[a^{p-1}\equiv 1(mod\space p) \]

那么:

\[a\times a^{p-2}\equiv 1(mod\space p)\\ 故a的逆元a^{-1}=a^{p-2} \]

然后用快速冪求出逆元即可，

拓展：快速冪

// a是底數，b是指數
LL pow(LL a, LL b, LL mod) {
    LL ans = 1;
    while(b) {
        if(b&1) ans = ans * a % mod;
        b >>= 1, a = a * a % mod;
    }
    return ans;
}

拓展：慢速乘，防止快速冪程序中乘法運算爆long long，只需把快速冪中的乘法替換為mul（）即可，一般情況用不上，

// a為因數1，b為因數2
LL mul(LL a, LL b, LL mod) {
    LL ans = 0;
    while(b) {
        if(b&1) ans = (ans + a) % mod;
        b >>= 1, a = (a + a) % mod;
    }
    return ans;
}

歐拉定理

費馬小定理其實是歐拉定理的特殊情況，在模數p不為質數時，有：

\[a^{\phi(p)}\equiv 1(其中\phi()是歐拉函式) \]

同理可以得到：

\[a^{-1}\equiv a^{\phi(p)-1} \]

線性逆元表

有時我們需要快速得出1~p-1的所有逆元，這個時候我們需要一種O（n）的方法，這就是線性逆元表，

首先，1關于任意模p的逆元的逆元都是1，設p = k * i + r（r < i ，1 < i < p），那么有：

\[k\times i + r \equiv0(mod\space p) \]

利用逆元性質，兩邊同時乘上i的逆元以及r的逆元，得到：

\[k\times r^{-1}+i^{-1}\equiv 0(mod\space p) \]

移項得：

\[i^{-1}\equiv-k \times r^{-1}\equiv -?{p\over i}?\times (p\%i)^{-1}(mod\space p) \]

顯然p%i是小于i的，那么我們就可以利用之前的結果在O（1）時間算出單獨的逆元，

用代碼表示即：

// inv[i]對應i的逆元
LL inv[maxn];
void CalInv() {
	// 0沒有逆元，故不初始化
    inv[1] = 1;
    for (int i = 2; i < maxn; i++)
        inv[i] = (mod - mod / i) * inv[mod % i] % mod;
}

利用

\[(a\times c)\% mod = (a\%mod\times c\%mod)\%mod；\\ 由a\times a^{-1}\equiv1(mod \space p)得，\\ inv_{fac(i)}\equiv inv_{fac(i)}\times (i+1)^{-1}\times (i+1)\equiv inv_{fac(i+1)}\times (i+1)(mod\space p) \]

我們還可以在線性時間求出1~min(n，p)的階乘的逆元，代碼如下：

LL fac[maxn], inv[maxn];
void CalFacInv() {
    // 0的階乘是1
    fac[0] = fac[1] = 1;
    for (int i = 2; i < maxn; i++)
        fac[i] = fac[i - 1] * i % mod;
    inv[maxn - 1] = QPow(fac[maxn - 1], mod - 2);
    // 注意邊界是≥
    for (int i = maxn - 2; i >= 0; i--)
        inv[i] = inv[i + 1] * (i + 1) % mod;
}

例題：hdu1576

題意即把除法轉化為逆元乘法，3種方法的代碼：

#include<iostream>
#include<algorithm>
using namespace std;

typedef long long LL;
LL n, B;
const LL mod = 9973;

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

LL QPow(LL x, LL n)
{
    LL res = 1;
    while(n)
    {
        if(n & 1) res = res * x % mod;
        x = x * x % mod;
        n >>= 1;
    }
    return res;
}

LL inv[10000];
void CalInv() {
    inv[1] = 1;
    for (LL i = 2; i < 10000; i++) {
        inv[i] = (mod * mod - mod / i * inv[mod % i]) % mod;
    }
}


int main() {
    CalInv();
    int t;
    cin >> t;
    while (t--) {
        cin >> n >> B;
        // exgcd
        // LL x, y;
        // exgcd(B, mod, x, y);
        // cout << (x + mod) * n % mod << endl;
        
        // QPow
        // cout << QPow(B, mod - 2) * n % mod << endl;
        
        // 逆元打表
        cout << inv[B % mod] * n % mod << endl;
    }
    return 0;
}

轉載請註明出處，本文鏈接：https://www.uj5u.com/qita/122728.html

標籤：其他

上一篇：小伙伴們，你們誰有劉汝佳演算法競賽入門經典這本書的習題答案嘞?

下一篇：java進階_高級學習

最大公約數、擴展歐幾里得演算法、逆元