700字范文,内容丰富有趣,生活中的好帮手!
700字范文 > 蒙特卡洛法求概率密度函数matlab matlab的蒙特卡洛算法-----用概率算法求pi值

蒙特卡洛法求概率密度函数matlab matlab的蒙特卡洛算法-----用概率算法求pi值

时间:2022-10-23 09:20:15

相关推荐

蒙特卡洛法求概率密度函数matlab matlab的蒙特卡洛算法-----用概率算法求pi值

问题补充:

能提供一个例子看看吗?我就不懂MATLAB,想知道具体代码。因为现在写论文急,也没时间仔细看书了

蒙特卡罗模拟

就是随机数相关的东西,你只要知道随机数是怎么得到。其它的事就要好办了。

rand(m,n)产生m*n均匀随机数。

ex:

用概率方法求pi

N=100000;

x=rand(N,1);

y=rand(N,1);

count=0;

for i=1:N

if (x(i)^2 y(i)^2<=1)

count=count 1;

end

end

PI=4*count/N

-------------------------------

有关使用matlab进行蒙特卡罗模拟的程序问题

检举|-10-1

15:37

提问者:轩fang|浏览次数:1425次

已经分布是均匀分布(连续),区间为(12,62),请问各位大侠,如何用matlab编程实现此蒙特卡罗模拟,我想模拟2000次,得到概率密度图与累积概率密度图,程序应该如何编,麻烦大家指教,期待帮助,谢谢。小妹现在没有分数,还是希望大家能帮我,谢谢1

问题补充:

我的邮箱是512608029@,谢谢大家的帮忙。n=2000;

%随机点数(可增加点数) x=12 (62-12)*rand(1,n); %产生2000个12到62的随机数 xx=12:2:62;

%画概率密度图的区间 nx=histc(x,xx); %计算x在xx每个小区间内的点数。 px=nx/n;

sumpx=cumsum(px); subplot(1,2,1) bar(xx(1:end-1),px(1:end-1));

title('概率密度') subplot(1,2,2) plot(xx(1:end-1),sumpx(1:end-1));

title('累积概率密度')-----------------------------------

我这儿正好有份程序,希望有所帮助。

代码:

一份蒙特卡洛程序

count=input('input the count:');%输入模拟粒子数

sigmaedata=[28370,13845,6908,2555,1223,2602,1925,905.3,479.7,164.5,74.24,23.86,66.60,36.62,22.29,9.978,5.298,1.668,0.7378,0.2361,0.1099,0.06211,0.03939,0.02030];

sigmacdata=[0.0220,0.0393,0.0568,0.0904,0.1209,0.1479,0.1722,0.2136,0.2480,0.3092,0.3486,0.3932,0.4153,0.4268,0.4319,0.4291,0.4215,0.3969,0.3691,0.3269,0.2944,0.2709,0.2512,0.2209];

Edata=[1,1.5,2,3,4,5,6,8,10,15,20,30,40,50,60,80,100,150,200,300,400,500,600,800];%截面数据

channel=zeros(1,ceil(662/5) 10);%多道数组

nget=0;%探测到的总计数

ntotal=0;%进入探测器的总计数

for ii=1:count%count个粒子循环

collidetime=0;%当前粒子碰撞次数

%粒子状态初始化

E0=622;

E=E0;

z=-2;

r=0;

theta=2*pi*rand(1);% 源抽样,z,r,theta坐标

miu=2*rand(1)-1;

fai=2*pi*rand(1);%方向角抽样

if miu

continue;else

z=0;

r=2*sqrt(1-miu^2)/miu;

theta=fai;

end

while E>1%一个粒子在闪烁体中的输运过程

sigmae=interp1(Edata,sigmaedata,E,'linear');

sigmac=interp1(Edata,sigmacdata,E,'linear');

sigmat=sigmae sigmac;%线性插值得到截面数据

L=-log(rand(1))/sigmat;%下次碰撞的距离

%计算下次碰撞位置坐标

rnew=sqrt(r^2 L^2*(1-miu^2)

2*r*L*sqrt(1-miu^2)*cos(fai-theta));

z=z L*miu;

cdth=(rnew^2 r^2-L^2*(1-miu^2))/2/r/rnew;

sdth=L*sqrt(1-miu^2)*sin(fai-theta)/rnew;

dtheta=asin(sdth);

if

cdth<0

dtheta=pi-dtheta;

end

theta=theta dtheta;

r=rnew;

if(r>2)|(z>=4)|(z<0)%判断是否在闪烁体内

break;

else

collidetime=collidetime

1;

end

if rand(1)

E=0;else%康普顿散射

alpha=E/511;

flag=0;

while

flag==0

if

rand(1)<=27/(4*alpha 29)

x=(1

2*alpha)/(1 2*alpha*rand(1));

if

rand(1)<=0.5*((alpha 1-x/alpha)^2 1)

flag=1;

end

else

x=1

2*alpha*rand(1);

if

rand(1)<=27/4*((x-1)^2)/x^3

flag=1;

endend

end

本内容不代表本网观点和政治立场,如有侵犯你的权益请联系我们处理。
网友评论
网友评论仅供其表达个人看法,并不表明网站立场。