グラフの彩色数求値 $O(2^n n)$ や $O(2^n)$ を定数倍高速化したもの

この記事の第 1 部は、競プロ Advent Calender 2 日目として公開されています。

• 第 1 部 2021-12-02 投稿
• 第 2 部 2023-12-14 投稿

第 1 部：問題

Library Checker – Chromatic Number

Library Checker

Library Checker

judge.yosupo.jp

問題文中の $N$ は、ここでは代わりに $n$ と表します。頂点全体の集合を $V$ とします。 $\vert V \vert = n$ です。 タイトルのように、この問題を解く計算量 $O(2^n n)$ や $O(2^n)$ のアルゴリズムを実装し、その定数倍高速化をしました。

グラフの彩色数を用いた解法がある問題の例： (1-1:cp-unspoiler)

また、同じアルゴリズムを用いて和集合に関するある問題を高速に解けます。例： (2-1)

第 1 部：そもそも、そのアルゴリズムは何

頂点集合の部分集合 $s$ について、 $k$ 色で塗ることができるかどうかの真偽値を $\mathrm{dp}_k[s]$ とします。実際は添え字は非負整数を用いて、ビット列で集合を表現します。 $\mathrm{dp}_1$ を求めておきます。

bitwise OR convolution を用いる $O(2^n n^2)$

テーブル $X$ をゼータ変換して得られるテーブルを $\hat{X}$ と表します。

頂点集合 $s$ を $k$-彩色でき、頂点集合 $t$ を $l$-彩色できるとき、 $s \cup t$ は $(k+l)$-彩色可能です。よって、 $2$ つのテーブル $\mathrm{dp}_k$ と $\mathrm{dp}_l$ から、 bitwise OR convolution を用いて $\mathrm{dp}_{k+l}$ を求められます。 ここで、 $\mathrm{dp}_k[V]$ ( $V$ は頂点全体の集合、 $1 \leq k \leq n$ ) に注目すると、 $k$-彩色が判定でき、よって彩色数が求まります。この時点では計算量は $O(2^n n^2)$ 、ダブリングと二分探索を用いると $O(2^n n \log n)$ となります。

$1$ 色塗るごとに計算量が指数的に減少 $O(2^n n)$

彩色後のグラフについて、任意の $k$ 頂点を選ぶとそれらは $k$ 色以下で塗られていることが明らかです。ここで、頂点 $n-k,n-k+1, \ldots ,n-1$ を始めの $k$ 色で塗ることにします。色に色 $0$ 、色 $1$ 、……と番号を振るなら、頂点 $n-k$ の色を色 $x$ として $x \leq k-1$ が成り立ちます。すると、色 $0,1, \ldots ,k$ を塗る頂点の集合を決めた後は $k$ 個の頂点を無視できることになります。したがって、テーブル $\mathrm{dp}_k$ のサイズは色を塗るごとに半分にしてもよく、つまり $k$ の増加とともに指数的に小さくなります。これで計算量 $O(2^n n)$ が達成できます。

ちなみに、この方法では次の紹介する方法に比べて具体的な彩色の構築が容易です。

適当な値で割った余りを用いる $O(2^n n)$

以上で説明したアルゴリズムでは、 bitwise OR convolution の結果を真偽値に変換するため、 bitwise OR convolution の累乗を高速に求める方法(計算量 $O(2^n n)$ ) を適用できません。しかし、真偽値に変換せずに適当な値 $P$ で割ったあまりを用いることで、出力が正しい保証と引き換えに高速な累乗を可能にします。

bitwise OR convolution は以下の手順で行われます。

1. 高速ゼータ変換
2. 各点積
3. 高速メビウス変換

まず、 $\hat{\mathrm{dp}_1}$ の各点を $k$ 乗します。興味があるのは $\mathrm{dp}_k[V]$ のみなので、高速メビウス変換の代わりに $O(2^n)$ で値を求められます。 $\mathrm{dp}_k[V] \neq 0$ ならグラフは $k$ 彩色可能です。 $k=1,2, \ldots ,n$ と順に求めると (各点の $k$ 乗は順に求まるので) $O(2^n n)$ 時間になります。このアルゴリズムが(Library Checker の提出で)よく用いられています。

