ゼファーネットのロゴ

高速フーリエ変換:多点評価のスケーリング

日付:

高速フーリエ変換(FFT)は、大きな数の乗算や多項式の乗算など、多くのアルゴリズムの重要な構成要素です。 フーリエ変換は、信号処理、量子力学、その他の分野でも重要な用途があり、世界経済の重要な部分を実現するのに役立ちます。 XNUMXつの異なる操作について説明します。多点多項式評価(次数の評価)とその逆多項式補間(次数の評価が与えられた場合)です。

画像

トリガー警告:特殊な数学的トピック

フィードバックを提供してくれたKarlFloerschに特に感謝します

数論でより興味深いアルゴリズムのXNUMXつは、高速フーリエ変換(FFT)です。 FFTは、次のような多くのアルゴリズムの重要な構成要素です。 大きな数の非常に高速な乗算、多項式の乗算、およびの非常に高速な生成と回復 イレイジャーコード.

特にイレイジャーコードは非常に用途が広いです。 フォールトトレラントなデータの保存と回復における基本的なユースケースに加えて、イレイジャーコードには次のようなより高度なユースケースもあります。 スケーラブルなブロックチェーンでデータの可用性を確保する & スターク。 この記事では、高速フーリエ変換とは何か、およびそれらを計算するためのより単純なアルゴリズムのいくつかがどのように機能するかについて説明します。

経歴

画像
画像
画像

わかりました。フーリエ変換は、信号処理、量子力学、その他の分野でも非常に重要な用途があり、世界経済の重要な部分を実現するのに役立ちます。 しかし、さあ、象は涼しいです。

オリジナル フーリエ変換 は、「周波数領域」と「時間領域」の間でデータを変換することとしてよく説明される数学演算です。 これがより正確に意味するのは、データがある場合、アルゴリズムを実行すると、周波数と振幅が異なる正弦波のコレクションが生成され、それらを合計すると元のデータに近似するということです。 フーリエ変換は、次のような素晴らしいことに使用できます。 従円と周転円を介して正方形の軌道を表現する & 象を描くことができる方程式のセットを導出する:

この投稿で説明するフーリエ変換の種類は、 連続的な フーリエ変換 実数または複素数、 それは 離散フーリエ変換 が 有限体 (「モジュラー数学の幕間」セクションを参照してください。 ここ 有限体とは何かについての復習のために)。

「周波数領域」と「時間領域」の間の変換について説明する代わりに、ここではXNUMXつの異なる操作について説明します。 多点多項式評価 (学位の評価 <N での多項式 N 異なるポイント)とその逆、 多項式補間 (学位の評価を与えられた <N での多項式 N 異なる点、多項式の回復)。 たとえば、モジュラス5の素数体で操作している場合、多項式は y =x²+ 3 (便宜上、係数を昇順で書くことができます: 【3,0,1])ポイントで評価 【0,1,2] 値を与える 【3,4,2] (会員登録はお済みでしょうか? 【3,4,7] なぜなら、数値が5)でラップアラウンドする有限体で動作しているため、実際に評価を行うことができるからです。 【3,4,2] そしてそれらが評価された座標(【0,1,2])元の多項式を復元する 【3,0,1].

マルチポイント評価と補間の両方に、どちらの操作も実行できるアルゴリズムがあります。 O(N²) 時間。 マルチポイント評価は単純です。各ポイントで多項式を個別に評価するだけです。 これを行うためのPythonコードは次のとおりです。

def eval_poly_at(self, poly, x, modulus): y = 0 power_of_x = 1 for coefficient in poly: y += power_of_x * coefficient power_of_x *= x return y % modulus

アルゴリズムは、すべての係数を通過するループを実行し、係数ごとにXNUMXつのことを実行するため、次のように実行されます。 オン) 時間。 マルチポイント評価では、この評価を次の場所で行います。 N 異なるポイントなので、合計実行時間は O(N²).

ラグランジュ補間はより複雑です(「ラグランジュ補間」を検索してください) ここ より詳細な説明については)。 基本戦略の重要な構成要素は、あらゆるドメインの構成要素です。 D そしてポイント x、を返す多項式を作成できます 1 for x & 0 の任意の値 D 以外の x。 例えば、もし D=【1,2,3,4] & x=1、多項式は次のとおりです。

