蒙特卡洛估算定积分的平均值法的Python实现

蒙特卡洛估算定积分的第二种方法

平均值法

前几天，我们利用蒙特卡洛的随机投点法实现了y=x^2在到1上的定积分的估算（传送门），今天，我们介绍另一种蒙特卡洛估算定积分的方法：平均值法

原理：辛钦大数定律

这里的大数定律，指的是伯努利大数定律，即

绝对值里面的前者为平均观测值，后者为概率（由二点分布的的公式得到），这就是说，当n很大时，后者（概率）可以由前者（平均观测值）近似代替。而后者又可以写成期望的形式，连续型随机变量求期望就是求对应的定积分，我们由此来估算定积分

具体的推导证明不再叙述（因为我也没看太懂），原理有点复杂，但是我们用代码来实现起来还是挺容易的。

直接上代码，还是以y=x^2为例

# -*- coding: cp936 -*-

import time

import numpy as np

def f(x):

return x**2

def n(N):

start_time=time.time()

a=1/3.0

x=np.random.uniform(0,1,N)#随机生成N个0-1之间的一维点

c=f(x)#把x的值代入f（x）计算

s=0

for i in c:

s=s+i

j=s*1.00/N

print"蒙特卡洛估计值为： ",j

print"与真实值之间的误差为：",abs(a-j)

end_time=time.time()

clap=end_time-start_time

result=round(clap,3)

print"程序运行耗时： ",result,"秒"

下面进行效果测试：

当N取100000000时，误差已经非常小了，只是计算时间变成了80多秒，有点耗时，但，毕竟，欲速则不达，对吧。

从开学到现在感觉还没有睡醒，头脑昏昏沉沉的，错误之处，还望指出，睡觉了，晚安。