出力が正しい保証がないと言いましたが、高々 $n$ 個の値のどれかが $P$ の倍数かつ正になって初めて失敗しうるので、出力が正しくない可能性はかなり低いです。

調べた限りでは、この方法で定数倍が小さいものが Library Checker によく提出されています。

以上をまとめて $O(2^n)$

以上の $2$ つの考えを合わせて計算量を $O(2^n)$ にできます。つまり、 $\hat{\mathrm{dp}_k}$ は、 $P$ で割った余りを計算しながら、色を塗るごとにサイズを半分にします。

$\mathrm{dp}_1$ を高速ゼータ変換したテーブルは、場合分けを用いて計算量 $O(2^n)$ で直接計算できます（ Nyaan さんの提出 を参考にしました ）。定数倍高速化の項で改めて解説します。

$P \leftarrow 10^9+9$
$T[i] \leftarrow 1$ $(i=0,1,\ldots 2^n-1)$
$A \leftarrow \hat{\mathrm{dp}_1}$
$T[i] \leftarrow T[i] \cdot (-1)^{\text{popcount}(i) \bmod 2}\ (i=0,1,\ldots 2^n-1)$
$\text{for }k \leftarrow 1,2, \ldots ,n$
$\hspace{20px}z \leftarrow 2^{n-k}$
$\hspace{20px}T[i] \leftarrow (T[i]A[i]) \bmod P\ (i=0,1,\ldots 2z-1)$
$\hspace{20px}T[i] \leftarrow (T[i]+T[i+z]) \bmod P\ (i=0,1,\ldots z-1)$
$\hspace{20px}\text{if }(\sum_{i=0}^{z-1}T[i]) \bmod P =0$
$\hspace{40px}\text{return }k$

いったい何を数え上げているんだ…?

適当な値で割った余りを用いる $O(2^n n)$ では、

• 各頂点を複数の色で塗ってよいとしたときの塗り方の個数

を数え上げることになります。

以上をまとめて $O(2^n)$ では、

• 各頂点を複数の色で塗ってよいが、頂点 $n-1-i$ は番号が $i$ 以下の色しか塗ってはいけないときの塗り方の個数

を数え上げることになります。

第 1 部：これを定数倍高速化します

概要

出力が $1$ の場合は自明なので、別で判定し、テーブルの最大サイズを $2^n$ から $2^{n-1}$ に減らします。つまり、各テーブルは初めから頂点 $n-1$ を塗る前提で構築します。

詳細

まず変数を定義します。

• 非負整数 $E_i$ の $j$ ビット目が立つのは、頂点 $i$ に頂点 $j$ が隣接しているときである。
• $n$ 以下の正整数 $k$ 、 $\{0,1,2,\ldots ,n-1-k\}$ の部分集合 $d$ について、 $\mathrm{dp2}_k[d]$ が $0$ でないのは、 $d$ で表される頂点集合に頂点 $n-k,n-k+1, \ldots n-1$ を足した頂点集合は、 $k$ 色で塗れるとき。ある $k$ についてこのテーブルをまとめて $\mathrm{dp2}_k$ と呼ぶ。
  • 実際は $1$ つの非負整数を bitset として $d$ に用いるため、この意味で非負整数と頂点集合を同一視する場合がある。
  • $P$ で割った余りを用いて高速化するときはこの値が間違っている可能性がある。
• $\hat{\mathrm{dp2}_k}$ は $\mathrm{dp2}_k$ をゼータ変換したもの。
  • 上と同様、この値が間違っている可能性がある。

$\hat{\mathrm{dp}_1}$ の計算

$\hat{\mathrm{dp}_1}[d]$ の値は、

• $d$ の部分集合であって $1$ 色で塗れるものの個数

