IEEE754な浮動小数点の任意区間でのカウントと列挙

Question

今、
配列を重複なく乱数で埋めるには - Qiita
http://qiita.com/yumetodo/items/48d77f5d554df84f66f7
C++で効率よく重複のない乱数列を生成する = Qiita
http://qiita.com/hmito/items/9f4bdc8442b6f6b3c7bc
の浮動小数点版を作ろうと思っています。

ここで、make_rand_array_select(リンク先参照)のアルゴリズムで生成しようと思っています

template<typename type> 
std::vector<type> make_rand_array_select(const size_t size, type rand_min, type rand_max) ;

1. 生成個数(第1引数)と生成範囲(第2, 3引数)を見た時重複なく生成可能かを調べ、不可能な場合は例外を投げたい
2. 生成範囲(rand_min, rand_max自身も含む)のすべての数が重複なく列挙されたvectorを作りたい

と思っています。

Twitterで質問を投げた際には

• std::nextafterで全部辿ってしまう(@ytomino)
• ゼロ近傍を跨いでいるとすごくいろいろ面倒ですが（あと端のほうも無限大とか変な奴がいる）、それを別にすればIEEE 754方式の浮動小数点型の中身は、基本的に整数として比較できるようにできてるのでそのようにして(@ksmakoto)
• https://gist.github.com/mskashi/bd59f51d934d3c7ed42cこんな感じでいけると思います。全てのdoubleに連番をふる感じ。ただし、+0と-0は同一視しています。列挙はそれこそnextafterですかね。(@mkashi)

などという意見を頂きましたが、正直Twitterはこういうのを聞くのに向かないな、と思いまとめました。(なにより私が混乱してきた)

環境等

• 浮動小数点型はfloat, double, (long double)のことを指す
• 浮動小数点型の内部仕様はIEEE754に従う
• (size_t)((int)(rand_max) - (int)(rand_min)) < size
• sizeof(long double) == sizeof(uint64_t)

追記

浮動小数点表現での一定範囲内の全列挙…だと?! | okkyの日記 | スラド
にて、全列挙のコストの高さを突きつけられたので別のアルゴリズムを模索中です。

https://twitter.com/fjs_kyousosama/status/704681468535009281

@yumetodo ヒント1: 64bitのビット列を完全に出力させて、それをdouble形式として読み込み、分布を調べると 0.0近辺に集中していることが判る。「それはランダムな浮動小数点分布」として欲しい分布なのか?

自力で検証していないけど0.0付近に集中したらアカンなぁ・・・
完全に整数版で作ったbit列を浮動小数点に変換する方向に考えていたので。

区間を生成個数で分割して各区間ごとに乱数生成してvectorに放り込んでshuffleとかしなきゃいけないんだろうか

Answer 1

こんにちは。

難しそうなテーマですね。
リンク先はざっとしか見れなかったので、主に質問の文面から回答してみます。
（頭で考えただけで実証してないので間違ってたらごめんなさい。）

1.生成個数(第1引数)と生成範囲(第2, 3引数)を見た時重複なく生成可能かを調べ、不可能な場合は例外を投げたい

rand_min, rand_maxの間に存在するIEEE754フォーマットの浮動小数点の数を計算できれば良い筈ですね。指数部と仮数部に分けて考えるのが良さそうに思います。

まず、指数部を固定した時に表現できる浮動小数点値の個数は仮数部が取りうる値の数と同じですね。
この原理にもとづいて下記３つの数を計算できると思います。
①rand_minの指数値についてrand_minの仮数部の値以上となる仮数部の個数は計算できます。
②rand_min～rand_maxの間(それぞれを含まない)指数部の値の個数も計算できます。
③rand_maxの仮数部の値の個数について①と同じ原理で仮数部の個数を計算できます。
この３つを足して、第一パラメータと比較すれば判定できるのではないでしょうか？

2.生成範囲(rand_min, rand_max自身も含む)のすべての数が重複なく列挙されたvectorを作りたい