画像

あなたは精神的にプラグインすることができます 123 & 4 上記の式に変換し、それが返されることを確認します 1 for x=1 & 0 他のXNUMXつの場合。

これらの多項式を乗算およ​​び加算することにより、特定のドメインで任意の出力セットを提供する多項式を復元できます。 上記の多項式を呼び出すと P1、および同等のもの X = 2X = 3X = 4P2P3 & P4、次に戻る多項式 【3,1,4,1] ドメイン上 【1,2,3,4] 単に 3⋅P1+P2+4⋅P3+P4。 円周率多項式の計算には オン²) 時間(最初に、ドメイン全体で0に戻る多項式を作成します。これには時間がかかります。 O(N²) 時間、それからそれを別々に割る (x-xi) それぞれ xi)、および線形結合の計算には別の時間がかかります O(N²)時間なので、 O(N²) 実行時の合計。

高速フーリエ変換

このはるかに高速なアルゴリズムを使用するために支払わなければならない代償があります。つまり、任意のフィールドと任意のドメインを選択することはできません。 一方、ラグランジュ補間では、任意のx座標とy座標、および任意のフィールドを選択でき(単純な古い実数に対しても実行できます)、FFTを使用してそれらを通過する多項式を取得できます。 、有限体を使用する必要があり、ドメインは 乗法部分群 フィールドの(つまり、「ジェネレータ」値の累乗のリスト)。

たとえば、モジュロを法とする整数の有限体を使用できます。 337、およびドメインでの使用 【1,85,148,111,336,252,189,226] (それがの力です 85 フィールドで、例えば。 85 ^ 3%337 = 111; で止まります 226 の次の力のため 85 1)に戻ります。 さらに、乗法部分群にはサイズが必要です 2n (フォームの数に対して機能させる方法があります 2 ^m⋅3^ n そしておそらくわずかに高い素数冪ですが、それははるかに複雑で非効率的になります)。 モジュロを法とする整数体の有限体 59たとえば、位数の乗法的サブグループしかないため、機能しません。 2、29 & 58; 2 小さすぎて面白くない、そして要因 29 大きすぎてFFTに対応できません。 サイズの乗法群から生じる対称性 2 ^ n はるかに少ない作業量から必要な結果を非常に巧妙に計算する再帰的アルゴリズムを作成できます。

アルゴリズムと、実行時間が短い理由を理解するには、再帰の一般的な概念を理解することが重要です。 再帰的アルゴリズムは、アルゴリズムへの入力が十分に小さく、出力を直接与えることができる「基本ケース」と、必要な計算が何らかの「グルー計算」で構成される「再帰的ケース」のXNUMXつのケースがあるアルゴリズムです。さらに、同じアルゴリズムをXNUMXつ以上使用して、入力を小さくします。 たとえば、リストの並べ替えに再帰的なアルゴリズムが使用されているのを見たことがあるかもしれません。 リストがある場合(例: 【1,8,7,4,5,6,3,2,9])、次の手順を使用して並べ替えることができます。

  • 入力にXNUMXつの要素がある場合、それはすでに「ソート」されているため、入力を返すだけで済みます。
  • 入力に複数の要素がある場合は、リストの前半と後半を別々に並べ替えてから、並べ替えられたXNUMXつのサブリストをマージします(それらを呼び出します)。 A & B) 次のように。 XNUMXつのカウンターを維持し、 アポス bpos、両方ともゼロから始まり、空で始まる出力リストを維持します。 どちらかまで アポス or bpos 対応するリストの最後にある場合は、 A[アポス] or B[bpos] 小さいです。 どちらか小さい方の値を出力リストの最後に追加し、そのカウンターを1増やします。これが完了したら、完全に処理されていない残りのリストを出力リストの最後に追加し、出力を返します。リスト。

XNUMX番目の手順の「接着剤」にはランタイムがあることに注意してください オン):XNUMXつのサブリストのそれぞれに N 要素の場合、各リストのすべてのアイテムをXNUMX回実行する必要があるため、次のようになります。 オン) 計算合計。 したがって、アルゴリズムは全体として、サイズの問題をとることによって機能します N、そしてそれをサイズのXNUMXつの問題に分割します N2、プラス オン) 「接着剤」の実行の。