と解釈できます。頂点 $h$ を塗るかどうかで場合分けすると、 $2^h \leq d \lt 2^{h+1}$ のとき

• 頂点 $h$ を塗らない場合の個数は $\hat{\mathrm{dp}_1}[d \cap \overline{\{h\}}]$
• 頂点 $h$ を塗る場合、 $h$ に隣接する頂点は塗れないので、求める個数は $\hat{\mathrm{dp}_1}[d \cap \overline{\{h\}} \cap \overline{E_h}]$

として、より小さい $d$ についての結果から計算でき、 $\hat{\mathrm{dp}_1}[0]=1$ と合わせると動的計画法でテーブル全体を計算できます。

$\hat{\mathrm{dp2}_1}$ の計算

同様に場合分けし、動的計画法を導きます。

• $\hat{\mathrm{dp2}_1}[0]=1$
• 頂点 $h$ を塗らない場合の個数は $\hat{\mathrm{dp}_1}[d \cap \overline{\{h\}}]$
• 頂点 $h$ を塗る場合、
  • 頂点 $h$ が頂点 $n-1$ に隣接する場合、求める個数は $0$ 。
  • 隣接しない場合、 $h$ に隣接する頂点は塗れないので、求める個数は $\hat{\mathrm{dp}_1}[d \cap \overline{\{h\}} \cap \overline{E_h}]$

オーバーフローの吟味

(2023/01/20 更新)

実は、 $32$ ビット整数で計算すると、ゼータ変換したテーブルどうしで各要素を掛け合わせたときに最大で $2n$ ビット程度の値が現れ、オーバーフローします。ここで、 C++ の符号なし $32$ ビット整数型は加減乗算を $\bmod2^{32}$ で合同な値を用いて計算することを用いて正しく計算します。

(以降、 2023/01/20 追記) この部分は当時の勘違いにもとづいて書かれており、結果的には誤りではなかったものの誤解を招いた可能性がありました。申し訳ございません。

$n=20$ の場合、 bitwize OR convolution の計算結果に表れる値の上界に関して $3^{20}\leq 3.5\times 10^9\lt 2^{32}$ が成り立つので、正しい結果が得られることが示せますが、例えば $n=26$ などとすると間違う可能性を否定できません。 $64$ ビット整数型、より一般に $2n$ ビット以上の表現で計算すれば安全になります。

全体をもう一度

除数 $P$ を用いない計算量 $O(2^n n)$ のアルゴリズムは次のようになります。

$A \leftarrow \hat{\mathrm{dp}_1}$
$T \leftarrow \hat{\mathrm{dp2}_1}$
$\text{if }E_i=\emptyset\ (i=0,1,\ldots n-1)$
$\hspace{20px}\text{return }1$
$\text{for }k \leftarrow 2,3, \ldots ,n$
$\hspace{20px}z \leftarrow 2^{n-k}$
$\hspace{20px}T[i] \leftarrow T[i]A[i]\ (i=0,1,\ldots 2z-1)$
$\hspace{20px}T[i] \leftarrow -T[i]+T[i+z]\ (i=0,1,\ldots z-1)$
$\hspace{20px}T \leftarrow \text{first half of }T$
$\hspace{20px}T \leftarrow \text{Möbius transformed }T$
$\hspace{20px}\text{if }(\sum_{i=0}^{z-1}T[i]) \bmod P =0$
$\hspace{40px}\text{return }k$
$\hspace{20px}T \leftarrow \text{zeta transformed }T$

除数 $P$ を用いる計算量 $O(2^n)$ のほうは次のようになります。

