ぎるばーとの日記

もっともっと遠くへ行きたい 空が広く見える場所まで

log2(x)

 JavaのMathクラスにはlog2(x)関数が見当たりません。
 他のプログラミング言語でもlog2(x)がない場合には役立つ話かも……?

 さて、下のやり方で何か問題があるでしょうか?

Math.log(x) / Math.log(2)

 対数の底の変換は、
\log_a x = \frac{\log x}{\log a} = (\log_a e) \log x
なので数学的には自然というか素朴な方法ですね。

 しかし、浮動小数点演算では予想外の結果になってしまうことがあります。

    System.out.println(Math.log(1 << 29) / Math.log(2));
29.000000000000004

 xがぴったり2のn乗なら、正確にnを返してほしいと思いませんか?(思わなければ続きは読まなくて大丈夫です。)
 かといって2のn乗を特別扱いすると、単調性(a \le bならばf(a) \le f(b))が保たれるか分からないので、それはそれでうーんな感じです。
 ということで、なるべく誤差を発生させない計算方法を考えます。

 まず、
x = s \times 2^k
と表します。(kは整数で、sの範囲を制限すると決まる。)
\log_2 x = \log_2 s + k
 ここで、\log_2 sの絶対値が小さいとき、浮動小数点数の性質(0に近いほど細かく値の区別ができる)から誤差も小さく計算できます。

 ここではsの範囲を、
\cases{0.5 \le s < 1 & \text{($x < 1$)}\\ 1 \le s < 2 & \text{($x \ge 1$)}}
とします。(絶対にこうすべきというわけでもない。)
 あとは非正規化数に注意しつつコード化すると次のようになります。

public static double log2(double x) {
    // 特殊な値
    if (Double.isNaN(x) || x < 0.0) return Double.NaN;
    if (x == Double.POSITIVE_INFINITY) return Double.POSITIVE_INFINITY;
    if (x == 0.0) return Double.NEGATIVE_INFINITY;
    // ここから
    int k = Math.getExponent(x);
    if (k < Double.MIN_EXPONENT) {
        // 非正規化数は取扱い注意!
        k = Math.getExponent(x * 0x1.0p52) - 52;
    }
    if (k < 0) {
        k++;
    }
    double s = Math.scalb(x, -k);
    final double LOG2_E = 1.4426950408889634;
    return k + LOG2_E * Math.log(s);
}

 比較テスト。(x = 2^n;\ n = -1074, \ldots, 1023

    int j = 0, k = 0;
    for (int i = -1074; i <= 1023; i++) {
        double x = Math.scalb(1.0, i);
        if (log2(x) != i) j++;
        if (Math.log(x) / Math.log(2) != i) k++;
    }
    System.out.println("log2(x) function errors: " + j);
    System.out.println("Math.log(x) / Math.log(2) errors: " + k);
log2(x) function errors: 0
Math.log(x) / Math.log(2) errors: 441

 ライブラリでも作るのでなければこれで十分かと……。