Numpy向量化嵌套的’for’循环

huangapple go评论57阅读模式
英文:

Numpy Vectorization for Nested 'for' loop

问题

我尝试编写一个程序,用于绘制任何给定函数的等值线。

```python
rmin = -5.0
rmax = 5.0
c = 4.0
x = np.arange(rmin, rmax, 0.1)
y = np.arange(rmin, rmax, 0.1)
x, y = np.meshgrid(x, y)
f = lambda x, y: y**2.0 - 4*x
realplots = []
for i in range(x.shape[0]):
    for j in range(x.shape[1]):
        if abs(f(x[i, j], y[i, j]) - c) < 1e-4:
            realplots.append([x[i, j], y[i, j]])

但由于它是一个嵌套的for循环,所以需要很长时间。非常感谢任何帮助,可以将上述代码向量化或采用新的等值线绘制方法。(注意:在运行时,函数'f'将发生变化。因此,必须在不考虑函数属性的情况下进行向量化)

我尝试通过以下方式向量化:

ans = np.where(abs(f(x, y) - c) < 1e-4, np.array([x, y]), [0, 0])

但它给出了错误消息“operands could not be broadcast together with shapes (100,100) (2,100,100) (2,)”,我在np.where中添加了[0,0]作为对else条件的一种逃避,这是错误的。


<details>
<summary>英文:</summary>

I was trying to write a program which plots level set for any given function.
```lang-python
rmin = -5.0
rmax = 5.0
c = 4.0
x = np.arange(rmin,rmax,0.1)
y = np.arange(rmin,rmax,0.1)
x,y = np.meshgrid(x,y)
f = lambda x,y: y**2.0 - 4*x
realplots = []
for i in range(x.shape[0]):
    for j in range(x.shape[1]):
        if abs(f(x[i,j],y[i,j])-c)&lt; 1e-4:
            realplots.append([x[i,j],y[i,j]])`

But it being a nested for loop, is taking lot of time. Any help in vectorizing the above code/new method of plotting level set is highly appreciated.(Note: The function 'f' will be changed at the time of running.So, the vectorization must be done without considering the function's properties)

I tried vectorizing through
ans = np.where(abs(f(x,y)-c)&lt;1e-4,np.array([x,y]),[0,0])
but it was giving me operands could not be broadcast together with shapes (100,100) (2,100,100) (2,)
I was adding [0,0] as an escape from else condition in np.where which is indeed wrong.

答案1

得分: 2

你可以直接使用掩码来索引 xy,参考“布尔数组索引”部分的文档。

def vectorized(x, y, c, f, threshold):
    mask = np.abs(f(x, y) - c) < threshold
    x, y = x[mask], y[mask]
    return np.stack([x, y], axis=-1)

您的参考函数:

def op(x, y, c, f, threshold):
    res = []
    for i in range(x.shape[0]):
        for j in range(x.shape[1]):
            if abs(f(x[i, j], y[i, j]) - c) < threshold:
                res.append([x[i, j], y[i, j]])
    return res

测试:

rmin, rmax = -5.0, +5.0
c = 4.0
threshold = 1e-4

x = np.arange(rmin, rmax, 0.1)
y = np.arange(rmin, rmax, 0.1)
x, y = np.meshgrid(x, y)

f = lambda x, y: y**2 - 4 * x

res_op = op(x, y, c, f, threshold)
res_vec = vectorized(x, y, c, f, threshold)

assert np.allclose(res_op, res_vec)
英文:

Since you get the values rather than the indexes, you don't really need np.where.

You can directly use the mask to index x and y, look at the "Boolean array indexing" section of the documentation.

It is straightforward:

def vectorized(x, y, c, f, threshold):
    mask = np.abs(f(x, y) - c) &lt; threshold
    x, y = x[mask], y[mask]
    return np.stack([x, y], axis=-1)

Your function for reference:

def op(x, y, c, f, threshold):
    res = []
    for i in range(x.shape[0]):
        for j in range(x.shape[1]):
            if abs(f(x[i, j], y[i, j]) - c) &lt; threshold:
                res.append([x[i, j], y[i, j]])
    return res

Tests:

rmin, rmax = -5.0, +5.0
c = 4.0
threshold = 1e-4

x = np.arange(rmin, rmax, 0.1)
y = np.arange(rmin, rmax, 0.1)
x, y = np.meshgrid(x, y)

f = lambda x, y: y**2 - 4 * x

res_op = op(x, y, c, f, threshold)
res_vec = vectorized(x, y, c, f, threshold)

assert np.allclose(res_op, res_vec)

huangapple
  • 本文由 发表于 2023年2月10日 12:02:21
  • 转载请务必保留本文链接:https://go.coder-hub.com/75406815.html
匿名

发表评论

匿名网友

:?: :razz: :sad: :evil: :!: :smile: :oops: :grin: :eek: :shock: :???: :cool: :lol: :mad: :twisted: :roll: :wink: :idea: :arrow: :neutral: :cry: :mrgreen:

确定