$P \leftarrow 10^9+9$
$A \leftarrow \hat{\mathrm{dp}_1}$
$T \leftarrow \hat{\mathrm{dp2}_1}$
$T[i] \leftarrow T[i] \cdot (-1)^{\text{popcount}(i) \bmod 2}\ (i=0,1,\ldots 2^n-1)$
$\text{if }E_i=\emptyset\ (i=0,1,\ldots n-1)$
$\hspace{20px}\text{return }1$
$\text{for }k \leftarrow 2,3, \ldots ,n$
$\hspace{20px}z \leftarrow 2^{n-k}$
$\hspace{20px}T[i] \leftarrow (T[i]A[i]) \bmod P\ (i=0,1,\ldots 2z-1)$
$\hspace{20px}T[i] \leftarrow T[i]+T[i+z]\ (i=0,1,\ldots z-1)$
$\hspace{20px}\text{if }(\sum_{i=0}^{z-1}T[i]) \bmod P =0$
$\hspace{40px}\text{return }k$

第 1 部：完成したものがこちら

修正 (20 ms)https://t.co/TTOwLimcMo
運ゲー (6 ms fastest)https://t.co/Xt159wBURL

— Nachia (@NachiaVivias) November 18, 2021

(2024-01-15 コンパイルエラー修正)

• $O(2^n n)$ , 19 ms : https://judge.yosupo.jp/submission/183424
• $O(2^n)$ , 6 ms : https://judge.yosupo.jp/submission/183425

第 2 部：決定性 $O(2^n)$ 時間へ

アルゴリズム

「以上をまとめて $O(2^n)$ 」のアルゴリズムにおいて、 ${}\bmod P$ を取らずに多倍長整数で計算することを考えます。実はこれが $O(2^n)$ 時間の deterministic アルゴリズムになります。（！？）

高速ゼータ変換をした後の配列に含まれる値は $2^n$ 以下です。よって、 bitwise OR convolution で $k$ 乗した配列を高速ゼータ変換した列に含まれる値は $(2^n)^k$ 以下です。これを $k$ ワードで管理します。

$k$ 色目を塗る処理にかかる計算量は $O(k2^{n-k})$ 時間と表せるので、全体の計算量は

$$\sum_{k=1}^{n}k2^{n-k}\leq 2^n \sum_{k=1}^{\infty}k2^{-k} =2^n\sum_{s=1}^{\infty}2^{-s+1}=2^n 2^2$$

と評価して $O(2^n)$ です。

実装

実装しました。 (2024-01-15 修正)

• (9 ms) https://judge.yosupo.jp/submission/183422

第 1 部での制作と比較すると、オーバーフローや除数による不正確さといった問題点が取り除かれたうえで、処理速度を維持することができました。

この実装では、正しく動作する範囲として $1\leq n\leq 31$ を想定しています。

実装のポイント

メモリの戦略

第 1 部では一貫して、番号の大きい頂点から取り除くことで、アクセスする配列の範囲を狭めることを優先しました。今回は多倍長整数を扱うので $1$ 要素あたりのメモリが一定ではありません。

このため、頂点を取り除くときは番号の小さい頂点から順に取り除くことにして、配列を間引くように空間を減らすようにします。これにより、 $k$ 頂点取り除いた後は $1$ 要素あたり $2^k$ ワードを使えるようになり、必要な $k+1$ ワードを確保できます。

多倍長整数の演算

ゼータ変換の長さを半分にするステップでは、

• $1$ ワードの非負整数値 $a,b$ と、 $k$ ワードの非負整数値 $x,y$ が与えられたとき、
• $k+1$ ワードの非負整数値 $ax-by$ を計算する

という計算が繰り返されます。 $ax$ と $by$ でそれぞれ繰り上がりを計算してから単に引き算することにより、引き算の繰り上がりの処理が簡潔になるようにしました。

また、彩色可能か判定するためにテーブルの値を足し引きするときは、符号・桁ごとに足し合わせてから繰り上がりを処理し、正の成分と負の成分それぞれの総和が一致するかどうか判定するように実装しました。

おわりに

第 1 部では地道な高速化で差を付けることにしていましたが、第 2 部では競プロでよく用いられる方法と比べて確かに高速な解法プログラムを作ることができました。ソースコード長もそんなに大きくありません。 Library Checker の問題設定において、競プロ用ライブラリとしてはかなり良いところに落ち着いているのではないでしょうか。