Learning Algorithms

アルゴリズムの勉強メモ

彩色数を求めるアルゴリズム

はじめに

彩色数を求めるアルゴリズムについて書きました。

English version is available here.

彩色数とは

彩色数(chromatic number)とは、無向グラフにおいて、辺で繋がれた頂点同士が、互いに異なる色でなければいけないという制約のもとで、すべての頂点に彩色をするために最低必要な色の数のことです。例えば、辺が無いグラフだと彩色数は $1$ で、完全グラフでは彩色数は頂点数に等しくなります。

さらに、彩色多項式(chromatic polynomial)とは、あるグラフに対して、その頂点をちょうど $K$ 色で塗り分ける方法の総数を $K$ の関数として求めた多項式のことです。

ここで例として木の彩色多項式 $P(K)$ を簡単に求めてみます。まず、木には少なくとも $1$ つの葉が存在するので、それを任意の色で塗ります。するとその隣り合う頂点には今塗った色以外の任意の色が塗れます。木の構造から、一つ前に塗った色以外を順に塗っていけるので、結局木の彩色多項式は $K \geq 2$ のとき以下のようになります。ただし $n$ は頂点数です。

$P(K) = K(K - 1)^{n - 1}$

$K = 1$ のときは明らかに塗り分ける方法が存在しないことから $P(K) = 0$ となり、任意の $K (\geq 1)$ について上の式が成り立つので、これが木の彩色多項式になります。

その他のいくつかの特殊なグラフの彩色多項式についてはWikipediaに書かれているので、興味がある方はご覧ください。

さて、実際に彩色数を求めるときには、次のようにします。

1. 彩色多項式 $P(K)$ を求める
2. $K = 1$ から順に $P(K)$ を評価していく
3. $P(K)$ が初めて $0$ より大きくなったとき、その $K$ が彩色数になる

シンプルですね。

再び例として、木の彩色数を求めてみると、$K(1) = 0$ であり、$K(2) = 2 > 0$ となるので、彩色数は $2$ であるとわかります(計算結果の $2$ は木を $2$ 色で塗る方法の総数になっていて、最初にどちらの色を塗るかの $2$ 通りです)。

ここからは、任意の無向グラフの彩色数を求めていきます。

解説

さて、任意の無向グラフの彩色数を求めるには、以下のようなことをすればいいです。

1. 次のスライドを見る
指数時間アルゴリズム入門

2. 実装する

シンプルですね。

補足

少しだけ下手な補足をします。

まず $K$ 色で彩色可能であるという条件を、$K$ 個の独立集合に分割可能、と言い換えます。独立集合は必ず $1$ 色で彩色可能であることからこのような言い換えができます。こうすると、かなり見通しが良くなります。

そこで、グラフの頂点集合を $V$ として、任意の頂点集合 $U \subseteq V$ について、「$U$ の部分集合であって、独立集合をなすものの個数」を求め、これを $ind(U)$ と書くことにします。

ここから $K$ 個選ぶ方法の総数は $ind(U)^K$ なので、包除原理より次の彩色多項式を得ることができます。

$P(K) = \sum (-1)^{|V| - |U|} ind(U)^K$

これをそのまま実装します。$ind(U)$ の計算はふつうの $bitDP$ をしています。また、今回は彩色多項式の値自体には興味がないので、計算途中の値が大きくなってしまっても大丈夫なように、十分大きい素数による剰余をとっています(ちゃんとランダムにとったほうがよいです)。

実装
int ChromaticNumber(const vector<int> &g) {
        int n = g.size();
        if (n == 0) return 0;
        //randomly choose a large prime
        const int modulo = 1077563119;
        int all = 1 << n;
        vector<int> ind(all), s(all);
        for (int i = 0; i < all; i ++) s[i] = ((n - __builtin_popcount(i)) & 1 ? -1 : 1);
        ind[0] = 1;
        for (int i = 1; i < all; i ++) {
                int ctz = __builtin_ctz(i);
                ind[i] = ind[i - (1 << ctz)] + ind[(i - (1 << ctz)) & ~g[ctz]];
                if (ind[i] >= modulo) ind[i] -= modulo;
        }
        //compute the chromatic number (= \sum (-1)^{n - |i|} * ind(i))^k)
        for (int k = 1; k < n; k ++) {
                long long sum = 0;
                for (int i = 0; i < all; i ++) {
                        long long cur = ((s[i] * (long long) ind[i]) % modulo);
                        s[i] = cur;
                        sum += cur;
                }
                if (sum % modulo != 0) return k;
        }
        return n;
}

これで彩色数が求められました。計算量は $O(2^n\ n)$ です。

