2017-02-05 5 views
1

私は外側と内側のforループをベクトル化する方法をチェックしています。これらにはいくつかの計算があり、それらの中にも削除があります - それはそれほど単純ではありません。これらが計算と削除を含むときに、外側と内側のループをベクトル化する

これはどのようにして最適化されますか?

import numpy as np 

flattenedArray = np.ndarray.tolist(someNumpyArray) 

#flattenedArray is a python list of lists. 
c = flattenedArray[:] 
for a in range (len(flattenedArray)): 
    for b in range(a+1, len(flattenedArray)): 
     if a == b: 
      continue 
     i0 = flattenedArray[a][0] 
     j0 = flattenedArray[a][1] 
     z0 = flattenedArray[a][2] 
     i1 = flattenedArray[b][0] 
     i2 = flattenedArray[b][1] 
     z1 = flattenedArray[b][2] 

     if ((np.square(z0-z1)) <= (np.square(i0-i1) + (np.square(j0-j2)))): 
      if (np.square(i0-i1) + (np.square(j0-j1))) <= (np.square(z0+z1)): 
         c.remove(flattenedArray[b]) 

答えて

0

あなたはマスクで作業している場合、内側のループが問題のあまりないですベクトル化:

import numpy as np 

# I'm using a random array 
flattenedArray = np.random.randint(0, 100, (10, 3)) 
mask = np.zeros(flattenedArray.shape[0], bool) 

for idx, row in enumerate(flattenedArray): 
    # Calculate the broadcasted elementwise addition/subtraction of this row 
    # with all following 
    added_squared = np.square(row[None, :] + flattenedArray[idx+1:]) 
    subtracted_squared = np.square(row[None, :] - flattenedArray[idx+1:]) 
    # Check the conditions 
    col1_col2_added = subtracted_squared[:, 0] + subtracted_squared[:, 1] 
    cond1 = subtracted_squared[:, 2] <= col1_col2_added 
    cond2 = col1_col2_added <= added_squared[:, 2] 
    # Update the mask 
    mask[idx+1:] |= cond1 & cond2 

# Apply the mask 
flattenedArray[mask] 

あなたも1放送でそれを行う必要があり、外側のループをベクトル化したい場合は、そのただし、O(n)の代わりに多くのメモリO(n**2)を使用します。クリティカルな内部ループが既にベクトル化されていることを考えると、外部ループをベクトル化することによって多くの高速化はありません。

+0

ありがとうございます!申し訳ありませんが、私はこれらの2つの式で間違いを犯したことに気付きました。私はあなたがやったことで遊んでみましたが、それを得ることはできません。私は元の投稿の2つの公式を修正しました。もう一度それをチェックしてもよろしいですか? – user1812626

+0

@ user1812626私はそれが正しいと思う、見てください:) – MSeifert

1

@MSeifertはもちろん、しばしば正しいです。したがって、次の完全なベクトル化は、「完了した方法」を示すことだけです。

import numpy as np 

N = 4 
data = np.random.random((N, 3)) 

# vectorised code 
j, i = np.tril_indices(N, -1) # chose tril over triu to have contiguous columns 
           # useful later 
sqsum = np.square(data[i,0]-data[j,0]) + np.square(data[i,1]-data[j,1]) 
cond = np.square(data[i, 2] + data[j, 2]) >= sqsum 
cond &= np.square(data[i, 2] - data[j, 2]) <= sqsum 
# because equal 'b's are grouped together we can use reduceat: 
cond = np.r_[False, np.logical_or.reduceat(
    cond, np.add.accumulate(np.arange(N-1)))] 
left = data[~cond, :] 


# original code (modified to make it run) 
flattenedArray = np.ndarray.tolist(data) 

#flattenedArray is a python list of lists. 
c = flattenedArray[:] 
for a in range (len(flattenedArray)): 
    for b in range(a+1, len(flattenedArray)): 
     if a == b: 
      continue 
     i0 = flattenedArray[a][0] 
     j0 = flattenedArray[a][1] 
     z0 = flattenedArray[a][2] 
     i1 = flattenedArray[b][0] 
     j1 = flattenedArray[b][1] 
     z1 = flattenedArray[b][2] 

     if ((np.square(z0-z1)) <= (np.square(i0-i1) + (np.square(j0-j1)))): 
      if (np.square(i0-i1) + (np.square(j0-j1))) <= (np.square(z0+z1)): 
       try: 
         c.remove(flattenedArray[b]) 
       except: 
         pass 

# check they are the same 
print(np.alltrue(c == left)) 
+0

それは賢いです(+1)!しかし、今私は小さなアレイではずっと速く、大きなアレイではそれよりずっと遅いことがわかります。 :-) – MSeifert

+0

私を打つ。それは空想的なインデックスのコストかもしれませんか?あるいは、大きすぎる配列は効率的ではありませんか?なぜそれが小さい、ないアイデアの方が速いのですか? –

関連する問題