[精讚] [會員登入]
625

[PERL] 常態分配亂數產生函數

一般而言我們使用程式內建的 rand 函數取得的亂數是平均分配的,但如果需要常態分配的亂數產生器該怎麼做?

此文完整連結 http://n.sfs.tw/10536

複製連結 [PERL] 常態分配亂數產生函數@新精讚
(文章歡迎轉載,務必尊重版權註明連結來源)
2017-01-09 11:53:50 最後編修
2017-01-03 01:43:56 By 張○○
 

一般而言我們使用程式內建的 rand 函數取得的亂數是平均分配的,但如果需要常態分配的亂數產生器該怎麼做?

因此我採用Perl改寫亂數產生函式,使其產生1,000,000筆亂數,在我的SERVER上跑的效能極佳,只需要大約4.6秒即可以完成(HP DL380G3 Xeon3.0x2)。

以下是測試的結果,均值為6,標準差為1,我把每隔1統計一筆資料:
以下是程式產生亂數後再分群統計數量,由結果看的確是符合常態分配的函數:

Program Start....
1 3
2 235
3 5919
4 60437
5 242095
6 382825
7 241478
8 60621
9 6177
10 209
11 1
Elapsed time=00:00:04.671079


至於實測10,000,000筆的資料,在 linux OS和INTEL雙CPU(非雙核喔)表現上(沒有寫多執行緒,也許沒法最佳化),所費的時間毫無意外是1,000,000筆的10倍。這裡如果拿 IA64 來測的話,不知效果如何… 。

Program Start....
1 31
2 2231
3 59565
4 605552
5 2417906
6 3831267 
<==均值
7 2415974
8 605805
9 59397
10 2245
11 27
Elapsed time=00:00:46.395085[


以下是 Perl的程式碼片段,分群組只花了三行:

#!/usr/bin/perl -w

use Time::Elapse;
print "Program Start....\n";
$avg = 6; # 均值
$sd =1; # 標準差
$num = 10000000; # 產生組數
$dif = 1; # 統計柱狀圖區間幾個 $sd
# ====== 以上是參數的設定部分

my %h=(); # 分組的HASH
Time::Elapse->lapse(my $now);
for(1..$num){
    do{
        $v1 = (2 * rand 1) - 1;
        $v2 = (2 * rand 1) - 1;
        $r = $v1 **2 + $v2 **2;
    }while($r >= 1);
    $fac = sqrt(-2 * log($r)/ log(exp(1)) / $r);
    $gauss = $v2 * $fac;
    $vr= $gauss * $sd + $avg;     # $vr 就是得到的亂數
    # 計算群組中間值
    $sddif = $sd*$dif;
    $hid= int( ( $vr + $sddif /2 ) / $sddif ) * $sddif;
    # 置入hash
    $h{$hid}++;
}

foreach $key (sort {$a<=>$b} keys %h) {
    print "$key \t\t$h{$key}\n";
}
print "Elapsed time=". $now;

如果只要一個「常態分配」的亂數產生函數而不要像上面的統計分析,我把他抽出來變副程式:

#!/usr/bin/perl -w
use strict;

#&mean_dis_fun(均值, 標準差);
print &mean_dis_fun(5, 1);

sub mean_dis_fun{
    my $avg =shift;
    my $sd = shift;
    my ($r, $v1, $v2);
    do{
        $v1 = (2 * rand 1) - 1;
        $v2 = (2 * rand 1) - 1;
        $r = $v1 **2 + $v2 **2;
    }while($r >= 1);
    my $fac = sqrt(-2 * log($r)/ log(exp(1)) / $r);
    my $gauss = $v2 * $fac;
    return $gauss * $sd + $avg;
}


原文撰於 2008-10-11,重新整理於 2009-11-15

你可能感興趣的文章

[PERL] 21- 使用資料庫 Perl 連Mysql

[PERL] 08-陣列 #2 --操作 更多perl陣列的操作,如拆開、黏合、取出、加入等等

[PERL] 17-參照 PERL的參照,就是指標

[PERL] 03-條件式判斷 perl的條件式判斷

[PERL] 位元運算 Bitwise operation PERL的位元運算整理

[PERL] 自寫 trim, ltrim, rtrim 函數 常用的 perl 函數中沒有trim, ltrim, rtrim這三個函數,那就自己來寫,非常容易

自行撰寫syslog server建立資訊安全控管中心#3 -- 程式下載和安裝 利用PERL將syslog收攏到資料庫的程式,程式下載和安裝

[PERL] 22-日期和時間 Perl 取得日期時間的方法整理,使用 localtime等多種函數

關於Perl,你得知道這幾點 Perl 有不少奇異的地方,值得看一看

[PERL] 05-運算子 #1 PERL的運算子介紹,總共有21種

我有話要說


限制:留言最高字數1000字,超過部分會被截掉。請注意:留言不可帶有網址,會被濾掉。 限制:未登入訪客,每則留言間隔需超過10分鐘,每日最多5則留言。

訪客留言

[無留言]

隨機好文

為什麼要買長達二十年的保單? 為什麼要買長達二十年的保單?找一個可以說服我買二十年保單的理由。

PHP for sphinx 函式庫安裝 PECL/sphinx PHP>= 5.2.2 已經能原生支援 sphinx,可是預設的沒有裝,我們得自己裝才能用

維修冰箱 維修冰箱

魔球中小女孩唱的歌 The show 魔球中小女孩唱的歌 The show

[AS3] 我做的唯一一個Flash As3遊戲UFO INVADSION [AS3] 我做的唯一一個Flash As3遊戲,是第一個也是最後一個,後來就沒再寫as3,不過as3還滿好玩的。