Re: 浮動小数点数の二分探索

rsk0315.hatenablog.com

こういう記事を過去に書いたんですが、どうにもなんだかな〜という気がします。

導入

というのも、64 bit の浮動小数点数を使うとして 90 回もチェックするのはうれしくなさそうです。 （IEEE 754 準拠というのを仮定するとして）double のビット表現は、整数として見ても順序がある程度保たれます。

であれば、整数として見てにぶたんしちゃいたいな〜となります。そうすれば 64 回程度で済むので。

そこで、double としての値とビット表現を並べてみます。

値	ビット表現
-NaN	1111111111111000000000000000000000000000000000000000000000000000
-inf	1111111111110000000000000000000000000000000000000000000000000000
-1.7976931348623157e308	1111111111101111111111111111111111111111111111111111111111111111
-1.0	1011111111110000000000000000000000000000000000000000000000000000
-2.2250738585072014e-308	1000000000010000000000000000000000000000000000000000000000000000
-5e-324	1000000000000000000000000000000000000000000000000000000000000001
-0.0	1000000000000000000000000000000000000000000000000000000000000000
0.0	0000000000000000000000000000000000000000000000000000000000000000
5e-324	0000000000000000000000000000000000000000000000000000000000000001
2.2250738585072014e-308	0000000000010000000000000000000000000000000000000000000000000000
1.0	0011111111110000000000000000000000000000000000000000000000000000
1.7976931348623157e308	0111111111101111111111111111111111111111111111111111111111111111
inf	0111111111110000000000000000000000000000000000000000000000000000
NaN	0111111111111000000000000000000000000000000000000000000000000000

値が負 (-0.0 含む) であるところに関しては、ビット表現を整数として読んだ順と逆、正 (+0.0 含む) に関しては同じになっていることがわかります。 最上位ビットをよしなにすることも考え、次のような変換を考えます。Rust です。

fn f2u(f: f64) -> u64 {
    let u = f.to_bits();
    if u >> 63 == 1 { !u } else { u | 1 << 63 }
}

値は次のようになります。

f	f2u(f) のビット表現
-NaN	0000000000000111111111111111111111111111111111111111111111111111
-inf	0000000000001111111111111111111111111111111111111111111111111111
-1.7976931348623157e308	0000000000010000000000000000000000000000000000000000000000000000
-1.0	0100000000001111111111111111111111111111111111111111111111111111
-2.2250738585072014e-308	0111111111101111111111111111111111111111111111111111111111111111
-5e-324	0111111111111111111111111111111111111111111111111111111111111110
-0.0	0111111111111111111111111111111111111111111111111111111111111111
0.0	1000000000000000000000000000000000000000000000000000000000000000
5e-324	1000000000000000000000000000000000000000000000000000000000000001
2.2250738585072014e-308	1000000000010000000000000000000000000000000000000000000000000000
1.0	1011111111110000000000000000000000000000000000000000000000000000
1.7976931348623157e308	1111111111101111111111111111111111111111111111111111111111111111
inf	1111111111110000000000000000000000000000000000000000000000000000
NaN	1111111111111000000000000000000000000000000000000000000000000000

単調増加になっていることがわかります。それから、NaN が inf より外側にいるのも都合がいいです。

というわけで、これに関してにぶたんすればよいです。

実装例

ABC 026 D を使います。

f2u の逆関数 u2f を用意しておきます。

use std::f64::consts::PI;

use proconio::input;

fn main() {
    input! {
        a: f64,
        b: f64,
        c: f64,
    }

    let f = |x: f64| a * x + b * (c * x * PI).sin();

    let res_u = {
        let mut lb = f2u(0.0);
        let mut ub = f2u(1.0e9);
        while ub - lb > 1 {
            let mid_u = lb + (ub - lb) / 2;
            let mid = u2f(mid_u);
            *(if f(mid) < 100.0 { &mut lb } else { &mut ub }) = mid_u;
        }
        ub
    };
    let res = u2f(res_u);

    println!("{}", res);
}

fn f2u(f: f64) -> u64 {
    let u = f.to_bits();
    if u >> 63 == 1 { !u } else { u ^ 1 << 63 }
}

fn u2f(u: u64) -> f64 {
    f64::from_bits(if u >> 63 == 1 { u ^ 1 << 63 } else { !u })
}

atcoder.jp

最上位ビットを持ってきてがちゃがちゃしたいときは、算術シフトを使えば if がいらないんですね。

fn mask(u: u64) -> u64 { ((u as i64 >> 63) as u64 >> 1) | 1 << 63 }

fn f2u(f: f64) -> u64 {
    let u = f.to_bits();
    u ^ mask(u)
}

fn u2f(u: u64) -> f64 { f64::from_bits(u ^ mask(!u)) }

atcoder.jp

人はなぜ

C++ での実装例をよこせと言われそうな気がしますが、C++ は標準で直接 double のビット表現を得る関数がないんですね。bit_cast というのを C++20 まで待つ必要があるらしいです。それから雑にやろうとすると未定義動作にもなります。嫌だねえ。

というので、ちょっと面倒ですが、以下を参考にしつつ書きます。

en.cppreference.com

提出コードのリンクだけ載せておきます。

atcoder.jp

おきもち

定数回ループの方が直感的に思う人も多そうだし、好きな方を使えばいいんじゃないかなとも思います。

それはそれとして、\( (-\infty, +\infty)\) の探索をしたいとき、整数に変換してやれば上から bit を決めていけばよくなって楽かもみたいな気持ちもあります。NaN には注意する必要があるかも。

おわり

ねこにゃんこ。