[BZOJ2004][HNOI2010]Bus 公交线路
题目描述
小Z所在的城市有\(N\)个公交车站，排列在一条长(\(N-1\))km的直线上，从左到右依次编号为\(1\)到\(N\)，相邻公交车站间的距离均为\(1\)km。 作为公交车线路的规划者，小Z调查了市民的需求，决定按下述规则设计线路：
1.设共\(K\)辆公交车，则\(1\)到\(K\)号站作为始发站，\(N-K+1\)到N号台作为终点站。
2.每个车站必须被一辆且仅一辆公交车经过（始发站和终点站也算被经过）。
3.公交车只能从编号较小的站台驶往编号较大的站台。
4.一辆公交车经过的相邻两个站台间距离不得超过\(P\)km。
在最终设计线路之前，小Z想知道有多少种满足要求的方案。由于答案可能很大，你只需求出答案对\(30031\)取模的结果。
Input
仅一行包含三个正整数\(N, K ,P\)，分别表示公交车站数，公交车数，相邻站台的距离限制。
\(N<=10^9，1<P<=10，K<N，1<K<=P\)
Output
仅包含一个整数，表示满足要求的方案数对\(30031\)取模的结果。
Sample Input
10 3 3
5 2 3
10 2 4
Sample Output
1
3
81
看到\(k,p\)如此之小自然就是状压\(dp\)了。
我们可以十分想到\(O((n-k)*p*2^p)\)的\(dp\)，即\(dp[i][j]\)表示第\(i\)千米，第\(i-p+1\)到\(i\)千米是否有车的状态为\(j\)的方案数。
转移是很显然的枚举一个\(1\)把它放到后面去，由于直接转移的话会算重复，所以我们要限制一位。
我们只要限制最低位一定要有一个\(1\)在，这要的话就可以保证状态的转移不会重复。
\(dp\)部分见下：
dp[k][(1 << k) - 1] = 1;
ret(i, k, n) drep(K, (1 << P) - 1, 1) if (dp[i][K]) {
if (!(K & 1))continue;
ret(p, 0, P) if (K & 1 << p) {
int Now = (K ^ (1 << p)) % (1 << P - 1);
Now <<= 1, Now |= 1;
dp[i + 1][Now] += dp[i][K], Mod(dp[i + 1][Now]);
}
}
但是！！！这是过不了的。
我们发现每次dp的转移都是一样的
每次都是一个状态转移到另一个固定的状态，和\(i\)没有关系。
于是我们想到了矩阵加速，利用矩阵加速这个递推。
还没完！！！
如果把所有状态都算上去的话，矩阵的大小为\(1024*1024\)，每次矩阵乘法的复杂度为\(O(1024*1024*1024)\)无法承受。
但是我们有一个限制条件，即状态\(i\)必须保证\(i\)的最低位为\(1\)，这样的话我们优化到了\(O(512*512*512)\)的复杂度。
好像还是有一点难过掉，怎么办呢？
我们一开的的状态\(1\)的个数为\(k\)，而转移时并不会改变\(i\)的\(1\)的个数，所以我们只要把所有二进制位中\(1\)的个数为\(k\)的状态利用矩阵转移即可。
那么矩阵的大小为多少呢，为\(C(p-1,k)\)，最大为\(C(9,4)=126\)，这样的话一次矩阵乘法的复杂度为\(O(126*126*126)\)，我们就可以过掉这道题了。
#include <bits/stdc++.h>
using namespace std;
#define LL long long
#define reg register
#define clr(a,b) memset(a,b,sizeof a)
#define Mod(x) (x>=mod)&&(x-=mod)
#define abs(a) ((a)<0?-(a):(a))
#define debug(x) cerr<<#x<<"="<<x<<endl;
#define debug2(x,y) cerr<<#x<<"="<<x<<" "<<#y<<"="<<y<<endl;
#define debug3(x,y,z) cerr<<#x<<"="<<x<<" "<<#y<<"="<<y<<" "<<#z<<"="<<z<<endl;
#define rep(a,b,c) for(reg int a=(b),a##_end_=(c); a<=a##_end_; ++a)
#define ret(a,b,c) for(reg int a=(b),a##_end_=(c); a<a##_end_; ++a)
#define drep(a,b,c) for(reg int a=(b),a##_end_=(c); a>=a##_end_; --a)
#define erep(i,G,x) for(int i=(G).Head[x]; i; i=(G).Nxt[i])
#pragma GCC optimize(2)
#pragma GCC optimize(3)
#pragma GCC optimize(3,"Ofast","inline")
inline int Read(void) {
int res = 0, f = 1;
char c;
while (c = getchar(), c < 48 || c > 57)if (c == '-')f = 0;
do res = (res << 3) + (res << 1) + (c ^ 48);
while (c = getchar(), c >= 48 && c <= 57);
return f ? res : -res;
}
template<class T>inline bool Min(T &a, T const&b) {return a > b ? a = b, 1 : 0;}
template<class T>inline bool Max(T &a, T const&b) {return a < b ? a = b, 1 : 0;}
const int N = 3e2 + 5, M = 1e5 + 5, K = 10, mod = 30031;
bool MOP1;
int n, k, P, cnt, Sz[1 << K], B[130], C[1 << K];
struct Matrix {
int Num[N][N];
inline void clear(void) {clr(Num, 0);}
inline void Init(void) {rep(i, 1, cnt)Num[i][i] = 1;}
inline Matrix operator*(Matrix _)const {
Matrix Ans;
Ans.clear();
rep(i, 1, cnt)rep(j, 1, cnt)rep(k, 1, cnt)Ans.Num[i][j] = (Ans.Num[i][j] + 1ll * Num[i][k] * _.Num[k][j]) % mod;
return Ans;
}
} us;
inline Matrix qpow(Matrix A, int k) {
Matrix res;
res.clear(), res.Init();
while (k) {
if (k & 1)res = res * A;
A = A * A, k >>= 1;
}
return res;
}
bool MOP2;
void _main(void) {
n = Read(), k = Read(), P = Read();
ret(i, 1, 1 << P) {
Sz[i] = Sz[i ^ (i & -i)] + 1;
if (!(i & 1))continue;
if (Sz[i] == k)B[++cnt] = i, C[i] = cnt;
}
Matrix res;
res.clear();
rep(i, 1, cnt) {
int K = B[i];
ret(p, 0, P) if (K & 1 << p) {
if (k - P > p)continue;
int Now = ((K ^ (1 << p)) & ((1 << P - 1) - 1)) << 1 | 1;
res.Num[i][C[Now]] = 1;
}
}
Matrix Ans = qpow(res, n - k);
printf("%lld\n", Ans.Num[C[(1 << k) - 1]][C[(1 << k) - 1]]);
}
signed main() {
_main();
return 0;
}
