JZOJ.5335【NOIP2017模拟8.24】早苗

 

Description

 

Input

Output

 

Sample Input

3 3

Sample Output

21
 

Data Constraint

 

Hint

考虑DP，我们设g[i]表示第i-m+1天~第i天为第一次出现非法方案的方案数,f[i]为第i天之前合法方案的方案数，h[i]为第i天前非法方案的方案数。

根据定义容易得出：h[i]=h[i-1]*m+g[i]（即第i-1天前的方案数乘以第i天可召唤的种类加上g[i]）

f[i]=mi-h[i]

g[i]=m!*f[i-m]-$\sum _{k=1}^{m-1}k!\ast g\left[ i-k\right]$=m!*mi-m-m!*h[i-m]-$\sum _{k=1}^{m-1}k!\ast g\left[ i-k\right]$

对于g[i]，在m!*f[i-m]个方案中会存在一些不符合g[i]定义的方案，即不是第i-m+1天~第i天才开始出现第一次的非法方案。

如：对于序列.....5 4 1 2 3 4 5......（红色部分为g[i]的范围）

我们发现这个实际上在g[i-2]就出现非法方案了，不是g[i]的非法方案，我们需要减去这部分的方案。

我们发现，在上面这个序列，第i-2-5天到第i-2天是一个排列，第i-5天到第i天也是一个排列。

所以在第i-2-5天到第i-2天的某一个排列里，在第i-2天到第i天必有2!个方案多算的。

所以进一步整理我们可以得出g[i]=m!*f-[i-m]-1!*g[i-1]-2!*g[i-2]-…-(m-1)!*g[i-m+1],也就是上面那个例子。

所以我们得到了递推式：

g[i]=m!*mi-m-m!*h[i-m]-$\sum _limits{k=1}^{m-1}k!\ast g\left[ i-k\right]$

h[i]=h[i-1]*m+g[i]

不过这会超时，但是因为这是递推式，我们可以采用矩阵优化。

这是一个关于i的递推式，我们考虑列出矩阵出来。

mi-m
g[i-m+1]
g[i-m+2]
……
g[i-1]
hi-m

 

 

 

 

 

考虑它要转移到i+1的矩阵

mi-m+1
g[i-m+2]
g[i-m+3]
……
g[i]
hi-m+1

 

 

 

 

 

对照以上递推式我们得出转移的系数：

m           * mi-m
    1       g[i-m+1]
      ……     g[i-m+2]
        1   ……
m! -(m-1)! -(m-2)! …… -1! -m! g[i-1]
  1       m hi-m

 

 

 

 

 

 

（例：由mi-m到mi-m+1是乘以一个m

因为hi-m+1=hi-m*m+g[i-m+1],所以g[i-m+1]的系数为1，h[i-m]的系数为m，其余同理）

所以用个矩阵快速幂后就可以可以得出答案为mn-hn,即幂运算后的系数矩阵第一列的第一个数减去最后一次数即可。

 1 #include <cstring>
 2 #include <cstdio>
 3 #include <iostream>
 4 #include <algorithm>
 5 #define mo 1000000007
 6 #define N 109
 7 using namespace std;
 8 long long ans,n,m;
 9 struct data{
10     long long ju[N][N];
11     data operator * (const data &a) const{
12         data b;
13         for (long long i=1;i<=m+1;++i)
14          for (long long j=1;j<=m+1;++j){
15              b.ju[i][j]=0;
16              for (long long k=1;k<=m+1;++k)
17                  b.ju[i][j]=(b.ju[i][j]+ju[i][k]*a.ju[k][j]%mo+mo)%mo;
18          }
19         return b;
20     }
21 }qwq,tmp;
22 long long jie(long long x){
23     long long s=1;
24     for (long long i=2;i<=x;i++)
25         s=s*i%mo;
26     return s;
27 }
28 void kuai(){
29     long long b=n;
30     for (long long i=1;i<=m+1;++i) tmp.ju[i][i]=1;
31     while (b){
32         if (b&1) tmp=tmp*qwq;
33         qwq=qwq*qwq;
34         b>>=1;
35     }
36 }
37 void init(){
38     scanf("%lld%lld",&n,&m);
39     qwq.ju[1][1]=m;
40     for (long long i=2;i<=m-1;++i)
41         qwq.ju[i][i+1]=1;
42     for (long long i=2;i<=m;++i)
43         qwq.ju[m][i]=-jie(m-i+1);
44     qwq.ju[m][1]=jie(m);
45     qwq.ju[m][m+1]=-qwq.ju[m][1];
46     qwq.ju[m+1][2]=1;
47     qwq.ju[m+1][m+1]=m;
48 }
49 void solve(){
50     ans=tmp.ju[1][1]-tmp.ju[m+1][1];
51     printf("%lld\n",(ans+mo)%mo);
52 }
53 int main(){
54     init();
55     kuai();
56     solve();
57     return 0;
58 }
神奇的代码

（人生第一次写矩阵快速幂QAQ）

posted @ 2017-08-24 22:45  ~Lanly~  阅读(3591)  评论(0)  编辑  收藏  举报