2024-04-07

2次元遅延セグ木の見果てぬ夢

$2$ 次元遅延セグ木は (かなり制限した設定でも) 実現不能と考えられている (以下の資料を参照せよ)。

勉強会をしました。今日のテーマは「2次元セグメント木」です。スライドを公開します。https://t.co/tL8fMsdBv3
— こたつがめ (@kotatsugame_t) 2024年3月18日

本記事では、 $2$ 次元遅延セグ木になれそうでなれないデータ構造たちを紹介する。

可換モノイドの長方形加算長方形和

$(M, +, 0)$ を可換モノイドとし、 $M$ の $2$ 次元列 $a \in M ^ {H \times W}$ を考える。以下のクエリを $\Theta ( \log (H) \log (W) )$ の時間計算量で処理する。

$\mathtt { add }$ : $l, r, d, u \in \mathbb{Z} , x \in M$ が指定され、 $(i , j) \in \lbrack l, r) \times \lbrack d, u)$ について $a _ {i , j} \gets a _ {i , j} + x$ と更新する。
$\mathtt { sum }$ : $l, r, d, u \in \mathbb{Z}$ が指定され、 $\displaystyle \sum _ { (i , j) \in \lbrack l, r) \times \lbrack d, u) } a _ { i , j }$ を出力する。

基本となる知識は区間mul 区間積 O(log N) - よすぽの日記この記事で説明されているデータ構造である。まず、上述のセグ木を $H$ 本用意し、 $i$ 本目のセグ木は $a _ { i , \ast }$ を管理するものとする。セグ木のノードを配列に展開すれば ( $2$ 種類の値を管理するが、適当に並べておく) このセグ木群は $2$ 次元配列と見做すことができ、 $\mathrm { seg } _ {i , k}$ が $i$ 本目のセグ木の $k$ 番目のノードを表す。

$\mathtt { add }$ クエリの処理を考えると、 $i \in \lbrack l , r )$ について、 $\mathrm{seg} _ i$ に $\lbrack d , u )$ への $x$ の区間加算を行うことになる。さらに区間加算の内部の処理を考えると、「 $k, c$ が与えられて、 $i \in \lbrack l , r )$ について $\mathrm{seg} _ {i , k} \gets \mathrm{seg} _ {i , k} + c$ 」という操作 $\Theta ( \log (W) )$ 回で表されることが分かる *1。すると、 $\mathrm { seg }$ を転置して $b _ k := ( \mathrm{seg} _ {i , k} ) _ {i = 0} ^ {H - 1}$ を $1$ つの列として管理したとき、 $b _ k$ への区間加算で操作が表されることになる。これは $b _ k$ たちを区間加算セグ木で管理すれば $\Theta ( \log (H) \log (W) )$ で実行できる。

$\mathtt { sum }$ クエリも、和の順序を適切に入れ替えれば $b _ k$ に対する区間和クエリを用いてシミュレートできることが分かる。

あり得そうな使いどころとしては、 $(M , + , 0) = (\mathbb{Z} , + , 0)$ や $(M , + , 0) = (\mathbb{Z} , \min , \infty)$ が挙げられる。その場合は冪乗が $\Theta ( 1 )$ で求まるため、上記のアルゴリズムのうち複数の部分が簡略化、効率化されるだろう。

kd-tree

可換モノイド $M$ と作用 $F$ が遅延セグ木の要件を満たすとする。 $2$ 次元平面に $N$ 個の点からなる点群 $P$ があり、各点 $p \in P$ が値 $a _ p \in M$ を持っているとき、以下のクエリを $\Theta ( \sqrt N )$ の時間計算量で処理する。

$\mathtt { apply }$ : $l, r, d, u \in \mathbb{Z} , f \in F$ が指定され、 $p \in \lbrack l, r) \times \lbrack d, u)$ について $a _ p \gets f (a _ p)$ と更新する。
$\mathtt { prod }$ : $l, r, d, u \in \mathbb{Z}$ が指定され、 $\displaystyle \prod _ { p \in \lbrack l, r) \times \lbrack d, u) } a _ p$ を出力する。

kd-tree については検索すれば様々な資料が見つかるはずである。注意点として、kd-tree はユークリッド距離での最近点の検索によく使われるが、そちらの計算量を保証するには点群のランダム性が必要となる。

kd-tree は $P$ を $x$ 座標の大小で半分に分け、それらを $y$ 座標の大小で半分に分け、... と軸を交互に切り替えながら点群を分割していく過程を $2$ 分木にしたものである。各ノードが対応する点群についての $a _ p$ の総積などを持てば、( $1$ 次元の) 遅延セグ木と全く同様にしてクエリを処理することができる。すなわち、クエリの長方形 $\lbrack l , r ) \times \lbrack d , u )$ とノードの点群の bounding box の関係を見て、共通部分がないかクエリに包含されていれば直ちに return し、いずれでもなければ子に対して再帰的に処理を行う。

問題は、上記のように再帰を行ったときに訪れるノードの個数を抑えることである。観察すると、あるノードについての処理が子に再帰するとき、クエリの長方形の $4$ 辺のうち少なくとも $1$ 本がノードの bounding box と交差していることが分かる。そこで、軸に平行な直線 $\ell$ を任意にとり、 $\ell$ と bounding box が交差するようなノードの個数を考える。根から深さを $1$ ずつ増やしていくと、 $\ell$ と交差するノードの個数は深さが $2$ 増える毎に $2$ 倍になることが分かる。 kd-tree の深さは $\log _ 2 (N)$ であるから、 $\ell$ と交差するノードは全体で $\Theta ( \sqrt N )$ 個あることが分かる。長方形クエリの場合は $4$ 辺があるため、この $4$ 倍で抑えられることになる。

