BSGS及其拓展

BSGS

介绍

这是一个求解 a x ≡ b ( m o d p ) a ^ {x} \equiv b \pmod p axb(modp),的方法。并且 p p p是质数, a , p a, p a,p互质,费马小定理可知,这个式子有周期性,

我们一般取 m = s q r t ( p ) m = sqrt(p) m=sqrt(p),假设 x = i ∗ m + j , 0 < = i , j < = m x = i * m + j, 0 <= i, j <= m x=im+j0<=i,j<=m,则有

a i ∗ m + j ≡ b ( m o d p ) a ^ {i * m + j} \equiv b \pmod p aim+jb(modp)

a i + m ≡ b a j ( m o d p ) a ^ {i + m} \equiv \frac {b} {a ^ j} \pmod p ai+majb(modp)

为了方便我们设 x = i ∗ m − j x = i * m - j x=imj,则有

a i ∗ m ≡ b ∗ a j ( m o d p ) a ^ {i * m} \equiv b * a ^ {j} \pmod p aimbaj(modp)

所以我们只要通过一次枚举 j j j,记录下 b ∗ a j b * a ^ j baj,再一次枚举 i i i去刚才枚举过的 j j j中寻找有没有符合要求的答案即可,整体复杂度是 p \sqrt {p} p

P2485 [SDOI2011]计算器

模板题

/*
  Author : lifehappy
*/
#pragma GCC optimize(2)
#pragma GCC optimize(3)
#include <bits/stdc++.h>

// #include <cstdio>
// #include <iostream>
// #include <stdlib.h>
// #include <algorithm>
// #include <cmath>


#define mp make_pair
#define pb push_back
#define endl '\n'

using namespace std;

typedef long long ll;
typedef unsigned long long ull;
typedef pair<int, int> pii;

const double pi = acos(-1.0);
const double eps = 1e-7;
const int inf = 0x3f3f3f3f;

inline ll read() {
    ll f = 1, x = 0;
    char c = getchar();
    while(c < '0' || c > '9') {
        if(c == '-')    f = -1;
        c = getchar();
    }
    while(c >= '0' && c <= '9') {
        x = (x << 1) + (x << 3) + (c ^ 48);
        c = getchar();
    }
    return f * x;
}

void print(ll x) {
    if(x < 10) {
        putchar(x + 48);
        return ;
    }
    print(x / 10);
    putchar(x % 10 + 48);
}

ll quick_pow(ll a, ll n, ll mod) {
    ll ans = 1;
    while(n) {
        if(n & 1) ans = (ans * a) % mod;
        n >>= 1;
        a = (a * a) % mod;
    }
    return ans;
}

ll quick_mult(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;
}

ll ex_gcd(ll a, ll b, ll & x, ll & y) {
    if(!b) {
        x = 1, y = 0;
        return a;
    }
    ll gcd = ex_gcd(b, a % b, x, y);
    ll temp = x;
    x = y;
    y = temp - a / b * y;
    return gcd;
}


void BSGC(ll a, ll b, ll p) {
    map<ll, ll> MP;
    int m = sqrt(p) + 1;
    ll x = b, nex = quick_pow(a, m, p);
    for(int i = 0; i <= m; i++) {
        MP[x] = i;
        x = (x * a) % p;
    }
    x = 1;
    if(nex == 0) {
        if(b == 0) {
            puts("0");
        }
        else {
            puts("Orz, I cannot find x!");
        }
        return ;
    }
    for(int i = 0; i <= m; i++) {
        if(MP.count(x)) {
            printf("%d\n", i * m - MP[x]);
            return ;
        }
        x = (x * nex) % p;
    }
    puts("Orz, I cannot find x!");
}

int main() {
    // freopen("in.txt", "r", stdin);
    // freopen("out.txt", "w", stdout);
    // ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
    int n = read(), k = read();
    for(int i = 1; i <= n; i++) {
        ll a = read(), b = read(), p = read();
        if(k == 1) {
            printf("%lld\n", quick_pow(a, b, p));
        }
        else if(k == 2) {
            ll x, y;
            ll gcd = ex_gcd(a, p, x, y);
            if(b % gcd) {
                puts("Orz, I cannot find x!");
            }
            else {
                b /= gcd;
                x = (((x % p + p) % p) * b) % p;
                printf("%lld\n", x);
            }
        }   
        else {
            BSGC(a, b % p, p);
        }
    }
	return 0;
}

2019牛客暑期多校训练营(第五场)C generator 2

思路

x 0 = x 0 x_0 = x_0 x0=x0

x 1 = a ∗ x 0 ∗ b x_1 = a * x_0 * b x1=ax0b

x 2 = a ∗ x 1 + b = a 2 ∗ x 0 + a ∗ b + b x_2 = a * x_1 + b = a ^{2} * x_0 + a * b + b x2=ax1+b=a2x0+ab+b

容易发现后项是一个等比数列求和

x n = a n x 0 + b ( 1 − a n ) 1 − a x_n = a ^ {n} x_0 + \frac {b (1 - a ^ n)} {1 - a} xn=anx0+1ab(1an)

