这是一个小代码,演示了我遇到的错误。
import numpy as np
r=4.0
L=20.0
ratio = 4*np.pi / 3.0 * (r/L)**3
for i in range(5, 16):
n = 10**i
m = int(ratio * n)
print i,n,m
ip = np.random.random_integers(100, size=(n,3))
jp = np.random.random_integers(100, size=(m,3))
a = np.expand_dims(ip, -1) == jp.T
b = np.where( a.all(axis=1).any(axis=1) )[0]
我得到以下输出:
5 100000 3351
6 1000000 33510
Traceback (most recent call last):
File "example.py", line 16, in <module>
b = np.where( a.all(axis=1).any(axis=1) )[0]
AttributeError: 'bool' object has no attribute 'all'
有人知道这是怎么回事吗?
另外,在ip中索引jp的元素位置的合理快速方法也可以使用。我可能会从这里开始采用第二种解决方案
您正在广播ip
反对jp
创建非常大的数组。当i==6
您拥有100GB阵列时。
一种解决方案是遍历数组:
for i in range(2,6):
t=time.time()
n = 10**i+1
m = int(ratio * n)
print i,n,m
ip = np.random.random_integers(10, size=(n,3))
jp = np.random.random_integers(10, size=(m,3))
chunksize=10000
if chunksize>ip.shape[0]:
chunksize=ip.shape[0]
span=ip.shape[0]/chunksize
remainder=(ip.shape[0]-span*chunksize)
out=[]
start=0
for n in xrange(span):
end=start+chunksize
a = np.expand_dims(ip[start:end], -1) == jp.T
b = np.where( a.all(axis=1).any(axis=1) )[0]
out.append(b+start)
start+=chunksize
if remainder!=0:
a = np.expand_dims(ip[-remainder:], -1) == jp.T
b = np.where( a.all(axis=1).any(axis=1) )[0]
out.append(b+end)
end=np.sort(np.concatenate(out))
print time.time()-t,end.shape
时间大约是10秒,i==6
因此i==7
大约需要20分钟。
本文收集自互联网,转载请注明来源。
如有侵权,请联系 [email protected] 删除。
我来说两句