蒙特卡洛算法实验
这么看来蒙特卡洛方法的理论支撑其实是概率论或统计学中的大数定律。基本原理简单描述是先大量模拟,然后计算一
个事件发生的次数,再通过这个发生次数除以总模拟次数,得到想要的结果。下面我们以三个经典的小实验来学习下蒙
特卡洛算法思想。
1.计算圆周率 pi(π)值
实验原理:在正方形内部有一个相切的圆,圆面积/正方形面积之比是(PixRxR)/(2Rx2R)= Pi/4。在这个正方形内随机产
生 n 个点,假设点落在圆内的概率为 P,那么 P=圆面积/正方形面积,则 P= Pi/4。如何计算点落在圆内的概率 P?可以
计算点与中心点的距离,判断是否落在圆的内部,若这些点均匀分布,用 M 表示落到圆内投点数 , N 表示总的投点数,
则圆周率 Pi=4P=4xM/N。
实验步骤:
(1)将圆心设在原点(0,0),以 R 为半径形成圆,则圆面积为 PixRxR
(2)将该圆外接正方形, 坐标为(-R,-R)(R,-R)(R, R)(-R,R),则该正方形面积为 R*R
(3)随即取点(X,Y),使得-R <=X<=R 并且-R <=Y<=R,即点在正方形内
(4)通过公式 XxX+YxY<= RxR 判断点是否在圆周内(直角三角形边长公式)。
(5)设所有点(也就是实验次数)的个数为 N,落在圆内的点(满足步骤 4 的点)的个数为 M,则 P=M/N,于是
Pi=4xM/N。
(6)运行结果为 3.143052
def cal_pai_mc(n=1000000):
r = 1.0
a, b = (0.0, 0.0)
x_neg, x_pos = a - r, a + r
y_neg, y_pos = b - r, b + r
m = 0
for i in range(0, n+1):
x = random.uniform(x_neg, x_pos)
y = random.uniform(y_neg, y_pos)
if x**2 + y**2 <= 1.0:
m += 1
return (m / float(n)) * 4
2.计算函数定积分值