洛谷 P3306 [SDOI2013]随机数生成器 题解【BSGS】【数列】

点击量：140

数据比较强的 BSGS 题，很多细节都卡到了√

目录

题目描述

小 W 喜欢读书，尤其喜欢读《约翰克里斯朵夫》。最近小 W 准备读一本新书，这本书一共有 $P$ 页，页码范围为 $0 \cdots P-1$。

小 W 很忙，所以每天只能读一页书。为了使事情有趣一些，他打算使用 NOI2012 上学习的线性同余法生成一个序列，来决定每天具体读哪一页。

我们用 $X_i$ 来表示通过这种方法生成出来的第 $i$ 个数，也即小 W 第 $i$ 天会读哪一页。这个方法需要设置 $3$ 个参数 $a,b,X_1$，满足 $0\le a,b,X_1\leq p-1$ ，且 $a,b,X_1$ 都是整数。按照下面的公式生成出来一系列的整数：$X_{i+1}\equiv aX_i+b \pmod p$ 其中 $\text{mod}$ 表示取余操作。

但是这种方法可能导致某两天读的页码一样。

小 W 要读这本书的第 $t$ 页，所以他想知道最早在哪一天能读到第 $t$ 页，或者指出他永远不会读到第 $t$ 页。

输入输出格式

输入格式：

输入含有多组数据，第一行一个正整数 $T$，表示这个测试点内的数据组数。

接下来 $T$ 行，每行有五个整数 $p$，$a$，$b$，$X_1$，$t$，表示一组数据。保证 $X_1$ 和 $t$ 都是合法的页码。 注意：$P$ 一定为质数。

输出格式：

共 $T$ 行，每行一个整数表示他最早读到第 $t$ 页是哪一天。如果他永远不会读到第 $t$ 页，输出 $-1$。

输入输出样例

输入样例#1：

3
7 1 1 3 3
7 2 2 2 0
7 2 2 2 1

输出样例#1：

1 
3 
-1

说明

$0\le a\le P-1,0\le b\le P-1,2\le P\le 10^9$。

题解：

一开始感觉要处理多项式，因为要迭代很多层 $AX_i+B$。

不过这个式子还是很好化的。

令 $x=X_1$，则 $X_1=x$，$X_2=Ax+B$，$X_3=AX_2+B=A(Ax+B)+B=A^2x+AB+B$。

可以发现，在这个迭代中，$x$ 永远只有一个一次项。继续往下化可以发现更多规律。

$X_4=AX_3+B=A(A^2x+AB+B)+B=A^3x+A^2B+AB+B$。

此时断言 $X_n=A^{n-1}x+A^{n-2}B+A^{n-3}B+\cdots+B$。

具体证明可以参考数学归纳法。

此时还没有什么优秀的表达。把后面的 $B$ 提出来，得到 $X_n=A^{n-1}x+B(A^{n-2}+A^{n-3}+\cdots+1)$。

后面的括号中是一个等比数列，用等比求和可以构造出一个 $A^{n-1}$，这样就能和前面的合并了。
$$
\begin{aligned}
X_n&=A^{n-1}x+B\cdot\frac{1-A^{n-1}}{1-A}\\
&=A^{n-1}x+\frac{B}{1-A}-\frac{A^{n-1}B}{1-A}\\
&=A^{n-1}(x-\frac{B}{1-A})+\frac{B}{1-A}
\end{aligned}
$$
现在我们需要解的是 $n$，式子中也恰好只有 $n$ 这一个未知数。

移项之后得到 $X_n+\frac{B}{A-1}=(x+\frac{B}{A-1})A^{n-1}$。

但是这个题有一些细节需要注意。

在上面的式子中，分母是 $A-1$，在实数意义下计算时分母不能为 $0$；同样地，在模意义下计算时如果分母为 $0$ 也会产生很多错误。

同样，当 $A=1$ 时，等比数列的公比是 $0$，也不能用这个公式。

所以在读入时特判 $A=1$ 的情况，$X_n=x+(n-1)B$。

直接解出 $(n-1)$ 的值然后再 $+1$ 即可。这里的 $n-1$ 是在 BSGS 的解的范围 $[0,p-1]$ 内的，但是 $+1$ 后不能模 $p$，因为是从第一天开始算的，所以答案应该 $\ge 1$。

时间复杂度 $O(\sqrt p)$。

Code：

#include<cmath>
#include<cstdio>
#include<cstring>
#define gc() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++)
char buf[2100000],*p1=buf,*p2=buf;
int read()
{
    int x=0;
    char ch=gc();
    while(ch<'0'||ch>'9')
        ch=gc();
    while(ch>='0'&&ch<='9')
    {
        x=x*10+(ch&15);
        ch=gc();
    }
    return x;
}
#define ll long long
ll qpow(ll x,ll y,ll p)
{
    ll ans=1;
    while(y)
    {
        if(y&1)
            ans=ans*x%p;
        x=x*x%p;
        y>>=1;
    }
    return ans;
}
struct edge
{
    int x,nxt,v;
}e[100000];
int head[10007],ecnt=-1;
void add(int x,int v)
{
    e[++ecnt]=(edge){x,head[x%10007],v};
    head[x%10007]=ecnt;
}
int Find(int x)
{
    for(int i=head[x%10007];~i;i=e[i].nxt)
        if(e[i].x==x)
            return e[i].v;
    return 0;
}
int main()
{
#ifdef wjyyy
    freopen("a.in","r",stdin);
#endif
    int T=read(),p,a,b,X,t;
    while(T--)
    {
        p=read();
        a=read();
        b=read();
        X=read();
        t=read();
        if(t==X)
        {
            puts("1");
            continue;
        }
        if(a==0)
        {
            if(t==b)
                puts("2");
            else
                puts("-1");
            continue;
        }
        else if(a==1)
        {
            if(b==0)
                puts("-1");
            else
                printf("%lld\n",qpow(b,p-2,p)*(t-X+p)%p+1);
            continue;
        }
        for(int i=0;i<10007;++i)
            head[i]=-1;
        ecnt=-1;
        int inv=qpow(a-1,p-2,p)*b%p,sq=sqrt(p),flag=0;
        t=(t+inv)%p;
        X=(X+inv)%p;
        ll f=a;
        for(int i=1;i<=sq;++i,f=f*a%p)
            add(f*t%p,i);
        f=f*qpow(a,p-2,p)%p;
        ll g=f*X%p,k;
        for(int q=sq;q-sq<p;q+=sq,g=g*f%p)
            if(k=Find(g))
            {
                printf("%lld\n",q-k+1);
                flag=1;
                break;
            }
        if(!flag)
            puts("-1");
    }
    return 0;
}