なお、kd-tree はオーダーが $\Theta ( \sqrt N )$ である上に定数倍があまり良くないため、思ったほど振り回せる道具ではない。

$4$ 分木

(※ 四分木 (Quadtree) と言うと kd-tree のように点群に対するデータ構造を指すようなので、この節のように行列に対するものには別の名前が付いているかもしれない。)

可換モノイド $M$ と作用 $F$ が遅延セグ木の要件を満たすとし、 $M$ の $2$ 次元列 $a \in M ^ {N \times N}$ を考える。以下のクエリを $\Theta ( N )$ の時間計算量で処理する。

$\mathtt { apply }$ : $l, r, d, u \in \mathbb{Z} , f \in F$ が指定され、 $(i, j) \in \lbrack l, r) \times \lbrack d, u)$ について $a _ { i , j} \gets f (a _ { i , j} )$ と更新する。
$\mathtt { prod }$ : $l, r, d, u \in \mathbb{Z}$ が指定され、 $\displaystyle \prod _ { (i , j) \in \lbrack l, r) \times \lbrack d, u) } a _ { i , j}$ を出力する。

正方形領域を縦横に $4$ 分割することを繰り返すと、その過程は $4$ 分木として表される。この木に対しては、遅延セグ木と同じようにクエリを行うことができる。アクセスするノードの個数は $\Theta ( N )$ になってしまうが、これは冒頭で挙げた資料にある通り限界なのでどうしようもない。

この木は、 $N \times N$ の $N ^ 2$ 個の点に対して kd-tree を構築したと見ても本質的には同じである。ただし、一般の点群と比べると bounding box が簡単に計算できるなど特殊な事情があるため、定数倍は良くなるかもしれない。

*1:一般の遅延セグ木クエリだと $c$ が $i$ に依存してしまうので、この議論はできない

2024-03-26

文字列の連結と辞書順の関係を冪級数で証明

概要

cp-unspoiler

この問題の公式解説にある補題群を冪級数を用いて証明する。

定義

この記事内では以下のように定義する。

$\displaystyle \Sigma := \lbrace 1 , 2 , \dots , \sigma \rbrace \\ \displaystyle \Sigma ^ { + } := \Sigma \text{ 上の空でない有限列全体 } \\ \displaystyle \Sigma ^ { \ast } := \Sigma \text{ 上の有限または無限列全体 }$

また、列の大小関係を辞書順で定める。

辞書順関連

$\varepsilon \gt 0$ を十分小さい正の実数とし、 $f \colon \Sigma ^ { \ast } \to \mathbb R$ を $\displaystyle f ( S ) := \sum _ { i = 0 } ^ { \infty } S _ i \varepsilon ^ i$ で定める。

このとき $S , T \in \Sigma ^ { \ast }$ について $S \lt T \iff f ( S ) \lt f ( T )$ となる。

定理 $1$ : $S , T \in \Sigma ^ { + }$ について、 $S ^ { \infty } \lt T ^ { \infty } \iff S T \lt T S$

証明:

\begin{align} & & S ^ { \infty } & \lt T ^ { \infty } \cr & \iff & f ( S ^ { \infty } ) & \lt f ( T ^ { \infty } ) \cr & \iff & \frac { f ( S ) } { 1 - \varepsilon ^ { \lvert S \rvert } } & \lt \frac { f ( T ) } { 1 - \varepsilon ^ { \lvert T \rvert } } \cr & \iff & f ( S ) \left ( 1 - \varepsilon ^ { \lvert T \rvert } \right ) & \lt f ( T ) \left ( 1 - \varepsilon ^ { \lvert S \rvert } \right ) \cr & \iff & f ( S ) + \varepsilon ^ { \lvert S \rvert } f ( T ) & \lt f ( T ) + \varepsilon ^ { \lvert T \rvert } f ( S ) \cr & \iff & f ( S T ) & \lt f ( T S ) \cr & \iff & S T & \lt T S \quad \blacksquare \end{align}

以下の系は、文字列を並び替えてその連結の辞書順を最小化あるいは最大化したい場合に用いる事実として有名である。

系 $2$ : $\Sigma ^ { + }$ 上の順序 $\preceq$ を $S \preceq T \iff S T \leq T S$ で定めると弱順序となる。

定理 $3$ : $S \in \Sigma ^ { + } , T \in \Sigma ^ { \ast }$ について、 $S ^ { \infty } \lt T \iff S T \lt T$

証明:

\begin{align} & & S ^ { \infty } & \lt T \cr & \iff & f ( S ^ { \infty } ) & \lt f ( T ) \cr & \iff & \frac { f ( S ) } { 1 - \varepsilon ^ { \lvert S \rvert } } & \lt f ( T ) \cr & \iff & f ( S ) & \lt f ( T ) \left ( 1 - \varepsilon ^ { \lvert S \rvert } \right ) \cr & \iff & f ( S ) + \varepsilon ^ { \lvert S \rvert } f ( T ) & \lt f ( T ) \cr & \iff & f ( S T ) & \lt f ( T ) \cr & \iff & S T & \lt T \quad \blacksquare \end{align}

LCP関連

