CF1117D Magic Gems 矩阵快速幂 DP

原题链接:CF1117D Magic Gems

题目大意

一开始有无穷多个牛逼石头,每个牛逼石头可以分解成\(m\)个普通石头,而普通石头不能再次分解出普通石头.你有一个大小为\(n\)的背包,牛逼石头和普通石头的大小都是1个单位,问恰好填满整个背包的方案数有多少.答案对\(10^9+7\)取模

此处讨论的方案数包含牛逼石头位置不同,但是数量相同的方案,即不认为数量相同的位置不同的属于同构的方案.

数据范围:

\(1 \leq n \leq 10^{18}\)

\(2 \leq m \leq 100\)

思路

看到\(n\)大的离谱而\(m\)很小,支持\(O(m^3logn)\)很可能就是一个矩阵快速幂了,顺着可以先把递推表达式搞出来.

这里有一个比较浅显的思路就是直接\(dp\)以\(f[i][j]\)表示当前一共有\(i\)个石头\(j\)个位置的方案数,不过这样去做有很多组合数上的问题不太好搞,有个更好的做法是:\(f[i]\)表示占用\(i\)个单位的方案数.那么转移也比较显然,就是枚举这个状态的最后一部分石头是怎么来的:要么从前一个状态加了一个不分解的牛逼石头,也就是\(f[i] += f[i - 1]\),要么把一个牛逼石头分解了再加进来,也就是\(f[i] += f[i - m]\)合起来就是\(f[i] = f[i - 1] + f[i - m]\).这一步比较显然,不是这个题的重点.

考虑优化,直接暴力\(dp\)显然是\(O(nm)\)的,需要矩阵快速幂:

• 状态矩阵\(A = [f[0],f[1],f[2],...f[m-1]]\),状态矩阵是一个\(1*m\)的行向量.
• 转移矩阵\(B\).我们希望能有这样的一个结果:\([f[0],f[1],f[2],...f[m-1]] *B = [f[1],f[2],f[3]...f[m]]\).那么可以注意到除了最后一列元素\(f[m-1]\)比较怪之外,其他的元素其实只是后面的往前面挪动了一位,对此可以先构造转移矩阵\(B\)的前\(m-1\)列,对于这些列,根据矩乘的式子可以找到就是让第\(i\)列的\(i+1\)行之外的所有元素为\(0\)而他系数是\(1\)就可以了.对于最后一列元素直接单独对他做就可以了:\(f[m]=f[0] + f[m-1]\),只需要这一列本来的元素以及最左侧的\(f[0]\)就可以了,系数就是\(B[1][m] = 1\)以及\(B[m][m]=1\).
• 考虑转移:首先\(A*B\)的结果会使得\(A[1][1]\)的位置变成\(f[1]\),所以乘上\(B^i\)之后,\(A[1][1]=f[i]\).那么结果需要的是\(f[n]\)直接让\(A*B^n\)即可.后者可以矩阵快速幂解决.

最后答案取\(A[1][1]\)输出即可.

(由于矩乘部分也是copy别人的代码来的,不知道会不会有锅,有锅别找我(逃))

代码

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
#define forn(i,x,n) for(int i = x;i <= n;++i)	

const int N = 105,MOD = 1e9+7;
ll n,m;
struct Matrix
{
	ll a[N+5][N+5];
	ll* operator[](int x){return a[x];}
	Matrix operator*(Matrix& rhs){
		Matrix res;
		for(int i=1;i<=m;++i){
			for(int j=1;j<=m;++j){
				for(int k=1;k<=m;++k)
				{
					res[i][j] = (res[i][j] + (ll)a[i][k] * rhs[k][j] % MOD) % MOD;
				}
			}
		}
		return res;
	}
	Matrix()	{memset(a,0,sizeof(a));}
};
Matrix mat_pow(Matrix x,ll i)
{
	Matrix y;
	for(int j=1;j<=m;++j)y[j][j]=1;
	while(i)
	{
		if(i&1LL)y=y*x;
		x=x*x;
		i>>=1;
	}
	return y;
}

int main()
{
	scanf("%lld%lld",&n,&m);
	Matrix A,B;
	forn(i,1,m)	A.a[1][i] = 1;
	forn(i,1,m - 1)	B.a[i + 1][i] = 1;
	B.a[1][m] = 1;B[m][m] = 1;
	
	B = mat_pow(B,n);
	A = A * B;
	printf("%lld",A[1][1]);
	
    return 0;
}