状压 DP

枚举子集

有时，对于一个二进制状态 \(i\)，我们只希望从它的子集转移，因为这样能有效的减少时间复杂度。

Y555 最优组队

题目大意

给定 \(n\) 个人（\(n\le 16\)），你需要把他们分成若干组，满足所有组的权值和最大。题目给定所有可能出现的组对应的权值。

我们注意到 \(n\) 很小，可以指数复杂度通过，因此考虑状压。设二进制状态 \(i\)，表示对应集合中的人都已组好队，最大的权值和。很容易写出状态转移方程：

\[ f_i=\max_{j\subset i}\{f_j+w_{i-j}\} \]

对于状态 \(i\)，暴力枚举所有状态 \(j\)，再判断是否满足 \(j\subset i\)，决定是否转移。时间复杂度为 \(O(4^n)\)，“不能”通过本题。

for(int i = 1; i <= (1 << n) - 1; i++) {
    for(int j = 0; j < i; j++) {
        if((i | j) == i){
            f[i] = max(f[i], f[j] + a[i ^ j]);
        }
    }
}

cout << f[(1 << n) - 1] << endl;

我也不知道这段暴力为什么就水过了，可能是带一个小常数。

我们注意到，大量的时间被浪费在了枚举那些不是 \(i\) 的子集的 \(j\) 上面。在此我们提供一种可以只遍历到所有真子集的方法：

for(int i = 1; i <= (1 << n) - 1; i++) {
    int j = i;
    while(j) {
        j = (j - 1) & i;
        f[i] = max(f[i], f[j] + a[i ^ j]);
    }
}

通过此种优化，时间复杂度降为 \(O(3^n)\)，可以在更短时间内通过本题。

时间复杂度证明

设一个状态有 \(i\) 个元素，其子集一共有 \(2^i\) 个。有 \(i\) 个元素的状态一共有 \(C_n^i\) 种方案，总转移时间复杂度为：

\[ \sum_{i=1}^n{C_n^i\times 2^i}=3^n \]

此种方法比暴力转移要更优。

简单状压

Y556 最短路径

题目大意

给定一个 \(n\) 个点 \(m\) 条边的有向图，有 \(k\) 个标记点，要求从规定的起点按任意顺序经过所有标记点到达规定的终点，问最短的距离是多少。

\(n\le 5\times 10^4,\ m\le 10^5,\ k \le 10\)。

我们可以在图上跑一遍 Dijkstra，同时通过状压记录经过了哪些标记点。时间复杂度 \(O((n\times 2^k + m\times 2 ^ k)\log (m\times 2^k))\)，不能通过本题，考虑优化。

通过仔细分析该算法，我们发现各个标记点之间的最短路被重复计算了 \(2^k\) 次，这是完全不必要的，考虑优化这些重复计算。

我们可以先求出这 \(k\) 个点之间两两的最短路，再在这 \(k\) 个点的简化图上跑状压 Dijkstra。这样，时间复杂度就被优化为 \(O\big(k(n+m)\log m+k^2\log k^2\big)\)，可以通过本题。

代码

合理利用 namespace 区分名称相似的函数可以避免很多不必要的麻烦。

#include<iostream>
#include<cstring>
#include<queue>
#define int long long

using namespace std;
const int N = 5E4 + 10;
const int M = 1E5 + 10;
const int K = 15;

struct Edge {
    int v, w;
    int next;
} pool[2 * M];
int nn, head[N];

void addEdge(int u, int v, int w) {
    pool[++nn] = {v, w, head[u]};
    head[u] = nn;
}

int n, m, k, s, t;
int tag[N];
int b[K], dis[K][K];

namespace Dij1 {

    struct Node {
        int u, d;
        inline bool operator<(const Node &other) const {
            return d > other.d;
        }
    };

    int vis[N], d[N];
    priority_queue<Node> que;

    void Dijkstra(int x) {

        while(!que.empty()) que.pop();
        memset(d, 0x3f, sizeof(d));
        memset(vis, 0, sizeof(vis));
        d[b[x]] = 0;
        que.push({b[x], 0});
        while(!que.empty()) {
            int u = que.top().u;
            que.pop();
            if(vis[u]) continue;
            vis[u] = 1;
            for(int i = head[u]; i; i = pool[i].next) {
                int v = pool[i].v, w = pool[i].w;
                if(d[u] + w < d[v]) {
                    d[v] = d[u] + w;
                    que.push({v, d[v]});
                }
            }
        }
        for(int i = 1; i <= k + 2; i++) {
            dis[x][i] = d[b[i]];
        }

    }

}

namespace Dij2 {

