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

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

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

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

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

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

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

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

$(x_{i+1}+K)=a\cdot (x_{i}+K) (mod\;p)$

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

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

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

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

$a==0$的情况特判

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

  1 #include <cmath>
  2 #include <cstdio>
  3 #include <cstring>
  4 #include <algorithm>
  5 #define N1 400010
  6 #define M1 40010
  7 #define ll long long
  8 #define dd double
  9 #define cint const int 
 10 #define cll const long long 
 11 #define inf 23333333333333333ll
 12 using namespace std;
 13 
 14 int gcd(int a,int b){ if(!b) return a; return gcd(b,a%b); }
 15 void exgcd(ll a,ll b,ll &x,ll &y)
 16 {
 17     if(!b){ x=1; y=0; return; }
 18     exgcd(b,a%b,x,y); ll t=x; x=y; y=t-a/b*y;
 19 }
 20 int P,A,B,X1,T,K;
 21 namespace S0{
 22 int solve()
 23 {
 24     int ans,g; ll inv,invy; T=(T-X1+P)%P;
 25     g=gcd(B,P); if(T%g!=0) return -1;
 26     T/=g; B/=g; P/=g; 
 27     exgcd(B,P,inv,invy); inv=(inv%P+P)%P;
 28     ans=inv*T%P; P*=g; ans=1ll*ans*g%P;
 29     return ans+1; 
 30 }
 31 };
 32 namespace S1{
 33 #define maxn 400000
 34 struct Hsh{
 35 int head[N1],to[M1],nxt[M1],val[M1],cte;
 36 void ins(int x,int w)
 37 {
 38     int u=x%maxn,j,v;
 39     for(j=head[u];j;j=nxt[j])
 40     {
 41         v=to[j];
 42         if(v==x) return;
 43     }
 44     cte++; to[cte]=x; nxt[cte]=head[u];
 45     head[u]=cte; val[cte]=w;
 46 }
 47 int find(int x)
 48 {
 49     int u=x%maxn,j,v;
 50     for(j=head[u];j;j=nxt[j])
 51     {
 52         v=to[j];
 53         if(v==x) return val[j];
 54     }
 55     return -1;
 56 }
 57 }h;
 58 ll qpow(ll x,ll y,cll &P)
 59 {
 60     ll ans=1;
 61     while(y){
 62         if(y&1) ans=ans*x%P;
 63         x=x*x%P; y>>=1;
 64     }return ans;
 65 }
 66 #undef maxn
 67 int solve(int p)
 68 {
 69     int i,j,sq;ll K,pw,now,ans,inv,invy,Y,tmp; cll P=p;
 70     exgcd(A-1,P,inv,invy); inv=(inv%P+P)%P; K=inv*B%P; 
 71     T=(T+K)%P; X1=(X1+K)%P;
 72     exgcd(X1,P,inv,invy); inv=(inv%P+P)%P; T=T*inv%P;
 73     sq=sqrt(P); 
 74     for(pw=qpow(A,sq,P),now=1,i=1;(i-1)*sq<P;i++)
 75     {
 76         now=1ll*now*pw%P;
 77         h.ins(now,i);
 78     }
 79     for(ans=inf,now=T,i=0;i<sq;i++)
 80     {
 81         tmp=h.find(now);
 82         if(tmp!=-1) ans=min(ans,tmp*sq-i);
 83         now=now*A%P;
 84     }
 85     memset(&h,0,sizeof(h));
 86     if(ans==inf) return -1;
 87     return ans+1;
 88 }
 89 };
 90 int Case;
 91 
 92 int main()
 93 {
 94     scanf("%d",&Case);
 95     while(Case--)
 96     {
 97         scanf("%d%d%d%d%d",&P,&A,&B,&X1,&T);
 98         if(X1==T){ puts("1"); continue; }
 99         if(A==0)
100         {
101             if(T==B) puts("2");
102             else puts("-1");
103             continue;
104         }
105         if(A==1) printf("%d\n",S0::solve());
106         else printf("%d\n",S1::solve(P));
107     }
108     return 0;
109 }

 

posted @ 2019-01-10 18:58  guapisolo  阅读(204)  评论(1编辑  收藏  举报