Python-蒙特卡罗计算圆周率

蒙特卡洛法

        蒙特卡罗方法于20世纪40年代美国在第二次世界大战中研制原子弹的“曼哈顿计划”计划的成员S.M.乌拉姆和J.冯·诺伊曼首先提出。数学家冯·诺伊曼用驰名世界的赌城—摩纳哥的Monte Carlo—来命名这种方法，为它蒙上了一层神秘色彩。在这之前，蒙特卡罗方法就已经存在。1777年，法国数学家布丰（Georges Louis Leclere de Buffon，1707—1788）提出用投针实验的方法求圆周率π。这被认为是蒙特卡罗方法的起源。

简单来说：就是一个正方形框里有一个内切圆（如下图），往框里随机投放大量的点，利用落在圆里的点的数目来求出圆周率

 我们知道，正方形和内切圆面积之比的4：π（不信你可以算一下），所以圆周率=落在圆内的随机点数*4/投点次数

                                                   

 

不多说，我们直接上代码！！！

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.patches import  Circle

#投点次数
n=10000

#圆的基本信息
r=1.0
a,b=(0.,0.)

#随机点的正方形区域边界
x_min,x_max=a-r,a+r
y_min,y_max=b-r,b+r

#开始投放随机点
x=np.random.uniform(x_min,x_max,n)  #均匀分布，左闭右开
y=np.random.uniform(y_min,y_max,n)

#计算随机点到圆心的距离
d=np.sqrt((x-a)**2+(y-b)**2)

#统计落在单位圆的随机点数目
count=sum(np.where(d<r,1,0))

#蒙特卡罗求圆周率
pi=4*count/n

print(pi)

#图形表示
fig=plt.figure()
axes=fig.add_subplot(111)
axes.plot(x,y,'ro',markersize=1)
plt.axis('equal')   #防止图形变形

circle=Circle(xy=(a,b),radius=r,alpha=0.5)
axes.add_patch(circle)

plt.show()

结果显示：

                       

可能有些语句函数会看不懂，下面我来简单的说一下，仅供参考：

1）random.uniform(low,high,size)：在区间[low,high)内随机生成size个样本数目，需要注意的是uniform不能单独使用，需要通过random库调用

2）np.where(condition,x,y)：根据条件生成新的数组，condition是条件，符合条件的就是x，不符合条件的就等于y，所以np.where(d<r,1,1)意思就是落在圆内的随机点记为1，反之为0，再用sum()就可以直接算出落在圆内的随机点的数目了

3）add_subplot(111)：意思是把画布分割成一行一列，图形在第一个位置，创建完figure对象后不能直接画图

4）plt.axis('equal')：图形的x轴和y轴的范围相等

5）add_patch()：在之前的111位置上添加一个东西，这里添加的是圆形

小结一下画图的步骤：

（1）创建figure对象

（2）创建一个或多个subplot子图

（3）在子图上绘画图形

（4）展示图形

 

总结：蒙特卡罗法就是通过多次投放随机点，模拟概率，

参考来源：https://www.cnblogs.com/hhh5460/p/6713287.html