    struct Node {
        int u, sta;
        int d;
        inline bool operator<(const Node &other) const {
            return d > other.d;
        }
    };

    int vis[K][1 << 12], d[K][1 << 12];
    priority_queue<Node> que;

    void Dijkstra() {

        if(s == t && k == 0) {
            cout << 0 << endl;
            exit(0);
        }

        memset(d, 0x3f, sizeof(d));
        d[k + 1][tag[s] ? (1 << tag[s] - 1) : 0] = 0;
        que.push({k + 1, tag[s] ? (1 << tag[s] - 1) : 0, 0});
        while(!que.empty()) {
            int u = que.top().u;
            int sta = que.top().sta;
            que.pop();
            if(vis[u][sta]) continue;
            vis[u][sta] = 1;
            if(u == k + 2 && (sta & ((1 << k) - 1)) == (1 << k) - 1) {
                cout << d[u][sta] << endl;
                exit(0);
            }
            for(int v = 1; v <= k + 2; v++) {
                if(u == v) continue;
                int tmp = sta | (1 << v - 1);
                if(d[u][sta] + dis[u][v] < d[v][tmp]) {
                    d[v][tmp] = d[u][sta] + dis[u][v];
                    que.push({v, tmp, d[v][tmp]});
                }
            }
        }

    }

}

signed main() {

    cin >> n >> m >> k >> s >> t;
    for(int i = 1; i <= m; i++) {
        int u, v, w;
        cin >> u >> v >> w;
        addEdge(u, v, w);
    }
    for(int i = 1; i <= k; i++) {
        int x;
        cin >> x;
        tag[x] = i;
        b[i] = x;
    }

    b[k + 1] = s;
    b[k + 2] = t;
    for(int i = 1; i <= k + 2; i++) {
        Dij1::Dijkstra(i);
    }

    Dij2::Dijkstra();

    cout << -1 << endl;

    return 0;
}

连通性状压 & 最小表示法

P2109 [NOI2007] 生成树计数

题意

给定排成一行的 \(n\) 个节点，相邻节点之间的距离都为 \(1\)。每个点和所有与其距离不超过 \(k\) 的节点连边。问这个图的生成树的数量对 \(65521\) 取模的结果。

\(n\le 10^{15},k\le 5\)

\(n\le 10^{15}\)，但不易找到通项公式，因此考虑矩阵快速幂优化 dp。

考虑第 \(i\) 个节点，枚举其向 \([i-k,i-1]\) 中的哪些节点有连边。注意到 \(k\le 5\)，因此当前节点最多只能与前面的 \(5\) 个节点连边，dp 也只需要记录最近 \(5\) 个节点的状态。然而生成树要求图上不能出现环，因此我们需要保证第 \(i\) 个节点连接的 \([i-k,i-1]\) 里的点原本是不连通的。这就需要我们记录最后 \(k\) 个点的连通性。

要想记录所有可能出现的连通性状态，可以先暴力搜索出所有合法的连边方案，再把连通性相同的状态合并到一起，重新编号（这叫做最小表示法）。经过实测，\(k=5\) 时只有 \(52\) 种连通性状态。

int id[55560], st[60], cnt[60], scnt;
int fa[K];

int find(int x) {
    if(fa[x] == x) return x;
    return fa[x] = find(fa[x]);
}

// 参数 s 表示边集 E 的二进制状压形式。
void check(int s) {
    for(int i = 1; i <= k; i++) fa[i] = i;
    for(int i = 1, c = 0; i <= k; i++) {
        for(int j = 1; j < i; j++, c++) {
            // 判环
            if(s & (1 << c)) {
                if(find(i) == find(j)) return;
                fa[max(find(i), find(j))] = min(find(i), find(j));
            }
        }
    }
    int sta = 0;
    for(int i = 1; i <= k; i++) {
        sta = sta * 10 + find(i);
    }
    if(!id[sta]) {
        id[sta] = ++scnt;
        st[scnt] = sta;
    }
    cnt[id[sta]]++;
}

int main() {

    // something...

    for(int i = 0; i < (1 << (k * (k - 1) / 2)); i++) {
        check(i);
    }

    // something...

}

这样，我们就得到了所有连通性状态的最小表示法和出现次数。接下来，我们只需要枚举第 \(i\) 个节点向 \([i-k,i-1]\) 的连边，判断是否有环，再找出新的连通性状态的最小表示法，就可以知道 \([i-k-1,i-2]\) 的每一个连通性状态 \(s_1\) 对 \([i-k,i-1]\) 的连通性状态 \(s_2\) 的贡献（转移系数）。