我们要求 x n = v x_n = v xn=v,化简

a n x 0 + b ( 1 − a n ) 1 − a = v a ^ {n} x_0 + \frac {b (1 - a ^ n)} {1 - a} = v anx0+1ab(1an)=v

a n x 0 ( 1 − a ) + b ( 1 − a n ) = v ( 1 − a ) a ^{n} x_0(1 - a) + b (1 - a ^ n) = v(1 - a) anx0(1a)+b(1an)=v(1a)

a n ( x 0 − a x 0 − b ) = v ( 1 − a ) − b a ^ n (x_0 - ax_0 - b) = v(1 - a) - b an(x0ax0b)=v(1a)b

a n = v ( 1 − a ) − b x 0 − a x 0 − b a^n =\frac {v (1 - a) - b} {x_0 - a x_0 - b} an=x0ax0bv(1a)b

上面式子都是 m o d    p \mod p modp下的同余等式,为了方便写了 = = =

看到这里就简单了,我们要求解的是 n n n,显然右边这一坨都是已知的,我们假定为 B B B,求解 a n = B a ^ n = B an=B,这不就是个裸题了吗。

这道题目还要稍加分类讨论一下:

  • a == 0

x 0 = x 0 , x i = b ( i > = 1 ) x_0 = x_0, x_i = b(i >= 1) x0=x0,xi=b(i>=1)

  • a = = 1 a == 1 a==1

因为这种情况上面不能直接相除,所以我们需要特殊考虑 a i = b + i ∗ a a_i = b + i * a ai=b+ia

也就是求解 i ∗ a = v − b i * a = v - b ia=vb,这个时候只要左右两边同时乘以 a a a的逆元即可得到我们要的 i i i

代码

/*
  Author : lifehappy
*/
#pragma GCC optimize(2)
#pragma GCC optimize(3)
#include <bits/stdc++.h>
#define mp make_pair
#define pb push_back
#define endl '\n'

using namespace std;

typedef long long ll;
typedef unsigned long long ull;
typedef pair<int, int> pii;

const double pi = acos(-1.0);
const double eps = 1e-7;
const int inf = 0x3f3f3f3f;

inline ll read() {
    ll f = 1, x = 0;
    char c = getchar();
    while(c < '0' || c > '9') {
        if(c == '-') f = -1;
        c = getchar();
    }
    while(c >= '0' && c <= '9') {
        x = (x << 1) + (x << 3) + (c ^ 48);
        c = getchar();
    }
    return f * x;
}

void print(ll x) {
    if(x < 10) {
        putchar(x + 48);
        return ;
    }
    print(x / 10);
    putchar(x % 10 + 48);
}

ll quick_pow(ll a, ll n, ll mod) {
    ll ans = 1;
    while(n) {
        if(n & 1) ans = (ans * a) % mod;
        a = (a * a) % mod;
        n >>= 1;
    }
    return ans;
}

int main() {
    // freopen("in.txt", "r", stdin);
    // freopen("out.txt", "w", stdout);
    // ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
    int T = read();
    while(T--) {
        ll n = read(), x0 = read(), a = read(), b = read(), p = read();
        ll x = 1, Unit = quick_pow(a, 1000, p);
        unordered_map<ll, int> MP;
        for(int i = 1; i <= 1000000; i++) {
            x = (x * Unit) % p;
            if(!MP.count(x)) {
                MP[x] = i * 1000;
            }
        }
        ll inv = quick_pow(((x0 - a * x0 - b) % p + p) % p, p - 2, p);
        int t = read();
        while(t--) {
            ll v = read();
            if(a == 0) {
                if(v % p == x0 % p) {
                    puts("0");
                }
                else if(v % p == b % p) {
                    puts("1");
                }
                else {
                    puts("-1");
                }
                continue;
            }
            if(a == 1) {
                ll ans = (((((v - x0) % p + p) % p) * quick_pow(b, p - 2, p))) % p;
                if(ans < n) {
                    printf("%lld\n", ans);
                }
                else {
                    puts("-1");
                }
                continue;
            }
            v = (((v * (1 - a) - b) % p + p) % p * inv) % p;
            if(Unit == 0) {
                if(v == 0) {
                    puts("0");
                }
                else {
                    puts("-1");
                }
                continue;
            }
            x = v;
            int ans = p + 1;
            for(int i = 0; i <= 1000; i++) {
                if(MP.count(x)) {
                    ans = min(ans, MP[x] - i);
                }
                x = (x * a) % p;
            }
            if(ans != p + 1 && ans < n) {
                printf("%d\n", ans);
            }
            else {
                puts("-1");
            }
        }
    }
	return 0;
}

BSGS拓展

介绍

求解 a x ≡ b ( m o d p ) a ^ x \equiv b \pmod p axb(modp),但是 a , b a, b a,b不互质,这里就要用到我们的 E X B S G S EXBSGS EXBSGS了。

假定 a x + p y = b a ^ x + py = b ax+py=b

