【NOI2019模拟2019.7.4】朝夕相处 (动态规划+BM)

Description:


【NOI2019模拟2019.7.4】朝夕相处 (动态规划+BM)

题解:


这种东西肯定是burnside引理:
\(\sum置换后不动点数 \over |置换数|\)

一般来说,是枚举置换\(i\),则\(对所有x,满足a[x+i]=a[i]\),然后a还要满足题目条件,但是仔细想一想,设\(d=gcd(i,n)\),只要a[0..d-1]满足就好了,所以:
\(Ans=\sum_{d|n}f(d)*\phi(n/d),f(d)表示\)不考虑循环同构时的答案。

然后考虑dp:

枚举0这一列的块是什么,然后设\(dp[i][1..4][1..4][1..3][1..3]\),表示确定了后\(i\)列,第一行的颜色,和长度,第二行的颜色和长度,1、2表示第0列的颜色, 3、4表示不同于1、2列的颜色。

然后暴写dp,我死了……

n<=3的时候要特判,那么n更大的时候要矩阵乘法。

事实上有矩阵乘法的做法可以改成BM,最短递推式大概不超过35项

Code:

#include<bits/stdc++.h>
#define fo(i, x, y) for(int i = x, B = y; i <= B; i ++)
#define ff(i, x, y) for(int i = x, B = y; i <  B; i ++)
#define fd(i, x, y) for(int i = x, B = y; i >= B; i --)
#define ll long long
#define pp printf
#define hh pp("\n")
using namespace std;

const int mo = 1e9 + 7;

ll ksm(ll x, ll y) {
    ll s = 1;
    for(; y; y /= 2, x = x * x % mo)
        if(y & 1) s = s * x % mo;
    return s;
}

namespace cz {
    ll a[7] = {0, 0, 0, 24, 312, 1320, 3720};
    ll cz(ll x) {
        if(x <= 6) return a[x];
        ll s = 0;
        fo(i, 0, 6) {
            ll xs = 1;
            fo(j, 0, 6) if(i != j) xs = xs * (x - j) % mo * ksm(i - j, mo - 2) % mo;
            s += xs * a[i];
        }
        return s % mo;
    }
}


const int N = 200;

ll num, n, k, a[N];
ll f[N][5][5][4][4];

ll dp(int n, int cnt, int a0, int a1, int b0, int b1) {
    ll m = k - cnt;
    memset(f, 0, sizeof f);
    f[n][1][cnt][a0][b0] = 1;
    fd(i, n, 1) {
        int c0 = 0, c1 = 0, c3 = 0, c4 = 0;
        if(i - 1 < a0 || i - 1 >= n - a1) c0 = 1;
        if(i - 1 < b0 || i - 1 >= n - b1) c1 = cnt;
        if(i - 1 == a0 || i - 1 == n - a1 - 1) c3 = 1;
        if(i - 1 == b0 || i - 1 == n - b1 - 1) c4 = cnt;
        fo(x, 1, 4) fo(y, 1, 4) fo(lx, 1, 3) fo(ly, 1, 3) if(f[i][x][y][lx][ly]) {
            ll F = f[i][x][y][lx][ly] % mo;
            fo(u, 1, 3) {
                ll xs = 1;
                if(c0 && u != c0) xs = 0;
                if(c1 && u != c1) xs = 0;
                if(c3 && u == c3) xs = 0;
                if(c4 && u == c4) xs = 0;
                if(u == 2 && cnt == 1) xs = 0;
                if(u <= 2 && (x == u || y == u)) xs = 0;
                if(!xs) continue;
                if(u == 3) {
                    xs *= m;
                    if(x >= 3) xs --;
                    if(y >= 3) xs --;
                    if(x == 3 && y == 3) xs ++;
                }
                xs = max(xs, 0ll);
                f[i - 1][u][u][1][1] += xs * F % mo;
            }
            fo(u, 1, 3) fo(v, 1, 4) if(u != v && v != 3) {
                ll xs = 1;
                if(c0 && u != c0) xs = 0;
                if(c1 && v != c1) xs = 0;
                if(c3 && u == c3) xs = 0;
                if(c4 && v == c4) xs = 0;
                if((u == 2 || v == 2) && cnt == 1) xs = 0;
                if(!xs) continue;
                fo(r1, 0, 1) fo(r2, 0, 1) {
                    if(r1 && (u != x || lx == 3)) continue;
                    if(r2 && (v != y || ly == 3)) continue;
                    if(!r1 && u <= 2 && u == x) continue;
                    if(!r2 && v <= 2 && v == y) continue;
                    if(!r1 && !r2 && x != y) continue;
                    if(x == y && (r1 || r2)) continue;
                    int nlx = r1 ? lx + 1 : 1, nly = r2 ? ly + 1 : 1;
                    xs = 1;
                    if(u <= 2 || v <= 2 || r1 || r2) {
                        if(r1 && r2) {
                            xs = x != y;
                        } else
                        if(!r1 && !r2) {
                            xs *= (u <= 2 || r1) ? 1 : (m - (x >= 3));
                            xs *= (v <= 2 || r2) ? 1 : (m - (y >= 3));
                        } else {
                            if(r2) {
                                xs *= (u <= 2) ? 1 : (m - (x >= 3) - (y >= 3 && v >= 3));
                            } else {
                                xs *= (v <= 2) ? 1 : (m - (x >= 3 && u >= 3) - (y >= 3));
                            }
                        }
                    } else {
                        if(x <= 2 && y <= 2) xs *= m * (m - 1); else
                        if(x > 2 && y > 2) {
                            if(x == y) xs *= (m > 2) ? (m - 1) * (m - 2) : 0; else
                            xs *= (m > 3) ? (m - 2) * (m - 3) : 0;
                        } else {
                            xs *= (m > 1) ? (m - 1) * (m - 1) : 0;
                        }
                    }
                    xs = max(xs, 0ll);
                    xs %= mo;
                    f[i - 1][u][v][nlx][nly] += xs * F % mo;
                }
            }
        }
    }
    ll ans = 0;
    fo(x, 1, 4) fo(y, 1, 4) fo(lx, 1, 3) fo(ly, 1, 3) if(f[0][x][y][lx][ly]) {
        ans = (ans + f[0][x][y][lx][ly]) % mo;
    }
    return ans;
}

