円周率.jp (http://xn--w6q13e505b.jp/method/fft/dwt.html)

離散荷重変換

DFT を行う前後で、データ列と同じ長さの適当な係数列 $\{a_i\}$ を付ける変換を考えてみよう。

\[ X_k = \sum_{j=0}^{n-1}a_jx_j\omega^{jk}, \quad x_j = a_j^{-1} n^{-1}\sum_{k=0}^{n-1}X_k\omega^{-jk} \]

この変換を畳み込み乗算のルーチンに組み込むと

\[ \begin{eqnarray} z_l &=& a_l^{-1}n^{-1} \sum_{k=0}^{n-1} \left(\sum_{j=0}^{n-1} a_j x_j \omega^{jk} \right) \left(\sum_{j'=0}^{n-1} a_{j'} y_{j'} \omega^{j'k} \right) \omega^{-kl}\\ &=& \sum_{(j+j') \bmod n = l} (a_l^{-1} a_j a_{j'}) (x_j y_{j'}) \end{eqnarray} \]

という形になり、$\{a_j\}$ の定義次第で計算のバリエーションを増やすことができる。 この計算方法を離散荷重変換 (Discrete Weighted Transform) という。 以下に使い勝手の良い例をいくつか挙げる。

負巡回畳み込み乗算

$\{a_j\} = \{\omega^{j/2}\}$ とすると

\[ \begin{eqnarray} z_l = \sum a_{j+j'-l} (x_j y_{j'}) &=& \sum_{j+j'=l} a_0 x_j y_{j'} + \omega^{n/2} \sum_{j+j'=l+n} x_j y_{j'}\\ &=& \sum_{j+j'=l} x_j y_{j'} - \sum_{j+j'=l+n} x_j y_{j'} \end{eqnarray} \]

となり、負巡回畳み込み乗算形になる。

複素巡回畳み込み乗算

計算が複素数体でなされている時、$\{a_j\} = \{\omega^{j/4}\}$ とすると $\omega^{n/4} = i$ となるので

\[ \begin{eqnarray} z_l = \sum a_{j+j'-l} (x_j y_{j'}) &=& \sum_{j+j'=l} x_j y_{j'} + \omega^{n/4} \sum_{j+j'=l+n} x_j y_{j'}\\ &=& \sum_{j+j'=l} x_j y_{j'} + i \sum_{j+j'=l+n} x_j y_{j'} \end{eqnarray} \]

となり、巡回部分が虚数部へ巡回する畳み込み乗算になる。 この DWT は実数 DFT のバリエーションの 1 つでもある。

無理数基底離散荷重変換 (IBDWT)

多倍長整数での剰余計算をするとき、通常なら基底 $R$ を定めて

\[ X = \sum_{j=0}^{n-1} x_jR^j \]

というように数列 $\{a_j\}$ に分割定義して計算する。($0\leq x_j \lt R$)

このとき $R$ は $2^{64}$ や $10^{19}$ の様に解釈 しやすい定数を用いるのが普通であるが、 これを少し変えて $R$ にバリエーションを持たせてみる。 具体的には $r$ を 2 や 10 の小さな単位での基数として、

\[ 0\leq x_j \lt R_j = r^{d_j} \]

のように各要素に詰める桁数をバラけさせる。 さらに

\[ D = \sum_{j=0}^{n-1} d_j \]

とすれば多倍長整数 $X$ は $r$ 進数で $D$ 桁の数と解釈できる。 もしこの計算空間で巡回畳み込み乗算ができれば $\pmod{r^D-1}$ の計算が容易になるが、それを実現する方法がここで紹介する 無理数基底離散荷重変換 (IBDWT ; Irrational Base DWT) である。 ちなみに長々と一般的な話で紹介しているが、 現在、メルセンヌ素数($r=2$、$D\in P$)の確認以外に使いみちは無い。

計算方法

$\{d_j\}$、$\{a_j\}$ を以下のように定義する。

\[ d_j = \left \lceil \frac{D(j+1)}{n} \right \rceil - \left \lceil \frac{Dj}{n} \right \rceil \] \[ a_j = r^{\left\lceil \frac{Dj}{n} \right\rceil - \frac{Dj}{n}} \Leftrightarrow \log_r a_j = \left\lceil \frac{Dj}{n} \right\rceil - \frac{Dj}{n} \]

わかりにくいが、 $a_j$ は $\dfrac{Dj}{n}$ の小数部分 $f_j$ を使って $r$ の $(1-f_j)$ 乗を取った数である。厳密には $\dfrac{Dj}{n}$ が整数になる場合は $a_j=r^0=1$ である。 この定義を使って DWT を行うと

\[ z_l = \sum a_{j}a_{j'}a_l^{-1} \cdot x_j y_{j'} = \sum r^{E(j,j',l)} (x_j y_{j'}) \] \[ \left( E(j,j',l) = \left\lceil \dfrac{Dj}{n} \right\rceil + \left\lceil \dfrac{Dj'}{n} \right\rceil - \left\lceil \dfrac{Dl}{n} \right\rceil - \dfrac{D(j+j'-l)}{n} \right) \]

となる。 ちなみに $j+j'-l=\{0,n\}$ ということを考えると $E(j,j',l)=\{0,1\}$ となり、計算結果が整数になることが保証される。