P1357 花园

P1357 花园 - 洛谷 | 计算机科学教育新生态 (luogu.com.cn)

为了方便理解，数学公式不是完全的数学。

m 很小，相邻 m 个花圃的情况也就 \(2^m = 32\) 种，可以状态压缩（0为p，1为c），对于一个状态 S = 1010，可以看做是由 X101 加上0（X1010）后，去掉 X 推过来的，可以认为整个花圃就是这样一个接一个推出来的，那么就要用 DP 统计，链。设 \(f[i][s]\) 为有 i 个花圃且最后 m 个状态为 s，那么很容易得到递推式子：

\[f[i][s] = f[i - 1][s >> 1] + f[i - 1][s >> 1 | (1 << m - 1)] \]

其中 \(s >> 1 | (1 << m - 1)\) 与 \(s >> 1\) 分别是 X 为 1 和为 0 状态，自己容易求出。

而花圃是环，要破环为链，正常复制原序列不好统计，于是转换思路（这里是变环为链），一个有n个点的环与一个有n+m个点且前m个点与后m个点完全相同的链的种类数是一样的。这 m 个点就可以用我们的状态 S。而这些状态 S 不能同时作为初始状态即：\(f[m][s] = 1\) 否则最后得到的 \(f[n + m][s]\) 无法分出前后 m 个是否相同。所以要单独统计，而前 m 个数包含所有合法 S，所以每个都要统计，算入总和即为答案。

分别枚举每个初始合法状态 S，以 \(f[m][s] = 1\) 为递推初值，其余（包括其他 \(f[m][s]\)） \(f\) 为 \(0\)，那么最后得出 \(f[n + m][s]\) 就是前面 m 个点和后面 m 个点相同的序列的种类数的一部分。

因为 n 很大，为了加速要使用矩阵快速幂，加快递推（也因为递推方程简单）。建议用 dfs 写转移矩阵。

#include <iostream>
#include <cstring>
#include <cstdio>
#include <algorithm>

#define x first
#define y second

using namespace std;

typedef long long LL; // 记得开longlong
typedef pair<int, int> PII;

const int M = (1 << 5) + 2, mod = 1e9 + 7;

LL n, m, k;
PII id[M];
int sum, cnt;
int st[M], g[M];

struct Mat
{
    int l, r;
    int a[M][M];
}A, B; // 初始矩阵，转移矩阵

void print(Mat &A) // 用来输出矩阵，方便调试
{
    for (int i = 1; i <= A.l; i ++ )
    {
        for (int j = 1 ; j <= A.r; j ++ )
            printf("%d ", A.a[i][j]);
        puts("");
    }
    puts("");
}


int get(int x)
{
    int cnt = 0;
    while (x)
    {
        if (x & 1) cnt ++ ;
        x >>= 1;
    }
    return cnt;
}

void dfs(int x)
{
    if (x > cnt) return;
    
    int i = g[x];
    int a = id[i].x, b = id[i].y;
    
    if (st[a]) B.a[st[a]][x] = 1;
    if (st[b]) B.a[st[b]][x] = 1;
    dfs(x + 1);
}

Mat Mat_mul(Mat &A, Mat &B, int p) // 矩阵乘
{
    Mat C;
    C = {A.r, A.r};

    for (int i = 1; i <= A.l; i ++ )
        for (int j = 1; j <= B.r; j ++ )
            for (int k = 1; k <= A.r; k ++ )
                C.a[i][j] = (C.a[i][j] + 1ll * A.a[i][k] * B.a[k][j] % mod) % mod;
    // print(C);   
    return C;
}

Mat Mat_qmi(Mat A, LL k, int p) // 矩阵快速幂
{
    Mat res = {A.l, A.r};
    for (int i = 1; i <= A.l; i ++ ) res.a[i][i] = 1;
    // print(res);
    while (k)
    {
        if (k & 1) res = Mat_mul(res, A, p);
        k >>= 1;
        A = Mat_mul(A, A, p);
        // print(res);
    }
    
    return res;
}


int main()
{
    cin >> n >> m >> k;
    
    
    for (int i = 0; i < (1 << m); i ++ )
    {
        id[i] = {i >> 1, i >> 1 | (1 << m - 1)};
        if (get(i) <= k) 
        {
            st[i] = ++ cnt;
            g[cnt] = i;
            // printf("%d %d %d %d\n", id[i].x, id[i].y, i, cnt);
        }
    }
    A.l = 1;
    A.r = cnt;
    
    B.l = B.r = cnt;
    dfs(1); 
    B = Mat_qmi(B, n, mod);
    int ans = 0;
    for (int i = 0; i < (1 << m); i ++ )
    {
        if (st[i] == 0) continue;
        Mat C = A;
        C.a[1][st[i]] = 1;
        
        Mat x = Mat_mul(C, B, mod);
        ans = (ans + x.a[1][st[i]]) % mod;
    }
    
    
    cout << ans << endl;
    
    return 0;
}