$g \colon \Sigma ^ { \ast } \to \mathbb R \lbrack \lbrack x \rbrack \rbrack$ を $\displaystyle g ( S ) := \sum _ { i = 0 } ^ { \infty } S _ i x ^ i$ で定める。また、 $s \in \mathbb R \lbrack \lbrack x \rbrack \rbrack$ について $\mathrm { ord } ( s ) := \min \lbrace k \mid \lbrack x ^ k \rbrack s \neq 0 \rbrace$ とする。

このとき $S , T \in \Sigma ^ { \ast }$ について $\mathrm { len } ( \mathrm { LCP } ( S , T ) ) = \mathrm { ord } ( g ( S ) - g ( T ) )$ となる。

補題 $4$ : $s , t \in \mathbb R \lbrack \lbrack x \rbrack \rbrack$ が $\lbrack x ^ 0 \rbrack t \neq 0$ を満たすなら $\mathrm { ord } ( s ) = \mathrm { ord } ( s t )$

証明: 何かを乗じることで $\mathrm { ord }$ が減少することはない。仮定より $t$ は乗法逆元を持つので、 $\mathrm { ord } ( s ) \leq \mathrm { ord } ( s t ) \leq \mathrm { ord } ( s t t ^ { - 1 } ) = \mathrm { ord } ( s )$ 。 $\blacksquare$

定理 $4$ : $S \in \Sigma ^ { + } , T \in \Sigma ^ { \ast }$ について、 $\mathrm { LCP } ( S ^ { \infty } , T ) = \mathrm { LCP } ( S T , T )$

証明:

\begin{align} \mathrm { len } ( \mathrm { LCP } ( S ^ { \infty } , T ) ) & = \mathrm { ord } ( g ( S ^ { \infty } ) - g ( T ) ) \cr & = \mathrm { ord } \left ( \frac { g ( S ) } { 1 - x ^ { \lvert S \rvert } } - g ( T ) \right ) \cr & = \mathrm { ord } \left ( g ( S ) - g ( T ) \left ( 1 - x ^ { \lvert S \rvert } \right ) \right ) \cr & = \mathrm { ord } \left ( \left ( g ( S ) + x ^ { \lvert S \rvert } g ( T ) \right ) - g ( T ) \right ) \cr & = \mathrm { ord } ( g ( S T ) - g ( T ) ) \cr & = \mathrm { len } ( \mathrm { LCP } ( S T , T ) ) \end{align}

LCP の長さが等しく、 $T$ が共通しているから、 $\mathrm { LCP } ( S ^ { \infty } , T ) = \mathrm { LCP } ( S T , T )$ である。 $\blacksquare$

定理 $5$ : $S , T \in \Sigma ^ { + }$ について、 $S T \neq T S$ ならば $\mathrm { LCP } ( S ^ { \infty } , T ^ { \infty } ) = \mathrm { LCP } ( S T , T S )$

証明:

定理 $4$ を繰り返し用いて、 $\mathrm { LCP } ( S ^ { \infty } , T ^ { \infty } ) = \mathrm { LCP } ( S T ^ { \infty } , T ^ { \infty } ) = \mathrm { LCP } ( S T ^ { \infty } , T S T ^ { \infty } )$ 。 $S T \neq T S$ より、 $\mathrm { LCP } ( S T ^ { \infty } , T S T ^ { \infty } ) = \mathrm { LCP } ( S T , T S )$ 。 $\blacksquare$

定理 $4$ から定理 $5$ を導いたように、定理 $3$ から定理 $1$ を導くこともできる。定理 $5$ を直接示そうとすると、LCP の長さが等しいことは分かるが、LCP 自体が等しいことは上手く示せそうにない。

追記

LCP の一致と辞書順の一致から、LCP の次の文字の組み合わせも一致しているのではないかという推測ができ、これは正しい。文字集合を不定元によって $\Sigma := \lbrace x _ 1, x _ 2, \dots, x _ { \sigma } \rbrace$ と定め、 $\mathbb Z \lbrack x _ 1 , x _ 2 , \dots , x _ { \sigma } \rbrack \lbrack \lbrack x \rbrack \rbrack$ 上で LCP の節と同様の議論を行う。その際、 $\mathrm { ord }$ だけでなく $\lbrack x ^ { \mathrm { ord } ( s ) } \rbrack s$ も $g ( S ^ { \infty } ) - g ( T )$ と $g ( S T ) - g ( T )$ で一致していることが示され、証明を得る。

定理 $6$ : $S \in \Sigma ^ { + } , T \in \Sigma ^ { \ast }$ について、 $k := \mathrm { len } ( \mathrm { LCP } ( S ^ { \infty } , T ) ) = \mathrm { len } ( \mathrm { LCP } ( S T , T ) )$ とする。 $k$ が有限なら $S ^ { \infty } _ k = ( S T ) _ k$

定理 $7$ : $S , T \in \Sigma ^ { + }$ について、 $k := \mathrm { len } ( \mathrm { LCP } ( S ^ { \infty } , T ^ { \infty } ) ) = \mathrm { len } ( \mathrm { LCP } ( S T , T S ) )$ とする。 $k$ が有限なら $S ^ { \infty } _ k = ( S T ) _ k$ かつ $T ^ { \infty } _ k = ( T S ) _ k$

2024-03-16

FPS の合成と逆関数、冪乗の係数列挙 Θ(n (log(n))^2)

概要

温めておいたのですが、ABC で似た問題が出題されてしまったので (想定解は違いましたが) 起源が問題になる前に公開しようと思います。

問題