と呼ばれる定理があります マスター定理 これにより、このようなアルゴリズムの合計実行時間を計算できます。 多くのサブケースがありますが、サイズNの実行をサイズのk個のサブケースに分割する場合 N / k   オン) 接着剤(ここの場合のように)、結果は実行に時間がかかるということです O(N・log(N)).

画像

FFTも同じように機能します。 サイズの問題を取ります N、サイズのXNUMXつの問題に分割します N2、 そして、やります オン) 接着剤は、小さなソリューションを大きなソリューションに結合するために機能します。 O(N・log(N)) 実行時の合計– はるかに高速 より O(N2)。 これが私たちのやり方です。 最初に、マルチポイント評価(つまり、一部のドメイン)にFFTを使用する方法について説明します。 D および多項式 P、計算 P(x) すべてのための x in D)、そして、わずかな調整で補間に同じアルゴリズムを使用できることがわかりました。

与えられた定義域がの累乗であるFFTがあると仮定します x ある分野では、 x ^ 2 ^ k = 1 (例えば、上で紹介した場合、ドメインはの力です 85 モジュロ 33785 ^ 2 ^ 3 = 1)。 いくつかの多項式があります。 y=6×7+2×6+9×5+5×4+x3+4×2+x+3 (私たちはそれを次のように書きます p = [3,1,4,1,5,9,2,6])。 ドメイン内の各ポイントでこの多項式を評価したいと思います。 のXNUMXつの力のそれぞれで 85。 これが私たちの仕事です。

まず、多項式をXNUMXつの部分に分割します。これを次のように呼びます。 偶数 & オッズ:偶数=【3,4,5,2] & オッズ 【1,1,9,6] (または 偶数=2×3+5×2+4x+3 & オッズ=6×3+9×2+x+1; はい、これは偶数次の係数と奇数次の係数を使用しているだけです)。 ここで、数学的観察に注目します。 p(x)= evens(x2)+x⋅odds(x2) & p(−x)= evens(x2)−x⋅odds(x2) (これについて自分で考え、先に進む前に理解していることを確認してください)。

「接着剤」は比較的簡単です(そして オン) 実行時):偶数とオッズの評価を次のように受け取ります サイズ-N / 2 リストなので、単純に 各インデックスiのp [i] = evens_result [i] + domain [i]⋅odds_result[i]およびp [N2 + i] = evens_result [i] -domain [i]⋅odds_result[i]。

完全なコードは次のとおりです。

def fft(vals, modulus, domain): if len(vals) == 1: return vals L = fft(vals[::2], modulus, domain[::2]) R = fft(vals[1::2], modulus, domain[::2]) o = [0 for i in vals] for i, (x, y) in enumerate(zip(L, R)): y_times_root = y*domain[i] o[i] = (x+y_times_root) % modulus o[i+len(L)] = (x-y_times_root) % modulus return o

実行してみることができます:

>>> fft([3,1,4,1,5,9,2,6], 337, [1, 85, 148, 111, 336, 252, 189, 226])
[31, 70, 109, 74, 334, 181, 232, 4]

そして、結果を確認できます。 位置での多項式の評価 85たとえば、実際には結果が得られます 70。 これは、ドメインが「正しい」場合にのみ機能することに注意してください。 それは形である必要があります [range(n)内のiのx ^ i%モジュラス] コラボレー x ^ n = 1.

逆FFTは驚くほど単純です。

def inverse_fft(vals, modulus, domain): vals = fft(vals, modulus, domain) return [x * modular_inverse(len(vals), modulus) % modulus for x in [vals[0]] + vals[1:][::-1]]

基本的に、FFTを再度実行しますが、結果を逆にして(最初の項目が所定の位置にとどまる場合を除く)、すべての値をリストの長さで除算します。

>>> domain = [1, 85, 148, 111, 336, 252, 189, 226]
>>> def modular_inverse(x, n): return pow(x, n - 2, n)
>>> values = fft([3,1,4,1,5,9,2,6], 337, domain)
>>> values
[31, 70, 109, 74, 334, 181, 232, 4]
>>> inverse_fft(values, 337, domain)
[3, 1, 4, 1, 5, 9, 2, 6]