我们将转移系数写成一个矩阵，并将 \([i-k,i-1]\) 的每种连通性状态对应的方案数写在一个 \(1\times n\) 的向量里，将向量和矩阵相乘，即可得到后一个位置 \(i+1\) 的状态。

最后使用矩阵快速幂计算出 \(op^{n-k}\)，再以 \(cnt\) 为初始状态，得到最终状态。输出 \(k\) 个节点全部连通的状态对应的方案数即可。

代码

#include<iostream>
#include<cstring>
#include<cmath>
#include<vector>
#define cint const int&
#define int long long
using namespace std;
const int K = 7;
const int MOD = 65521;
const int F[] = {-1, -1, 11, 111, 1111, 11111};

struct Matrix {
    int n;
    int a[60][60];
    inline Matrix(int _n) { n = _n; memset(a, 0, sizeof(a)); }
    inline int* operator [] (int index) { return a[index]; }
    inline const int* operator [] (int index) const { return a[index]; }
    inline Matrix operator * (const Matrix &b) const {
        Matrix res(n);
        for(int i = 1; i <= n; i++) {
            for(int j = 1; j <= n; j++) {
                for(int k = 1; k <= n; k++) {
                    res[i][j] = (res[i][j] + a[i][k] * b[k][j] % MOD) % MOD;
                }
            }
        }
        return res;
    }
};

inline vector<int> operator*(const vector<int> &a, const Matrix &b) {
    vector<int> res(b.n + 1);
    for(int i = 1; i <= b.n; i++) {
        for(int j = 1; j <= b.n; j++) {
            res[i] = (res[i] + a[j] * b[j][i]) % MOD;
        }
    }
    return res;
}

int n, k, lim;

int id[55560], st[60], cnt[60], scnt;
vector<int> beg, ed;

int p[K], fa[K];

int find(int x) {
    if(fa[x] == x) return x;
    return fa[x] = find(fa[x]);
}

void check(int s) {
    for(int i = 1; i <= k; i++) fa[i] = i;
    for(int i = 1, c = 0; i <= k; i++) {
        for(int j = 1; j < i; j++, c++) {
            if(s & (1 << c)) {
                if(find(i) == find(j)) return;
                fa[max(find(i), find(j))] = min(find(i), find(j));
            }
        }
    }
    int sta = 0;
    for(int i = 1; i <= k; i++) {
        sta = sta * 10 + find(i);
    }
    if(!id[sta]) {
        id[sta] = ++scnt;
        st[scnt] = sta;
    }
    cnt[id[sta]]++;
}

inline void add(Matrix &op, cint s1, cint s2) {
    op[id[s1]][id[s2]]++;
}

Matrix qpow(Matrix a, int k) {
    Matrix res(a.n);
    for(int i = 1; i <= a.n; i++) res[i][i] = 1;
    while(k) {
        if(k & 1) res = res * a;
        a = a * a;
        k >>= 1;
    }
    return res;
}

signed main() {

    cin >> k >> n;

    for(int i = 0; i < (1 << (k * (k - 1) / 2)); i++) {
        check(i);
    }

    beg.resize(scnt + 1);
    for(int i = 1; i <= scnt; i++) beg[i] = cnt[i];

    Matrix op(scnt);
    for(int s = 1; s <= scnt; s++) {
        int sta = st[s];
        for(int i = k; i >= 1; i--) p[i] = sta % 10, sta /= 10;
        for(int i = 0; i < (1 << k); i++) {
            for(int j = 1; j <= k; j++) fa[j] = p[j];
            fa[k + 1] = k + 1;
            bool flag = true;
            for(int j = 1; j <= k; j++) {
                if(i & (1 << j - 1)) {
                    if(find(k + 1) == find(j)) {
                        flag = false;
                        break;
                    }
                    fa[max(find(k + 1), find(j))] = min(find(k + 1), find(j));
                }
            }
            if(!flag) continue;
            int sta2 = 0, tag = 0;
            for(int i = 2; i <= k + 1; i++) if(find(i) == 1) { tag = i; break; }
            if(!tag) continue;
            for(int i = 2; i <= k + 1; i++) {
                sta2 = sta2 * 10 + (find(i) == 1 ? tag : find(i)) - 1;
            }
            add(op, st[s], sta2);
        }
    }

    op = qpow(op, n - k);
    ed = beg * op;

    cout << ed[id[F[k]]] << endl;

    return 0;
}

C250420 D1 最小生成树

题意

给定一张 \(n\) 个点 \(m\) 条边的 无向连通图，图中的第 \(i\) 条边有 \(a_i,b_i\)​ 两个权重。

