4

最初にデータを補間してサーフェスを作成することにより、二重積分を実行しようとしています。numba を使用してこのプロセスを高速化しようとしていますが、時間がかかりすぎます。

これが私のコードで、コードを実行するために必要な画像はhereとhereにあります。

4

1 に答える 1

8

あなたのコードには 4 重にネストされた一連の for ループがあることに注意して、内側のペアの最適化に焦点を当てました。古いコードは次のとおりです。

for i in xrange(K.shape[0]):
    for j in xrange(K.shape[1]):

        print(i,j)
        '''create an r vector '''
        r=(i*distX,j*distY,z)

        for x in xrange(img.shape[0]):
            for y in xrange(img.shape[1]):
                '''create an ksi vector, then calculate
                   it's norm, and the dot product of r and ksi'''
                ksi=(x*distX,y*distY,z)
                ksiNorm=np.linalg.norm(ksi)
                ksiDotR=float(np.dot(ksi,r))

                '''calculate the integrand'''
                temp[x,y]=img[x,y]*np.exp(1j*k*ksiDotR/ksiNorm)

        '''interpolate so that we can do the integral and take the integral'''
        temp2=rbs(a,b,temp.real)
        K[i,j]=temp2.integral(0,n,0,m)

K と img はそれぞれ約 2000x2000 であるため、最も内側のステートメントを 16 兆回実行する必要があります。これは Python を使用して単純に実用的ではありませんが、NumPy を使用してベクトル化する C や Fortran に作業を移すことができます。結果が一致することを確認するために、この手順を一度に 1 つずつ慎重に実行しました。これが私が最終的に得たものです：

'''create all r vectors'''
R = np.empty((K.shape[0], K.shape[1], 3))
R[:,:,0] = np.repeat(np.arange(K.shape[0]), K.shape[1]).reshape(K.shape) * distX
R[:,:,1] = np.arange(K.shape[1]) * distY
R[:,:,2] = z

'''create all ksi vectors'''
KSI = np.empty((img.shape[0], img.shape[1], 3))
KSI[:,:,0] = np.repeat(np.arange(img.shape[0]), img.shape[1]).reshape(img.shape) * distX
KSI[:,:,1] = np.arange(img.shape[1]) * distY
KSI[:,:,2] = z

# vectorized 2-norm; see http://stackoverflow.com/a/7741976/4323                                                    
KSInorm = np.sum(np.abs(KSI)**2,axis=-1)**(1./2)

# loop over entire K, which is same shape as img, rows first                                                        
# this loop populates K, one pixel at a time (so can be parallelized)                                               
for i in xrange(K.shape[0]):                                                                                    
    for j in xrange(K.shape[1]):                                                                                

        print(i, j)

        KSIdotR = np.dot(KSI, R[i,j])
        temp = img * np.exp(1j * k * KSIdotR / KSInorm)

        '''interpolate so that we can do the integral and take the integral'''
        temp2 = rbs(a, b, temp.real)
        K[i,j] = temp2.integral(0, n, 0, m)

ループの内側のペアは完全になくなり、前もって実行されたベクトル化された演算に置き換えられました (スペース コストは入力のサイズに比例します)。

これにより、 Numba を使用せずに、Macbook Air 1.6 GHz i5 で外側の 2 つのループの反復あたりの時間が 340 秒から 1.3 秒に短縮されます。反復あたり 1.3 秒のうち、0.68 秒がrbs関数に費やされscipy.interpolate.RectBivariateSplineます。おそらく、さらに最適化する余地があります。いくつかのアイデアを次に示します。

  1. Numba を再度有効にします。私のシステムにはありません。この時点では大きな違いはないかもしれませんが、簡単にテストできます。
  2. 実行中の基本的な計算を単純化するなど、よりドメイン固有の最適化を行います。私の最適化は無損失であることを意図しており、問題のドメインがわからないため、あなたができるほど深く最適化することはできません.
  3. 残りのループをベクトル化してみてください。scipy RBS 関数を呼び出しごとに複数の計算をサポートするものに置き換えようとしない限り、これは難しいかもしれません。
  4. より高速な CPU を入手してください。私のはかなり遅いです。私の小さなラップトップよりも優れたコンピューターを使用するだけで、おそらく少なくとも 2 倍のスピードアップを得ることができます。
  5. データをダウンサンプリングします。テスト画像は 2000x2000 ピクセルですが、詳細はほとんど含まれていません。直線の寸法を 2 ～ 10 倍にすると、大幅なスピードアップが得られます。

今のところは以上です。これはあなたをどこに置きますか？わずかに優れたコンピューターで、それ以上の最適化作業がないと仮定すると、最適化されたコードでさえ、テスト画像を処理するのに約 1 か月かかります。これを一度だけ行う必要がある場合は、おそらくそれで問題ありません。より頻繁に実行する必要がある場合、またはさまざまなことを試すときにコードを反復する必要がある場合は、最適化を続ける必要があります。現在半分以上の時間を消費している RBS 関数から始めます。

kおまけのヒント:とのようなほぼ同一の変数名がなく、変数名としても複素数接尾辞 ( ) としてKも使用されていない場合、コードははるかに扱いやすくなります。j0j

于 2013-07-13T04:12:51.540 に答える