蒲丰氏投针问题的模拟过程

蒲丰氏投针问题的模拟过程

ID:13252923

大小:23.50 KB

页数:5页

时间:2018-07-21

蒲丰氏投针问题的模拟过程_第1页
蒲丰氏投针问题的模拟过程_第2页
蒲丰氏投针问题的模拟过程_第3页
蒲丰氏投针问题的模拟过程_第4页
蒲丰氏投针问题的模拟过程_第5页
资源描述:

《蒲丰氏投针问题的模拟过程》由会员上传分享,免费在线阅读,更多相关内容在行业资料-天天文库

1、蒲丰氏投针问题的模拟过程,随机数发生器也是自编的,以供大家参考和提出建议。谢谢。(seed1和seed2最好选择3和5,为了使投针次数达到1000000,CVF进行如下设置Project->settings->link->output,将stackallocationsreserve:设为1000000000)programgetpi      implicitnone    real,parameter::a=5,L=4,pi=3.14159      integer::n1,i,counter=0      real,allocatable::

2、R1(:),R2(:)    real::theta,x,pi1    write(*,*)'inputthesizeofthearray:'      read(*,*)n1    allocate(R1(n1))    allocate(R2(n1))      callrandom(n1,R1,R2)    doi=1,n1    x=a*(2*R1(i)-1)    theta=pi*R2(i)      if(abs(x)

3、unter*1.0)*(2.0*L/a)    write(*,*)'thisisPI:'      write(*,*)pi    write(*,"('thisisratioofcountertototalnumber',F10.6)")1.0    &*counter/n1    stop      endprogram        subroutinerandom(n,R1,R2)    implicitnone      integern,seed1,seed2,i,little,m    real::R1(n),R2(n)      i

4、nteger::temp1(n),temp2(n)    write(*,*)'inputseed1'      read(*,*)seed1      write(*,*)'inputseed2'    read(*,*)seed2      m=2**30      m=m*2-1      temp1(1)=397204094      temp2(1)=temp1(1)      R1(1)=(1.0*temp1(1))/(1.0*m)      R2(1)=(1.0*temp2(1))/(1.0*m)    little=0    if(R

5、1(1)<0.5)little=little+1      doi=1,n-1    temp1(i+1)=mod(seed1*temp1(i),m)      R1(i+1)=(1.0*temp1(i+1))/(1.0*m)    temp2(i+1)=mod(seed2*temp2(i),m)      R2(i+1)=(1.0*temp2(i+1))/(1.0*m)    if(R1(i+1)<0.5)little=little+1.      enddo;    write(*,*)'ratioofnumberwhichislittletha

6、n0.5'      write(*,*)1.0*little/n    return    endsubroutine拟蒙特卡洛方法(Quasi-MonteCarlo)积分实例%使用Matlab提供的函数求积分,exp(-1/2*x^2)在(0,1)间积分formatlong;symsxa=sym(1/2);f=exp(-a*x^2);ezplot(f)disp(int(f,-1,1));fprintf('integralresult:%1.18f.',double(int(f,0,1)));%disp(double(int(f,0,1)));

7、复制代码%使用拟蒙特卡洛方法积分%得到拟蒙特卡洛序列,即低偏差序列,halton法%如果有相关的工具箱的话,可以用Matlab里面的haltonset,faureset,sobolset函数实现,x=halton(10000,2,5577);n=length(x);mju=0;fori=1:n  mju=mju+exp(-0.5*x(i)^2);endmju=mju/n;fprintf('Quasi-MonteCarloresult:%1.18f.',mju);%disp(mju);%使用蒙特卡洛方法积分%得到Uniform序列,x=random

8、('unif',0,1,10000,1);n=length(x);mju=0;fori=1:n  mju=m

当前文档最多预览五页,下载文档查看全文

此文档下载收益归作者所有

当前文档最多预览五页,下载文档查看全文
温馨提示:
1. 部分包含数学公式或PPT动画的文件,查看预览时可能会显示错乱或异常,文件下载后无此问题,请放心下载。
2. 本文档由用户上传,版权归属用户,天天文库负责整理代发布。如果您对本文档版权有争议请及时联系客服。
3. 下载前请仔细阅读文档内容,确认文档内容符合您的需求后进行下载,若出现内容与标题不符可向本站投诉处理。
4. 下载文档时可能由于网络波动等原因无法下载或下载错误,付费完成后未能成功下载的用户请联系客服处理。