对于图中的每条边，你可以选择 \(a_i,b_i\)​ 其中之一作为该条边的边权。

对于所有满足 \(0\le k\le m\) 的整数 \(k\)，你需要求出，若选择 恰好 \(k\) 个 \(a_i\)​ 作为对应边的边权，恰好 \(m−k\) 个 \(b_i\)​ 作为对应边的边权，该图的最小生成树的边权和 最大 是多少。

\(n\le 9,\ m\le 100\)

\(n\le 9\) 可以状压 dp。我们设 \(f[i][s][j]\) 表示考虑了前 \(i\) 条边，包含了 \(j\) 条 \(a\) 边，连通性为 \(s\)，最大的 mst 边权和是多少。

为了保证求得的是最小生成树，我们应该像 kk 算法一样，按边权从小到大考虑。对于一条边 \((u,v,a,b)\)，我们将它拆成两条边 \((u,v,a)\) 和 \((u,v,b)\)，然后将所有边按边权升序排序。

我们现在尝试向森林 \(s\) 中加入一条边 \((u,v,w)\)，分讨几种情况：

• 森林中 \(u,v\) 未连通：
• 是 \(a\) 边：
  • \(a\le b\)，可选可不选；
  • \(a>b\)，说明排在前面的 \(b\) 边没有选，由于 \(a,b\) 边必须要选一条，因此这条边必须选；
• 是 \(b\) 边：
  • \(a\le b\)，说明排在前面的 \(a\) 边没有选，由于 \(a,b\) 边必须要选一条，因此这条边必须选；
  • \(a>b\)，可选可不选；
• 森林中 \(u,v\) 已经连通：是否选择这条边不会对森林产生影响；

考虑这些情况如何进行转移。

