洛谷P5110 块速递推 题解

题目链接：P5110 块速递推

题意：给定一个数列 $a$ 满足递推式
$$\begin{gathered} a_0=0,a_1=1 \\ {a}_n=233a_{n-1}+666a_{n-2} \end{gathered}$$
求 $a_n\mathrm{m}\mathrm{o}\mathrm{d}(10^9+7)$

多组询问

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

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

也就是，如果数列 $a_n$ 在模 $M$ 意义下存在循环节 $p$ ，则有 $a_n\equiv a_{n\mathrm{m}\mathrm{o}\mathrm{d}p}\mathrm{m}\mathrm{o}\mathrm{d}M$

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

具体方法如下

对于二阶线性递推数列
$$\begin{gathered} a_0=A,a_1=B \\ {a}_n=p{a}_{n-1}+q{a}_{n-2},n\ge 2 \end{gathered}$$
考虑使用特征方程求解

$a_n=p{a}_{n-1}+q{a}_{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{array}{l}\alpha +\beta =A\\ \alpha x_1+\beta x_2=B\end{array}$$
解出 $\alpha ,\beta$ 即可

本题的通项公式为
$$a_n=\frac{1}{\sqrt{56953}}\left({\left(\frac{233+\sqrt{56953}}{2}\right)}^n-{\left(\frac{233-\sqrt{56953}}{2}\right)}^n\right)$$
注意到这里有个 $\sqrt{56953}$ ，而我们要求它模意义下的值

也就是求出所有的 $x$
$$x^2\equiv 56953\mathrm{m}\mathrm{o}\mathrm{d}(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) 

转载请说明出处