【jzoj5335】早苗
题目描述
Sanae准备对大结界进行一次连续n天的风祭。Sanae每天可以召唤一种神风。她一共有m种神风可以召唤，但是如果有连续m天刮了m种不同的神风，环境就会遭到破坏。现在，Sanae想知道这n天她有多少种召唤神风的方案。
由于答案可能很大，所以你只需要告诉她方案数对1000000007取模的结果即可。
输入
两个正整数n，m，分别代表风祭的天数和Sanae能召唤的风的种数。
输出
一个整数，为答案对1000000007取模的结果。
样例输入
3 3
样例输出
21
数据范围
对于8%的数据，m=2
对于另16%的数据，n<=10,m<=4
对于48%的数据，n<=100000,m<=10
对于80%的数据，n<=100000,2<=m<=100
对于100%的数据，2<=m<=100,m<=n<=10^16
思路
题目大意：求用1~m组成有n个数的序列，任意连续m个数不能互不相同的方案数
？？？任意连续m个数不能互不相同？如果我第i个放了数字j，那我前面必须放哪些数？i和i-1有什么关系？任意连续m个数不能互不相同的情况太多了吧……
好像不管怎么想都得用减法啊……
多试试逆向
设f[i]是用了i个数的合法方案数，h[i]是用了i个数的不合法方案数。
f[i]+h[i]=mn
求不出f[i]，那我们来求h[i]吧~
emmmmmmm，考虑h[i]和h[i-1]的关系，可以看出h[i]必定包含h[i-1]*m（意思是最后一位随便放，比f好多了），也就是h[i]=h[i-1]*m+g[i],其中g[i]为加上第i数后新出现的不合法方案数。显然，g[i]的最后m个数是个全排列，也就是g[i]包含了m!*f[i-m]。
然而，我们这样算会算重。
“g[i]为加上第i数后新出现的不合法方案数。”，新出现，也就是只有i-m+1~i是全排列，那么i-m~i-1，i-m-1~i-2……都不能是全排列。m!*f[i-m]不能保证这点。
例如：m=5,你想算g[9],5 3 2 3 1 5 4 2 3中,f[4]确实包含5 3 2 3这个方案，但g[9]中不包含5 3 2 3 1 5 4 2 3，因为5 3 2 3 1 5 4早就已经不合法了。
但我们可以去重啊。
对于m!*f[i-m]中的每一个算重的序列，我们找到这个序列第一次出现全排列的位置假设为i-j-m+1~i-j(1<=j<=m-1),那么i-j-m+1~i-j是全排列而i-m+1~i也是全排列，那么i-j+1~i也是有j个元素的全排列（数字不一定正好是1~j，但肯定是j个元素）。那么我们就要减掉j!*g[i-j]个方案数.
总递推式:g[i]=m!*f[i-m]-$\sum _limits{j=1}^{m-1}j!\ast g\left[ i-j\right]$
=m!*(mi-m-h[i-m])-$\sum _limits{j=1}^{m-1}j!\ast g\left[ i-j\right]$
其中h[i-m]=h[i-m-1]*m+g[i-m].
目前的复杂度为O(n*m)，80分。
优化
递推式的主体可以看作是n，而n的范围n<=10^16.
妥妥的矩阵乘法优化。
根据递推式里，与i有关的有mi-m,h[i-m],g[i-j]，所以我们把矩阵转移方向定为:
总转移：
转移矩阵就让大家自己思考啦。
代码：
1 #include <iostream> 2 #include <cstdio> 3 #include <cstdlib> 4 #include <cmath> 5 #include <algorithm> 6 #include <cstring> 7 #include <string> 8 #include <queue> 9 #define fo(p,q,r) for(p=q;p<=r;++p) 10 #define fow(p,q,r) for(p=q;p>=r;--p) 11 #define arclr(p,q) memset(p,q,sizeof(p)) 12 using namespace std; 13 typedef long long LL; 14 15 const LL mo=1000000007; 16 17 LL mr[105][105][2],temp[105][105],st[105],en[105],ans; 18 19 LL jc[105],n,m,i,j,cnt; 20 21 void mul(LL p,LL q,LL r) 22 { 23 LL i1,j1,k1; 24 arclr(temp,0); 25 fo(i1,1,m+1) 26 fo(j1,1,m+1) 27 fo(k1,1,m+1) 28 temp[i1][j1]=(temp[i1][j1]+((mr[i1][k1][p]*mr[k1][j1][q])%mo))%mo; 29 fo(i1,1,m+1) 30 fo(j1,1,m+1) 31 mr[i1][j1][r]=temp[i1][j1]; 32 } 33 34 void build() 35 { 36 mr[1][1][0]=m; 37 fo(i,2,m-1) 38 mr[i][i+1][0]=1; 39 mr[m][m+1][0]=(2*mo-jc[m])%mo; 40 mr[m][1][0]=jc[m]; 41 fo(i,2,m) mr[m][i][0]=(2*mo-jc[m-i+1])%mo; 42 mr[m+1][2][0]=1; mr[m+1][m+1][0]=m; 43 44 fo(i,1,m+1) 45 mr[i][i][1]=1; 46 } 47 48 void mi(LL p) 49 { 50 LL i1=p; 51 while (i1>0) 52 { 53 if (i1&1) mul(1,0,1); 54 mul(0,0,0); 55 i1>>=1; 56 } 57 } 58 59 int main() 60 { 61 scanf("%lld%lld",&n,&m); 62 63 jc[0]=1; 64 fo(i,1,m) jc[i]=(jc[i-1]*i)%mo; 65 66 st[1]=1; 67 fo(i,2,m+1) st[i]=0; 68 69 build(); 70 71 fo(i,1,m+1) 72 { 73 fo(j,1,m+1) printf("%lld ",mr[i][j][0]); 74 printf("\n"); 75 } 76 77 mi(n); 78 79 arclr(en,0); 80 fo(i,1,m+1) 81 fo(j,1,m+1) 82 en[i]=(en[i]+((mr[i][j][1]*st[j])%mo))%mo; 83 84 ans=(en[1]-en[m+1]+mo)%mo; 85 printf("%lld\n",ans); 86 }
