(ex)BSGS/（擴展）大步小步演算法 學習筆記

在即將暫時退役之際殺掉了P4195的毒瘤模板題，于是來寫篇學習筆記，

謹此為我初中三年擺爛的OI生涯畫上一個句號，（距離中考還有20天！）

BSGS

link

求\(a^x\equiv b\pmod p\)的非負整數解，其中\(a, p\)互質，

演算法思路

我們不妨令\(t=\lceil{\sqrt{p}\rceil}\)，\(j\lt t\)，\(i\leq t\)

原式轉化為\(a^{it-j} \equiv b \pmod p\)

即 \(\left(a^t\right)^i\equiv b\cdot a^j \pmod p\)

于是我們可以這么在\(\Theta\left(\sqrt{n}\right)\)（不考慮hash表）內求出解：

• 列舉\(j \in [0, t)\)，求出所有的\(b\cdot a^j \mod p\)，用hash表記錄下來

• 列舉\(i \in [0, t]\)，求出所有的\(\left(a^t\right)^i \mod p\)，在hash表中查找是否有對應的\(j\)值

需要注意的是，當\(a^t \mod p=0\)時，若\(b=0\)，則解為\(x=1\)；
否則無解，

正確性討論

由Euler定理，我們有\(a^{\varphi\left(p\right)}\equiv 1\pmod p\)，其中\(\varphi\left(x\right)\)為Euler函式，

也就是說，\(\mod p\)意義下，\(a^x\mod p\)在\(x=n\cdot\varphi\left(p\right)\)（\(n\)遍歷非負整數）后回圈，

我們知道對于任意整數\(p \gt 1\)，\(p>\varphi\left(p\right)\)，而我們取的\(it-j\)可以遍歷\([0, x]\)，因此能夠取到\(a^p \mod p\)的一切情況，故一定不會漏解，

代碼實作

#include <bits/stdc++.h>

using namespace std;

#define int long long

int qpow(int a, int n, int p) {
    a%=p;
    int res=1;
    while (n) {
        if (n&1) res=res*a%p;
        a=a*a%p;
        n>>=1;
    }
    return res;
}

signed bsgs(int a, int b, int p) {
    b%=p;
    int t=sqrt(p)+1;
    unordered_map<int, int> hash; hash.clear();
    for (int j=0; j<t; ++j) {
        int power=b*qpow(a,j,p)%p;
        hash[power]=j;
    }
    a=qpow(a,t,p);
    if (a==0) return b==0?1:-1;
    for (int i=0; i<=t; ++i) {
        int val=qpow(a,i,p);
        int j=hash.find(val)==hash.end()?-1:hash[val];
        if (j>=0 && i*t-j>=0) return i*t-j;
    }
    return -1;
}

signed main() {
    int a, b, p;
    cin>>p>>a>>b;
    int res=bsgs(a,b,p);
    if (res==-1) puts("no solution");
    else cout<<res<<endl;
    return 0;
}

這份代碼是可以通過P3846的，我們可以考慮對它進行優化，

我們發現列舉\(a^j\)和\(\left(a^t\right)^i\)的時候，\(i, j\)是遞增的，于是我們可以直接用一個變數來維護\(a^j\)和\(\left(a^t\right)^i\)，

另外，用unordered_map來實作hash表似憾訓快一些，

inline int bsgs(int a, int b, int p) {
    int t=(int)(sqrt(p))+1;
    unordered_map<int,int> h; h.clear();
    int powa=1;
    for (reg int j=0; j<t; ++j) {
        int val=powa*b%p;
        h[val]=j;
        powa=powa*a%p;
    }
    a=qpow(a,t,p);
    powa=1;
    for (reg int i=0; i<=t; ++i) {
        int val=powa%p;
        int j=h.find(val)==h.end()?-1:h[val];
        if (j>=0 && i*t-j>=0) return i*t-j;
        powa=powa*a%p;
    }
    return -1;
}

為exBSGS進行修改

我們現在來考慮\(ka^x\equiv b\pmod p\)（\(a, p\)互質，\(k\)為正整數）的式子，我們可以同樣地將它們變形為

于是稍微修改一下上述代碼就可以了，

inline int bsgs(int a, int b, int k, int p) {
    int t=(int)(sqrt(p))+1;
    unordered_map<int,int> h; h.clear();
    int powa=1;
    for (reg int j=0; j<t; ++j) {
        int val=powa*b%p;
        h[val]=j;
        powa=powa*a%p;
    }
    a=qpow(a,t,p);
    powa=1;
    for (reg int i=0; i<=t; ++i) {
        int val=k*powa%p;
        int j=h.find(val)==h.end()?-1:h[val];
        if (j>=0 && i*t-j>=0) return i*t-j;
        powa=powa*a%p;
    }
    return -1;
}

exBSGS