d 1 = g c d ( a , p ) , 如 果 有 解 则 一 定 : g c d ( a , p ) ∣ b d_1 = gcd(a, p),如果有解则一定:gcd(a, p) \mid b d1=gcd(a,p)gcd(a,p)b

得到 a x − 1 a d 1 + p d 1 y = b d 1 a ^{x - 1}\frac{a} {d_1} + \frac{p}{d_1}y = \frac {b} {d_1} ax1d1a+d1py=d1b

如果 g c d ( a , p d 1 ) ! = 1 gcd(a, \frac{p} {d_1}) != 1 gcd(a,d1p)!=1,继续化简

d 2 = g c d ( a , p d 1 ) − > a x − 2 a 2 d 1 d 2 + p d 1 d 2 y = b d 1 d 2 d_2 = gcd(a, \frac{p} {d_1})->a ^{x - 2} \frac{a ^ 2} {d_1d_2} + \frac{p} {d_1d_2}y = \frac{b}{d_1d_2} d2=gcd(a,d1p)>ax2d1d2a2+d1d2py=d1d2b

如果 g c d ( a , p 2 d 1 d 2 ) ! = 1 gcd(a, \frac{p ^ 2} {d_1d_2}) != 1 gcd(a,d1d2p2)!=1

重复上面操作,最后得到式子

a x − n a n d 1 d 2 … … d n + p d 1 d 2 … … d n y = b d 1 d 2 … … d n a ^{x - n} \frac {a ^ n} {d_1d_2……d_n} + \frac{p}{d_1d_2……d_n}y = \frac{b}{d_1d_2……d_n} axnd1d2dnan+d1d2dnpy=d1d2dnb

A = a x − n , A ′ = a n d 1 d 2 … … d n , P = p d 1 d 2 … … d n , B = b d 1 d 2 … … d n A = a {x - n}, A' = \frac{a^n}{d_1d_2……d_n}, P = \frac{p}{d_1d_2……d_n}, B = \frac{b} {d_1d_2……d_n} A=axn,A=d1d2dnan,P=d1d2dnp,B=d1d2dnb

则变成里求 A x − n ≡ B A ′ − 1 ( m o d P ) A^{x - n} \equiv B A'^{-1} \pmod{P} AxnBA1(modP)

记录进行了多少次 g c d gcd gcd求解,通过求 A ′ A' A的逆元化简式子,再通过一次 B S G S BSGS BSGS求得 x − n x - n xn,即可得到我们得答案 x x x

P4195 【模板】扩展BSGS

一定注意求逆元不能用费马小定理,我就入了这个坑。

/*
    Author : lifehappy
*/
#include <bits/stdc++.h>

using namespace std;

typedef long long ll;

ll gcd(ll a, ll b) {
    return b ? gcd(b, a % b) : a;
}

ll quick_pow(ll a, ll n, ll mod) {
    ll ans = 1;
    while(n) {
        if(n & 1) ans = (ans * a) % mod;
        a = (a * a) % mod;
        n >>= 1;
    }
    return ans;
}

int exgcd(int a, int b, int & x, int & y) {
    if(!b) {
        x = 1, y = 0;
        return a;
    }
    int gcd = exgcd(b, a % b, x, y);
    int temp = x;
    x = y;
    y = temp - a / b * y;
    return gcd;
}

int inv(int a, int b) {
    int x, y;
    exgcd(a, b, x, y);
    return (x % b + b) % b;
}

int BSGS(ll a, ll b, ll p) {
    int m = sqrt(p) + 1;
    unordered_map<int, int> mp;
    int x = b;
    for(int i = 0; i <= m; i++) {
        mp[x] = i;
        x = (1ll * x * a) % p;
    }
    x = 1;
    int Unit = quick_pow(a, m, p);
    if(Unit == 0) {
        if(b == 0) {
            return 0;
        }
        else {
            return -1;
        }
    }
    for(int i = 1; i <= m; i++) {
        x = (1ll * x * Unit) % p;
        if(mp.count(x)) {
            return i * m - mp[x];
        }
    }
    return -1;
}

void EXBSGS(ll a, ll b, ll p) {
    a %= p, b %= p;
    if(b == 1 || p == 1) {
        puts("0");
        return ;
    }
    int cnt = 0, d, ad = 1;
    while((d = gcd(a, p)) != 1) {
        if(b % d) {
            puts("No Solution");
            return ;
        }
        cnt++;
        b /= d, p /= d;
        ad = (1ll * ad * a / d) % p;
        if(ad == b) {
            printf("%d\n", cnt);
            return ;
        }
    }
    int ans = BSGS(a % p, (1ll * b * inv(ad, p)) % p, p);
    if(ans == -1) {
        puts("No Solution");
    }
    else {
        printf("%d\n", ans + cnt);
    }
}

int main() {
    // freopen("in.txt", "r", stdin);
    ll a, b, p;
    while(scanf("%lld %lld %lld", &a, &p, &b) && (a || b || p)) {
        EXBSGS(a, b, p);
    }
    return 0;
}
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值