如果不考虑作图,这里的两个例子可以改写成下面的样子:
求圆周率
import random
'''
蒙特卡罗模拟
投点法计算圆周率
'''
# 投点游戏
def play_game():
# 圆
r = 1.0 # 半径
a, b = (0., 0.) # 圆心
# 正方形区域边界
x_min, x_max = a-r, a+r
y_min, y_max = b-r, b+r
# 在 正方形 区域内随机投点
x = random.uniform(x_min, x_max) # 均匀分布
y = random.uniform(y_min, y_max)
# 计算点到圆心距离
d = (x-a)**2 + (y-b)**2
# 根据落在圆内与否,返回1,0(为方便计数)
return [0, 1][d<r]
# 游戏次数
n = 100000
# 蒙特卡罗方法,模拟 n 次游戏
pi = 4.0 * sum((play_game() for _ in range(n))) / n # 4 倍频率近似于圆周率
print('pi: ', pi)
求定积分
import random
'''
蒙特卡罗模拟
投点法计算函数 y=x^2在[0,1]内的定积分
'''
# 函数 y=x^2
def f(x):
return x**2
# 投点游戏
def play_game():
# 矩形区域边界
x_min, x_max = 0, 1
y_min, y_max = 0, 1
# 在 矩形 区域内随机投点
x = random.uniform(x_min, x_max) # 均匀分布
y = random.uniform(y_min, y_max)
# 根据点落在函数 y=x^2图像下方与否,返回1,0(为方便计数)
return [0, 1][y<f(x)]
# 游戏次数
n = 100000
# 蒙特卡罗方法,模拟 n 次游戏
integral = sum((play_game() for _ in range(n))) / n # 频率近似于定积分
print('integral: ', integral)