さて、これを何に使うことができますか? これがXNUMXつの楽しいユースケースです。FFTを使用して数値を非常にすばやく乗算できます。 掛けたいとしましょう 1253 by 1895。 これが私たちがすることです。 まず、問題を少し簡単な問題に変換します。 多項式 【3,5,2,1] by 【5,9,8,1] (これは、昇順でXNUMXつの数字の数字だけです)次に、XNUMX回のパスを実行して数十桁を超えることにより、回答を数字に変換し直します。

多項式をFFTで乗算すると、多項式を次のように変換するとすぐに乗算できます。 評価表 (すなわち。 f(x) すべてのための x 一部のドメインでは D)の場合、評価を乗算するだけでXNUMXつの多項式を乗算できます。 だから私たちがすることは私たちのXNUMXつの数を表す多項式を取ることです 係数形式、FFTを使用してそれらを評価フォームに変換し、点ごとに乗算して、元に戻します。

>>> p1 = [3,5,2,1,0,0,0,0]
>>> p2 = [5,9,8,1,0,0,0,0]
>>> x1 = fft(p1, 337, domain)
>>> x1
[11, 161, 256, 10, 336, 100, 83, 78]
>>> x2 = fft(p2, 337, domain)
>>> x2
[23, 43, 170, 242, 3, 313, 161, 96]
>>> x3 = [(v1 * v2) % 337 for v1, v2 in zip(x1, x2)]
>>> x3
[253, 183, 47, 61, 334, 296, 220, 74]
>>> inverse_fft(x3, 337, domain)
[15, 52, 79, 66, 30, 10, 1, 0]

これにはXNUMXつのFFTが必要です(それぞれ O(N・log(N)) 時間)とXNUMXつの点ごとの乗算(オン) 時間)、それでそれはかかります O(N・log(N))全体の時間(技術的には少し長い O(N・log(N))、非常に大きな数の場合は置き換える必要があるため 337 モジュラスが大きくなると、乗算が難しくなりますが、十分に近くなります)。 これは はるかに高速 教科書の掛け算よりも O(N2) 時間:

 3 5 2 1 ------------
5 | 15 25 10 5
9 | 27 45 18 9
8 | 24 40 16 8
1 | 3 5 2 1 --------------------- 15 52 79 66 30 10 1

したがって、結果を取得して、数十桁を引き継ぐだけです(これは、「リストをXNUMX回ウォークスルーし、各ポイントでXNUMXつのことを実行する」アルゴリズムであるため、次のようになります。 オン) 時間):

[15, 52, 79, 66, 30, 10, 1, 0]
[ 5, 53, 79, 66, 30, 10, 1, 0]
[ 5, 3, 84, 66, 30, 10, 1, 0]
[ 5, 3, 4, 74, 30, 10, 1, 0]
[ 5, 3, 4, 4, 37, 10, 1, 0]
[ 5, 3, 4, 4, 7, 13, 1, 0]
[ 5, 3, 4, 4, 7, 3, 2, 0]

そして、数字を上から下に読むと、 2374435。 答えを確認しましょう…。

>>> 1253 * 1895
2374435

わーい! 出来た。 実際には、そのような小さな入力では、 O(N・log(N)) & O(N ^ 2) ありません それ アルゴリズムが単純であるという理由だけで、教科書の乗算はこのFFTベースの乗算プロセスよりも高速ですが、入力が大きい場合は非常に大きな違いがあります。

ただし、FFTは数値の乗算だけでなく便利です。 前述のように、多項式の乗算と多点評価は、イレイジャーコーディングを実装する上で非常に重要な操作です。これは、多くの種類の冗長フォールトトレラントシステムを構築するための非常に重要な手法です。 フォールトトレランスが好きで効率が好きなら、FFTはあなたの友達です。

FFTとバイナリフィールド

そこにある有限体の種類は素数体だけではありません。 別の種類の有限体(実際には、より一般的な概念の特殊なケース 拡張フィールド、複素数に相当する有限体のようなもの)はXNUMX進体です。 バイナリフィールドでは、各要素はすべてのエントリが次の多項式として表されます。 0 or 1、例えば。 x ^ 3 + x + 1.

