高速フーリエ変換の超丁寧な解説! Atcoder社の解説を詳解!
はじめに
学生最強コンのAが高速フーリエ変換そのものであると聞き、高速フーリエ変換を勉強しました。ネット上のフーリエ級数、変換の資料をあさりつつ、tsutajさん、Atcoder社の公式解説をいろいろ眺めて理解に努めようとしたが、頭が弱いため完全理解には6時間かかりました。
けんちょんさんがネットの集合知をまとめてくれることを敬虔に祈念しつつ、自分の勉強も兼ねて、これから高速フーリエ変換を勉強する方々の助けになれるようにこの記事を書きます。
この記事はAtcoder社の公式の高速フーリエ変換のスライドを見ながら、ページごとの詳解をするスタンスを取っています。スライドを同時に見ることを *強く* 推奨します!
また、スライドの後ろの発展的な話題については触れません。初心者向けの記事である上に僕はまだ理解できてません。
まだ始まってもないのですが、ソースコードを大変参考になったtsutajさん、この解説のベースとなったAtcoderの公式解説さん、理解が進まない頭の悪い僕に一生懸命説明したhotman78さんに強い感謝の気持ちをここで示します。
最低限の数学的知識
これさえわかれば、問題ありません。
また、この記事ではアルゴリズムと実装の理屈が分かって使えればよいので、数学的な厳密性は全く保証しません。
解説
これから、スライドの(Xページの)Y行目という表現を多用します。Y行目は、数式、タイトルも一行として数えます。
例えば、スライドの3ページでは、計7行あり、1行目にタイトル、4行目に数式があります。
1ページ
がんばるぞい!
ところでTypical Contestの難易度差激しすぎませんか?
2ぺージ
この問題の類題は2019年の学生最強コン決勝Aにも出ました。高速フーリエ変換を貼るだけで終わるみたいですが、もっと効率のいい解法が存在します。それについてもいずれかきます。(To do)
3ページ
主菜(A)の価格が分かれば、副菜(B)の価格もわかります。とが同じ時、Aは種類だけ、Bは種類だけ、それぞれ自由に選んでよいので、スライドのとおりの式になります。
後は、スライドの式のとおり、を1からまでループを回せば集計できます。
4, 5, 6ページ
唐突に多項式の係数にしてますが、スライドのとおり、二つの多項式の掛け算にすると、掛けてできた多項式の係数は、なんと先ほどの3ページのごとの計算や集計の結果となるのだ!
ちなみに、この際とする。こうすると、主菜や副菜のうち一品だけでになるのが答えに入らない。なぜなら、係数の掛け算で定数項が0なので、何を掛けても定数項は増えないし、ほかの項に影響を与えない。
つまり、定石みたいなものである!
ならば、
となる。スライド3の方法と一致することを確認できる。
7ページ
スライドにある通りです。
8ページ
スライド7ページのとおりに、うまく個以上の点を選び、そこからうまく補間すれば、掛けた式の係数が分かりそう!
9, 10ぺージ
急に複素数が出てきます!なぜ1のn乗根を持つのかというと、
10ページに書いてある通り
- 乗()間では異なるからは異なるが得られる。(3行目)
- 四行目のあるような式を計算すると、キレイサッパリします!なぜその式を計算するのかは、後々わかります。
つまり、簡単に言うと、をn個用意して、多項式に代入する形を取れば、計算上都合がいいのです!
11ページ
9ページに続いて、急によくわからないものを定義し始めました。スライドに書いてある通り、今まで考えてたのの計点を代入すると、個のになる。
番目のに代入した結果がの係数になるようなまた別の式を作っています。
急にそんなことをやりだす理由としては、こうやるとうまくいく!と今考えたほうがいいです。
ちなみに、この離散フーリエ変換(DFT)の形は僕の見てた参考資料と少し違います。競プロでFFTを扱う分には、フーリエ変換のほかの兄弟も知る必要はおそらくないので、気になる人は参考資料のサイトをご覧下さい。このアルゴリズムを理解する分には、こうやっておくとなんか計算がすごく楽!という理解で構いません。
12ページ
先ほど(11ページで)作ったばっかりのの式は、の足し合わせで表してました。
ここで、大もとで求めたかったはn次の多項式であることを思い出す。
多項式ならば、11ページでののように、と、係数を掛けたものを足し合わせた形で書けます(それはそう)
11ページのように、もその足し算で書いてみました。この際、それぞれのにかかってる係数はわからない(というかこれが知りたい)ので、ひとまずとスライドの中ではおいてます。
そして、12ページでの式変形の結果、どうやらもともと4行目で奥の()の中の足し算に入っていたを、外にくくりだせたようです。
次のソースコードの中の場所のビフォーアフターと同じです。
for(int i = 0; i < A; i++){ //今はココ for(int j = 0; j < B; j++){ //もともとココ } }
数学では、一つ外のループに出すことで式変形がうまくいく場合があります。今回もそのようです。
次のページで、なぜ1のn乗根をわざわざ代入したのか、なぜくくりだしたのか、解決します!
13, 14ページ
非実数であるを満たすは、
を満たすので、が成り立ちます。この時のは実数でない1の乗根のいずれでも成り立つ。
3行目での時、[\zeta_n^{j - k}]は必ず1のn乗根のうちのどれかになります。
そして、1の非実数のn乗根は先ほど述べた累乗を足した式を必ず満たします。3行目では、奥に入ってるΣではを、つまりからまで足し合わせるので、先ほどのn乗根の満たすべき式の右辺そのものになって、0になります!
の時は、0乗なので1になり、それが回足し算されるので、最終的にはスライドのとおり、と、係数の倍になります。つまり、で割ればの係数がわかる!この変換を離散フーリエ逆変換といいます(inverse DFT)。
いままでの変な式や代入するを選んだのは、うまく選ぶとこうなるからです。
15ページ
の定義で、とするところを、で書き直して、としてます。
スライドのとおり、どうやらフーリエ変換はを先に変換しても、掛ければの変換になります!
この性質は、次頁アルゴリズムで重要な役割を果たします。
16. 17, 18ページ
流れのまとめです。
ここから離散フーリエ変換と逆変換が高速でできれば幸せだなぁーとわかります。そこが時間計算量ボトルネックですので。
逆フーリエ変換は、最終的にフーリエ変換のソースコードのを-1倍して、最後にで割るだけなので、本質的にはフーリエ変換を頭よくやります。
19, 20ページ
を、のの偶奇によって、スライド通りに分けると
と となります。
をでくくりだすと、二つの式はで置換をしてみると
と
とかけます。
スライドの6行目の式のように、の値を知るには、どうやらとの値が分かればいいとなります。
この時点ではまだ意味不明のようだが、次の2,3枚スライドで意味が分かります。
そして、をすべて通り代入すると、スライド20ページにある値が必要となります。
21ページ
なんということでしょう!2行目のように式をいじくると、乗根のおかげで、になります!実をいうと、三角関数ので表すとまったく同じものであるともわかります!
これのおかげで、20ページで求めたかった2つの数列は7~8行目のように変形できます。この変形のうまあじは次のスライドで伝わります。
22ページ
乗根の周期性で、周期ががになったので、なんと20ページで求めたかった数列の前半と後半は実は同じだったのです!
例を挙げると
は、
と同じで、
それは
と同じである。そしてこれらは、2つで周期になってる。
23, 24ページ
22ページの性質を利用すると、
「個の値を代入した値を計算するためには、2つのの値の数列がわかればで計算できます。」
例えば、8個の点を代入した値の列がほしいのなら、4個の点を代入した列を2つ、求めればよい。そうすると、次図のように完全二分木みたいな感じの再帰関数になります。(次図はの離散フーリエ変換の例)
詳しい図の解説は次回の記事でします。
25ページ
具体的な実装となります。僕もこれを実質写経しました。
26~30ページ
見なかったことにしてください。おそらく初心者向きでは全くありません。
最後に
ここまで読んでいただきお疲れさまでした。ここまでで5555字だそうです。
上記の解説だけで、ピンとくる人は少ないと思うので、続編として実装例と、具体的な値を入れた時の挙動の解説をします。
何か誤字脱字、論理上のミスなどあったらご指摘お願いします。