如果满足另一个数组中值的条件,则对 numpy 数组中的值求和
Sum values from numpy array if condition on value in another array is met
我在对函数进行矢量化以使其有效应用于 numpy 数组时遇到问题。
我的节目条目:
- A pos_part Nb_particles 行的二维数组,3 列(基本上是 x、y、z坐标,只有z与困扰我的部分有关)Nb_particles可以达到几十万。
- 具有 Nb_particles 值的 prop_part 一维数组。这部分我已经介绍过了,创建是用一些很好的 numpy 函数完成的;我只是在这里放了一个与实际值相符的基本分布。
- A z_distances 一维数组,z=0 和 z=z_max.
之间的简单 np.arange
然后是需要时间的计算,因为我无法找到一种方法来仅通过数组的 numpy 操作来正确地做事。我想要做的是:
- 对于z_distances中的所有距离z_i,对[中的所有值求和=81=] 如果相应的粒子坐标 z_particle < z_i。这将 return 与 z_distances.
长度相同的一维数组
我目前的想法:
- 版本 0,for 循环,enumerate 和 np.where 检索我需要求和的值的索引。明明挺长的。
- 版本 1,在新数组上使用掩码(z 坐标和粒子属性的组合),并对掩码数组求和。似乎比 v0
- 版本 2,另一个掩码和 np.vectorize,但我知道它效率不高,因为矢量化基本上是 for环形。似乎仍然比 v0
- 第 3 版,我正在尝试在可以直接应用于 z_distances 的函数上使用掩码,但它目前无法正常工作。
所以,我来了。可能与排序和累计总和有关,但我不知道该怎么做,所以任何帮助将不胜感激。请在下面找到代码以使事情更清楚
提前致谢。
import numpy as np
import time
import matplotlib.pyplot as plt
# Creation of particles' positions
Nb_part = 150_000
pos_part = 10*np.random.rand(Nb_part,3)
pos_part[:,0] = pos_part[:,1] = 0
#usefull property creation
beta = 1/1.5
prop_part = (1/beta)*np.exp(-pos_part[:,2]/beta)
z_distances = np.arange(0,10,0.1)
#my version 0
t0=time.time()
result = np.empty(len(z_distances))
for index_dist, val_dist in enumerate(z_distances):
positions = np.where(pos_part[:,2]<val_dist)[0]
result[index_dist] = sum(prop_part[i] for i in positions)
print("v0 :",time.time()-t0)
#A graph to help understand
plt.figure()
plt.plot(z_distances,result, c="red")
plt.ylabel("Sum of particles' usefull property for particles with z-pos<d")
plt.xlabel("d")
#version 1 ??
t1=time.time()
combi = np.column_stack((pos_part[:,2],prop_part))
result2 = np.empty(len(z_distances))
for index_dist, val_dist in enumerate(z_distances):
mask = (combi[:,0]<val_dist)
result2[index_dist]=sum(combi[:,1][mask])
print("v1 :",time.time()-t1)
plt.plot(z_distances,result2, c="blue")
#version 2
t2=time.time()
def themask(a):
mask = (combi[:,0]<a)
return sum(combi[:,1][mask])
thefunc = np.vectorize(themask)
result3 = thefunc(z_distances)
print("v2 :",time.time()-t2)
plt.plot(z_distances,result3, c="green")
### This does not work so far
# version 3
# =============================
# t3=time.time()
# def thesum(a):
# mask = combi[combi[:,0]<a]
# return sum(mask[:,1])
# result4 = thesum(z_distances)
# print("v3 :",time.time()-t3)
# =============================
通过完全用 numpy 编写第一个版本,您可以获得更多的性能。将 pythons sum
替换为 np.sum
。而不是 for i in positions
列表理解,只需传递您正在创建的 positions
掩码。
事实上,np.where
不是必需的,我最好的版本是这样的:
#my version 0
t0=time.time()
result = np.empty(len(z_distances))
for index_dist, val_dist in enumerate(z_distances):
positions = pos_part[:, 2] < val_dist
result[index_dist] = np.sum(prop_part[positions])
print("v0 :",time.time()-t0)
# out: v0 : 0.06322097778320312
如果 z_distances 很长,您可以使用 numba
来加快速度。
运行 calc
第一次通常会产生一些开销,我们可以通过 运行 一些小的 `z_distances 函数来消除这些开销。
下面的代码在我的笔记本电脑上实现了比纯 numpy 大约两倍的加速。
import numba as nb
@nb.njit(parallel=True)
def calc(result, z_distances):
n = z_distances.shape[0]
for ii in nb.prange(n):
pos = pos_part[:, 2] < z_distances[ii]
result[ii] = np.sum(prop_part[pos])
return result
result4 = np.zeros_like(result)
# _t = time.time()
# calc(result4, z_distances[:10])
# print(time.time()-_t)
t3 = time.time()
result4 = calc(result4, z_distances)
print("v3 :", time.time()-t3)
plt.plot(z_distances, result4)
我在对函数进行矢量化以使其有效应用于 numpy 数组时遇到问题。
我的节目条目:
- A pos_part Nb_particles 行的二维数组,3 列(基本上是 x、y、z坐标,只有z与困扰我的部分有关)Nb_particles可以达到几十万。
- 具有 Nb_particles 值的 prop_part 一维数组。这部分我已经介绍过了,创建是用一些很好的 numpy 函数完成的;我只是在这里放了一个与实际值相符的基本分布。
- A z_distances 一维数组,z=0 和 z=z_max. 之间的简单 np.arange
然后是需要时间的计算,因为我无法找到一种方法来仅通过数组的 numpy 操作来正确地做事。我想要做的是:
- 对于z_distances中的所有距离z_i,对[中的所有值求和=81=] 如果相应的粒子坐标 z_particle < z_i。这将 return 与 z_distances. 长度相同的一维数组
我目前的想法:
- 版本 0,for 循环,enumerate 和 np.where 检索我需要求和的值的索引。明明挺长的。
- 版本 1,在新数组上使用掩码(z 坐标和粒子属性的组合),并对掩码数组求和。似乎比 v0
- 版本 2,另一个掩码和 np.vectorize,但我知道它效率不高,因为矢量化基本上是 for环形。似乎仍然比 v0
- 第 3 版,我正在尝试在可以直接应用于 z_distances 的函数上使用掩码,但它目前无法正常工作。
所以,我来了。可能与排序和累计总和有关,但我不知道该怎么做,所以任何帮助将不胜感激。请在下面找到代码以使事情更清楚
提前致谢。
import numpy as np
import time
import matplotlib.pyplot as plt
# Creation of particles' positions
Nb_part = 150_000
pos_part = 10*np.random.rand(Nb_part,3)
pos_part[:,0] = pos_part[:,1] = 0
#usefull property creation
beta = 1/1.5
prop_part = (1/beta)*np.exp(-pos_part[:,2]/beta)
z_distances = np.arange(0,10,0.1)
#my version 0
t0=time.time()
result = np.empty(len(z_distances))
for index_dist, val_dist in enumerate(z_distances):
positions = np.where(pos_part[:,2]<val_dist)[0]
result[index_dist] = sum(prop_part[i] for i in positions)
print("v0 :",time.time()-t0)
#A graph to help understand
plt.figure()
plt.plot(z_distances,result, c="red")
plt.ylabel("Sum of particles' usefull property for particles with z-pos<d")
plt.xlabel("d")
#version 1 ??
t1=time.time()
combi = np.column_stack((pos_part[:,2],prop_part))
result2 = np.empty(len(z_distances))
for index_dist, val_dist in enumerate(z_distances):
mask = (combi[:,0]<val_dist)
result2[index_dist]=sum(combi[:,1][mask])
print("v1 :",time.time()-t1)
plt.plot(z_distances,result2, c="blue")
#version 2
t2=time.time()
def themask(a):
mask = (combi[:,0]<a)
return sum(combi[:,1][mask])
thefunc = np.vectorize(themask)
result3 = thefunc(z_distances)
print("v2 :",time.time()-t2)
plt.plot(z_distances,result3, c="green")
### This does not work so far
# version 3
# =============================
# t3=time.time()
# def thesum(a):
# mask = combi[combi[:,0]<a]
# return sum(mask[:,1])
# result4 = thesum(z_distances)
# print("v3 :",time.time()-t3)
# =============================
通过完全用 numpy 编写第一个版本,您可以获得更多的性能。将 pythons sum
替换为 np.sum
。而不是 for i in positions
列表理解,只需传递您正在创建的 positions
掩码。
事实上,np.where
不是必需的,我最好的版本是这样的:
#my version 0
t0=time.time()
result = np.empty(len(z_distances))
for index_dist, val_dist in enumerate(z_distances):
positions = pos_part[:, 2] < val_dist
result[index_dist] = np.sum(prop_part[positions])
print("v0 :",time.time()-t0)
# out: v0 : 0.06322097778320312
如果 z_distances 很长,您可以使用 numba
来加快速度。
运行 calc
第一次通常会产生一些开销,我们可以通过 运行 一些小的 `z_distances 函数来消除这些开销。
下面的代码在我的笔记本电脑上实现了比纯 numpy 大约两倍的加速。
import numba as nb
@nb.njit(parallel=True)
def calc(result, z_distances):
n = z_distances.shape[0]
for ii in nb.prange(n):
pos = pos_part[:, 2] < z_distances[ii]
result[ii] = np.sum(prop_part[pos])
return result
result4 = np.zeros_like(result)
# _t = time.time()
# calc(result4, z_distances[:10])
# print(time.time()-_t)
t3 = time.time()
result4 = calc(result4, z_distances)
print("v3 :", time.time()-t3)
plt.plot(z_distances, result4)