Matlab模拟蒲丰投针实验计算Π值

1.前言

学习了概率论，决定自己动手编程实现蒲丰投针实验来计算$\pi$的值，其中本文代码部分参考了这篇文章.

2.实验目的

通过仿真模拟实现蒲丰投针实验来计算 $\pi$ 值.

3.实验原理

(1)问题描述

平面上画着许多间距为$a$的平行线，随机地向此平面任投掷一长度为$l(l<a)$的针$N$次,观察针与直线相交的次数记为$n$，计算此针与任意直线相交的概率。

(2)问题求解

以$x$表示针的中点到最近的一条平行线的距离，$\varphi$表示针与平行线的夹角。可知$0\le x\le \frac{a}{2},0\le \varphi \le \pi ,$以$S$表示长宽分别为$\frac{a}{2}$和$\pi$的长方形，要使针与平行线相交，必须保证$x\le \frac{l}{2}sin\varphi$，满足这个关系式的区域$G$为下图(图中$b$即为$l$)中黄色部分，则所求概率即表示为$p=\frac{G的面积}{S的面积}=\frac{\frac{1}{2}{\int }_0^{\pi }l\sin \varphi d\varphi }{\frac{1}{2}a\pi }=\frac{2l}{\pi a}$.用频率$\frac{n}{N}$代替概率$p$即可得$\pi =\frac{2lN}{an}$,并以此来计算$\pi$值.

(图来源网络，侵删)

4.实验过程

根据实验原理用Matlab编写程序，通过蒙特卡洛法生成规定量的随机数模拟蒲丰投针实验的结果，通过不断增加试验次数来得到更加精确结果.

5.实验结果

	$a$	$n$	$Pi$
第一次	1	10000	3.148779847809
第二次	1	10000	3.147128245476
第三次	1	100000	3.149440974227
第四次	1	100000	3.140539125883
第五次	1	1000000	3.144044393907
第六次	1	1000000	3.145742237881
第七次	1	10000000	3.139098147566
第八次	1	10000000	3.142210146825
第九次	1	20000000	3.140213681074
第十次	1	30000000	3.141194329202

6.实验结论

由实验结果可以看出，可以看出当$n$足够大时，频率接近概率，蒲丰投针模拟实验可以得到较为精确的$\pi$值。

7.Matlab代码

l = 0.6; %针长度
a = 1;   %平行线宽度,l和a均可修改，只要满足l<a即可
n = 100000; %做n次投针试验
y = unifrnd(0, a/2, 1, n); %在[0, a/2]内服从均匀分布随机产生n个数
theta = unifrnd(0, pi, 1, n); %在[0, pi]内服从均匀分布随机产生n个数
k=sum(y<(l/2)*sin(theta));%记录针与平行线相交的次数
pival = (2*l*n)/(a*k);
z=sprintf('%.12f',pival);
disp(z)