实验数据处理方法
第二部分,Monte Carlo模拟
第七章
均匀分布随机数的产生
? 随机数的定义和特性
? 随机数的产生
? 线性乘同余方法
第七章
均匀分布随机数的产生
7.1 随机数的定义和特性
7.1 随机数的定义和特性
什么是随机数？
?单个的数字不是随机数
?是指一个数列，其中的每一个体称为随机数，其值与数列中
的其它数无关；
?在一个均匀分布的随机数中，每一个体出现的概率是均等的；
?例如：在 [0,1]区间上均匀分布的随机数序列中，
0.00001与 0.5出现的机会均等
7.1 随机数的定义和特性
随机数应具有的基本特性
?考虑一个对高能粒子反应过程的模拟：需用随机数确定,
?出射粒子的属性：能量、方向,…
?粒子与介质的相互作用
?对这一过程的模拟应满足以下要求（相空间产生过程）,
? 出射粒子的属性应是互不相关的，即每一粒子的属
性的确定独立于其它的粒子的属性的确定；
? 粒子的属性的分布应满足物理所要求的理论分布；
?所模拟的物理过程要求随机数应具有下列特性,
?随机数序列应是独立的、互不相关的 (uncorrelated),
即序列中的任一子序列应与其它的子序列无关；
7.1 随机数的定义和特性
?长的周期 (long period),
实际应用中，随机数都是用数学方法计算出来的，这些
算法具有周期性，即当序列达到一定长度后会重复；
?均匀分布的随机数应满足均匀性 (Uniformity),
随机数序列应是均匀的、无偏的，即：如果两个子区间
的“面积”相等，则落于这两个子区间内的随机数的个
数应相等。
例如：对 [0,1)区间均匀分布的随机数，如果产生
了足够多的随机数，而有一半的随机数落
于区间 [0,0.1]?不满足均匀性
如果均匀性不满足，则会出现序列中的多组随机数相
关的情况 ?均匀性与互不相关的特性是有联系的
7.1 随机数的定义和特性
?有效性（ Efficiency),
模拟结果可靠
?模拟产生的样本容量大
?所需的随机数的数量大
?随机数的产生必须快速、有效，最好能
够进行并行计算。
第七章
均匀分布随机数的产生
7.2 随机数的产生
7.2 随机数的产生
? [0,1]区间上均匀分布的随机数是 Monte Carlo模拟的基础,
?[0,1]均匀分布的随机数的产生方法,
?利用一些具有内在的随机性的过程,
?放射性衰变过程（ radioactive decay)；
?热噪声 (thermal noise);
?宇宙线的到达时间（ cosmic ray arrival);
?…
?缺点：模拟的结果不可再现，使得模拟程序的找错困难
?利用事先制订好的随机数表,
?缺点：表的容量有限，不适合需要大量随机数的应用
?服从任意分布的随机数序列可以用 [0,1]区间均匀分布的
随机数序列作适当的变换或舍选后求得
7.2 随机数的产生
),...,,( 11 ???? ? knnnkn rrrTr
?利用数学递推公式在计算机中产生随机数
其中,T为某个函数，给定初值 r1,r2,…,r k,可按上式确
定 rn+1,n=1,2,… ?随机数序列,
算法：产生 [0,M]区间上的整数 In，然后利用公式
rn=In/M返回 [0,1]区间上的实数
优点,
–占用计算机的内存少；
–产生速度快；
–可以重复前次的模拟结果，便于程序的找错；
7.2 随机数的产生
缺点,
? 不满足随机数之间相互独立的要求：公式和初值确定
后，序列就唯一地确定了；
?伪随机数（ Pseudo-Random Number)
? 不满足均匀性：计算机能表示的 [0,1]区间内的数是有限
的（由字长确定）
?递推到一定次数后，出现周期性的重复现象
第七章
均匀分布随机数的产生
7.3 线性乘同余方法
（ Linear Congruential Method)
7.3 线性乘同余方法
（ Linear Congruential Method)
mcaII nn m o d)(1 ???
caIm
ca
,,
0,
0?
?
1948年由 Lehmer提出的一种产生伪随机数的方法，是最常用的方法。
1、递推公式,
其中,
I0,初始值（种子 seed)
a,乘法器 （ multiplier)
c,增值（ additive constant)
m,模数（ modulus)
mod：取模运算,(aIn+c)除以 m后的余数
a,c和 m皆为整数
?产生整型的随机数序列,随机性来源于取模运算
如果 c=0 ? 乘同余法：速度更快，也可产生长的随机数序列
7.3 线性乘同余方法
（ Linear Congruential Method)
]1,0[
)1(
)1,0[
)(
?
?
?
??
mf l o a t
I
r
mf l o a t
I
r
n
n
n
n
1??
?
mI
mI
n
n
2、实型随机数序列,
3、特点,
1）最大容量为 m,mI
n ??0
2）独立性和均匀性取决于参数 a和 c的选择
例,a=c=I0=7,m=10 ? 7,6,9,0,7,6,9,0,…
7.3 线性乘同余方法
（ Linear Congruential Method)
4、模数 m的选择,
? m 应尽可能地大，因为序列的周期不可能大于 m;
? 通常将 m取为计算机所能表示的最大的整型量，在 32位计算
机上,m=231=2x109
5、乘数因子 a的选择,
1961年,M,Greenberger证明：用线性乘同余方法产生的随机
数序列具有周期 m的条件是,
1,c和 m为互质数；
2,a-1是质数 p的倍数，其中 p是 a-1和 m的共约数；
3,如果 m是 4的倍数,a-1也是 4的倍数。
例,a=5,c=1,m=16,I0=1 ?周期 =m=16
1,6,15,12,13,2,11,8,9,14,7,4,5,10,3,0,1,6,15,12,13,2,.,
7.3 线性乘同余方法
（ Linear Congruential Method)
RANDU随机数产生器,
311 2m o d)6 5 5 3 9( nn II ???
1961年由 IBM提出
unsigned long seed = 9;
float randu() {
const unsigned long a = 65539;
const unsigned long m = pow(2,31);
unsigned long i1;
i1 = (a * seed) % m;
seed = i1;
return (float) i1/float(m);
}
void SetSeed(unsigned long i) {
seed = i;
}
7.3 线性乘同余方法
（ Linear Congruential Method)
存在严重的问题,
Marsaglia效用，存在于所有乘同余方法的产生器
void test() {
c1 = new TCanvas("c1",“Test of random number generator",200,10,700,900);
pad1 = new TPad("pad1",“one ",0.03,0.62,0.50,0.92,21);
pad2 = new TPad("pad2",“one vs one",0.51,0.62,0.98,0.92,21);
pad3 = new TPad("pad3",“one vs one vs one",0.03,0.02,0.97,0.57,21);
pad1->Draw();
pad2->Draw();
pad3->Draw();
TH1F * h1 = new TH1F("h1","h1",100,0.0,1.0);
TH2F * h2 = new TH2F("h2","h2",100,0.0,1.0,100,0.0,1.0);
TH3F * h3 = new TH3F("h3","h3",100,0.0,1.0,100,0.0,1.0,100,0.0,1.0);
7.3 线性乘同余方法
（ Linear Congruential Method)
for(int i=0; i < 10000; i++) {
float x = randu();
float y = randu();
float z = randu();
h1->Fill(x);
h2->Fill(x,y);
h3->Fill(x,y,z);
}
pad1->cd(); h1->Draw();
pad2->cd(); h2->Draw();
pad3->cd(); h3->Draw();
}
7.3 线性乘同余方法
（ Linear Congruential Method)
7.3 线性乘同余方法
（ Linear Congruential Method)
如果取 a=69069,将极大地改善结果
7.3 线性乘同余方法
（ Linear Congruential Method)
mIbIaI nnn m o d)( 21 ?? ????
1968年,Marsaglia对这一问题进行了研究，认为,
任何的乘同余产生器都存在这一问题：在三维和三维以
上的空间中，所产生的随机数总是集聚在一些超平面上
?随机数序列是关联的
对于 32位的计算机，在 d-维空间中超平面的最大数目为
d=3 2953
d=4 566
d=6 120
d=10 41
改进措施：将递推公式修改为
特点,1）需要两个初始值（种子）；
2）周期可大于 m;
7.3 线性乘同余方法
（ Linear Congruential Method)
#include <math.h>
unsigned long seed0 = 9;
unsigned long seed1 = 11;
float randac() {
const unsigned long a = 65539;
const unsigned long b = 65539;
unsigned long i2;
unsigned long m = pow(2,31);
i2 = (a * seed1 + b * seed0 ) % m;
seed0 = seed1; seed1 = i2;
return (float) i1/float(m);
}
void SetSeed(unsigned long i0,
unsigned long i1) {
seed0 = i0;
seed1 = i1;
}
7.3 线性乘同余方法
（ Linear Congruential Method)
a=b=65539,seed0=9,seed1=11
7.3 线性乘同余方法
（ Linear Congruential Method)
如何获取 [0,1]区间均匀分布的随机数产生器,
1,每一个 Monte Carlo模拟程序软件包都有自带的产生器,
? Jetset(LUND Monte Carlo模拟系列）：利用 Marsaglia等
所提出的算法，周期可达 1043
函数用法,r=rlu(idummy)
? Geant3(探测器模拟程序,FORTRAN),周期 =1018
Call grndm(vec*,len)
…,
2,利用 CERN程序库,
? Y=rndm(x),周期,5x108
? Y=rn32(dummy):乘同余法,a=69069,i0=65539
? Call ranmar(vec,len),周期,1043
? Call ranecu(vec,len,isq)
7.3 线性乘同余方法
（ Linear Congruential Method)
CLHEP(Class Library for High
Energy Physics)中的随机数产
生器
3,利用 CLHEP中的随机数产生器软件包,
http://hepg.sdu.edu.cn/~zhangxy/clhep/html/index.html
7.3 线性乘同余方法
（ Linear Congruential Method)
FORTRAN中使用随机数产生器应注意的问题,
在 FORTRAN中，如果随机数产生器是带 dummy变量的函数,
其中变量 idum在函数中不使用，应注意以下问题,
X=RAND(idum)
FORTRAN编译器在对程序进行优化时,
X=RAND(IDUM)+RAND(IDUM) ? X=2.0*RAND(IDUM)
DO I=1,10
X=RAND(IDUM)
…
END DO
X=RAND(IDUM)
DO I=1,10
…,
END DO
?
解决办法,DO I=1,10
IDUM = IDUM +1
X=RAND(IDUM)
…
END DO
第二部分,Monte Carlo模拟
第七章
均匀分布随机数的产生
? 随机数的定义和特性
? 随机数的产生
? 线性乘同余方法
第七章
均匀分布随机数的产生
7.1 随机数的定义和特性
7.1 随机数的定义和特性
什么是随机数？
?单个的数字不是随机数
?是指一个数列，其中的每一个体称为随机数，其值与数列中
的其它数无关；
?在一个均匀分布的随机数中，每一个体出现的概率是均等的；
?例如：在 [0,1]区间上均匀分布的随机数序列中，
0.00001与 0.5出现的机会均等
7.1 随机数的定义和特性
随机数应具有的基本特性
?考虑一个对高能粒子反应过程的模拟：需用随机数确定,
?出射粒子的属性：能量、方向,…
?粒子与介质的相互作用
?对这一过程的模拟应满足以下要求（相空间产生过程）,
? 出射粒子的属性应是互不相关的，即每一粒子的属
性的确定独立于其它的粒子的属性的确定；
? 粒子的属性的分布应满足物理所要求的理论分布；
?所模拟的物理过程要求随机数应具有下列特性,
?随机数序列应是独立的、互不相关的 (uncorrelated),
即序列中的任一子序列应与其它的子序列无关；
7.1 随机数的定义和特性
?长的周期 (long period),
实际应用中，随机数都是用数学方法计算出来的，这些
算法具有周期性，即当序列达到一定长度后会重复；
?均匀分布的随机数应满足均匀性 (Uniformity),
随机数序列应是均匀的、无偏的，即：如果两个子区间
的“面积”相等，则落于这两个子区间内的随机数的个
数应相等。
例如：对 [0,1)区间均匀分布的随机数，如果产生
了足够多的随机数，而有一半的随机数落
于区间 [0,0.1]?不满足均匀性
如果均匀性不满足，则会出现序列中的多组随机数相
关的情况 ?均匀性与互不相关的特性是有联系的
7.1 随机数的定义和特性
?有效性（ Efficiency),
模拟结果可靠
?模拟产生的样本容量大
?所需的随机数的数量大
?随机数的产生必须快速、有效，最好能
够进行并行计算。
第七章
均匀分布随机数的产生
7.2 随机数的产生
7.2 随机数的产生
? [0,1]区间上均匀分布的随机数是 Monte Carlo模拟的基础,
?[0,1]均匀分布的随机数的产生方法,
?利用一些具有内在的随机性的过程,
?放射性衰变过程（ radioactive decay)；
?热噪声 (thermal noise);
?宇宙线的到达时间（ cosmic ray arrival);
?…
?缺点：模拟的结果不可再现，使得模拟程序的找错困难
?利用事先制订好的随机数表,
?缺点：表的容量有限，不适合需要大量随机数的应用
?服从任意分布的随机数序列可以用 [0,1]区间均匀分布的
随机数序列作适当的变换或舍选后求得
7.2 随机数的产生
),...,,( 11 ???? ? knnnkn rrrTr
?利用数学递推公式在计算机中产生随机数
其中,T为某个函数，给定初值 r1,r2,…,r k,可按上式确
定 rn+1,n=1,2,… ?随机数序列,
算法：产生 [0,M]区间上的整数 In，然后利用公式
rn=In/M返回 [0,1]区间上的实数
优点,
–占用计算机的内存少；
–产生速度快；
–可以重复前次的模拟结果，便于程序的找错；
7.2 随机数的产生
缺点,
? 不满足随机数之间相互独立的要求：公式和初值确定
后，序列就唯一地确定了；
?伪随机数（ Pseudo-Random Number)
? 不满足均匀性：计算机能表示的 [0,1]区间内的数是有限
的（由字长确定）
?递推到一定次数后，出现周期性的重复现象
第七章
均匀分布随机数的产生
7.3 线性乘同余方法
（ Linear Congruential Method)
7.3 线性乘同余方法
（ Linear Congruential Method)
mcaII nn m o d)(1 ???
caIm
ca
,,
0,
0?
?
1948年由 Lehmer提出的一种产生伪随机数的方法，是最常用的方法。
1、递推公式,
其中,
I0,初始值（种子 seed)
a,乘法器 （ multiplier)
c,增值（ additive constant)
m,模数（ modulus)
mod：取模运算,(aIn+c)除以 m后的余数
a,c和 m皆为整数
?产生整型的随机数序列,随机性来源于取模运算
如果 c=0 ? 乘同余法：速度更快，也可产生长的随机数序列
7.3 线性乘同余方法
（ Linear Congruential Method)
]1,0[
)1(
)1,0[
)(
?
?
?
??
mf l o a t
I
r
mf l o a t
I
r
n
n
n
n
1??
?
mI
mI
n
n
2、实型随机数序列,
3、特点,
1）最大容量为 m,mI
n ??0
2）独立性和均匀性取决于参数 a和 c的选择
例,a=c=I0=7,m=10 ? 7,6,9,0,7,6,9,0,…
7.3 线性乘同余方法
（ Linear Congruential Method)
4、模数 m的选择,
? m 应尽可能地大，因为序列的周期不可能大于 m;
? 通常将 m取为计算机所能表示的最大的整型量，在 32位计算
机上,m=231=2x109
5、乘数因子 a的选择,
1961年,M,Greenberger证明：用线性乘同余方法产生的随机
数序列具有周期 m的条件是,
1,c和 m为互质数；
2,a-1是质数 p的倍数，其中 p是 a-1和 m的共约数；
3,如果 m是 4的倍数,a-1也是 4的倍数。
例,a=5,c=1,m=16,I0=1 ?周期 =m=16
1,6,15,12,13,2,11,8,9,14,7,4,5,10,3,0,1,6,15,12,13,2,.,
7.3 线性乘同余方法
（ Linear Congruential Method)
RANDU随机数产生器,
311 2m o d)6 5 5 3 9( nn II ???
1961年由 IBM提出
unsigned long seed = 9;
float randu() {
const unsigned long a = 65539;
const unsigned long m = pow(2,31);
unsigned long i1;
i1 = (a * seed) % m;
seed = i1;
return (float) i1/float(m);
}
void SetSeed(unsigned long i) {
seed = i;
}
7.3 线性乘同余方法
（ Linear Congruential Method)
存在严重的问题,
Marsaglia效用，存在于所有乘同余方法的产生器
void test() {
c1 = new TCanvas("c1",“Test of random number generator",200,10,700,900);
pad1 = new TPad("pad1",“one ",0.03,0.62,0.50,0.92,21);
pad2 = new TPad("pad2",“one vs one",0.51,0.62,0.98,0.92,21);
pad3 = new TPad("pad3",“one vs one vs one",0.03,0.02,0.97,0.57,21);
pad1->Draw();
pad2->Draw();
pad3->Draw();
TH1F * h1 = new TH1F("h1","h1",100,0.0,1.0);
TH2F * h2 = new TH2F("h2","h2",100,0.0,1.0,100,0.0,1.0);
TH3F * h3 = new TH3F("h3","h3",100,0.0,1.0,100,0.0,1.0,100,0.0,1.0);
7.3 线性乘同余方法
（ Linear Congruential Method)
for(int i=0; i < 10000; i++) {
float x = randu();
float y = randu();
float z = randu();
h1->Fill(x);
h2->Fill(x,y);
h3->Fill(x,y,z);
}
pad1->cd(); h1->Draw();
pad2->cd(); h2->Draw();
pad3->cd(); h3->Draw();
}
7.3 线性乘同余方法
（ Linear Congruential Method)
7.3 线性乘同余方法
（ Linear Congruential Method)
如果取 a=69069,将极大地改善结果
7.3 线性乘同余方法
（ Linear Congruential Method)
mIbIaI nnn m o d)( 21 ?? ????
1968年,Marsaglia对这一问题进行了研究，认为,
任何的乘同余产生器都存在这一问题：在三维和三维以
上的空间中，所产生的随机数总是集聚在一些超平面上
?随机数序列是关联的
对于 32位的计算机，在 d-维空间中超平面的最大数目为
d=3 2953
d=4 566
d=6 120
d=10 41
改进措施：将递推公式修改为
特点,1）需要两个初始值（种子）；
2）周期可大于 m;
7.3 线性乘同余方法
（ Linear Congruential Method)
#include <math.h>
unsigned long seed0 = 9;
unsigned long seed1 = 11;
float randac() {
const unsigned long a = 65539;
const unsigned long b = 65539;
unsigned long i2;
unsigned long m = pow(2,31);
i2 = (a * seed1 + b * seed0 ) % m;
seed0 = seed1; seed1 = i2;
return (float) i1/float(m);
}
void SetSeed(unsigned long i0,
unsigned long i1) {
seed0 = i0;
seed1 = i1;
}
7.3 线性乘同余方法
（ Linear Congruential Method)
a=b=65539,seed0=9,seed1=11
7.3 线性乘同余方法
（ Linear Congruential Method)
如何获取 [0,1]区间均匀分布的随机数产生器,
1,每一个 Monte Carlo模拟程序软件包都有自带的产生器,
? Jetset(LUND Monte Carlo模拟系列）：利用 Marsaglia等
所提出的算法，周期可达 1043
函数用法,r=rlu(idummy)
? Geant3(探测器模拟程序,FORTRAN),周期 =1018
Call grndm(vec*,len)
…,
2,利用 CERN程序库,
? Y=rndm(x),周期,5x108
? Y=rn32(dummy):乘同余法,a=69069,i0=65539
? Call ranmar(vec,len),周期,1043
? Call ranecu(vec,len,isq)
7.3 线性乘同余方法
（ Linear Congruential Method)
CLHEP(Class Library for High
Energy Physics)中的随机数产
生器
3,利用 CLHEP中的随机数产生器软件包,
http://hepg.sdu.edu.cn/~zhangxy/clhep/html/index.html
7.3 线性乘同余方法
（ Linear Congruential Method)
FORTRAN中使用随机数产生器应注意的问题,
在 FORTRAN中，如果随机数产生器是带 dummy变量的函数,
其中变量 idum在函数中不使用，应注意以下问题,
X=RAND(idum)
FORTRAN编译器在对程序进行优化时,
X=RAND(IDUM)+RAND(IDUM) ? X=2.0*RAND(IDUM)
DO I=1,10
X=RAND(IDUM)
…
END DO
X=RAND(IDUM)
DO I=1,10
…,
END DO
?
解决办法,DO I=1,10
IDUM = IDUM +1
X=RAND(IDUM)
…
END DO
