產生指數隨機數的兩種方法的比較…

800)  r--->0  r末尾溢出
ln(r)  -inf 
所以連乘, n不易過大.....
改進NewRandom_exp2.m
對於任意n, 分段處理,可分m+1段
m=[n/100](取整)    尾段=mod(n,100) 模數
100對於連乘適中..  對於誤差適中...
test.m測試程式...(測試環境 P-m 1.6G  512m DDR400)
產生100000個。。。。（S)
>> [Time1,Time2]=test(100000)
Time1 is the runtime of random_exp1 Time1 =   78.0422m Time2 =    1.0715
                                                                                         made by aris  2006 3 20

 

%exponent distribution1........
function [ExpRandom,RunTime]=random_exp1(n);
t=cputime;
ExpTemp=[];
for i=1:n
    ExpTemp(i)=-log(rand);
end
RunTime=cputime-t;
ExpRandom=ExpTemp;

%exponent distribution2........
function [ExpRandom,RunTime]=random_exp2(n)
t=cputime;
ExpTemp=[];
RandTemp1=rand(1,n);
RandTemp2=rand(1,n-1);
r=1;
for i=1:n
    r=r*RandTemp1(i);
end
r=log(r);
RandTemp2=[0,RandTemp2,1];
RandTemp2=sort(RandTemp2);
for i=1:n
    ExpTemp(i)=(RandTemp2(i)-RandTemp2(i+1))*r;
end
RunTime=cputime-t;
ExpRandom=ExpTemp;

 

%exponent distribution........
function [ExpRandom,RunTime]=NewRandom_exp2(n)
t=cputime;
m=int16(n/100);
l=mod(n,100);
RandTemp=[];
for i=1:m
    RandTemp(i,:)=random_exp2(100);
end
RandTemp=vec(RandTemp);
RandTempl=random_exp2(l);
ExpRandom=[RandTemp';RandTempl];
RunTime=cputime-t;

test the [ExpRandom,RunTime]=random_exp1(n)
%         [ExpRandom,RunTime]=random_exp2(n)
%         [ExpRandom,RunTime]=NewRandom_exp2(n)
function [Time1,Time2]=test(n);
[ExpRandom,RunTime1]=random_exp1(n);
[ExpRandom,RunTime2]=NewRandom_exp2(n);
disp('Time1 is the runtime of random_exp1');
Time1=RunTime1;
Time2=RunTime2;