詳細はこちらをどうぞ→Mersenne Twister Home Page
少なくとも,ABで使っている乱数よりは性能がいいそうです。
また,不必要に乗算を使用しないためかABの乱数よりも少しだけ速いです。
<処理時間比較>
MT法 ABのrand (時間:msec)
100,000回 71 76
500,000回 371 394
1,000,000回 759 878
5,000,000回 3855 4276
大学の授業に出てきたのですが,興味があったので原版のCプログラムをABに移植してみました。
※関数リファレンスは省略しました。使用法ドキュメントをお読み下さい。
コードを表示する [ここをクリックすると内容が表示されます]
関数定義部使用例
コード: 全て選択
' ############################################################
' ## ##
' ## メルセンヌツイスターによる乱数発生ライブラリ AB版 ##
' ## Sinryow http://www.sinryow.net/ ##
' ## ##
' ############################################################
' 【注記】
' 「/* ~ */」 の形のコメントの部分は,原プログラム
' http://www.math.sci.hiroshima-u.ac.jp/~ ... t19937ar.c
' に書かれていたコメントそのものです。私が追加したコメントは全
' て「' ~」の形になっています。
'
' また英語の訳の正確性についての保証は致しかねます。もし間違い
' がございましたら,私のホームページ,あるいは私のメールアドレ
' ス mail (at-mark) sinryow.netにご連絡下さい。
/*
A C-program for MT19937, with initialization improved 2002/1/26.
Coded by Takuji Nishimura and Makoto Matsumoto.
Before using, initialize the state by using init_genrand(seed)
or init_by_array(init_key, key_length).
Copyright (C) 1997 - 2002, Makoto Matsumoto and Takuji Nishimura,
All rights reserved.
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions
are met:
1. Redistributions of source code must retain the above copyright
notice, this list of conditions and the following disclaimer.
2. Redistributions in binary form must reproduce the above copyright
notice, this list of conditions and the following disclaimer in the
documentation and/or other materials provided with the distribution.
3. The names of its contributors may not be used to endorse or promote
products derived from this software without specific prior written
permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
Any feedback is very welcome.
http://www.math.sci.hiroshima-u.ac.jp/~ ... T/emt.html
email: m-mat @ math.sci.hiroshima-u.ac.jp (remove space)
*/
' 【訳】
' MT19937(メルセンヌツイスター)のCのプログラムです。初期化の
' 方法を2002年1月26日に改良しました。
'
' 使用する前に,init_genrand(seed) あるいは init_by_array(init_key, key_length).
' を用いて乱数の状態を初期化して下さい。
'
' Copyright (C) 1997 - 2002, Makoto Matsumoto and Takuji Nishimura,
' All rights reserved.
'
' 以下の条件にあてはまる限り,ソースコードでもバイナリ形式でも,
' 改変があってもなくても再配布・利用が可能です。(訳注1)
'
' <利用条件>
' 1. ソースコードの再配布においては,上記の著作権表示,この利
' 用条件および下記の免責条項を残さなければならない。
' 2. バイナリ形式での再配布においては,それと共に配布されるド
' キュメント等に,上記の著作権表示,この利用条件および下記
' の免責条項を複製して載せなければならない。
' 3. このプログラムの作成に関わった者は,このプログラムを利用
' して作られた別のプログラムについては,事前に明記された許
' 可がない限り,そのプログラムを発展させた者の名前として用
' いてはならない。
'
' <免責条項>
' 本ソフトウェアは著作権者によって,「現に提供されている状態の
' ままで」提供されるものとする。本ソフトウェアについては、明示
' あるいは黙示を問わず,商用品として通常そなえるべき品質をそな
' えているとの保証も,特定の目的に適合するとの保証を含め,何の
' 保証もないものとする。事由のいかんを問わず,損害発生の原因い
' かんを問わず,且つ,責任の根拠が契約であるか厳格責任であるか
' (過失その他の)不法行為であるかを問わず,著作権者も寄与者も,
' 仮にそのような損害が発生する可能性を知らされていたとしても,
' 本ソフトウェアの使用から発生した直接的な・間接的な・偶発的な・
' 特別な・懲罰的な・結果的な損害のいずれに対しても(代替品また
' はサービスの提供,使用機会,データまたは利益の損失の補償,ま
' たは,業務の中断に対する補償を含め)責任を一切負わないものと
' する。
'
' ご意見・ご感想はこちらへお願いします。(訳注2)
' http://www.math.sci.hiroshima-u.ac.jp/~ ... T/emt.html
' メール: m-mat @ math.sci.hiroshima-u.ac.jp
' (スパム防止のためスペースを入れてあります)
'
' 【訳注1】
' 日本語訳の部分は法的には有効ではありませんので,実際にこのプ
' ログラムの再配布を行う場合は,「英語の部分を」再配布するもの
' に付加する必要があります。
'
' 【訳注2】
' AB版のプログラムについてのお問い合わせは,私(Sinryow)にお
' 願いします。松本氏(C言語版作者)には問い合わせないようお願
' いします。
Const __MT_N = 624
Const __MT_M = 397
Const __MT_MATRIX_A = &H9908b0df /* constant vector a */
Const __MT_UPPER_MASK = &H80000000 /* most significant w-r bits */
Const __MT_LOWER_MASK = &H7fffffff /* least significant r bits */
Dim __mt_mt[ELM(__MT_N)] As DWord /* the array for the state vector */
Dim __mt_mti = (__MT_N)+1 As Long /* mti==N+1 means mt[N] is not initialized */
/* initializes mt[N] with a seed */
Sub init_genrand(s As DWord)
__mt_mt[0] = s
'__mt_mt[0] = s And &Hffffffff
For __mt_mti = 1 To (__MT_N)-1
__mt_mt[(__mt_mti)] = _
(1812433253 As DWord * (__mt_mt[(__mt_mti)-1] Xor (__mt_mt[(__mt_mti)-1] >> 30)) _
+ (__mt_mti))
/* See Knuth TAOCP Vol2. 3rd Ed. P.106 for multiplier. */
/* In the previous versions, MSBs of the seed affect */
/* only MSBs of the array mt[]. */
/* 2002/01/09 modified by Makoto Matsumoto */
'__mt_mt[__mt_mti] = __mt_mt[__mt_mti] And &Hffffffff
/* for >32 bit machines */
Next
End Sub
/* initialize by an array with array-length */
/* init_key is the array for initializing keys */
/* key_length is its length */
/* slight change for C++, 2004/2/26 */
Sub init_by_array(init_key As DWordPtr, key_length As Long)
Dim i As Long, j As Long, k As Long
init_genrand(19650218)
i=1: j=0
If (__MT_N)>key_length Then
k = (__MT_N)
Else
k = key_length
End If
While k
__mt_mt = _
(__mt_mt Xor ((__mt_mt[i-1] Xor (__mt_mt[i-1] >> 30)) * 1664525 As DWord)) _
+ init_key[j] + j /* non linear */
'__mt_mt = __mt_mt And &Hffffffff; /* for WORDSIZE > 32 machines */
i = i+1: j = j+1
If i>=(__MT_N) Then
__mt_mt[0] = __mt_mt[(__MT_N)-1]
i=1
End If
If j>=key_length Then j=0
k = k-1
Wend
For k=(__MT_N)-1 To 1 Step -1
__mt_mt = _
(__mt_mt Xor ((__mt_mt[i-1] Xor (__mt_mt[i-1] >> 30)) * 1566083941 As DWord)) _
- i /* non linear */
'__mt_mt = __mt_mt And &Hffffffff /* for WORDSIZE > 32 machines */
i = i+1
If i>=(__MT_N) Then
__mt_mt[0] = __mt_mt[(__MT_N)-1]
i=1
End If
Next
__mt_mt[0] = &H80000000; /* MSB is 1; assuring non-zero initial array */
End Sub
/* generates a random number on [0,0xffffffff]-interval */
Dim __mt_mag01[ELM(2)]=[0, __MT_MATRIX_A] As DWord
/* mag01[x] = x * MATRIX_A for x=0,1 */
Function genrand_int32() As DWord
Dim y As DWord
If __mt_mti >= (__MT_N) Then
/* generate N words at one time */
Dim kk As Long
If __mt_mti = (__MT_N)+1 Then init_genrand(5489)
/* if init_genrand() has not been called, */
/* a default initial seed is used */
For kk=0 To (__MT_N)-(__MT_M)-1
y = (__mt_mt[kk] And __MT_UPPER_MASK) Or (__mt_mt[kk+1] And __MT_LOWER_MASK)
__mt_mt[kk] = __mt_mt[kk+(__MT_M)] Xor (y >> 1) Xor __mt_mag01[y And &H1]
Next
While kk<(__MT_N)-1
y = (__mt_mt[kk] And __MT_UPPER_MASK) Or (__mt_mt[kk+1] And __MT_LOWER_MASK)
__mt_mt[kk] = __mt_mt[kk+(__MT_M)-(__MT_N)] Xor (y >> 1) Xor __mt_mag01[y And &H1]
kk = kk+1
Wend
y = (__mt_mt[(__MT_N)-1] And __MT_UPPER_MASK) Or (__mt_mt[0] And __MT_LOWER_MASK)
__mt_mt[(__MT_N)-1] = __mt_mt[(__MT_M)-1] Xor (y >> 1) Xor __mt_mag01[y And &H1]
__mt_mti = 0
End If
y = __mt_mt[__mt_mti]
__mt_mti = __mt_mti+1
/* Tempering */
y = y Xor (y >> 11)
y = y Xor ((y << 7) And &H9d2c5680 As DWord)
y = y Xor ((y << 15) And &Hefc60000 As DWord)
y = y Xor ((y >> 18))
genrand_int32 = y
End Function
/* generates a random number on [0,0x7fffffff]-interval */
Function genrand_int31() As Long
genrand_int31 = genrand_int32()>>1
End Function
/* generates a random number on [0,1]-real-interval */
Function genrand_real1() As Double
genrand_real1 = (genrand_int32() As Double)*(1.0 As Double/4294967295.0 As Double)
/* divided by 2^32-1 */
End Function
/* generates a random number on [0,1)-real-interval */
Function genrand_real2() As Double
genrand_real2 = (genrand_int32() As Double)*(1.0 As Double/4294967296.0 As Double)
/* divided by 2^32 */
End Function
/* generates a random number on (0,1)-real-interval */
Function genrand_real3() As Double
genrand_real3 = ((genrand_int32() As Double) + 0.5)*(1.0 As Double/4294967296.0 As Double)
/* divided by 2^32 */
End Function
/* generates a random number on [0,1) with 53-bit resolution*/
Function genrand_res53() As Double
Dim a = genrand_int32()>>5 As DWord, b = genrand_int32()>>6 As DWord
genrand_res53 =
((a As Double)*(67108864.0 As Double)+(b As Double))* _
(1.0 As Double/9007199254740992.0 As Double)
End Function
/* These real versions are due to Isaku Wada, 2002/01/09 added */
コード: 全て選択
#include "Mersenne.idx"
#console
Declare Function timeGetTime Lib "winmm.dll" () As Long
' ↓ ここからプログラムが実行されます
Dim i As Long
Dim init[ELM(4)]=[&H123, &H234, &H345, &H456] As DWord
Dim length=4 As Long
init_by_array(init, length)
' 以下,乱数発生テスト
/*
Print "1000 outputs of genrand_int32()"
For i=0 To 999
Print genrand_int32()
Next
Print
Print "1000 outputs of genrand_real2()"
For i=0 To 999
Print genrand_real2()
Next
*/
' 以下,速度測定用
Dim StartTime As Long
Dim EndTime As Long
Dim REPEATTIME As Long
Input "Repeat"; REPEATTIME
StartTime=timeGetTime()
For i=1 To REPEATTIME
genrand_int32()
Next
EndTime = timeGetTime()
Print "Mersenne", EndTime-StartTime
StartTime=timeGetTime()
For i=1 To REPEATTIME
rand()
Next
EndTime = timeGetTime()
Print "rand()",, EndTime-StartTime
Dim s As String
Input s
※補足
AB4用ですが,「/* ~ */」型のコメントを「' ~ 」にして,Dim文中で変数の初期化を行わなければ,AB2・3でも使えるかもしれません。