これは第一パラメータは無視する前提ですよね？
であれば、単純に指数部と仮数部を連結したものを整数型とみなしてインクリメントしていき、rand_maxと同じビットパターンになった時点で停止すれば良さそうに思います。
負の時は絶対値を取って処理する等、工夫が必要と思いますが、同じ原理で計算できる筈と思います。

もし、よさそうだったら、検証はお任せします。
って、どうやれば検証できるんでしょう。更に難しそうです。

---

【追記】
質問の追記へのコメントです。

コストの高さの件：
第一パラメータで指定した個数の要素を持つ配列を確保するのですよね？
当然、第一パラメータで許容する値は、現実的な範囲に制限するのでは？

分布の偏りの件：
浮動小数点数は0.0に近い程密に数値を表現できますね。
ということは、浮動小数点で表現できる「全て」の数値を生成したら、0.0に近い方が密度は高くなります。

なんだか当たり前のことを悩まれているようにも見えます。
課題の2の「全ての数」の部分が間違いですね？
「生成範囲(rand_min, rand_max自身も含む)の数が重複なく均等に分布するvectorを作りたい」でしょうか？

Comments

[yumetodo]

>ということは、浮動小数点で表現できる「全て」の数値を生成したら、0.0に近い方が密度は高くなります。

うーん、そうだよな、そう言われればそうだよな、指数だし・・・。
純粋にuniform_real_distribusionで一様生成した時もやはり0にかたよるんだろうか。

[Chironian]

http://cpprefjp.github.io/reference/random/uniform_real_distribution.html
に簡単なヒストグラムがありました。これを見る限り普通に「一様」ですね。

-1.0～1.0の間にIEEE754浮動小数点で表現できる値の丁度半分の個数の値が入りますね。
そして、binary64なら、-1x2の+1023乗～-1.0、1.0～1x2の+1023乗の範囲に残りの半分が入ります。
分布の偏り、天文学的と思いません？

[yumetodo]

さすが指数、おそるべし。
uniform_real_distributionは一体どうやって「一様」を得ているっていうんだ。

[takotakot]

uniform_real_distribusion は、rand()/RAND_MAX *(max - min) + min と同じような働きをしているだけ、得られる分布が確率的に一様というだけで、各ビットについてどうということはないのではないでしょうか。
各ビットについてのある程度のランダム性も期待するのであれば、2分探索っぽくなにかやるのかなぁ…

[yumetodo]

やはりそういうことなのかな・・・。
すると0.0~1.0の間で重複なく一様な4個の乱数を得るとしたら区間を[0.0-0.25],[0.25-0.5],[0.5-0.75],[0.75-1.0]の4つに分けてその中で一回生成、ということしかないのかな

[YsMana]

その場合、[0.0-0.25]の出現回数が0回とか2回以上という組み合わせが出現しなくなりますが、それは構わないと言うことでしょうか。
重複の定義を「0.25で割った時の整数部が等しい」と考えればこれで良いとも言えそうですが。

[Chironian]

binary64なら1+2048+52=2101ビット整数の一様な乱数を生成してbinary64へ写像すれば偏りは無くせると思います。
2101ビット整数≒263バイト整数とか、ぞっとしますし、ほとんど全てのビットが無駄になるむちゃくちゃ勿体無い使い方です。遥かに効率の良い方法があるはずです。

ああ、大きな方の桁から「順次」生成していけば一様に作れそうな気がします。
binary64なら、最初に53ビットの一様乱数を生成し、MSB(第52ビット)が1になるまで左シフトして、左シフトした分の乱数を再度生成して埋める感じです。そうしてできた53ビットの内の下位52ビットを仮数部とします。
指数部は最初1023で左へ1ビットシフトする度に1引きます。

最初に生成した乱数やそれに続いて生成した乱数が0だった場合は順次52ビット左シフトします。
符号ビットは別途生成するのが簡単かも。

[yumetodo]

うーん、わたしの理解力を超えてしまっていて・・・。
ちょっと手元でなんかコードを書いて考えてみます

Answer 2

ずいぶん放置してましたが、もともと浮動小数に対してやるような話じゃなかったという方向で質問を閉じます。