多項式の加法はモジュロで行われます 2、および減算は加算と同じです(として −1 = 1mod2)。 モジュラスとして既約多項式を選択します(例: x ^ 4 + x + 1; x ^ 4 + 1 なぜならうまくいかない x ^ 4 + 1 因数分解することができます (x2 + 1)⋅(x2 + 1) したがって、「既約」ではありません)。 乗算はそのモジュラスを法として行われます。 たとえば、バイナリフィールドmodで x ^ 4 + x + 1、乗算 x2 + 1 x 3 + 1 与えるだろう x5 + x3 + x2 + 1 掛け算をするだけなら x^5+x^3+x^2+1=(x^4+x+1)⋅x+(x^3+x+1)、結果は余りです x ^ 3 + x + 1.

この例は九九として表現できます。 最初に掛ける 【1,0,0,1] (すなわち。 x ^ 3 + 1)によって 【1,0,1] (すなわち。 x ^ 2 + 1):

 1 0 0 1 --------
1 | 1 0 0 1
0 | 0 0 0 0
1 | 1 0 0 1 ------------ 1 0 1 1 0 1

乗算結果には、 x ^ 5 減算できるように項 (x4 + x + 1)⋅x:

 1 0 1 1 0 1 - 1 1 0 0 1 [(x⁴ + x + 1) shifted right by one to reflect being multipled by x] ------------ 1 1 0 1 0 0 

そして、結果が得られます。 【1,1,0,1] (or x3 + x + 1).

画像

バイナリフィールドmodの加算および乗算テーブル x ^ 4 + x + 1。 フィールド要素は、バイナリから変換された整数として表されます(例: x ^ 3 + x ^ 2→1100→12)

バイナリフィールドは、XNUMXつの理由で興味深いものです。 まず、バイナリデータをイレイジャーコーディングする場合は、バイナリフィールドが非常に便利です。 N データのバイトは、バイナリフィールド要素として直接エンコードできます。また、データに対して計算を実行して生成するバイナリフィールド要素も次のようになります。 N バイト長。 素数フィールドのサイズは正確にXNUMXの累乗ではないため、素数フィールドでこれを行うことはできません。 たとえば、すべてをエンコードできます 2 からの数としてのバイト 0 ... 65536 素数体を法として 65537 (これは素数です)が、これらの値に対してFFTを実行すると、出力に次の情報が含まれる可能性があります。 65536、XNUMXバイトで表現することはできません。

第二に、足し算と引き算が同じ演算になるという事実、そして 1 + 1 = 0、いくつかの非常に興味深い結果につながるいくつかの「構造」を作成します。 バイナリフィールドの特に興味深く、有用な奇妙な点のXNUMXつは、「一年生の夢」定理: (x + y)^ 2 = x ^ 2 + y ^ 2 (そして指数についても同じです 4,8,16 ... 基本的にXNUMXの累乗)。

ただし、イレイジャーコーディングにバイナリフィールドを使用し、それを効率的に使用する場合は、バイナリフィールドに対して高速フーリエ変換を実行できる必要があります。 しかし、問題があります。バイナリフィールドでは、 (自明ではない)位数の乗法群はありません 2 ^ n。 これは、乗法群がすべて順序であるためです 2 ^ n-1。 たとえば、モジュラスのあるバイナリフィールドでは x ^ 4 + x + 1、の連続する累乗の計算を開始した場合 x + 1、サイクルバックします 1 After 15 ステップ–ではありません 16。 その理由は、フィールド内の要素の総数が 16、ただし、そのうちのXNUMXつはゼロであり、フィールドでゼロ以外の値をそれ自体で乗算してもゼロに達することはありません。 x + 1 ゼロ以外のすべての要素を循環するため、周期の長さは次のようになります。 15はなく、 16。 どうしようか?

以前に2n個の要素を持つ乗法群の「構造」を持つドメインが必要だった理由は、ドメイン内の各数値をXNUMX乗することにより、ドメインのサイズをXNUMX分のXNUMXに縮小する必要があったためです。 【1,85,148,111,336,252,189,226] に削減されます 【1,148,336,189] なぜなら 1 両方の二乗です 1 & 336、148 両方の二乗です 85 & 252、など。

