如何用C語言生成(0,1)之間的隨機浮點數？

02-02

不包括0和1，樓主剛自學編程不久，論文的演算法需要用到隨機浮點數，我試了好多辦法，生成的數還是包括0和1。請大家給點建議，謝謝。

Wikipedia: Rejection sampling

double r; do { r = someRNG(); // output [0, 1) } while (r == 0.0);

-----

8月30日補充

因為看到@北極的答案有問題，我在這裡再講講一般是如何生成隨機浮點數的。

我們知道一般浮點數是用n-bit尾數（mantissa，前面有hidden bit為1），數學表示為：

$1.m_1m_2...m_n imes2^e$

我們可以使用整數PRNG生成尾數，並把e設為0，那麼可以得到一個[1, 2)的浮點數。例如給定64位的整數PRNG，可以這樣生成[0, 1)的double：

double RandomDouble() { union { double d; uint64_t u; }x; x.u = (RandomUint64() &>&> 12) | 0x3FF0000000000000ULL; return x.d - 1.0; } // @北極使用的bit-field是有endianness問題，這種實現是endianness友好的。

實際上，e=0表達的是一個定點數（fixed point number），所以才能確保整數隨機數的均勻性質能保存下來。所以不應該把e也設置隨機數。

參考

[1] Saito, Mutsuo, and Makoto Matsumoto. "A PRNG specialized in double precision floating point numbers using an affine transition." Monte Carlo and Quasi-Monte Carlo Methods 2008. Springer Berlin Heidelberg, 2009. 589-602. http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/ARTICLES/dSFMT.pdf

-----

Rejection sampling的補充

Rejection sampling是有缺點的，例如不能確定生成的時間、SIMD計算中會出現分支等問題。

但在這個問題中，rejection sampling應該是可以接受的。對於

萬一連續N次都重試呢？萬一不是用lib C的標準庫做隨機而碰巧連續多次都重試怎麼辦？不能認為可能性低就完全忽視這種可能性吧。

這個疑問，我們可以用概率去了解。對於 52-bit 的尾數，出現一次0、連續兩次0、n次0的概率是多少呢？假設每個隨機數採樣是獨立事件，而採樣是均勻分布的，那麼

出現一次0（第二次不是0）的概率為：

$P( ext{0 occurs once}) = frac{1}{2^{52}} cdot (1 - frac{1}{2^{52}}) approx 2.22044605 imes10^{-16}$

連續出現兩0（第三次不是0）的概率為：

$P( ext{0 repeated twice})=left (frac{1}{2^{52}} ight )^2left(1 - frac{1}{2^{52}} ight) approx 4.93038066 imes 10^{-32}$

那麼n次也是這樣，可以看出這是一個等比數列。

然後我們看看實際的情況。假設我們只需要1納秒來生成1個[0, 1)的double隨機數，一天也只能產生 $10^9 imes 60 imes 60 imes 24 = 8.64 imes 10^{13}$ 個隨機數。

如果沒算錯的話，大約36天才有一半概率會出現一個0，大約10^14年才會有一半概率出現連續兩個0。宇宙的年齡才10^10年的數量級。

----

Trick

如果不介意對精度有極輕微的影響，只需把上述代碼修改一個字元，把最低位設置為1：

x.u = (RandomUint64() &>&> 12) | 0x3FF0000000000001ULL;

那麼就不用rejection sampling也可保證隨機數不為0，而又非常接近均勻分布。

贊同@vczh的說法，不贊同目前回答里的所有遇到0重試的方法，萬一連續N次都重試呢？萬一不是用lib C的標準庫做隨機而碰巧連續多次都重試怎麼辦？不能認為可能性低就完全忽視這種可能性吧。

@vczh在實踐中可能實現有點偏差，因為精度控制很麻煩，如果要完全符合精度，我倒是有個十分麻煩但肯定精準的方法，當然僅在x86平台，32位，C語言標準環境下可用：

根據浮點數的定義在內存中分別是：

符號位+指數位+有效數位

那麼根據這個規則，以單精度浮點數為例：

單精度浮點數在0-1之間的數字的範圍是：

符號位永遠為0

當有效數在1 ~ Max-1（Max=0x7FFFFF）之間的時候：指數取值是0到126

當有效數是0的時候：指數取值是1到126

當有效數是Max的時候：指數的取值是0到125

所以有代碼如下：

#include &


union FloatRand

{

    struct

    {

        unsigned long Frac:23;

        unsigned long Exp:8;

        unsigned long Signed:1;

    } BitArea;

    float Value;

    unsigned long Binary; /* for debug only */

};
float GetFloatRand()

{

    union FloatRand r;
    r.BitArea.Signed = 0;

    r.BitArea.Exp = 1;

    r.BitArea.Frac = (rand() * rand()) % 0x800000;

    if (r.BitArea.Frac == 0x7FFFFF)

        r.BitArea.Exp = rand() % 0x7E;

    else if (r.BitArea.Frac == 0)

        r.BitArea.Exp = rand() % 0x7E + 1;

    else

        r.BitArea.Exp = rand() % 0x7F;

return r.Value; }

驗證代碼

int main(int argc, char * argv[]) { int i; float r;

for (i = 1; i &< 1000000000; i++) { r = GetFloatRand(); printf("%G ", r); if (r == 0.0 || r == 1.0) { printf("failed 0x%08X ", r); return -1; } } return 0; }

同樣的double類型的也可以做出類似效果，肯定是隨機的，並且不需要重試。

補充：評論里有人提到會破壞均勻性，這確實是一個問題，但只要精心調整Exp的出現概率（具體計算就太複雜了）還是可以做到均勻分布的。

改良演算法是：

float GetFloatRand() { union FloatRand r;


    r.BitArea.Signed = 0;

    r.BitArea.Exp = 1;

    r.BitArea.Frac = (rand() * rand()) % 0x800000;

    if (r.BitArea.Frac == 0x7FFFFF)

        r.BitArea.Exp = 0x7D;

    else if (r.BitArea.Frac == 0)

        r.BitArea.Exp = 0x7E;

    else

        r.BitArea.Exp = 0x7E;

return r.Value; }

我再次聲明我的觀點連續出現0的概率確實非常小，但不能因為非常小就完全忽略，既然邏輯上有能力避免，那麼就應該避免，否則就是故意留下一個邏輯上的漏洞。大部分人都對隨機函數有信心，但可惜的是，隨機函數屬於這個代碼里的「外部函數」，不屬於自己設計的，那麼在一個高可靠性的系統里，不應該以「概率」的眼光去看待外部函數，萬一它被黑客攻擊了一直都產生0呢？

說點跟話題無關的內容：當年IPv4設計出來的時候，計算機很少，那時候從概率上來說是不是IP地址永遠都用不完？不應該把自己的代碼漏洞置於概率之上，也不應該把數學上小概率事件放到計算機系統里。另外Milo Yip說的，概率確實很低，但要是float的話，概率恐怕沒那麼低吧。

你生成的隨機數 * 0.999998 + 0.000001 就好了

do { r = random() } while (r === 0 || r === 1)

static int seed=1;

static unsigned int MaxRand=0xffffffff;

void MySrand(int s)

{

seed = s;

}

float MyRand()

{

seed = seed * 214013 + 2531011;

return (float )((seed&>&>16)(0x7FFF))/MaxRand;

}

#include&

...

double r;

{

r=rand()/RAND_MAX;

}while(r==0.0||r==1.0)

...

double foo() { double r; while(1) r = random(); if (0.000001 &< r r &< 0.999999) { return r; } } }

非常簡單：當結果等於0或1時，重新產生一個隨機數。

可以保證結果仍然是均勻分布。

支持重試，實踐中遇到連續n次0乃至影響性能的概率比遇到太陽黑子爆發宕機小几個數量級

(rand() + 1.0) / (RAND_MAX + 2.0); //這樣行不行？

寫得不嚴謹，包括兩個端點是閉區間。應該是［a，b］。至於怎麼生成random哈函數就是，但是不包括兩個端點，此時你可以定個精度，再判斷減速0或減去1是不是無限接近這個精度。是的話就取零或者1。

KR的C程序設計語言里講了，第八章。

(double)rand()/RAND_MAX

等於0或1的情況重新生成就可以了。