BZOJ 3122 [SDOI2013]随机数生成器 (BSGS)

题目大意：略 洛谷传送门 BZOJ传送门

观察式子$x_{i+1}=acdot x_{i}+b (mod;p)$

我们能用欧几里得算法求解$ax=c(mod;b)$的线性同余方程问题

还能用$BSGS$算法求解$a^{x}=c(mod;b)$的指数同余方程问题

而式子里的$b$这一项显得十分棘手

。

。

。

还记得等比数列吗，是数学老师教的

考虑构造一个"等比数列”方便我们递推

$(x_{i+1}+K)=acdot (x_{i}+K) (mod;p)$

解得$K=frac{b}{a-1}$

这样，原式化成了$A^{X}=B(mod;p)$的形式

好像求出指数$x$就行了？

于是乎这道题就变成了一个裸的$BSGS$问题？

$a==0$的情况特判

$a==1?$原式变成了等差数列，用$exgcd$求解线性同余方程即可

#include <cmath>
#include <cstdio>
#include <cstring>
#include <algorithm>
#define N1 400010
#define M1 40010
#define ll long long
#define dd double
#define cint const int 
#define cll const long long 
#define inf 23333333333333333ll
using namespace std;

int gcd(int a,int b){ if(!b) return a; return gcd(b,a%b); }
void exgcd(ll a,ll b,ll &x,ll &y)
{
    if(!b){ x=1; y=0; return; }
    exgcd(b,a%b,x,y); ll t=x; x=y; y=t-a/b*y;
}
int P,A,B,X1,T,K;
namespace S0{
int solve()
{
    int ans,g; ll inv,invy; T=(T-X1+P)%P;
    g=gcd(B,P); if(T%g!=0) return -1;
    T/=g; B/=g; P/=g; 
    exgcd(B,P,inv,invy); inv=(inv%P+P)%P;
    ans=inv*T%P; P*=g; ans=1ll*ans*g%P;
    return ans+1; 
}
};
namespace S1{
#define maxn 400000
struct Hsh{
int head[N1],to[M1],nxt[M1],val[M1],cte;
void ins(int x,int w)
{
    int u=x%maxn,j,v;
    for(j=head[u];j;j=nxt[j])
    {
        v=to[j];
        if(v==x) return;
    }
    cte++; to[cte]=x; nxt[cte]=head[u];
    head[u]=cte; val[cte]=w;
}
int find(int x)
{
    int u=x%maxn,j,v;
    for(j=head[u];j;j=nxt[j])
    {
        v=to[j];
        if(v==x) return val[j];
    }
    return -1;
}
}h;
ll qpow(ll x,ll y,cll &P)
{
    ll ans=1;
    while(y){
        if(y&1) ans=ans*x%P;
        x=x*x%P; y>>=1;
    }return ans;
}
#undef maxn
int solve(int p)
{
    int i,j,sq;ll K,pw,now,ans,inv,invy,Y,tmp; cll P=p;
    exgcd(A-1,P,inv,invy); inv=(inv%P+P)%P; K=inv*B%P; 
    T=(T+K)%P; X1=(X1+K)%P;
    exgcd(X1,P,inv,invy); inv=(inv%P+P)%P; T=T*inv%P;
    sq=sqrt(P); 
    for(pw=qpow(A,sq,P),now=1,i=1;(i-1)*sq<P;i++)
    {
        now=1ll*now*pw%P;
        h.ins(now,i);
    }
    for(ans=inf,now=T,i=0;i<sq;i++)
    {
        tmp=h.find(now);
        if(tmp!=-1) ans=min(ans,tmp*sq-i);
        now=now*A%P;
    }
    memset(&h,0,sizeof(h));
    if(ans==inf) return -1;
    return ans+1;
}
};
int Case;

int main()
{
    scanf("%d",&Case);
    while(Case--)
    {
        scanf("%d%d%d%d%d",&P,&A,&B,&X1,&T);
        if(X1==T){ puts("1"); continue; }
        if(A==0)
        {
            if(T==B) puts("2");
            else puts("-1");
            continue;
        }
        if(A==1) printf("%d
",S0::solve());
        else printf("%d
",S1::solve(P));
    }
    return 0;
}