蒙特.卡罗方法求解圆周率
背景
蒙特卡洛是一个地名,位于赌城摩纳哥,象征概率。蒙特卡洛(Monte Carlo)方法是由大名鼎鼎的数学家冯·诺伊曼提出的,诞生于上世纪40年代美国的“曼哈顿计划”。原理是通过大量随机样本,去了解一个系统,进而得到所要计算的值。
原理
对于某些不能精确求解的问题,蒙特.卡罗方法是一种非常巧妙的寻找近似解的方法。
以求解圆周率的问题为例,假设有一个单位圆及其外切正方形,我们往正方形内扔飞镖,当扔的次数足够多以后,“落在圆内的次数/落在正方形内的次数
”这个比值会无限接近“圆的面积/正方形的面积
”这个比值,也就是圆周率的四分之一
。模拟扔飞镖的次数越多,圆周率的近似结果越精确。
故:“落在圆内的次数/落在正方形内的次数
”这个比值会无限接近“圆的面积/正方形的面积
”这个比值,也就是圆周率的四分之一
。
思路:
当我们在(0,1)的范围内随机选择一个坐标(x, y)时,每个坐标点被选中的概率相等。则坐标落在直径为1的正方形中的圆的概率为:
#1
from random import random
from time import perf_counter
DARTS = 1000*1000
hits = 0.0 # 落入圆的次数
start = perf_counter()
for i in range(0, DARTS):
x, y = random(), random() #random()返回随机生成的一个实数,它在[0,1)范围内。
dist = x ** 2 + y ** 2
if dist <= 1.0:
hits = hits + 1
pi = 4 * (hits/DARTS)
print("圆周率值是: {}".format(pi))
print("运行时间是: {:.5f}s".format(perf_counter() - start))
# 2
import numpy as np
from time import perf_counter
def pi(n):
batch = 1000
t = 0
for i in range(n // batch):
# 随机获取坐标(batch个点)
# 以给定的形状创建一个数组(batch*2),并在数组中加入在[0,1)之间均匀分布的随机样本。
p = np.random.rand(batch, 2)
# 计算坐标平方和
p = (p * p).sum(axis=1) # 数组元素对位相乘,最外围axis=0
# 平方和小于1的即为落在圆中的点
t += (p <= 1).sum()
return 4 * t / n
start = perf_counter()
Pi = pi(10 ** 8)
print(Pi)
print("运行时间是: {:.5f}s".format(perf_counter() - start))