概率论与数理统计实验报告1
- 格式：doc
- 大小：61.00 KB
- 文档页数：9
概率论与数理统计实验报告
实验题目：蒙特卡洛算法计算积分
实验时间：2012.06.01
姓名：王文栋
学号：2110904023
班级：物理试验班12
实验报告
一．实验目的
1. 初步了解蒙特卡洛算法，以及用其计算一些高等数学中不能直接计算出的积分；
2. 计算出的真值与蒙特卡洛法得值的差值，比较其有效性。
二．实验原理
1. 蒙特卡洛法的思想简述
当我们所求解问题是某种随机事件出现的概率，或者是某个随机变量的期望值时，通过某种“实验”的方法，以这种事件出现的频率估计这一随机事件的概率，或者得到这个随机变量的某些数字特征，并将其作为问题的解。 有一个例子我们可以比较直观地了解蒙特卡洛方法：假设我们要计算一个不规则图形的面积，那么图形的不规则程度和分析性计算（比如，积分）的复杂程度是成正比的。蒙特卡洛方法是如下计算的：假想有一袋豆子，把豆子均匀地朝这个图形上撒，然后数这个图形之中有多少颗豆子，这个豆子的数目就是图形的面积。当豆子越小，撒的越多的时候，结果就越精确。在这里我们要假定豆子都在一个平面上，相互之间没有重叠。
2. 蒙特卡洛法与积分
通常蒙特卡洛方法通过构造符合一定规则的随机数来解决数学上的各种问题。对于那些由于计算过于复杂而难以得到解析解或者根本没有解析解的问题，蒙特卡洛方法是一种有效的求出数值解的方法。一般蒙特卡洛方法在数学中最常见的应用就是蒙特卡洛积分。非权重蒙特卡洛积分，也称确定性抽样，是对被积函数变量区间进行随机均匀抽样，然后对被抽样点的函数值求平均，从而可以得到函数积分的近似值。此种方法的正确性是基于概率论的中心极限定理。
3. 本实验原理简述
在本实验中，我们主要是计算积分值与误差比较。在计算积分时，我们要选择合适的变量分布，其中有均匀分布，有正态分布，要视情况而选择。在利用蒙特卡洛方法计算积分时，我们要分情况。
①对于积分为这种形式，我们可以转化为这种形式，然后利用其等于（b-a）E(x)的计算结果。E(x)可利用求随机变量的均值来得到。
②对于积分为这种形式的，我们依然可以利用上面的算法计算，将其化为即可得。E(x，y)可利用求随机变量的均值来得到。
对于有无穷的区间，我们应选择正态分布。
三．实验内容
1.计算下面的积分，并与真值比较
①20sinxxdx ②20xedx ③.22221xyxyedxdy
1.计算下面的积分，并与平均值比较计算方差
④.231xedx ⑤104011dxx ⑥4011dxx ⑦2222411xydxdyxy
计算过程及程序如下：
①20sinxxdx
n=10000;m=10;S=0;I=0;
r=unifrnd(0,pi/2,m,n)
for j=1:m
s=0;
for i=1:n
s=s+r(j,i)*sin(r(j,i));
end
S(j)=pi/2*s/n
end
D=0;d=0;
for j=1:m
D=D+(S(j)-1)^2;
end
d=D/(m-1)
结果S =1.009 1.0023 0.9898 0.9859 0.9948 0.9828
0.9963 1.0049 0.9943 0.9924
d =8.5885e-005
②20xedx
n=10000;m=10;S=0;sum=0；I=0;
r=normrnd(0,1,m,n);
for i=1:m
s=0;
for j=1:n
s=s+sqrt(pi/2)*exp(-0.5*r(i,j)^2);
end
S(i)=s/n;sum=sum+S(i);
end
I=sum/m；D=0;d=0;
for i=1:m
D=D+(S(i)-I)^2;
end
d=D/(m-1）
结果I =0.8851；d =3.6034e-006
③.22221xyxyedxdy
m=10000;sum=0;n=50;D=0;
X=unifrnd(-1,1,n,m);Y=unifrnd(-1,1,n,m);
for i=1:n
a=0;
for j=1:m
if(X(i,j)^2+Y(i,j)^2<=1)
Z(i,j)=exp(X(i,j)^2+Y(i,j)^2);
a=a+Z(i,j);
end
end
S(i)=a/m;sum=sum+S(i);
end
I=sum/n*4
for i=1:n
D=D+(S(i)*4-pi*(exp(1)-1))^2;
end
d=D/n
结果I =5.4089；d =0.0011
④.231xedx
I=0;S=0;sum=0;m=1000n=10000
r=unifrnd(1,3,m,n)
for j=1:m
s=0;
for i=1:n
s=s+2*exp(r(j,i)^2);
end
S(i)=s/n；sum=sum+S(i);
end
I=sum/m
for j=1:m
D=D+(S(i)-I)^2;
end
d=D/(m-1)
结果I =1.4433e+003；d = 53.7036
⑤104011dxx
n=10000;m=50;sum=0;I=0;S=0;
r=unifrnd(0,10,m,n);
for j=1:m
s=0;
for i=1:n
s=s+10/sqrt(1+r(j,i)^4);
end
S(j)=s/n； sum=sum+S(j);
end
I=sum/m；D=0;d=0;
for j=1:m
D=D+(I-S(j))^2;
end
d=D/(m-1)
结果I =1.7505 d =7.6373e-004
⑥4011dxx
n=10000;m=10;S=0;sum=0;I=0;
r=normrnd(0,1,m,n);
for j=1:m
s=0;
for i=1:n
s=s+sqrt(pi/(2*(1+r(j,i)^4)))*exp(0.5*r(j,i)^2);
end
S(j)=s/n；sum=sum+S(j);
end
I=sum/m；d=0;D=0;
for j=1:m
D=D+(I-S(j))^2;
end
d=D/(m-1)
I =1.6619；d = 0.0088
⑦2222411xydxdyxy
n=1000;m=100;sum=0;S=0;I=0;
x=unifrnd(-2,2,m,n);y=unifrnd(-2,2,m,n);
for j=1:m
s=0;
for i=1:n
if x(j,i)^2+y(j,i)^2<=4
s=s+16/sqrt(1+x(j,i)^4+y(j,i)^2);
end
end
S(j)=s/n;sum=sum+S(j);
end
I=sum/m；D=0;d=0;
for j=1:m
D=D+(I-S(j))^2;
end
d=D/(m-1)
结果I =7.4531；d =0.0188
四．实验心得
通过本次试验，我掌握了用蒙特卡洛方法计算积分。在此之外，还了解了蒙特卡洛方法的一些其他用法。
