多項式近似
三角関数や対数のような計算の難しい関数を、掛け算と足し算だけでできる多項式で近似する方法と、割り算や平方根を繰り返しの計算で求める方法について説明します。
<注意>
このページでは、数式を表示するためにMathJaxを使用しています。ブラウザの設定でスクリプト実行をブロックしている場合、ページ内の数式が正しく表示されません。
なぜ多項式で近似するのか
CPUが直接できるのは、足し算・引き算・掛け算(SNESでは乗算器と除算器を使う)と、時間のかかる割り算だけです。
\(\ \sin x\ \)や\(\ \log x\ \)、\(\ \sqrt{x}\ \)を直接求める命令はありません。
多項式\(\ c_0 + c_1x + c_2x^2 + \cdots\ \)なら、掛け算と足し算だけで計算できます。
そこで、求めたい関数に近い多項式を探して、その値で代用します。
ゲームでは、表を引く方法(時間計算量と空間計算量)がよく使われますが、表の大きさが足りないときや、表の間の値が必要なときに、多項式の近似を使ったり組み合わせたりします。
テイラー展開
なめらかな関数は、ある点の近くで、多項式の無限の和として表せます。
定理 (テイラー展開)
\[ \sin x = x - \frac{x^3}{6} + \frac{x^5}{120} - \cdots, \qquad \cos x = 1 - \frac{x^2}{2} + \frac{x^4}{24} - \cdots \]\[ e^x = 1 + x + \frac{x^2}{2} + \frac{x^3}{6} + \cdots, \qquad \log(1+x) = x - \frac{x^2}{2} + \frac{x^3}{3} - \cdots \]
\(\ x = 0\ \)の近くで、次の式が成り立つ。
途中までの項で打ち切れば、近似の多項式になります。
\(\ x = 0\ \)の近くではよく合いますが、0から離れるほど外れていきます。
\(\ 0 \le x \le \pi/2\ \)での最大の誤差は、次のとおりです。
| 打ち切った式 | 最大の誤差 |
|---|---|
| \(\ x\ \) | 0.571 |
| \(\ x - x^3/6\ \) | 0.075 |
| \(\ x - x^3/6 + x^5/120\ \) | 0.0045 |
| 7次の項まで | 0.00016 |
範囲を狭める
近似の式は、使う範囲が狭いほど、少ない項で精度が出ます。
そこで、関数の性質を使って、計算する範囲を狭めてから近似します。
- 三角関数: 対称性(\(\ \sin(\pi - x) = \sin x\ \)など、三角関数)を使えば、\(\ 0 \le x \le \pi/2\ \)の範囲だけ計算すれば済みます。256方位なら0～64です。
- 対数: \(\ x = 2^k \times m\ \)(\(\ 1 \le m < 2\ \))と分けると、\(\ \log_2 x = k + \log_2 m\ \)です。\(\ k\ \)は最上位のビットの位置で、\(\ m\ \)はそれより下のビットなので、近似が必要なのは\(\ 1 \le m < 2\ \)の範囲の\(\ \log_2 m\ \)だけです(対数とデシベル)。
たとえば、\(\ \log_2 m \approx m - 1\ \)とするだけでも、最大の誤差は約0.086で、ビットの位置\(\ k\ \)と合わせて大まかな対数が求まります。
2次式\(\ 1.357(m-1) - 0.367(m-1)^2\ \)にすると、誤差は約0.0097になります。
誤差をならす
テイラー展開は1点の近くで最もよく合うように作られた式なので、範囲の端で誤差が大きくなります。
範囲全体で誤差が小さくなるように係数を決め直すと、同じ次数でずっと精度が上がります。
\(\ 0 \le x \le \pi/2\ \)で、\(\ \sin x\ \)を\(\ ax + bx^3\ \)の形で近似した場合の誤差です。
| 式 | 最大の誤差 |
|---|---|
| \(\ x - 0.1667x^3\ \)(テイラー展開) | 0.075 |
| \(\ 0.9888x - 0.1451x^3\ \)(誤差の2乗の合計が最小になるように決めた) | 0.0090 |
係数を少し変えただけで、最大の誤差が約8分の1になります。
誤差の最大値そのものが最も小さくなるように係数を決めたもの(ミニマックス近似)は、さらに誤差が均等にならされます。
係数はパソコンで一度計算しておけばよいので、計算の手間はSNESの側にはかかりません。
多項式の計算の順番
多項式は、次のように\(\ x\ \)をくくり出して計算します(ホーナー法)。
\[ c_0 + c_1x + c_2x^2 + c_3x^3 = c_0 + x\,(c_1 + x\,(c_2 + x\,c_3)) \]\(\ n\ \)次式が、掛け算\(\ n\ \)回と足し算\(\ n\ \)回で計算できます。
\(\ x^2,\ x^3\ \)を別々に求めるより掛け算が少なく、固定小数点数で計算したときの丸めの誤差も小さくなります。
固定小数点数で計算するときは、各係数と途中の値が、選んだ形式の範囲に収まるかを確かめておきます(固定小数点数)。
表と組み合わせる
範囲をいくつかの区間に分け、区間ごとに別の短い多項式を使う方法もあります。
区間ごとの係数を表に持てば、表の大きさと計算の量を調整できます。
1次式(直線)を区間ごとに使うのが、補間法の線形補間です。
繰り返しで求める(ニュートン法)
平方根や逆数のように、多項式で近似しにくい関数は、答えに近づいていく計算を繰り返して求めます。
定理 (ニュートン法)
\[ x \leftarrow x - \frac{g(x)}{g'(x)} \]
\(\ g(x) = 0\ \)の解を求めるには、\(\ x\ \)をで繰り返し更新する(\(\ g'\ \)は\(\ g\ \)の微分)。答えに近い値から始めれば、1回ごとに正しい桁数がおよそ2倍になる。
平方根
\(\ g(x) = x^2 - a\ \)とすると、\(\ \sqrt{a}\ \)を求める次の式になります。
\[ x \leftarrow \frac{1}{2}\left(x + \frac{a}{x}\right) \]\(\ \sqrt{2} = 1.4142136\cdots\ \)を\(\ x = 1\ \)から求めると、1.5、1.41667、1.414216と進み、3回目で上から6桁が合います。
割り算が必要ですが、最初の値を表から引いておけば、1～2回で十分な精度になります。
ベクトルの長さを正確に求めるときに使えます。
割り算をしない逆数
\(\ g(x) = 1/x - a\ \)とすると、\(\ 1/a\ \)を求める次の式になります。
\[ x \leftarrow x\,(2 - a\,x) \]この式には割り算がなく、掛け算2回と引き算1回だけです。
\(\ 1/3\ \)を\(\ x = 0.3\ \)から求めると、0.33、0.3333、0.33333333と、正しい桁数が倍々に増えます。
最初の値は、\(\ a\ \)の上位のビットで小さな表を引いて決めます。
逆数が求まれば、割り算\(\ b/a\ \)は掛け算\(\ b \times (1/a)\ \)になります。