最小クリーク被覆について

ついでに、最小クリーク被覆(minimum clique cover)について触れておきます。グラフをいくつかのクリークに分割して、重複なく頂点集合を覆うとき、その分割する個数が最小になるようなもののことを最小クリーク被覆といいます。ところで、あるグラフにおけるクリークは、その補グラフにおいて独立集合をなすのでした。彩色数は、「最小独立集合被覆」とも言えるので、これらのことから次が言えます。

あるグラフの彩色数は、その補グラフにおける最小クリーク被覆のサイズに等しい

彩色数については以上です。

例題

以上で述べたことを使うと、以下の問題が解けるので、解いてみてください!

Problem - H - Codeforces

解説はまた別のページに書くかもしれませんが、実装だけ載せておきます。

#include <iostream>
#include <cstdio>
#include <algorithm>
#include <string>
#include <vector>
#include <map>
#include <set>
using namespace std;

unsigned long xor128() {
        static unsigned long x = 123456789, y = 362436069, z = 521288629, w = 88675123;
        unsigned long t = (x ^ (x << 11));
        x = y; y = z; z = w;
        return (w = (w ^ (w >> 19)) ^ (t ^ (t >> 8)));
}
int RandomPrime() {
        auto prime = [&](int n) {
                for (int i = 2; i * i <= n; i ++) if (n % i == 0) return false;
                return true;
        };
        int modulo = 1000000000;
        modulo += (int) (xor128() % 100000000);
        while (!prime(modulo)) modulo ++;
        return modulo;
}
int ChromaticNumber(const vector<vector<bool>> &tmpg) {
        int n = tmpg.size();
        if (n == 0) return 0;
        vector<int> g(n);
        for (int i = 0; i < n; i ++) {
                for (int j = 0; j < n; j ++) {
                        if (tmpg[i][j]) {
                                g[i] |= (1 << j);
                        }
                }
        }
        int all = 1 << n;
        vector<int> a(all), s(all);
        a[0] = 1;
        int modulo = RandomPrime();
        for (int i = 1; i < all; i ++) {
                int j = __builtin_ctz(i);
                a[i] = a[i - (1 << j)] + a[(i - (1 << j)) & ~g[j]];
                if (a[i] >= modulo) a[i] -= modulo;
        }
        for (int i = 0; i < all; i ++) {
                s[i] = ((n - __builtin_popcount(i)) & 1 ? -1 : 1);
        }
        for (int k = 1; k < n; k ++) {
                long long sum = 0;
                for (int i = 0; i < all; i ++) {
                        long long cur = ((s[i] * (long long) a[i]) % modulo);
                        s[i] = (int) cur;
                        sum += cur;
                }
                if (sum % modulo != 0) return k;
        }
        return n;
}

int main() {
        int n;
        scanf("%d", &n);
        vector<string> c(n);
        vector<vector<bool>> reach(n, vector<bool> (n));
        for (int i = 0; i < n; i ++) cin >> c[i];
        for (int i = 0; i < n; i ++) {
                for (int j = 0; j < n; j ++) {
                        reach[i][j] = (i == j || c[i][j] == 'A');
                }
        }
        for (int k = 0; k < n; k ++) {
                for (int i = 0; i < n; i ++) {
                        for (int j = 0; j < n; j ++) {
                                if (reach[i][k] && reach[k][j]) {
                                        reach[i][j] = true;
                                }
                        }
                }
        }
        for (int i = 0; i < n; i ++) {
                for (int j = 0; j < n; j ++) {
                        if (reach[i][j] && c[i][j] == 'X') {
                                return !printf("-1\n");
                        }
                }
        }
        vector<int> idx(n, -1);
        vector<bool> used(n, false);
        int N = 0;
        for (int i = 0; i < n; i ++) if (!used[i]) {
                vector<int> v;
                for (int j = 0; j < n; j ++) if (reach[i][j]) v.push_back(j);
                if (v.size() >= 2) {
                        for (int j = 0; j < v.size(); j ++) {
                                idx[v[j]] = N;
                        }
                        N ++;
                }
                for (int j = 0; j < v.size(); j ++) used[v[j]] = true; 
        }
        vector<vector<bool>> g(N, vector<bool> (N));
        for (int i = 0; i < n; i ++) {
                for (int j = 0; j < n; j ++) {
                        if (c[i][j] == 'X' && idx[i] != idx[j] && idx[i] >= 0 && idx[j] >= 0) g[idx[i]][idx[j]] = true;
                }
        }
        int ans = n - 1 + ChromaticNumber(g);
        printf("%d\n", ans);
        return 0;
}