[MtOI2019]幻想乡数学竞赛 解题报告
disangan233
2019-04-20 02:16:07
## 题意
存在一个数列$\{ a_n\} (n\in \{ 0,1,2,\cdots ,10^{18}\cdots \} )$。
已知$a_0=-3,a_1=-6,a_2=-12,a_n=3a_{n-1}+a_{n-2}-3a_{n-3}+3^n$。
* 现在给你给定的$n$，令$p=10^{9}+7$，请你求出$a_n \bmod p$。
## Solutions
### Sol 1
输出0，期望得分：$5~pts$。
### Sol 2
直接暴力递推，时间复杂度$O(tn)$，期望得分：$20~pts$。
### Sol 3
考虑将递推式转化为矩阵乘法：
$$\begin{bmatrix} a_{i-1}\\a_{i-2}\\a_{i-3}\\3^i \end{bmatrix} \times \begin{bmatrix} 3&1&-3&1\\1&0&0&0\\0&1&0&0\\0&0&0&3\end{bmatrix} = \begin{bmatrix} a_{i}\\a_{i-1}\\a_{i-2}\\3^{i+1} \end{bmatrix}$$
因为矩阵满足结合律，快速幂即可，时间复杂度$O\left(4^3 \times t\log {n}\right)$
* 加上对于数据点12的打表，期望得分：$20-60~pts$
### Sol 4
~~由题目名字得知这是一个数学题。~~
构造关于$x$的OGF：$f(x)=a_0x^0+a_1x^1+a_2x^2+a_3x^3+\cdots +a_nx^n+\cdots $，
那么可以得到：
$$(3x^3-x^2-3x+1)f(x)=a_0x^0+(a_1-3a_0)x^1+(a_2-3a_1-a_0)x^2+3^3x^3+\cdots+3^nx^n+\cdots$$
将$a_0=-3,a_1=-6,a_2=-12$代入后得：
$$(3x^3-x^2-3x+1)f(x)=-3+3^1x^1+3^2x^2+3^3x^3+\cdots+3^nx^n+\cdots$$
因式分解后得到：
$$(3x-1)(x+1)(x-1)f(x)=-3+3^1x^1+3^2x^2+3^3x^3+\cdots+3^nx^n+\cdots$$
用等比数列求和公式化简：
$$(3x-1)(x+1)(x-1)f(x)=-4+\dfrac{1}{1-3x}=\dfrac{12x-3}{1-3x} $$
$$f(x)=\dfrac{12x-3}{(1-3x)^2(1+x)(1-x)}$$
用待定系数法把这个式子拆开：
$$\dfrac{12x-3}{(1-3x)^2(1+x)(1-x)}=\dfrac{A}{(1-3x)^2}+\dfrac{B}{1+x}+\dfrac{C}{1-x}+\dfrac{D}{1-3x} $$
$$A=\dfrac{12x-3}{(1+x)(1-x)} \Bigg|_{x=\frac{1}{3}}=\dfrac{9}{8} $$
$$B=\dfrac{12x-3}{(1-3x)^2(1-x)} \Bigg|_{x=-1}=-\dfrac{15}{32} $$
$$C=\dfrac{12x-3}{(1-3x)^2(1+x)} \Bigg|_{x=1}=\dfrac{9}{8}$$
D单独求：
$$\because \dfrac{(12x-3)x}{(1-3x)^2(1+x)(1-x)}=\dfrac{Ax}{(1-3x)^2}+\dfrac{Bx}{1+x}+\dfrac{Cx}{1-x}+\dfrac{Dx}{1-3x} $$
$$\therefore \lim_{x \to \infty} (xf(x))=0 $$
$$\therefore \lim_{x \to \infty} (B-C-\dfrac{1}{3}D)=0 $$
$$\therefore D=3\times(B-C)=-\dfrac{153}{32}$$
将求得的$A$，$B$，$C$，$D$代入$f(x)$中：
$$\because A=\dfrac{9}{8},B=-\dfrac{15}{32},C=\dfrac{9}{8},D=-\dfrac{153}{32} $$
$$\therefore f(x)=\dfrac{9}{8(1-3x)^2}-\dfrac{15}{32(1+x)}+\dfrac{9}{8(1-x)}-\dfrac{153}{32(1-3x)} $$
$$\therefore f(x)=\dfrac{9}{8}\sum^{\infty}_{i=0}x^i-\dfrac{15}{32}\sum^{\infty}_{i=0}(-1)^ix^i-\dfrac{153}{32}\sum^{\infty}_{i=0}3^ix^i+\dfrac{9}{8}\sum^{\infty}_{i=0}3^ix^i(i+1) $$
$$\therefore f(x)=\dfrac{1}{32} \sum^{\infty}_{i=0} \left[ 3^{i+2}\times (4i-13)+36-15\times (-1)^i \right] x^i $$
$$\therefore a_n=\dfrac{3^{n+2}\times (4n-13)+36-15\times (-1)^n}{32}$$
快速幂即可，时间复杂度$O\left(t\log {n}\right)$，期望得分$30-60~pts$。
### Sol 5
分析算法后发现，**Sol 3**因为**常数巨大**而不够优秀，然而**Sol 4**已经**足够优秀**却仍然超时。
发现**Sol 4**中最烧时间的就是求$3^{n+2} \bmod p$的部分。
* 对，你想到了什么？欧拉定理！
**欧拉定理**
因为$a^b\equiv a^{b \bmod \varphi(p)} \pmod{p}$，我们把时间复杂度成功地优化到了$O\left(t\log {p}\right)$。
* 但是我们依然$\mathrm{TLE}$。
考虑继续优化求$3^p \bmod p$的部分。
**"光速"幂**
已知$p=10^9+7$，$\log_2{p}\leq 32$。所以我们预处理出$2^{16}$以内的所有快速幂情况，每一次快速幂就可以$O(1)$完成了！
分析后我们的算法时间复杂度为$O\left(2 \times 2^{16}+ 2t\right)$，空间复杂度为$O(2^{16})$，期望得分：$80-100~pts$。
仍然不够优秀，考虑常数优化，如：
* 合理使用`long long`来减少取模的常数。
* 使用循环展开来降低常数。
* 使用`inline`，`register`等常见的优化常数的技巧。
以上是std使用的优化方法，下面是一些神奇的优化方法（出题人没试过）：
* 修改`Mker`中常数大的地方来减少常数。
* 更多的编译优化~~&指令集？~~
* 在代码里注入`czakioi`等玄学的信仰优化
期望得分：$100~pts$
## Code
```cpp
#include<bits/stdc++.h>
using namespace std;
namespace Mker
{
// Powered By Kawashiro_Nitori
// Made In Gensokyo, Nihon
#include<climits>
#define ull unsigned long long
#define uint unsigned int
ull sd;int op;
inline void init() {scanf("%llu %d", &sd, &op);}
inline ull ull_rand()
{
sd ^= sd << 43;
sd ^= sd >> 29;
sd ^= sd << 34;
return sd;
}
inline ull rand()
{
if (op == 0) return ull_rand() % USHRT_MAX + 1;
if (op == 1) return ull_rand() % UINT_MAX + 1;
if (op == 2) return ull_rand();
}
}
#define re register int
#define ll long long
#define ull unsigned long long
#define in inline
const int base=16,lim=(1<<16)-1,inv32=281250002;
const ll p=1e9+7,pp=p*p;
int t,pre,ans,mul[2][lim+1];
ll m;
in ll multi(re a) {ll b=mul[0][a&lim];a>>=base;if(a)b=b*mul[1][a&lim]%p;return b;}
in void get_multi()
{
mul[0][0]=mul[1][0]=1;mul[0][1]=3;
for(re j=2;j<=lim;j++) mul[0][j]=1ll*mul[0][1]*mul[0][j-1]%p;
mul[1][1]=1ll*mul[0][lim]*mul[0][1]%p;
for(re j=2;j<=lim;j++) mul[1][j]=1ll*mul[1][1]*mul[1][j-1]%p;
}
int main()
{
scanf("%d",&t);Mker::init();get_multi();
while(t--)
{
ull n=Mker::rand();
pre=(n&1)?51:21;m=((n%p)<<2)-13+p;
ans^=(((multi((n+2)%(p-1))*m)%p+pre)*inv32+pp)%p;
}
printf("%d",ans);
return 0;
}
```