蒙特卡洛估算定积分的平均值法的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多秒,有点耗时,但,毕竟,欲速则不达,对吧。

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

  • 发表于:
  • 原文链接http://kuaibao.qq.com/s/20180307G1QVAB00?refer=cp_1026
  • 腾讯「腾讯云开发者社区」是腾讯内容开放平台帐号(企鹅号)传播渠道之一,根据《腾讯内容开放平台服务协议》转载发布内容。
  • 如有侵权,请联系 cloudcommunity@tencent.com 删除。

扫码关注腾讯云开发者

领取腾讯云代金券