#include <stdlib.h> #include <math.h> double gaussrand() { static double V1, V2, S; static int phase = 0; double X; if(phase == 0) { do { double U1 = (double)rand() / RAND_MAX; double U2 = (double)rand() / RAND_MAX; V1 = 2 * U1 - 1; V2 = 2 * U2 - 1; S = V1 * V1 + V2 * V2; } while(S >= 1 || S == 0); X = V1 * sqrt(-2 * log(S) / S); } else X = V2 * sqrt(-2 * log(S) / S); phase = 1 - phase; return X; }
其它的方法參見本文的擴展版本, 參見問題 20.36。
參考資料: [Knuth, Sec. 3.4.1 p. 117]; [Marsaglia&Bray]; [Press et al., Sec. 7.2 pp. 288-290]。
翻譯朱群英、孫雲, LaTeX2HTML 編譯 朱群英 (2005-06-23)