[精讚] [會員登入]
2170

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

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

分享此文連結 //n.sfs.tw/10536

分享連結 [PERL] 常態分配亂數產生函數@新精讚
(文章歡迎轉載,務必尊重版權註明連結來源)
2019-10-24 04:09:19 最後編修
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

END

你可能感興趣的文章

Freebsd11上安裝PERL模組 freebsd11+ perl module 在freebsd上安裝module還算簡單

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

[PERL] 18-套件及模組 套件和模組入門

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

Perl 的列舉寫法(map list) 使用列舉的寫法來設定陣列,非常的強大

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

我有話要說

>>

限制:留言最高字數1000字。 限制:未登入訪客,每則留言間隔需超過10分鐘,每日最多5則留言。

訪客留言

[無留言]

隨機好文

UTF-8 BOM (Byte Order Mark) 的問題 在 Michael Kaplan 那看到 Every character has a story #4: U+feff

[札記] 2016.7~12月札記 札記,只是札記

為什麼要重造輪子? 什麼輪子?造什麼輪子?我為什麼要重造輪子?

NETCRAFT發現你的網站及作業系統 NETCRAFT可以發現你的網站及作業系統

[PHP] 檢查檔案是否是圖檔 使用getimagesize函數檢查檔案是否是圖檔