洛谷P5110 块速递推 题解

题目链接：P5110 块速递推

题意：给定一个数列 \(a\) 满足递推式 \[ a_0=0,a_1=1 \\a_n = 233a_{n-1}+666a_{n-2} \] 求 \(a_n \bmod (10^9+7)\)

多组询问

这个题是有个循环节的，正好是 \(10^9+6\) ，据出题人说是凑好的

关于循环节怎么求的我也不太清楚，先留个坑，研究好了就补上来

upd.20220721 现在知道怎么求循环节了，详见本文最后

也就是，如果数列 \(a_n\) 在模 \(M\) 意义下存在循环节 \(p\) ，则有 \(a_n \equiv a_{n\,\bmod\, p} \bmod M\)

本题的解法就是手推通项公式

具体方法如下

对于二阶线性递推数列 \[ a_0=A,a_1=B \\a_n = pa_{n-1}+qa_{n-2} ,n\ge 2 \] 考虑使用特征方程求解

\(a_n=pa_{n-1}+qa_{n-2}\) 的特征方程为 \[ x^2=px+q \] 可以求出两个特解（不一定是实数） \(x_1,x_2\) ，则 \[ a_n=\alpha x_1^{n} + \beta x_2^{n} \]

注意，如果数列从 \(a_1\) 开始，这里就是 \(a_n=\alpha x_1^{n-1} + \beta x_2^{n-1}\)

然后将 \(a_0=A,a_1=B\) 代入可得 \[ \begin{cases} \alpha + \beta = A \\\alpha x_1+\beta x_2=B \end{cases} \] 解出 \(\alpha,\beta\) 即可

本题的通项公式为 \[ a_n=\dfrac{1}{\sqrt{56953}}\left(\left(\dfrac{233+\sqrt{56953}}{2}\right)^n-\left(\dfrac{233-\sqrt{56953}}{2}\right)^n\right) \] > 关于特征方程法的更多内容，详见 数列通项的求解

注意到这里有个 \(\sqrt{56953}\) ，而我们要求它模意义下的值

也就是求出所有的 \(x\) \[ x^2\equiv 56953 \bmod(10^9+7) \] 考虑二次剩余求解

什么？不会二次剩余？ 那么就打个暴力好了，最多几秒钟就跑出来了

#include <bits/stdc++.h>
using namespace std;
#define int long long
signed main()
{
    for(int i=1; i<=1000000000; i++)
        if(i*i%1000000007==56953)cout << i << endl;
    return 0;
}

然后有两个解 \(188305837,811694170\) 取个小点的代入就好了

则有 \[ a_n \equiv 233230706 \times(94153035^n-905847205^n) \] 注意到询问有 \(10^7\) 个，快速幂过不了

考虑光速幂，\(O(\sqrt{n})\) 预处理 \(f(x)=x^{65536t},g(x)=x^t\)

询问直接查询 \(f(x/65536)\times g(n\%65536)\) 即可

这里可以用位运算加速

代码：

#include <bits/stdc++.h>
using namespace std;
#define int long long
#define uint unsigned long long
#define INF 0x3f3f3f3f3f3f3f3f
namespace Mker
{
	unsigned long long SA,SB,SC;
	void init(){scanf("%llu%llu%llu",&SA,&SB,&SC);}
	unsigned long long rand()
	{
	    SA^=SA<<32,SA^=SA>>13,SA^=SA<<1;
	    unsigned long long t=SA;
		SA=SB,SB=SC,SC^=t^SA;return SC;
	}
}
#define N (int)(7e4+15)
const int p=(int)(1e9+7);
int Q;
uint ans;
uint pw1[2][N],pw2[2][N];
const int a=94153035;
const int b=905847205;
void init()
{
    pw1[0][0]=pw2[0][0]=1;
    for(int i=1; i<1<<16; i++)
        pw1[0][i]=pw1[0][i-1]*a%p;
    pw2[0][1]=pw1[0][(1<<16)-1]*a%p;
    for(int i=2; i<1<<16; i++)
        pw2[0][i]=pw2[0][i-1]*pw2[0][1]%p;

    pw1[1][0]=pw2[1][0]=1;
    for(int i=1; i<1<<16; i++)
        pw1[1][i]=pw1[1][i-1]*b%p;
    pw2[1][1]=pw1[1][(1<<16)-1]*b%p;
    for(int i=2; i<=1<<16; i++)
        pw2[1][i]=pw2[1][i-1]*pw2[1][1]%p;
}
int pow1(int n)
{
    return pw1[0][n&65535]%p*pw2[0][n>>16]%p;
}
int pow2(int n)
{
    return pw1[1][n&65535]%p*pw2[1][n>>16]%p;
}
uint solve(int n)
{
    return 233230706*(pow1(n)-pow2(n)+p)%p;
}
signed main()
{
    init();
    scanf("%lld",&Q);
    Mker::init();
    while(Q--)
        ans^=solve(Mker::rand()%(p-1));
    printf("%llu\n",ans);
    return 0;
}
//233230706×(94153035^n−905847205^n) 

如何求解循环节？

如果数列 \(a_n\) 在模 \(M\) 意义下存在循环节 \(p\) ，则有 \(a_n \equiv a_{n\,\bmod\, p} \bmod M\)

问了 @zx2017 老师 CCOrz

这个东西其实求起来非常简单，

在已知 \(a_n\) 递推式和 \(M\) 的情况下

直接暴力去枚举就好了

具体的，这道题 \[ a_0=0,a_1=1 \\a_n = 233a_{n-1}+666a_{n-2} \] 注意到 \(a_n\) 是从 \(a_{n-1}\) 和 \(a_{n-2}\) 推来的

那我们枚举 \(n\) ，直到出现了这样的情况 \[ \begin{cases} a_i\equiv a_2 \\\\a_{i-1}\equiv a_1 \\\\a_{i-2}\equiv a_0 \end{cases} \mod M \] 那么 \(i-2\) 就是循环节

是不是很神奇？我太菜了根本没想到 QAQ

当然可能会有更优秀的做法（比起暴力枚举）

不过起码知道怎么简单求就很好了 qwq