[BSGS]大步小步算法

問題

BSGS被用於求解離散對數,即同餘方程:
\[ A^x\equiv B\pmod{P} \]php

\(x\)的最小非負整數解。ios

保證\(A\perp P\)(互質)。算法

分析

首先,咱們根據費馬小定理,有
\[ A^{P-1}\equiv 1\pmod{P} \]ui

則顯然有
\[ A^{x-k(P-1)}\equiv A^x\pmod{P} \]spa


\[ A^{x\mod{P-1}}\equiv A^x\pmod{P} \]code

那麼顯然\(x<P-1\),咱們就獲得了一個\(O(P)\)的算法,然而太慢了。get

考慮分塊算法,對\(x\)\(m\)分一塊,則有
\[ A^{im-j}\equiv B\pmod{P} \]io

移項整理
\[ \left(A^m\right)^i\equiv A^j B\pmod{P} \]class

那麼咱們枚舉\(i\),就能夠求出\(A^j\)。再對於\(j\in[0,m-1]\)\(A^j\)存進哈希表/map,就能夠獲得\(x=im-j\)了。若是不考慮查詢哈希表/map的時間,則時間複雜度爲\(O(m+\frac{P}{m})\)stream

\(m\)應該取何值呢?求\(f(m)=m+\frac{P}{m}\)的駐點:
\[ \frac{\mathbb{d}f(m)}{\mathbb{d} m}=0 \]


\[ 1-\frac{P}{m^2}=0 \]

移項整理
\[ m^2=P \]

解得\(m=\sqrt{P}\)

那麼咱們令\(m=\lceil\sqrt{P}\rceil\),就獲得了一個\(O(\sqrt{P})\)的算法。

代碼

\(-1\)爲無解。

ll BSGS(ll a,ll b,ll p){
    if(!a)return b?-1:1;
    if(b==1)return 0;
    map<ll,ll>mp;
    ll m=ceil(sqrt(p)),ax=1;
    for(int i=0;i<m;i++){
        mp[ax]=i;
        ax=ax*a%p;
    }
    ll am=pow(a,m,p),aj=am*pow(b,p-2,p)%p;
    for(int i=1;i<=m;i++){
        if(mp.count(aj))return m*i-mp[aj];
        aj=aj*am%p;
    }
    return -1;
}

例題

[BZOJ2242][SDOI2011]計算器

#include<iostream>
#include<cstdio>
#include<cmath>
#include<map>
using namespace std;
typedef long long ll;
int t,k;
ll y,z,p;
ll pow(ll a,ll b,ll p){
    ll ans=1;
    while(b){
        if(b&1)ans=ans*a%p;
        a=a*a%p;
        b>>=1;
    }
    return ans;
}
ll BSGS(ll a,ll b,ll p){
    if(!a)return b?-1:1;
    if(b==1)return 0;
    map<ll,ll>mp;
    ll m=ceil(sqrt(p)),ax=1;
    for(int i=0;i<m;i++){
        mp[ax]=i;
        ax=ax*a%p;
    }
    ll am=pow(a,m,p),aj=am*pow(b,p-2,p)%p;
    for(int i=1;i<=m;i++){
        if(mp.count(aj))return m*i-mp[aj];
        aj=aj*am%p;
    }
    return -1;
}
int main(){
    scanf("%d%d",&t,&k);
    while(t--){
        scanf("%lld%lld%lld",&y,&z,&p);
        if(k==1)printf("%lld\n",pow(y,z,p));
        else if(k==2){
            if(y%p==0)printf("Orz, I cannot find x!\n");
            else printf("%lld\n",pow(y,p-2,p)*z%p);
        }else{
            ll ans=BSGS(y%p,z%p,p);
            if(~ans)printf("%lld\n",ans);
            else printf("Orz, I cannot find x!\n");
        }
    }
}
相關文章
相關標籤/搜索