0

数値積分は、予想よりも指数関数的に長くかかっています。メッシュ上で反復を実装する方法が要因になるかどうかを知りたいです。私のコードは次のようになります。

import numpy as np
import itertools as it

U = np.linspace(0, 2*np.pi)
V = np.linspace(0, np.pi)

for (u, v) in it.product(U,V):
    # values = computation on each grid point, does not call any outside functions
    # solution = sum(values)
return solution

計算は長いので省略しました。私の質問は、具体的には、パラメーター空間 (u, v) で計算を実装した方法に関するものです。次のような代替案を知っていnumpy.meshgridます。ただし、これらはすべて (非常に大きな) 行列のインスタンスを作成するようであり、それらをメモリに格納すると処理が遅くなると思います。

私のプログラムを高速化する代替手段はありit.productますか、それともボトルネックを他の場所で探す必要がありますか?

編集: 問題の for ループは次のとおりです (ベクトル化できるかどうかを確認するため)。

import random  
import numpy as np  
import itertools as it 

##########################################################################
# Initialize the inputs with random (to save space)
##########################################################################
mat1 = np.array([[random.random() for i in range(3)] for i in range(3)])
mat2 = np.array([[random.random() for i in range(3)] for i in range(3)]) 
a1, a2, a3 = np.array([random.random() for i in range(3)]) 
plane_normal = np.array([random.random() for i in range(3)])  
plane_point = np.array([random.random() for i in range(3)])  
d = np.dot(plane_normal, plane_point)  
truthval = True

##########################################################################
# Initialize the loop
##########################################################################
N = 100 
U = np.linspace(0, 2*np.pi, N + 1, endpoint = False) 
V = np.linspace(0, np.pi, N + 1, endpoint = False) 
U = U[1:N+1] V = V[1:N+1]

Vsum = 0
Usum = 0

##########################################################################
# The for loops starts here
##########################################################################   
for (u, v) in it.product(U,V):

    cart_point = np.array([a1*np.cos(u)*np.sin(v), 
                           a2*np.sin(u)*np.sin(v), 
                           a3*np.cos(v)])

    surf_normal = np.array(
            [2*x / a**2 for (x, a) in zip(cart_point, [a1,a2,a3])])


    differential_area = \
        np.sqrt((a1*a2*np.cos(v)*np.sin(v))**2 + \
        a3**2*np.sin(v)**4 * \
        ((a2*np.cos(u))**2 + (a1*np.sin(u))**2)) * \
        (np.pi**2 / (2*N**2)) 


    if (np.dot(plane_normal, cart_point) - d > 0) == truthval:
        perp_normal = plane_normal
        f = np.dot(np.dot(mat2, surf_normal), perp_normal)
        Vsum += f*differential_area
    else:
        perp_normal = - plane_normal
        f = np.dot(np.dot(mat2, surf_normal), perp_normal)
        Usum += f*differential_area

integral = abs(Vsum) + abs(Usum)
4

3 に答える 3

1

ipython では、50 x 50 グリッドスペースで単純な合計計算を行いました 

In [31]: sum(u*v for (u,v) in it.product(U,V))
Out[31]: 12337.005501361698

In [33]: UU,VV = np.meshgrid(U,V); sum(sum(UU*VV))
Out[33]: 12337.005501361693

In [34]: timeit UU,VV = np.meshgrid(U,V); sum(sum(UU*VV))
1000 loops, best of 3: 293 us per loop

In [35]: timeit sum(u*v for (u,v) in it.product(U,V)) 
100 loops, best of 3: 2.95 ms per loop

In [38]: timeit list(it.product(U,V))
1000 loops, best of 3: 213 us per loop

In [45]: timeit UU,VV = np.meshgrid(U,V); (UU*VV).sum().sum()
10000 loops, best of 3: 70.3 us per loop
# using numpy's own sum is even better

productproductそれ自体が遅いためではなく、ポイントごとの計算が原因で(10 倍) 遅くなります。2 (50,50) の配列を使用するように計算をベクトル化できる場合 (ループを一切行わずに)、全体の時間を短縮する必要があります。それが使用する主な理由ですnumpy。

于 2013-08-11T00:11:22.713 に答える
0

[k for k in it.product(U,V)]私にとっては2ミリ秒で実行され、itertoolパッケージは効率的に作られています。たとえば、最初に長い配列を作成しません（http://docs.python.org/2/library/itertools.html）。

犯人は、反復内のコード、またはリンスペースで多くのポイントを使用しているようです。

于 2013-08-10T20:37:07.333 に答える