matlab逆变换法产生随机数_生成伪随机数--Matlab实现

【来源】：上海交通大学BBS

随机数的产生

在利用计算机随机模拟方法(蒙特卡罗方法)时，需要产生随机数，为了避免一

些朋友们在这方面花费太多的时间，在这里对随机数的产生方法作一简单阐述：

一般说来，通过编程的方法产生随机数，是按照一定的计算方法产生的一列数，这

样产生的随机数并不是真正意义上的随机，称为伪随机数。但计算方法选的得当，它便

会近似于所需要的分布的随机数。

均匀分布的随机数是各类随机数产生的基础，在它的产生方法中最常用的为余数法

。

(由于插入符号较麻烦，我用中括号括起来的内容为其左边字符的脚标)

y[n+1]=w*y[n]+c(mod M)， (mod M)表示前面的表达式除以M后的余数，注意这

里是w*y[n]+c除以M。

初值y[0],乘子w,和模M取非负数，当c=0时，y[n+1]=w*y[n](mod M) ，

y[0]=a(a为奇数),x[0]=y[0]/M

得到的{X[n]}为需要的随机数，实际上是周期为L的伪随机数，L<=M.后面给出一个

基于MATLAB的M文件。

正态分布使用较频繁，在全赛的零件参数设计中采用蒙特卡罗方法就需要产生一系

列正态分布的随机数。产生N(0,1)正态分布的随机数中

较常用的有筛选法，变换法，近似法。这里介绍精度较高的筛选法，由于方法极其简单

，只有几行语句，MATLAB所用的语言几乎是和C一样的，感兴趣的朋友可阅读后面的程序

。

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%产生一组[0，1]分布的随机数，采用余数法

%从文献中得知下列参数组较为适用：

%y的初值为1，w=7,模M=10^10(伪随机数周期为5*10^7)

%y的初值为1，w=5^13,模M=2^36(伪随机数周期为2^34,约2*10^10)

%y的初值为1，w=5^17,模M=2^42(伪随机数周期为2^40,约10^12)

%--------------------------------------------------------

%function x=uniform(y,w,M,n)

%n为要产生的随机数个数

function x=uniform(y,w,M,n)

for i=1:n

x(i)=y./M;

y=rem(w*y,M);%rem()为求余函数

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%产生一个N(0,1)正态分布随机数

%采用筛选法，精度较高

%其他方法可参阅《现代应用数学手册--概率统计与随机过程卷》清华大学出版社马振华

主编

%function y=riddling()

function y=riddling()

sign=0;

while 1

x=rand(1,2);%产生两个[0，1]间均匀分布的随机数

v1=2*x(1)-1;v2=2*x(2)-1;

s=v1^2+v2^2;

if s<=1

a=sqrt(-2*log(s)/s);

y=v1*a;

sign=1;

end

if sign==1

break;

end

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%产生n个N(a,b)正态分布随机数

%其中a为均值，b为方差

%function x=normal(a,b,n)

function x=normal(a,b,n)

m=48;%m应尽量大，m取12时近似程度已较好

for i=1:n

r=rand(1,m);

x(i)=a+sqrt(b)*(sum(r)-m/2)/sqrt(m/12);

end