link

求\(a^x\equiv b\pmod p\)的非負整數解，其中\(a, p\)未必互質，

演算法思路

考慮利用同余式的約化性質，轉換成樸素的BSGS來求解，

我們有如下同余式的約化性質：
若\(a\equiv b\pmod p\)，\(d\mid a,d \mid b\)，則

\(\frac{a}{d}\equiv\frac{b}{d}\pmod {\frac{p}{\gcd(d,p)}}\)

我們回到\(a^x\equiv b\pmod p\)，令\(d_1=\gcd(a,p)\)，
當\(d_1 \mid b\)，我們可以變形為
\(\frac{a}{d_1}\cdot a^{x-1}\equiv \frac{b}{d_1} \pmod {\frac{p}{d_1}}\)
（若\(d_1 \nmid b\)，立刻得出無解）
若\(a,\frac{p}{d_1}\)仍不互質，我們可以繼續令\(d_2=\gcd(a,\frac{p}{d_1})\)，\(\cdots\)，直到\(a,\frac{p}{d_1d_2\cdots d_k}\)互質為止，

我們設一共進行了\(k\)次約化，\(D=d_1d_2\cdots d_k\)（此時\(a\)，\(\frac{p}{D}\)互質），原式可變形為

\(\frac{a^k}{D}\cdot a^{x-k} \equiv \frac{b}{D} \pmod {\frac{p}{D}}\)

我們令\(k=\frac{a^k}{D}, X=x-k, B=\frac{b}{D}, P=\frac{p}{D}\)，即

\(ka^X\equiv B \pmod P\)

于是可以利用修改后的BSGS演算法來求解，

注意求解之后得到\(X=x-k\)，故\(x=X+k\)，

Trick

\(\frac{a^k}{D}=\frac{a}{d_1}\frac{a}{d_2}\cdots\frac{a}{d_k}\)，于是可以在每一個回圈內單獨計算，

用cout<<'\n'似憾訓比用cout<<endl快一些，

可以預先將b%=p, a%=p，若取模過后\(b==1\)或者\(p==1\)，那么顯然\(x=0\)為解，

我們記\(D_k=\frac{a}{d1}\frac{a}{d2}\cdots \frac{a}{d_k}\)，當\(\frac{a^k}{D^k}==\frac{b}{D_k}\)時，有

即

于是\(x=k\)為解，其中\(k\)是正在進行的第\(k\)次約化，

代碼實作

#include <bits/stdc++.h>
using namespace std;

#define int long long
#define reg register
inline int qpow(int a, int n, int p) {
    a%=p; int res=1;
    while (n) {
        if (n&1) res=res*a%p;
        a=a*a%p; n>>=1;
    }
    return res;
}
inline int bsgs(int a, int b, int k, int p) {
    int t=(int)(sqrt(p))+1;
    unordered_map<int,int> h; h.clear();
    int powa=1;
    for (reg int j=0; j<t; ++j) {
        int val=powa*b%p;
        h[val]=j;
        powa=powa*a%p;
    }
    a=qpow(a,t,p);
    powa=1;
    for (reg int i=0; i<=t; ++i) {
        int val=k*powa%p;
        int j=h.find(val)==h.end()?-1:h[val];
        if (j>=0 && i*t-j>=0) return i*t-j;
        powa=powa*a%p;
    }
    return -1;
}
inline int exbsgs(int a, int b, int p) {
    a%=p, b%=p;
    if (b==1 || p==1) return 0;
    int A=a, NA=1, B=b, P=p, k=0, D=1;
    while (__gcd(a,P)>1) {
        int d=__gcd(a,P);
        if (B%d) return -1;
        k++; A/=d, B/=d, P/=d, NA=NA*(a/d)%P, D=D*d%P;
        if (NA==B) return k; // NA就是上文提到的(a^k)/(D^k)
    }
    int res=bsgs(a,B,NA,P);
    if (res==-1) return res;
    if ((qpow(a,res+k,p)-b)%p) return -1;
    return res+k;
}

signed main() {
    ios::sync_with_stdio(false);
    cin.tie(NULL); cout.tie(NULL);
    int a, b, p;
    while (cin>>a>>p>>b && a) {
        int res=exbsgs(a,b,p);
        if (res==-1) cout<<"No Solution\n";
        else cout<<res<<'\n';
    }
    return 0;
}

Record，\(2.46\rm{s}\)，可以通過本題（包括Hack資料），

后記

這道題算是比較毒瘤的，我是一共調了三天才過的（我太弱了）
還有\(20\)天就要中考了，而我還在這摸魚（悲）
就謹此紀念一下我的初中OI生涯罷，

順便在此祝福向宴醬中考順利！

標籤：其他

上一篇：讀改變未來的九大演算法筆記03_糾錯碼

下一篇：返回列表