標数の十分大きい体 $\mathbb{K}$ 上の形式的冪級数 $f ( x ) \in \mathbb{K} \lbrack \lbrack x \rbrack \rbrack$ の $k$ 次までの項が与えられる。 $\lbrack x ^ k \rbrack f ( x ) ^ 0 , \lbrack x ^ k \rbrack f ( x ) ^ 1 , \dots , \lbrack x ^ k \rbrack f ( x ) ^ {n - 1}$ を求めよ。ただし、 $\mathbb{K}$ を係数とする次数 $n$ の多項式乗算が $\mathrm{O} ( M ( n ) )$ で可能とする。

本記事では、この問題に対する時間計算量 $\mathrm{O} ( M ( n + k ) \log ( n + k ) )$ のアルゴリズムを提案する。

アルゴリズム

まず、 $\lbrack x ^ 0 \rbrack f ( x ) = 0$ の場合に帰着できる。何故なら、 $c := \lbrack x ^ 0 \rbrack f ( x )$ とすれば $\lbrack x ^ k \rbrack f ( x ) ^ i = \sum _ {j = 0} ^ {i} \binom{i}{j} c ^ j \lbrack x ^ k \rbrack ( f ( x ) - c ) ^ { i - j }$ であり、 $f ( x ) - c$ に対する答えから畳み込みで求まるためである。

$2$ 変数形式的冪級数 $\mathbb{K} \lbrack \lbrack y \rbrack \rbrack \lbrack \lbrack x \rbrack \rbrack$ 上で *1 、以下の値が $y$ についての多項式として $n - 1$ 次まで求められればいい。

\begin{align} [x ^ k] \sum _ {i = 0} ^ {\infty} y ^ i f ( x ) ^ i = [x ^ k] \frac{1}{1 - y f ( x )} \end{align}

$\lbrack x ^ 0 \rbrack f ( x ) = 0$ と仮定したため分母の $x ^ 0$ の係数は $1$ であり、分母が可逆であるからこの式変形は正当である。

これに Bostan-Mori 法 [1] を適用すると、分子と分母の $y$ についての次数がループごとに倍々になってしまう。しかし $x ^ k$ の係数を求めるためには分子と分母は $x$ について $k$ 次までしか必要ないため、 $x$ についての次数は半々になっていく。結局、 $t$ 回目のイテレーションにおいては $x$ について $k / 2 ^ t$ 次、 $y$ について $2 ^ t$ 次の $2$ 変数多項式の乗算を行うことになり、時間計算量は $\mathrm{O} (M (k) )$ となる。よって全体で時間計算量は $n + \mathrm{O} (M ( k ) \log (k) )$ であり、 $\lbrack x ^ 0 \rbrack f ( x )$ の場合に畳み込みが必要なことと併せても $\mathrm{O} ( M ( n + k ) \log ( n + k ) )$ で抑えられる。

応用

[Tutorial] Generating Functions in Competitive Programming (Part 2) - Codeforces この記事の Lagrange Inversion Formula の項で、 $f$ の逆関数 $g$ に対して $\left ( g ( x ) / x \right ) ^ k$ の計算ができれば冪乗の係数の列挙ができることが示されている。逆に前述のアルゴリズムにより冪乗の係数を列挙すれば、 $\left ( g ( x ) / x \right ) ^ k$ が分かるため $g$ を $\Theta ( n ( \log ( n ) ) ^ 2 )$ で求めることができる。

さらに、逆関数の計算と同じ時間計算量で合成の計算ができることが分かっているため [2]、合成も $\Theta ( n ( \log ( n ) ) ^ 2 )$ となる。

参考文献

[1] Bostan, A., & Mori, R. (2021). A simple and fast algorithm for computing the N-th term of a linearly recurrent sequence. In Symposium on Simplicity in Algorithms (SOSA) (pp. 118-132). Society for Industrial and Applied Mathematics.

[2] Brent, R. P., & Kung, H. T. (1978). Fast algorithms for manipulating formal power series. Journal of the ACM (JACM), 25(4), 581-595.

*1: 少し注意すると、 $\mathbb{K} \lbrack y \rbrack \lbrack \lbrack x \rbrack \rbrack$ で議論しても良い。

2024-02-06

ABC339-F Product Equality の確率を保証する

$\lbrack 10 ^ 9 , 2 \times 10 ^ 9 \rbrack$ から一様ランダムに素数 $P$ を選ぶ *1 。 $A _ i A _ j \neq A _ k$ となる $1$ 組の $( i , j , k )$ を固定し、誤判定する、すなわち $A _ i A _ j \equiv A _ k \pmod P$ となる確率を求める。これは $A _ i A _ j - A _ k$ が $P$ で割り切れることと同値である。

$\lbrack 10 ^ 9 , 2 \times 10 ^ 9 \rbrack$ の範囲に素数は $47374753$ 個存在する *2 。一方、 $0 \lt \lvert A _ i A _ j - A _ k \rvert \leq 10 ^ { 2000 }$ であるから、 $A _ i A _ j - A _ k$ は範囲内の素因数を $222$ 個以下しか持たない。よって一様ランダムに選んだとき、誤判定の確率は $222 / 47374753 \leq 10 ^ { - 5 }$ 以下である。

$4$ つの素数を一様ランダムかつ独立に選ぶと、全ての素数で同時に誤判定する確率は $\left ( 10 ^ { - 5 } \right ) ^ 4 = 10 ^ { - 20 }$ 以下である。独立に選ぶ部分が重要で、独立な事象が同時に発生する確率はその積となるという性質を用いて確率を抑えている。

