蒙特卡洛估算定积分的第二种方法
平均值法
前几天,我们利用蒙特卡洛的随机投点法实现了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多秒,有点耗时,但,毕竟,欲速则不达,对吧。
从开学到现在感觉还没有睡醒,头脑昏昏沉沉的,错误之处,还望指出,睡觉了,晚安。
领取专属 10元无门槛券
私享最新 技术干货