蒙特卡洛方法：Python实现

蒙特卡洛方法：Python实现

上期已介绍了蒙特卡洛方法，在这篇博客中，将使用Python实现蒙特卡洛算法，并逐步解释代码的作用。蒙特卡洛算法的一个经典应用是估算圆周率 π 的值。将通过在单位正方形内随机生成点，并统计落在单位圆内的点的比例来估算 π。

1. case 1

步骤1：导入必要的库

首先，我们需要导入一些必要的库，包括numpy用于数学计算和matplotlib用于可视化。

import numpy as np
import matplotlib.pyplot as plt

步骤2：生成随机点

我们在单位正方形内生成大量随机点，这些点的 x 和 y 坐标都在 [0, 1] 范围内。

# 生成随机点数量
num_points = 1000

# 生成随机点
points = np.random.rand(num_points, 2)

步骤3：判断点是否在单位圆内

通过判断点到原点的距离是否小于1，我们可以确定这个点是否在单位圆内。

# 计算每个点到原点的距离
distances = np.linalg.norm(points, axis=1)

# 判断点是否在单位圆内
inside_circle = distances < 1

步骤4：可视化生成的点

为了更直观地了解我们生成的随机点的分布情况，我们将它们可视化出来。

# 可视化生成的点
plt.scatter(points[:, 0], points[:, 1], c=inside_circle, cmap='viridis')
plt.title('Random Points in Unit Square')
plt.xlabel('X-axis')
plt.ylabel('Y-axis')
plt.show()

步骤5：估算圆周率 π

最后，我们将统计落在单位圆内的点的比例，并用这个比例估算圆周率 π。

# 统计落在单位圆内的点的数量
inside_circle_count = np.sum(inside_circle)

# 估算圆周率
estimated_pi = 4 * inside_circle_count / num_points

print(f'估算的圆周率 π 为: {estimated_pi}')

通过运行这个简单的蒙特卡洛模拟，可以得到一个对圆周率 π 的估算值。这个例子展示了蒙特卡洛算法的基本原理和实现过程。在实际应用中，可以根据问题的特点和需求进行更复杂的模拟和计算。

case-2

实现了蒙特卡洛算法估算圆周率和曲线下面积的完整代码：

# @copyright all reseved
# @author: Persist_Zhang
import numpy as np
import matplotlib.pyplot as plt

# 生成随机点数量
num_points = 10000

# 生成随机点
points = np.random.rand(num_points, 2)

# 计算每个点到原点的距离
distances = np.linalg.norm(points, axis=1)

# 判断点是否在单位圆内
inside_circle = distances < 1

# 可视化生成的点
plt.scatter(points[:, 0], points[:, 1], c=inside_circle, cmap='viridis')
plt.title('Random Points in Unit Square')
plt.xlabel('X-axis')
plt.ylabel('Y-axis')
plt.show()

# 统计落在单位圆内的点的数量
inside_circle_count = np.sum(inside_circle)

# 估算圆周率
estimated_pi = 4 * inside_circle_count / num_points

print(f'估算的圆周率 π 为: {estimated_pi}')

# 示例：使用蒙特卡洛方法估算曲线下面积

# 定义复杂曲线的函数
def complex_curve(x):
    return np.sin(x) + np.cos(2 * x)

# 生成随机点在曲线下方的比例
random_x = np.random.rand(num_points) * 10  # 假设曲线范围在 [0, 10]
random_y = np.random.rand(num_points) * 2 - 1  # 假设曲线高度在 [-1, 1]
below_curve = random_y < complex_curve(random_x)

# 估算曲线下面积的比例
estimated_area_ratio = np.sum(below_curve) / num_points

# 实际曲线下面积
actual_area = np.trapz(complex_curve(np.linspace(0, 10, 1000)), dx=0.01)

# 估算曲线下面积
estimated_area = 10 * 2 * estimated_area_ratio

print(f'实际曲线下面积: {actual_area}')
print(f'估算曲线下面积: {estimated_area}')

小结

蒙特卡洛算法是一种基于随机抽样的数值计算方法，通过在随机样本上进行统计分析，来估计数学问题的解。在本文中，我们通过两个实例演示了蒙特卡洛算法的应用：估算圆周率和估算曲线下面积。

对于圆周率的估算，我们在单位正方形内生成大量随机点，并通过统计落在单位圆内的点的比例来估算圆周率。这展示了蒙特卡洛算法在几何问题中的简单应用。

另外，我们还介绍了蒙特卡洛方法用于估算曲线下面积的例子。通过生成随机点并判断是否在曲线下方，可以估算复杂曲线的面积，而无需求解积分。

总的来说，蒙特卡洛算法具有通用性强、易于理解和实现的特点。它在处理高维、非线性或难以求解的问题时表现出色。然而，也需要注意到蒙特卡洛算法的收敛速度相对较慢，计算量可能较大，因此在实际应用中需要权衡选择。