$1$ つのケースで最大 $N ^ 3 \leq 10 ^ 9$ 組の $( i , j , k )$ に対して判定を行うので、 $1$ 回以上誤判定する確率は $10 ^ { - 20 } \times 10 ^ 9 = 10 ^ { - 11 }$ 以下である。この部分は union bound ( Boole's inequality - Wikipedia ) と呼ばれる直感的にはほとんど明らかな不等式を用いた。この $N ^ 3$ 組では同じ $P$ を使っているため、これらの事象は独立ではない。よって $1 - \left ( 1 - 10 ^ { - 20 } \right ) ^ {10 ^ 9}$ と計算するのは誤りである。

テストケース間では $P$ を取り直すのでそれぞれは独立となるが、独立性を気にせず union bound を使った方が簡単になる。全体で $1000$ ケースある場合、 $1$ ケース以上 WA となる確率は $10 ^ { - 11} \times 1000 = 10 ^ {- 8}$ 以下である。 $1000$ ケースという想定が甘いかどうかと、 $10 ^ {- 8 }$ という確率が許容可能かどうかは各自で判断するしかない。人生が終わるまでにジャッジできる量、宇宙線でアルゴリズムが破損するより低い確率、などで無理やり上限を付けることはできるかもしれない。

(おまけ) 乱択に工夫をするときに気を付けること

本問で独立に素数を選んだ部分で、既に選んだ素数を選ばないようにするという工夫が考えられる。同じ素数を選ぶことは得にならないのでこれは改善になっているが、十分多い候補から素数を選ばなくてはならない都合上、同じ素数を選ぶ確率自体がとてつもなく小さいため効果が無視できるほどしかない。

また、乱択繋がりで rolling hash の基数に原始根を選ぶという工夫がある。ランダムケースに対しては改善になっているかもしれないものの、もし作問者が参加者のそのような考えを見抜いて原始根を集中的に落とすケースを用意した場合は改悪となる *3 。

工夫によりアルゴリズムの挙動が複雑になると、確率の証明を付けることが難しくなりやすい。そもそも証明により安全性が保障されたシンプルなアルゴリズムがあるのだから、それをそのまま使えばよい。

*1: 範囲内から一様ランダムに整数を選び、素数判定することを繰り返す

*2: WolframAlpha に PrimePi[2*10⁹]-PrimePi[10⁹-1] と入力する

*3: 乱数の出目は作問者には想定のしようがないことに注意

2023-12-10

FFT の回数を削減するテクニック集

載ってないテクニックを募集しています。

前提

$\omega _ n := \exp \left ( 2 \pi i / n \right )$ と置く。以降 $i$ と書いたらそれは虚数単位ではなく添え字である。

数列 $\left ( a _ i \right ) _ { i = 0 } ^ { n - 1 }$ に対し、その DFT $\left ( \hat a _ j \right ) _ { j = 0 } ^ { n - 1 }$ を $\hat a _ j = \sum \limits _ { i = 0 } ^ { n - 1 } \omega _ n ^ { i j } a _ i$ で定義する。

基本事項

以下の事実は式変形で確認できる。

関係式 $\left ( 1 \right )$

$\left ( a _ i \right ) _ { i = 0 } ^ { 2 n - 1 }$ に対し、 $\left ( e _ i \right ) _ { i = 0 } ^ { n - 1 } , \left ( o _ i \right ) _ { i = 0 } ^ { n - 1 }$ を $e _ i = a _ i + a _ { n + i } , o _ i = \omega _ { 2 n } ^ i \left ( a _ i - a _ { n + i } \right )$ で定義する。このとき

\begin{align} \hat a _ { 2 j } &= \hat e _ j \cr \hat a _ { 2 j + 1 } &= \hat o _ j \end{align}

が成り立つ。すなわち、長さ $2 n$ の DFT の偶数部分と奇数部分は、それぞれが長さ $n$ のある列の DFT として与えられる。

参考: Cooley-Tukey Algorithm

関係式 $\left ( 2 \right )$

$\left ( a _ i \right ) _ { i = 0 } ^ { 2 n - 1 }$ に対し、 $\left ( e _ i \right ) _ { i = 0 } ^ { n - 1 } , \left ( o _ i \right ) _ { i = 0 } ^ { n - 1 }$ を $e _ i = a _ { 2 i } , o _ i = a _ { 2 i + 1 }$ で定義する。このとき

\begin{align} \hat a _ j &= \hat e _ { j \bmod n } + \omega _ { 2 n } ^ j \hat o _ { j \bmod n } \end{align}

が成り立つ。すなわち、長さ $2 n$ の DFT は、元の列の偶数部分と奇数部分それぞれの DFT を基に計算できる。

参考: https://hcpc-hokudai.github.io/archive/math_fft_001.pdf

DFT の使いまわし

もっとも基本的なテクニック。同じ列の DFT を複数回使うときは $1$ 度計算すれば十分である (それはそう)。やや気付きにくいケースとして、IDFT を行った列を後々もう一度 DFT し直している場合がある。多項式の $2$ 乗を求めるときや、有理式を通分して足し合わせるときなど、あらゆる場面で使う。

middle product

$\left ( a _ i \right ) _ { i = 0 } ^ { n + m - 1 }$ と $\left ( b _ i \right ) _ { i = 0 } ^ { n }$ を畳み込んだ列を $c$ として、 $\left ( c _ i \right ) _ { i = n } ^ { n + m - 1 }$ だけが必要となることがある。 DFT を用いた畳み込みは巡回畳み込みであるため普通は巡回しないよう十分な長さを取るが、この場合はいくらか巡回しても問題なく、長さが $n + m$ 以上の巡回畳み込みで足りる。特に $n = m = 2 ^ k$ の場合が頻出で、計算量が半分になる。ワイルドカードマッチング、転置原理を用いて導出されたアルゴリズム、オンライン畳み込みなどで用いる。

$2$ 倍の長さの DFT の計算

$\left ( a _ i \right ) _ { i = 0 } ^ { n - 1 }$ について $\hat a$ が求まっているとして、 $a$ を長さ $2 n$ に伸ばして $0$ 埋めした列 $\displaystyle \left ( b _ i \right ) _ { i = 0 } ^ { 2 n - 1 } , b _ i = \begin{cases} a _ i \quad \left ( i \lt n \right ) \cr 0 \quad \left ( i \geq n \right ) \end{cases}$ に対する $\hat b$ を計算したい場面がある。関係式 $\left ( 1 \right )$ より $\hat b _ { 2 j } = \hat a _ j$ であり、 $\left ( c _ i \right ) _ { i = 0 } ^ { n - 1 } , c _ i = \omega _ { 2 n } ^ i a _ i$ に対して $\hat b _ { 2 j + 1 } = \hat c _ j$ である。よって $\hat c$ を計算するだけでよい。 subproduct tree を作るときなどに用いる。

$f \left ( - x \right )$ の DFT の計算

$\left ( a _ i \right ) _ { i = 0 } ^ { 2 n - 1 }$ について $\hat a$ が求まっているとして、 $\displaystyle \left ( b _ i \right ) _ { i = 0 } ^ { 2 n - 1 } , b _ i = \left ( - 1 \right ) ^ i a _ i$ に対する $\hat b$ を計算したい場面がある。 $\hat b _ j = \sum \limits _ { i = 0 } ^ { 2 n - 1 } \omega _ { 2 n } ^ { i j } \left ( - 1 \right ) ^ i a _ i = \sum \limits _ { i = 0 } ^ { 2 n - 1 } \omega _ { 2 n } ^ { i \left ( j + n \right ) } a _ i = \hat a _ { \left ( j + n \right ) \bmod 2 n }$ であるから、 $\hat a$ から $\Theta \left ( n \right )$ で計算できる。 Bostan-Mori 法で用いる。

$f \left ( x ^ 2 \right )$ の DFT の計算

$\left ( a _ i \right ) _ { i = 0 } ^ { n - 1 }$ について $\hat a$ が求まっているとして、 $\displaystyle \left ( b _ i \right ) _ { i = 0 } ^ { 2 n - 1 } , b _ i = \begin{cases} a _ { i / 2 } \quad \left ( i \equiv 0 \pmod 2 \right ) \cr 0 \quad \left ( i \equiv 1 \pmod 2 \right ) \end{cases}$ に対する $\hat b$ を計算したい場面がある。関係式 $\left ( 2 \right )$ より $\hat b _ j = \hat a _ { j \bmod n }$ であるから、 $\hat a$ から $\Theta \left ( n \right )$ で計算できる。

前節の内容と組み合わせると以下の事実が分かる。次数 $2 n - 1$ の多項式 $g$ について $g \left ( x \right ) g \left ( - x \right ) = f \left ( x ^ 2 \right )$ が成り立っているとき、 $g$ の DFT の前半と後半を掛け合わせると $f$ の DFT が得られる。

$x f \left ( x \right )$ の DFT の計算

$\left ( a _ i \right ) _ { i = 0 } ^ { n - 1 }$ について $\hat a$ が求まっているとして、 $\displaystyle \left ( b _ i \right ) _ { i = 0 } ^ { n } , b _ i = a _ { \left ( i - 1 \right ) \bmod n }$ に対する $\hat b$ を計算したい場面がある。 $\hat b _ j = \sum \limits _ { i = 0 } ^ { n } \omega _ { n } ^ { i j } a _ { \left ( i - 1 \right ) \bmod n } = \sum \limits _ { i = 0 } ^ { n } \omega _ { n } ^ { \left ( i + 1 \right ) j } a _ i = \omega _ { n } ^ j \hat a _ j$ であるから、 $\hat a$ から $\Theta \left ( n \right )$ で計算できる。

偶数部分 (奇数部分) だけ IDFT

未知の $\left ( a _ i \right ) _ { i = 0 } ^ { 2 n - 1 }$ について $\hat a$ が分かっているとして、 $\left ( e _ i \right ) _ { i = 0 } ^ { n - 1 } , e _ i = a _ { 2 i }$ や $\hat e$ が求めたい場面がある。関係式 $\left ( 2 \right )$ と同様に、 $a$ の偶数部分と奇数部分をそれぞれ $e , o$ と定義する。関係式 $\left ( 2 \right )$ より、 $0 \leq j \lt n$ について $\hat a _ j = \hat e _ j + \omega _ { 2 n } ^ j \hat o _ j , \hat a _ { n + j } = \hat e _ j + \omega _ { 2 n } ^ { n + j } \hat o _ j = \hat e _ j - \omega _ { 2 n } ^ j \hat o _ j$ であるから、 $\hat e _ j = \left ( \hat a _ j + \hat a _ { n + j } \right ) / 2$ である。よって $\hat e$ は $\Theta \left ( n \right )$ で求まり、それを IDFT することで $e$ も求まる。 $o$ に関するものも同様に計算できる。

低次部分 (高次部分) の DFT の取り出し

未知の $\left ( a _ i \right ) _ { i = 0 } ^ { 2 n - 1 }$ について $\hat a$ が分かっているとして、 $\left ( \mathrm { lo } _ i \right ) _ { i = 0 } ^ { n - 1 } , \mathrm { lo } _ i = a _ i$ に対する $\hat { \mathrm { lo } }$ が求めたい場面がある。 $\left ( \mathrm { hi } _ i \right ) _ { i = 0 } ^ { n - 1 } , \mathrm { hi } _ i = a _ { n + i }$ と定義する。関係式 $\left ( 1 \right )$ より、 $\hat a$ の偶数部分は $\mathrm { lo } + \mathrm { hi }$ の DFT である。また、 $\hat a$ の奇数部分を IDFT することで $\omega _ { 2 n } ^ i \left ( \mathrm { lo } _ i - \mathrm { hi } _ i \right )$ が求まる。ここから $\mathrm { lo } - \mathrm { hi }$ が求まるため、これを DFT して $\hat a$ の偶数部分と足し合わせることで $2 \mathrm { lo }$ の DFT が求まる。 $\hat { \mathrm { hi } }$ も同様に計算できる。

長さ $2 n$ の IDFT ののち長さ $n$ の DFT を行うことと比べて $2 / 3$ 倍の時間計算量になる。

$f \left ( x ^ { - 1 } \right )$ や反転した列の DFT

$\left ( a _ i \right ) _ { i = 0 } ^ { n - 1 }$ について $\hat a$ が求まっているとして、 $\displaystyle \left ( b _ i \right ) _ { i = 0 } ^ { n } , b _ i = a _ { \left ( - i \right ) \bmod n }$ に対する $\hat b$ を計算したい場面がある。 $\hat b _ j = \sum \limits _ { i = 0 } ^ { n } \omega _ { n } ^ { i j } a _ { \left ( - i \right ) \bmod n } = \sum \limits _ { i = 0 } ^ { n } \omega _ { n } ^ { - i j } a _ i = \hat a _ { \left ( - j \right ) \bmod n }$ であるから、 $\hat a$ から $\Theta \left ( n \right )$ で計算できる。

IDFT の計算において回転子を $\omega _ n$ から $\omega _ n ^ { - 1 }$ に取り換えた DFT を計算するがこれも同じ式になるので、FFT の実装において IDFT 用に回転子の切り替えを可能にしておく必要はないことが分かる。もっとも、FFT に凝っている人は DFT の結果を bit 反転された状態で管理すると思われ、その場合 IDFT は別で実装することになるので回転子を一致させる旨味は少ないだろう。

$x f \left ( x \right )$ の DFT の計算と組み合わせると $a$ を反転した列の DFT も $\Theta \left ( n \right )$ で求まる。

2023-12-06

push_back / suffix の xor 基底クエリ Θ(log(A))

概要

最初空の数列 $A$ があり、以下のクエリを処理するアルゴリズムを説明する。

$x$ が与えられ、 $A$ の末尾に $x$ を追加する。
が与えられ、の xor に関する基底を出力する。
- ただし $x$ の highest set bit を $\operatorname { bsr } \left ( x \right ) := \max \left \lbrace k \mathrel { } \middle \vert \mathrel { } x _ k = 1 \right \rbrace$ と定義したとき $\operatorname { bsr } \left ( b _ j \right )$ は $j$ について単調減少である。

時間計算量はいずれも $\Theta \left ( \log \left ( \max A \right ) \right )$ である。

アルゴリズム

$0 \leq x \lt 2 ^ w$ のクエリだけが与えられるとして、各クエリを $\Theta \left ( w \right )$ で処理するアルゴリズムを説明する。多少工夫すれば $w$ を各クエリ時点での $\log \left ( \max A \right )$ に合わせることもできる。

基本となる操作は $l \lt r$ に対する $A _ l \leftarrow A _ l \oplus A _ r$ という操作である。 $A$ のいかなる suffix に着目してもこの操作は基本変形になっているから、この操作だけを行っている限りは $\operatorname { span } \left \lbrace A _ i , A _ { i + 1 } , \dots \right \rbrace$ は変化しない。

$A$ の添え字の列 $\left ( i _ 0 , i _ 1 , \dots , i _ { w - 1 } \right )$ を管理する。これらは以下の不変条件を維持する。

$\operatorname { bsr } \left ( A _ { i _ k } \right ) = k$
$i _ k$ として現れない $i$ について $A _ i = 0$

$A$ が最初からいくつかの要素を持っていてもクエリには影響しないため、最初 $A = \left ( 2 ^ 0 , 2 ^ 1 , \dots , 2 ^ { w - 1 } \right )$ として $i _ k = k$ とすれば条件を満たす。 $A = \left ( A _ 0 , A _ 1 , \dots , A _ { n - 1 } \right )$ の末尾に要素を追加した場合、以下の処理を実行し不変条件を回復する。

\begin{align} & j := n \cr & \mathtt { for } ~ k := w - 1 ~ \mathtt { to } ~ 0 \cr & \quad \quad \mathtt { if } ~ \left ( A _ j \right ) _ k = 0 \cr & \quad \quad \quad \quad \mathtt { continue } \cr & \quad \quad \mathtt { if } ~ j \gt i _ k \cr & \quad \quad \quad \quad \left ( j , i _ k \right ) \leftarrow \left ( i _ k , j \right ) \cr & \quad \quad A _ j \leftarrow A _ j \oplus A _ { i _ k } \end{align}

アルゴリズムの正当性は以下の事実に着目することで確認できる。

ループの始まりでは $\operatorname { bsr } \left ( A _ j \right ) \leq k$ であり、ループの終わりでは $\operatorname { bsr } \left ( A _ j \right ) \lt k$ である。
不変条件が、 $A _ j \neq 0$ となり得ることを除き常に満たされている。
$A$ を変更する操作は、 $l \lt r$ に対する $A _ l \leftarrow A _ l \oplus A _ r$ の形でしか行われていない。

処理の終わりでは $A _ j = 0$ となっているから、不変条件が完全に満たされる。 $\left \lbrace A _ i , A _ { i + 1 } , \dots \right \rbrace$ の基底は $\left \lbrace A _ { i _ k } \mathrel { } \middle \vert \mathrel { } i \leq i _ k \right \rbrace$ であるから、ただ出力すればよい。

実際には $A$ を陽には管理せず、 $\left ( i _ k , A _ { i _ k } \right )$ を管理すると空間計算量が抑えられる。

応用

区間 xor 基底クエリ

静的な区間 xor 基底クエリが $\Theta \left ( \left ( N + Q \right ) \log \left ( \max A \right ) \right )$ で処理できる。 $A _ 0 , A _ 1 , \dots , A _ i$ まで処理したときの $\left ( i _ k , A _ { i _ k } \right )$ を全て保持しておけば、オンラインクエリにも対応可能。

実装例

木上のパス xor 基底クエリ

静的な木上のパス xor 基底クエリが $\Theta \left ( N \log \left ( \max A \right ) + Q \left ( \log \left ( \max A \right ) \right ) ^ 2 \right )$ で処理できる。自由に根を決め、根から各頂点までのパスについての $\left ( i _ k , A _ { i _ k } \right )$ を計算しておく。これは各 $v$ について、根から $v$ の親へのパスに対する $\left ( i _ k , A _ { i _ k } \right )$ に $a _ v$ を push すればいいので $\Theta \left ( N \log \left ( \max A \right ) \right )$ で計算できる。 $uv$ パスに対するクエリにおいては、 $\operatorname { LCA } \left ( u , v \right )$ から $u , v$ それぞれへのパスの xor 基底をマージすればよい。それぞれは根から $u , v$ へのパスの suffix であるから本稿のアルゴリズムで求めることができる。

実装例 (LCA の計算に $\Theta \left ( \log \left ( N \right ) \right )$ のアルゴリズムを使っているので計算量は少し増えている)

拡張

アルゴリズムを観察すると、要素の追加が末尾である必要はないことが分かる。すなわち、最初 $S = \emptyset$ として、以下のクエリを処理できる。

$j , x$ が与えられ、 $S \leftarrow S \cup \left \lbrace \left ( i , x \right ) \right \rbrace$
が与えられ、の xor に関する基底を出力する。
- ただし $\operatorname { bsr } \left ( b _ k \right )$ は $k$ について単調減少である。

2023-12-05

木上の Disjoint Sparse Table

概要

$n$ 頂点の木の各頂点にモノイド $M$ の元 $a _ v$ が乗っているとして、指定されたパス上の値の総積を計算するクエリが、前計算 $\Theta \left ( n \log \left ( n \right ) \right )$ 、クエリ $\Theta \left ( 1 \right )$ で処理できる。

アルゴリズム

重心分解を行う。

木の重心を $c$ とし、各 $v \in V$ に対して、 $c v$ パスから $c$ を除いたものに対するクエリ結果を計算しておく *1 。 $u v$ パスが $c$ を含む場合、 $u v$ パスに対するクエリはこの前計算を用いて $\Theta \left ( 1 \right )$ で計算できる。 $c$ を含まない場合に対処するため、木から $c$ を取り除くことで得られる部分木たちに対し再帰的に同様の前計算を行う。

クエリとして $u , v$ が与えられた際に適切な $c$ を発見するためには、重心分解によって得られた木構造において LCA を計算すればよい。 (Disjoint) Sparse Table を用いた、前計算 $\Theta \left ( n \log \left ( n \right ) \right )$ 、クエリ $\Theta \left ( 1 \right )$ の LCA があるためそれを用いる。

良い性質

一般の $x , y \in M$ に対する $x y$ の計算はクエリの際に $1$ 度行われるだけで、前計算はすべて $x \in M, v \in V$ に対する $x a _ v , a _ v x$ の計算を用いている。よって、この種の計算が通常より高速に行える状況では計算量が改善する。例: 各頂点に非負整数が定められており、パス上の非負整数に対する xor 基底をクエリで求めたい場合。

これは Disjoint Sparse Table も持つ性質である。そもそも Disjoint Sparse Table はこのデータ構造をパスグラフに対して適用したものと一致する。したがってこのデータ構造は Disjoint Sparse Table の一般化であると見ることもできる。

LCA を計算する部分の対案

本稿のデータ構造でクエリを行う際に $c$ を計算する部分は、Disjoint Sparse Table においては整数の xor と bsr で行っていた。 1/3 重心分解などを用いて LCA を計算したい木を二分木にし、根から $v$ までのパスの左右を 01 でエンコードすれば、同様に xor と bsr で LCA を計算できる。正確には LCA ではなく LCA の深さが求まるが、そもそも前計算のデータを (頂点, 深さ) に対してデータを対応させるような形式で管理しておけばよい。

*1: モノイドが非可換なら双方向の向きについて計算しておく