Open3

"負のゼロ" のテストを書く

okuokuokuoku

言語処理系で実装される数には"負のゼロ"が存在することがある。

https://en.wikipedia.org/wiki/Signed_zero

この負のゼロは:

  • IEEE754で言うところの copySign https://cpprefjp.github.io/reference/cmath/copysign.html
  • -0.0 のようなソースコード上の記述を読みとったとき (Schemeの場合)
  • 演算の丸めの結果として生成される (-1 / DBL_MAXとか ※ これだとsubnormalのぶん足りない。後述。)

のようなケースで生成され、"正のゼロ" とは、

  • = の意味では一致する(Schemeで言うところの (= 0.0 -0.0) は真)
  • 区別を提供する処理系もある

SchemeではR6RS以降規定が割とクリアになっていて、処理系が負のゼロを扱えるなら (eqv? 0.0 -0.0) が偽、 (= 0.0 -0.0) が真になる。...が、処理系の仕様が常にバッチリ書かれているわけではないので、移植層を書く立場からするとテストが必要になる。

https://twitter.com/anohana/status/1430117898010337286

この辺の話はGauche devlogによくサマリされている。

http://blog.practical-scheme.net/gauche/20210826-negative-zero

イコールとは何ぞや

あと、これを拡張すると = の定義の問題になる( = が成立しなくなる最小のデルタの定義)。これはR6RSの事実上のリファレンス実装であるChezSchemeと、R7RSのリファレンス実装であるChibi-schemeの両方で超長いディスカッションになっている。

https://github.com/cisco/ChezScheme/issues/606

https://github.com/ashinn/chibi-scheme/issues/812

これらの処理系ではIEEE754の範囲を越えた正確数(exact numbers)の概念があり、ジェネリックなオペレータはそれらを混ぜて使えるため正確な定義が必要になる。

okuokuokuoku

Cで書く

Cで書いてみるとこうなる。つまり、文字列からの読み取りでも負のゼロが生成され、 -1 / DBL_MAX では微妙に負のゼロには足りないということになる。

#include <stdlib.h>
#include <stdio.h>
#include <float.h>
#include <math.h>

int
main(int ac, char** av){
    double z,p;
    z = 0;
    p = atof("-0.0");
    printf("%e\n", p);
    p = -1.0 / (DBL_MAX);
    printf("%e\n", p);
    p = copysign(0,-1);
    printf("%e\n", p);
    if(z == p){
        printf("0.0 == -0.0\n");
    }else{
        printf("0.0 != -0.0\n");
    }
    return 0;
}
-0.000000e+00
-5.562685e-309
-0.000000e+00
0.0 == -0.0

演算で負のゼロに至るか試すには割り続けるのが早そうだ。

#include <stdint.h>
#include <stdio.h>

int
main(int ac, char** av){
    union {
        double d;
        uint64_t i;
    }v;
    uint64_t d;
    d = 0;
    v.d = -1;

    do{
        v.d /= 2;
        d++;
        printf("d >> %d = %e (%lx)\n",
               d, v.d, v.i);
    } while (0.0 != v.d);

    return 0;
}

このプログラムは一見停止しなさそうに見えるが(-1をどんなに割ってもゼロにはならないから)、実際には 0.0 == -0.0 のルールにより、値が負または正のゼロに到達したときに停止することになる。

d >> 1 = -5.000000e-01 (bfe0000000000000)
d >> 2 = -2.500000e-01 (bfd0000000000000)
d >> 3 = -1.250000e-01 (bfc0000000000000)
d >> 4 = -6.250000e-02 (bfb0000000000000)
d >> 5 = -3.125000e-02 (bfa0000000000000)
d >> 6 = -1.562500e-02 (bf90000000000000)
/* snip */
d >> 1070 = -7.905050e-323 (8000000000000010)
d >> 1071 = -3.952525e-323 (8000000000000008)
d >> 1072 = -1.976263e-323 (8000000000000004)
d >> 1073 = -9.881313e-324 (8000000000000002)
d >> 1074 = -4.940656e-324 (8000000000000001)
d >> 1075 = -0.000000e+00 (8000000000000000)
okuokuokuoku

書いた

https://github.com/okuoku/yuni/issues/163

しかし、

  • そもそもソースコード上の -0.0 が読めるか
  • printしたときに -0.0 が出せるか
  • 演算の中で -0.00.0 を区別しているか

が全部別々の事項なので何とも纏めるのが難しいもんだな。。いやまぁIEEE754への準拠を強制すれば簡単ではあるけど。。