使用 NumPy 和 matplotlib 绘制的图表上未显示某些点答案

【问题标题】：Some points are not displayed on the graph plotted using NumPy and matplotlib使用 NumPy 和 matplotlib 绘制的图表上未显示某些点
【发布时间】：2021-04-07 11:23:05
【问题描述】：

对于下面的代码，其工作是对函数 f 执行蒙特卡罗积分，我想知道如果我将 f 定义为 y = sqrt(1-x^2) 会发生什么，这是一个单位季度的方程圆，并指定一个大于 1 的端点，因为我们知道 f 仅定义为 0

import numpy as np
import matplotlib.pyplot as plt

def definite_integral_show(f, x0, x1, N):
    """Approximate the definite integral of f(x)dx between x0 and x1 using
    N random points
    
    Arguments:
    f -- a function of one real variable, must be nonnegative on [x0, x1]
    N -- the number of random points to use
    
    
    """
    #First, let's compute fmax. We do that by evaluating f(x) on a grid
    #of points between x0 and x1
    #This assumes that f is generally smooth. If it's not, we're in trouble!
    x = np.arange(x0, x1, 0.01)
    
    y = f(x)
    print(y)
    f_max = max(y)
    
    
    #Now, let's generate the random points. The x's should be between
    #x0 and x1, so we first create points beterrm 0 and (x1-x0), and 
    #then add x0
    #The y's should be between 0 and fmax
    #
    #                  0...(x1-x0)
    x_rand = x0 + np.random.random(N)*(x1-x0)
    print(x_rand)
    
    y_rand = 0 +  np.random.random(N)*f_max
    
    
    
    #Now, let's find the indices of the poitns above and below
    #the curve. That is, for points below the curve, let's find
    #   i s.t. y_rand[i] < f(x_rand)[i]
    #And for points above the curve, find
    #   i s.t. y_rand[i] >= f(x_rand)[i]
    ind_below = np.where(y_rand < f(x_rand))
    ind_above = np.where(y_rand >= f(x_rand))
    
    
    #Finally, let's display the results
    plt.plot(x, y, color = "red")
    pts_below = plt.scatter(x_rand[ind_below[0]], y_rand[ind_below[0]], color = "green")
    pts_above = plt.scatter(x_rand[ind_above[0]], y_rand[ind_above[0]], color = "blue")
    plt.legend((pts_below, pts_above),
            ('Pts below the curve', 'Pts above the curve'),
            loc='lower left',
            ncol=3,
            fontsize=8)

def f1(x):
    return np.sqrt(1-x**2)

definite_integral_show(f1, 0, 6, 200)

令我惊讶的是，该程序仍然有效，并为我提供了以下图片。

我怀疑它可以工作，因为在 NumPy 中，数组中的 nan 在对数组执行操作时会被忽略。但是，我不明白为什么图片只包含 x 和 y 坐标都在 0 到 1 之间的点。不在此范围内但其值由计算的点在哪里

x_rand = x0 + np.random.random(N)*(x1-x0)
y_rand = 0 +  np.random.random(N)*f_max

【问题讨论】：

标签： python arrays numpy matplotlib

【解决方案1】：

您可以只打印出数组（例如通过仅生成一个随机点）并查看它们既不进入ind_below 也不进入ind_above...

这是因为所有涉及nan 的比较都返回False。（另见：What is the rationale for all comparisons returning false for IEEE754 NaN values?）。（所以y_rand < nan 和y_rand >= nan 都计算为False）

更改代码最简单的方法是

ind_below = np.where(y_rand < f(x_rand))
ind_above = np.where(~(y_rand < f(x_rand)))

（可选地只计算一次数组）

【讨论】：