一道矩阵乘法的神题 早上的时候我开挂了 想了2h想出来了。
关于这道题我推了很多矩阵 最终推出两个核心矩阵 发现这两个矩阵放在一起做快速幂就行了。
当k==1时 显然的矩阵乘法 多开一个位置维护前缀和即可。当然也可以 公式法:(f_1+f_2+...+f_n=f_{n+2}-1)
证明其 只需要数学归纳法即可。
当k==2时 不难发现 要求出((f_1+f_2+...f_n)+(f_2+...f_{n+1})+...(f_n+...f_{2n-1}))
把这个东西 画成图 可以发现是一个平行四边形。
考虑对这个东西求和 我开始想了一个比较麻烦的方法容斥。
可以发现 如果我把整个矩阵都给求出来 只需要减掉不合法方案即可(f_1)被多加了n-1次 (f_2)被多加了n-2次。
对这种形式求和怎么做？于是 我又构造了一个矩阵 多开两个位置 一个表示普通前缀和 一个表示当前前缀和。
发现这样就能把答案算出来了。
以上和正解毫无关系 只不过是我考试的时候的想法。
以下是正解：
还是考虑k==3 发现刚才的面形成了体 不过类似于平行四边形体？但是在题目角度来说其实是正方体。
对这个东西求和 发现很难 刚才的容斥不能做了。
不过此时可以发现每一面 由上一面加上上一面得到。
转回头看k==2 发现每一条线可以有上面和上上面两条线得到。
此时 我们就可以发现一个非常特殊的条件 每个维度也是由上一个维度给推出来的。
此时我们对于每个维度单独做 然后利用上个维度做当前维度。
这样我们得到了一个klogn的做法。注意此时的矩阵为(4cdot 4)的。
考虑优化。
可以发现求完一个维度之后 我们要求下一个维度。
这个过程还是一个重复的过程。
我们可以不手动调整 而是考虑构造一个矩阵帮我们进行调整。
这个矩阵也很容易构造。
剩下的就是先求出一个维度 然后利用维度生成维度是重复的。
矩阵快速幂来做即可。
总之 构造出来两个矩阵 可以发现 矩阵快速幂这两个矩阵的乘积即可。
非常巧妙。
const int MAXN=110;
int n,k,T,m;
ll f[5],w[5];
struct wy
{
ll a[5][5];
wy(){memset(a,0,sizeof(a));}
wy friend operator *(wy a,wy b)
{
wy c;
rep(1,m,i)rep(1,m,j)rep(1,m,k)
c.a[i][j]=(c.a[i][j]+a.a[i][k]*b.a[k][j])%mod;
return c;
}
wy friend operator ^(wy a,int p)
{
wy c;
rep(1,m,i)rep(1,m,j)c.a[i][j]=a.a[i][j];
while(p)
{
if(p&1)c=c*a;
a=a*a;p=p>>1;
}
return c;
}
}A,B,C;
int main()
{
freopen("fib.in","r",stdin);
freopen("fib.out","w",stdout);
get(T);
while(T--)
{
get(n);get(k);m=4;
if(n==1){puts("1");continue;}
A.a[1][1]=0;A.a[1][2]=1;A.a[1][3]=0;A.a[1][4]=1;
A.a[2][1]=1;A.a[2][2]=1;A.a[2][3]=1;A.a[2][4]=1;
A.a[3][1]=0;A.a[3][2]=0;A.a[3][3]=1;A.a[3][4]=0;
A.a[4][1]=0;A.a[4][2]=0;A.a[4][3]=0;A.a[4][4]=1;
B.a[3][1]=1;B.a[3][3]=1;B.a[4][2]=1;B.a[4][4]=1;
A=A^(n-2);C=A;
--k;
if(k)
{
A=A*B;
A=A^(k-1);
A=A*C;
}
f[1]=0;f[2]=1;f[3]=0;f[4]=1;
memset(w,0,sizeof(w));
rep(1,m,i)rep(1,m,j)w[i]=(w[i]+f[j]*A.a[j][i])%mod;
putl(w[m]);
}
return 0;
}