在 python 中使用 Monte Carlo 方法
Using Monte Carlo method in python
我正在处理下图中所写问题的第一个版本。我使用 rand
命令生成了一个随机点并测试了该点是否在圆圈内。我的代码是否接受 Monte Carlo 个函数作为输入值?我相信我选择了 N
,点数足够小,所以我不会 运行 内存不足。此外,当我 运行 这段代码时,它 运行 没有任何错误,但图表没有显示。在我可能出错的地方寻求帮助。
import numpy as np
import matplotlib.pyplot as plt
from random import random
xinside = []
yinside = []
xoutside = []
youtside = []
insidecircle = 0
totalpoints = 10**3
for _ in range(totalpoints):
x = random()
y = random()
if x**2+y**2 <= 1:
insidecircle += 1
xinside.append(x)
yinside.append(y)
else:
xoutside.append(x)
youtside.append(y)
fig, ax = plt.subplots()
ax.set_aspect('equal')
ax.scatter(xinside, yinside, color='g', marker='s')
ax.scatter(xoutside, youtside, color='r', marker='s')
fig.show()
没有显示的图表很神秘,也许试试plt.show()
。或者,您可以使用 savefig
保存绘图。这是代码第一部分的工作函数(只需修改问题中发布的代码)以及所需的输出图。
import numpy as np
import matplotlib.pyplot as plt
from random import random
def monte_carlo(n_points):
xin, yin, xout, yout = [[] for _ in range(4)] # Defining all 4 lists together
insidecircle = 0
for _ in range(n_points):
x = random()
y = random()
if x**2+y**2 <= 1:
insidecircle += 1
xin.append(x)
yin.append(y)
else:
xout.append(x)
yout.append(y)
print ("The estimated value of Pi for N = %d is %.4f" %(n_points, 4*insidecircle/n_points))
fig, ax = plt.subplots()
ax.set_aspect('equal')
ax.scatter(xin, yin, color='g', marker='o', s=4)
ax.scatter(xout, yout, color='r', marker='o', s=4)
plt.savefig('monte_carlo.png')
n_points = 10**4
monte_carlo(n_points)
> The estimated value of Pi for N = 10000 is 3.1380
矢量化方法 如果在函数中排除print
语句,则可以称为单行。时间分析留给你做作业
import numpy as np
import matplotlib.pyplot as plt
def monte_carlo(n_points, x, y):
pi_vec = 4*(x**2 + y**2 <= 1).sum()/n_points
print ("The estimated value of Pi for N = %d is %.4f" %(n_points, pi_vec))
# Generate points
n_points = 10**4
x = np.random.random(n_points)
y = np.random.random(n_points)
# x = [random() for _ in range(n_points)] # alternative way to define x
# y = [random() for _ in range(n_points)] # alternative way to define y
# Call function
monte_carlo(n_points, x, y)
> The estimated value of Pi for N = 10000 is 3.1284
或者,您可以通过使用单个 x 和 y 点数组来删除两个变量 x
和 y
,如下所示:
pts = np.random.uniform(0, 1, 2*n_points).reshape((n_points, 2))
monte_carlo(n_points, pts) # Function will be def monte_carlo(n_points, pts):
并在函数中使用
pi_vec = 4*(pts[:,0]**2 + pts[:,1]**2 <= 1).sum()/n_points
我正在处理下图中所写问题的第一个版本。我使用 rand
命令生成了一个随机点并测试了该点是否在圆圈内。我的代码是否接受 Monte Carlo 个函数作为输入值?我相信我选择了 N
,点数足够小,所以我不会 运行 内存不足。此外,当我 运行 这段代码时,它 运行 没有任何错误,但图表没有显示。在我可能出错的地方寻求帮助。
import numpy as np
import matplotlib.pyplot as plt
from random import random
xinside = []
yinside = []
xoutside = []
youtside = []
insidecircle = 0
totalpoints = 10**3
for _ in range(totalpoints):
x = random()
y = random()
if x**2+y**2 <= 1:
insidecircle += 1
xinside.append(x)
yinside.append(y)
else:
xoutside.append(x)
youtside.append(y)
fig, ax = plt.subplots()
ax.set_aspect('equal')
ax.scatter(xinside, yinside, color='g', marker='s')
ax.scatter(xoutside, youtside, color='r', marker='s')
fig.show()
没有显示的图表很神秘,也许试试plt.show()
。或者,您可以使用 savefig
保存绘图。这是代码第一部分的工作函数(只需修改问题中发布的代码)以及所需的输出图。
import numpy as np
import matplotlib.pyplot as plt
from random import random
def monte_carlo(n_points):
xin, yin, xout, yout = [[] for _ in range(4)] # Defining all 4 lists together
insidecircle = 0
for _ in range(n_points):
x = random()
y = random()
if x**2+y**2 <= 1:
insidecircle += 1
xin.append(x)
yin.append(y)
else:
xout.append(x)
yout.append(y)
print ("The estimated value of Pi for N = %d is %.4f" %(n_points, 4*insidecircle/n_points))
fig, ax = plt.subplots()
ax.set_aspect('equal')
ax.scatter(xin, yin, color='g', marker='o', s=4)
ax.scatter(xout, yout, color='r', marker='o', s=4)
plt.savefig('monte_carlo.png')
n_points = 10**4
monte_carlo(n_points)
> The estimated value of Pi for N = 10000 is 3.1380
矢量化方法 如果在函数中排除print
语句,则可以称为单行。时间分析留给你做作业
import numpy as np
import matplotlib.pyplot as plt
def monte_carlo(n_points, x, y):
pi_vec = 4*(x**2 + y**2 <= 1).sum()/n_points
print ("The estimated value of Pi for N = %d is %.4f" %(n_points, pi_vec))
# Generate points
n_points = 10**4
x = np.random.random(n_points)
y = np.random.random(n_points)
# x = [random() for _ in range(n_points)] # alternative way to define x
# y = [random() for _ in range(n_points)] # alternative way to define y
# Call function
monte_carlo(n_points, x, y)
> The estimated value of Pi for N = 10000 is 3.1284
或者,您可以通过使用单个 x 和 y 点数组来删除两个变量 x
和 y
,如下所示:
pts = np.random.uniform(0, 1, 2*n_points).reshape((n_points, 2))
monte_carlo(n_points, pts) # Function will be def monte_carlo(n_points, pts):
并在函数中使用
pi_vec = 4*(pts[:,0]**2 + pts[:,1]**2 <= 1).sum()/n_points