しかし、バイナリフィールドで、ドメインのサイズを半分にする別の方法があるとしたらどうでしょうか。 あることが判明しました:以下を含むドメインが与えられた 2 ^ k ゼロを含む値(技術的には、ドメインは 部分空間)、次のようにすることで、半分のサイズの新しいドメインD 'を構築できます。 x⋅(x + k) for x in D 特定のを使用して k in D。 元のドメインは部分空間であるため、 k ドメイン内にある、任意 x ドメイン内に対応する x + k ドメイン内でも、そして機能 f(x)=x⋅(x + k) と同じ値を返します x & x + k したがって、XNUMX乗で得られるのと同じ種類のXNUMX対XNUMXの対応が得られます。

画像

では、これに加えてFFTをどのように行うのでしょうか。 同じトリックを使用して、問題を変換します Nサイズの多項式と NサイズのドメインをXNUMXつの問題に分割し、それぞれに N / 2サイズの多項式と N / 2サイズのドメインですが、今回は異なる方程式を使用します。 多項式pをXNUMXつの多項式に変換します 偶数 & オッズ そのような p(x)= evens(x⋅(k-x))+x⋅odds(x⋅(k-x))。 私たちが見つけた偶数とオッズについては、 また 真実であること p(x + k)= evens(x⋅(k-x))+(x + k)⋅odds(x⋅(k-x))。 したがって、縮小されたドメインで偶数とオッズに対して再帰的にFFTを実行できます。 [Dのxのx⋅(k−x)]次に、これらXNUMXつの式を使用して、ドメインのXNUMXつの「半分」の答えを取得します。一方は他方からkだけオフセットされています。

pをに変換する 偶数 オッズ 上記のように、それ自体が重要であることが判明しました。 これを行うための「ナイーブな」アルゴリズムはそれ自体です O(N ^ 2)、しかし、バイナリフィールドでは、次の事実を使用できることがわかりました。 (x2−kx)2=x4−k2⋅x2、より一般的に (x2−kx)2i=x2i+1−k2i⋅x2i 、これを行うためのさらに別の再帰的アルゴリズムを作成するには O(N・log(N)) 時間。

そして、あなたがしたい場合 逆の FFT、補間を行うには、アルゴリズムのステップを逆の順序で実行する必要があります。 これを行うための完全なコードはここにあります: https://github.com/ethereum/research/tree/master/binary_fft、およびより最適なアルゴリズムの詳細が記載された論文はこちら: http://www.math.clemson.edu/~sgao/papers/GM10.pdf

では、この複雑さのすべてから何が得られるのでしょうか。 さて、「ナイーブ」の両方を備えた実装を実行してみることができます O(N ^ 2) マルチポイント評価と最適化されたFFTベースの評価、および両方の時間。 これが私の結果です:

>>> import binary_fft as b
>>> import time, random
>>> f = b.BinaryField(1033)
>>> poly = [random.randrange(1024) for i in range(1024)]
>>> a = time.time(); x1 = b._simple_ft(f, poly); time.time() - a
0.5752472877502441
>>> a = time.time(); x2 = b.fft(f, poly, list(range(1024))); time.time() - a
0.03820443153381348

そして、多項式のサイズが大きくなるにつれて、素朴な実装(_シンプル_フィート)FFTよりもはるかに速く遅くなります:

>>> f = b.BinaryField(2053)
>>> poly = [random.randrange(2048) for i in range(2048)]
>>> a = time.time(); x1 = b._simple_ft(f, poly); time.time() - a
2.2243144512176514
>>> a = time.time(); x2 = b.fft(f, poly, list(range(2048))); time.time() - a
0.07896280288696289

そして出来上がり、多項式を多点評価して補間するための効率的でスケーラブルな方法があります。 FFTを使用して、現在の場所でイレイジャーコーディングされたデータを回復する場合 行方不明 いくつかの部分、そしてこのためのアルゴリズム 存在するただし、単一のFFTを実行するよりも効率がやや劣ります。 楽しみ!

この記事はもともと「高速フーリエ変換

タグ

ハッカー正午に参加

無料のアカウントを作成して、カスタムの読書体験のロックを解除します。

PlatoAi。 Web3の再考。 増幅されたデータインテリジェンス。
アクセスするには、ここをクリックしてください。

ソース:https://hackernoon.com/fast-fourier-transform-scaling-multi-point-evaluation?source = rss

スポット画像

最新のインテリジェンス

スポット画像

私たちとチャット

やあ! どんな御用でしょうか?