2013-04-13 2 views
4

파이썬/싸이언트/낸피를 사용하여 많은 4x4 행렬을 매우 빠르게 곱할 수있는 방법을 찾고 있습니다. 아무도 제안 할 수 있습니까? 파이썬에서의 빠르고 작은, 반복적 인 행렬 곱셈

내 현재의 시도를 보여주기 위해, 내가 계산해야하는 알고리즘이

A_1 * A_2 * A_3 * ... * A_N 

모든

A_i != A_j 

파이썬에서 그것의 예 :

means = array([0.0, 0.0, 34.28, 0.0, 0.0, 3.4]) 
stds = array([ 4.839339, 4.839339, 4.092728, 0.141421, 0.141421, 0.141421]) 

def fn(): 
    steps = means+stds*numpy.random.normal(size=(60,6)) 
    A = identity(4) 
    for step in steps: 
     A = dot(A, transform_step_to_4by4(step)) 
%timeit fn() 

1000 loops, best of 3: 570 us per loop 

이 구현 Cython/Numpy의 알고리즘은 Eigen/C++을 사용하는 등가 코드보다 약 100 배 느립니다. ptimizations. 하지만 C++을 사용하고 싶지는 않습니다.

+0

C에서 플러그인을 코딩 하시겠습니까? 파이썬에는 꽤 쉬운 C API가있다. – Dave

+0

어떻게해야할까요? 'transform_step_to_4by4 (step) ' – askewchan

+0

@askewchan, 6D 벡터에서 4x4 매트릭스로의 좌표 변환이라고하면 충분합니다. – Mike

답변

2

나는 당신이 (4,4)(60,6) 배열을 설정하는 방법을 알고하지 않기 때문에 내가 당신의 방법에 속도를 비교할 수는 없지만,이 시퀀스의 점 걸릴 작동합니다

A = np.arange(16).reshape(4,4) 
B = np.arange(4,20).reshape(4,4) 
C = np.arange(8,24).reshape(4,4) 

arrs = [A, B, C] 

P = reduce(np.dot, arrs) 

을 그리고이 에 해당되어야하지만보다 빠른 :

P = np.identity(4, dtype = arrs[0].dtype) 
for arr in arrs: 
    P = np.dot(P, arr) 

타이밍 시험 :

In [52]: def byloop(arrs): 
    ....:  P = np.identity(4) 
    ....:  for arr in arrs: 
    ....:   P = np.dot(P, arr) 
    ....:  return P 
    ....: 

In [53]: def byreduce(arrs): 
    ....:  return reduce(np.dot, arrs) 
    ....: 

In [54]: byloop(arrs) 
Out[54]: 
array([[ 5104, 5460, 5816, 6172], 
     [ 15728, 16820, 17912, 19004], 
     [ 26352, 28180, 30008, 31836], 
     [ 36976, 39540, 42104, 44668]]) 

In [55]: byreduce(arrs) 
Out[55]: 
array([[ 5104, 5460, 5816, 6172], 
     [15728, 16820, 17912, 19004], 
     [26352, 28180, 30008, 31836], 
     [36976, 39540, 42104, 44668]]) 
In [56]: timeit byloop(arrs) 
1000 loops, best of 3: 1.26 ms per loop 

In [57]: timeit byreduce(arrs) 
1000 loops, best of 3: 656 us per loop 
+0

byloop는 3, 2 개의 점을 줄입니다. 더 빨라질 수도 있지만, 타이밍만큼이나 귀찮습니다. 또한 유형이 같지 않고 ID가 불필요한 시간을 소비합니다 (길이가 3 인 arrs는 많은 시간을 합산합니다). – seberg

+0

Thanks @seberg 긴 목록이있는 Timing이 지금 게시되었습니다. – askewchan

+1

여전히 dtype 문제 (그다지 작지는 않지만)는 순수 파이썬에서 제안 할 것이별로 없다는 것을 고려하면 ... – seberg

6

당신이 번식 할 행렬의 각을 생산하는 파이썬 함수 호출을 할 수있는 경우에, 당신은 기본적으로 성능 현명한 망했다 : len(arrs) = 1000. 당신이 transform_step_to_4by4 기능을 벡터화하고, 모양 (n, 4, 4)와 배열을 반환 할 수 있습니다하지만 당신은 matrix_multiply를 사용하여 시간을 절약 할 수 있습니다

import numpy as np 
from numpy.core.umath_tests import matrix_multiply 

matrices = np.random.rand(64, 4, 4) - 0.5 

def mat_loop_reduce(m): 
    ret = m[0] 
    for x in m[1:]: 
     ret = np.dot(ret, x) 
    return ret 

def mat_reduce(m): 
    while len(m) % 2 == 0: 
     m = matrix_multiply(m[::2], m[1::2]) 
    return mat_loop_reduce(m) 

In [2]: %timeit mat_reduce(matrices) 
1000 loops, best of 3: 287 us per loop 

In [3]: %timeit mat_loop_reduce(matrices) 
1000 loops, best of 3: 721 us per loop 

In [4]: np.allclose(mat_loop_reduce(matrices), mat_reduce(matrices)) 
Out[4]: True 

이제 로그를 (n)이 파이썬은 대신 N의 좋은 호출 2.5 배속으로 n = 1024인데 10 배에 가깝습니다. 분명히 matrix_multiply은 ufunc이므로 .reduce 메서드를 사용하면 코드가 파이썬에서 루프를 돌릴 수 없게됩니다. 나는 그것이 비록 실행하기 난해한 오류가 유지할 수 없었다 :이 작은 행렬에 대한

In [7]: matrix_multiply.reduce(matrices) 
------------------------------------------------------------ 
Traceback (most recent call last): 
    File "<ipython console>", line 1, in <module> 
RuntimeError: Reduction not defined on ufunc with signature 
1

, 모두 BLAS를 사용하지 않도록하는 것이 유리할 수 있습니다. 거기에는 작은 행렬 곱셈 라이브러리가 있습니다 (예 : libSMM) (더 있습니다).

f2py 또는 cython으로 감싸는 것이 가능하거나, cython 또는 fortran/f2py에서 직접 롤하는 것이 더 쉬울 수 있습니다.