• 森林中 \(u,v\) 未连通：
• 是 \(a\) 边：
  • \(a\le b\)，\(f[i-1][s][j]\to f[i][s'][j+1],\ f[i][s][j]\)
  • \(a>b\)，\(f[i-1][s][j]\to f[i][s'][j+1]\)
• 是 \(b\) 边：
  • \(a>b\)，\(f[i-1][s][j]\to f[i][s'][j],\ f[i][s][j]\)
  • \(a\le b\)，\(f[i-1][s][j]\to f[i][s'][j]\)

然而对于森林中 \(u,v\) 已经连通的情况，我们不知道连接 \(u,v\) 两个连通块的是 \((u,v)\) 对应的另一条边，还是一条不是 \((u,v)\) 的边。因此我们不知道对应的另一条边有没有被选择，也就不知道这一条边是否应该选择，不易处理 \(j\) 的转移。

考虑将 \(j\) 的贡献放在边权较小，先被考虑的一条边上面统计：

• 森林中 \(u,v\) 未连通：
• 是 \(b\) 边且 \(a>b\) 的情况：\(f[i-1][s][j]\to f[i][s'][j],\ f[i][s][j+1]\)
• 森林中 \(u,v\) 已连通：
• 是 \(a\) 边且 \(a\le b\)，\(f[i-1][s][j]\to f[i][s][j],\ f[i][s][j+1]\)
• 是 \(b\) 边且 \(a>b\)，\(f[i-1][s][j]\to f[i][s][j+1],\ f[i][s][j]\)
• 其余情况，\(f[i-1][s][j]\to f[i][s][j]\)，因为贡献已经被另一条边计算过了；

我们搜索出所有连通性状态 \(s\)，在 \(n=9\) 时有 \(Bell(9)=21147\) 种状态。时间复杂度为 \(O\big(m^2Bell(9)\big)\)。

由于 \(n\le 9\)，我们不能使用 \(2^{\binom n2}\) 的纯暴力枚举。考虑枚举当前节点向之前的哪些节点连边，再通过哈希使得每层的每一种连通性只会被处理一遍，可以将时间复杂度降低到 \(O\big(n2^nBell(n)\big)\)，可以通过本题，并且实际跑不满。

代码

#include<iostream>
#include<algorithm>
#include<unordered_set>
#include<unordered_map>
#include<cstring>
using namespace std;
const int N = 10;
const int M = 110;
const int INF = 0x3f3f3f3f;

#define gmax(x, y) x = max(x, y)

struct Edge {
    int u, v, w, t;
    inline bool operator<(const Edge &b) const {
        if(w != b.w) return w < b.w;
        return t > b.t;
    }
} edg[2 * M];
int ne;

int n, m;
int conn[22000][N], nn, ed;

int fa[N][N];
int find(int id, int x) { if(fa[id][x] == x) return x; return fa[id][x] = find(id, fa[id][x]); }
void merge(int id, int x, int y) { x = find(id, x); y = find(id, y); if(x < y) swap(x, y); fa[id][x] = y; }
unordered_set<int> his[N];
unordered_map<int, int> stat;

int hsh(int id) {
    int res = 0;
    for(int i = 0; i < n; i++) find(id, i);
    for(int i = 0; i < n; i++) res = res * n + fa[id][i];
    return res;
}

void dfs(int x) {
    int hs = hsh(x);
    if(his[x].count(hs)) return;
    his[x].insert(hs);
    if(x == n) {
        ++nn;
        for(int i = 0; i < n; i++) conn[nn][i] = fa[x][i];
        stat[hs] = nn;
        bool flag = true;
        for(int i = 0; i < n; i++) flag &= (fa[x][i] == 0);
        if(flag) ed = nn;
        return;
    }
    for(int i = 0; i < (1 << x); i++) {
        for(int j = 0; j < n; j++) fa[x + 1][j] = fa[x][j];
        for(int j = 0; j < x; j++) {
            if(i & (1 << j)) {
                merge(x + 1, x, j);
            }
        }
        dfs(x + 1);
    }
}

int f[22000][M], g[22000][M];

// #define FIO

int main() {

    #ifdef FIO
    freopen("mst.in", "r", stdin);
    freopen("mst.out", "w", stdout);
    #endif

    cin >> n >> m;
    for(int i = 1; i <= m; i++) {
        int u, v, a, b;
        cin >> u >> v >> a >> b;
        --u, --v;
        if(a <= b) {
            edg[++ne] = {u, v, a, 1};
            edg[++ne] = {u, v, b, 0};
        } else {
            edg[++ne] = {u, v, b, -1};
            edg[++ne] = {u, v, a, 0};
        }
    }

    sort(edg + 1, edg + 1 + ne);

    for(int i = 0; i < n; i++) fa[1][i] = i;
    dfs(1);

    for(int s = 0; s <= nn; s++) {
        for(int j = 0; j <= ne; j++) {
            g[s][j] = -INF;
        }
    }
    g[1][0] = 0;

    for(int i = 0; i < ne; i++) {
        int u = edg[i + 1].u, v = edg[i + 1].v, w = edg[i + 1].w, t = edg[i + 1].t;
        for(int s = 1; s <= nn; s++) {
            for(int j = 0; j <= m; j++) {
                f[s][j] = g[s][j];
                g[s][j] = -INF;
            }
        }
        for(int s = 1; s <= nn; s++) {
            for(int k = 0; k < n; k++) fa[0][k] = conn[s][k];
            merge(0, u, v);
            int s1 = stat[hsh(0)], flag = (s1 != s);
            for(int j = 0; j <= m; j++) {
                if(t == 1) {
                    gmax(g[s1][j + 1], f[s][j] + w * flag);
                    gmax(g[s][j], f[s][j]);
                } else if(t == -1) {
                    gmax(g[s1][j], f[s][j] + w * flag);
                    gmax(g[s][j + 1], f[s][j]);
                } else {
                    gmax(g[s1][j], f[s][j] + w * flag);
                }
            }
        }
    }

    for(int i = 0; i <= m; i++) {
        cout << g[ed][i] << '\n';
    }

    return 0;
}

多行状态

P2704 [NOI2001] 炮兵阵地

题意

YbtOj 链接

题面请参考洛谷链接。

我们发现第 \(i\) 行的状态受到第 \(i-1\) 行和 \(i-2\) 行的影响。状态应当同时包含 \(i\) 行和 \(i-1\) 行的信息。

我们设状态 \((i,j_1,j_2)\) 表示第 \(i\) 行的状态为 \(j_1\)，第 \(i-1\) 行的状态为 \(j_2\)，最多能摆放多少个大炮。

\[ f_{i,j_1,j_2}=\max_{j_3\cap j_1=\emptyset\text{ and }j_3\cap j_2=\emptyset}\{f_{i-1,j_2,j_3}+cnt(j_1)\} \]

再通过预处理每行的 valid 数组，找到最小表示法，能极大优化时间复杂度（因为每行的合法状态不超过 \(70\) 个）。

多维状压

有时，单个点的状态数可能大于 \(2\)，此时无法用二进制简单的表示出此状态。为此，对于一个状态 \(i\)，其第 \(j\) 个位置的状态可以表示为 \(\lfloor\frac{i}{k^j}\rfloor \bmod k\)。这样，对于 \(k\) 维的信息，我们也能使用状压处理了。只不过此种方法无法像二进制一样可以进行位运算。