void dp() {
    for(int n = 4; n <= 70; n ++) {
        a[n] = dp(n, 1, 1, 0, 1, 0) * k % mo;
        fo(i, 1, 3) fo(j, 1, i) fo(u, 1, 3) fo(v, 1, u)
            a[n] += dp(n, 2, j, i - j, v, u - v) * k % mo * (k - 1) % mo;
        a[n] %= mo;
    }
//  pp("%lld\n", a[5]);
    a[1] = 0;
    a[2] = (k * (k - 1) + k * (k - 1) % mo * (k - 2) * 2) % mo;
    a[3] = cz :: cz(k);
}

int m; ll b[N];

namespace bm {
    ll f[3][N], len[3], fail[3], delta[3];

    void cp(int x, int y) {
        fo(i, 1, len[x]) f[y][i] = f[x][i];
        len[y] = len[x]; fail[y] = fail[x]; delta[y] = delta[x];
    }

    void solve(ll *a, int n) {
        fo(i, 1, n) f[0][i] = 0; len[0] = 0;
        fo(i, 1, n) {
            ll tmp = a[i];
            fo(j, 1, len[0]) tmp -= f[0][j] * a[i - j] % mo;
            tmp = (tmp % mo + mo) % mo;
            delta[0] = tmp;
            if(!tmp) continue;
            fail[0] = i;
            if(!len[0]) {
                cp(0, 1); f[0][len[0] = 1] = 0;
                continue;
            }
            cp(0, 2);
            ll mul = delta[0] * ksm(delta[1], mo - 2) % mo;
            int st = i - fail[1]; len[0] = max(len[0], st + len[1]);
            f[0][st] = (f[0][st] + mul) % mo;
            fo(j, 1, len[1]) f[0][st + j] = (f[0][st + j] + (mo - mul) * f[1][j]) % mo;
            if(len[2] - fail[2] < len[1] - fail[1]) cp(2, 1);
        }
        m = len[0] + 1; b[m] = 1;
        fo(i, 1, len[0]) b[i] = -f[0][m - i];
    }
}

ll c[N * 2], d[N * 2], e[N * 2];

void qmo(ll *a) {
    fo(i, 0, 2 * m) a[i] %= mo;
    fd(i, 2 * m, m) if(a[i]) {
        ll v = a[i]; int c = i - m;
        fo(j, 1, m) a[c + j] = (a[c + j] - v * b[j]) % mo;
    }
}

ll solve(int n) {
    if(n <= 35) return a[n];
    fo(i, 0, 2 * m) c[i] = d[i] = 0;
    c[0] = 1; d[1] = 1;
    int y = n;
    for(; y; y /= 2) {
        if(y & 1) {
            fo(i, 0, 2 * m) e[i] = c[i], c[i] = 0;
            fo(i, 0, m) fo(j, 0, m) c[i + j] += e[i] * d[j] % mo;
            qmo(c);
        }
        fo(i, 0, 2 * m) e[i] = d[i], d[i] = 0;
        fo(i, 0, m) fo(j, 0, m) d[i + j] += e[i] * e[j] % mo;
        qmo(d);
    }
    ll ans = 0;
    fo(i, 1, m - 1) ans += c[i] * a[i] % mo;
    return ans % mo;
}

ll phi(int n) {
    ll s = n;
    for(int i = 2; i * i <= n; i ++) if(n % i == 0) {
        s = s / i * (i - 1);
        for(; n % i == 0; n /= i);
    }
    if(n > 1) s = s / n * (n - 1);
    return s;
}


int main() {
    freopen("experience.in", "r", stdin);
    freopen("experience.out", "w", stdout);
    scanf("%lld %lld %lld", &num, &n, &k);
    dp();
    bm :: solve(a, 70);
    ll ans = 0;
    for(int d = 1; d * d <= n; d ++) if(n % d == 0) {
        ans += solve(d) * phi(n / d) % mo;
        if(d != n / d) ans += solve(n / d) * phi(d) % mo;
    }
    ans = (ans % mo + mo) % mo * ksm(n, mo - 2) % mo;
    pp("%lld\n", ans);
}
上一篇:离散化


下一篇:6389. 【NOIP2019模拟2019.10.26】小w学图论