Tại sao einsum của numpy lại nhanh hơn các chức năng tích hợp của numpy?

Question

Cho phép bắt đầu với ba mảng dtype=np.double. Thời gian được thực hiện trên một CPU intel bằng cách sử dụng 1.7.1 xếp gọn được biên dịch với icc và được liên kết với số mkl của intel. Một cpu AMD với khối lượng 1.6.1 được biên dịch với gcc mà không cần mkl cũng được sử dụng để xác minh thời gian. Xin lưu ý timings quy mô gần tuyến tính với kích thước hệ thống và không phải do các chi phí nhỏ phát sinh trong các chức năng NumPy if báo cáo những sự khác biệt sẽ hiển thị trong micro không mili giây:

arr_1D=np.arange(500,dtype=np.double) 
large_arr_1D=np.arange(100000,dtype=np.double) 
arr_2D=np.arange(500**2,dtype=np.double).reshape(500,500) 
arr_3D=np.arange(500**3,dtype=np.double).reshape(500,500,500) 

Đầu tiên cho phép nhìn vào np.sum chức năng:

np.all(np.sum(arr_3D)==np.einsum('ijk->',arr_3D)) 
True 

%timeit np.sum(arr_3D) 
10 loops, best of 3: 142 ms per loop 

%timeit np.einsum('ijk->', arr_3D) 
10 loops, best of 3: 70.2 ms per loop 

Powers:

np.allclose(arr_3D*arr_3D*arr_3D,np.einsum('ijk,ijk,ijk->ijk',arr_3D,arr_3D,arr_3D)) 
True 

%timeit arr_3D*arr_3D*arr_3D 
1 loops, best of 3: 1.32 s per loop 

%timeit np.einsum('ijk,ijk,ijk->ijk', arr_3D, arr_3D, arr_3D) 
1 loops, best of 3: 694 ms per loop 

sản phẩm Outer:

np.all(np.outer(arr_1D,arr_1D)==np.einsum('i,k->ik',arr_1D,arr_1D)) 
True 

%timeit np.outer(arr_1D, arr_1D) 
1000 loops, best of 3: 411 us per loop 

%timeit np.einsum('i,k->ik', arr_1D, arr_1D) 
1000 loops, best of 3: 245 us per loop 

Tất cả những điều trên nhanh gấp hai lần với np.einsum. Đây nên là táo để so sánh táo như tất cả mọi thứ là cụ thể của dtype=np.double. Tôi mong chờ tốc độ lên trong một hoạt động như thế này:

np.allclose(np.sum(arr_2D*arr_3D),np.einsum('ij,oij->',arr_2D,arr_3D)) 
True 

%timeit np.sum(arr_2D*arr_3D) 
1 loops, best of 3: 813 ms per loop 

%timeit np.einsum('ij,oij->', arr_2D, arr_3D) 
10 loops, best of 3: 85.1 ms per loop 

Einsum có vẻ là ít nhất hai lần càng nhanh cho np.inner, np.outer, np.kron, và np.sum bất kể axes lựa chọn. Ngoại lệ chính là np.dot vì nó gọi DGEMM từ thư viện BLAS. Vì vậy, tại sao là np.einsum nhanh hơn các hàm numpy khác tương đương?

Trường hợp DGEMM cho đầy đủ:

np.allclose(np.dot(arr_2D,arr_2D),np.einsum('ij,jk',arr_2D,arr_2D)) 
True 

%timeit np.einsum('ij,jk',arr_2D,arr_2D) 
10 loops, best of 3: 56.1 ms per loop 

%timeit np.dot(arr_2D,arr_2D) 
100 loops, best of 3: 5.17 ms per loop 

---

Lý thuyết hàng đầu là từ @sebergs nhận xét rằng np.einsum có thể tận dụng SSE2, nhưng ufuncs NumPy sẽ phải đến NumPy 1.8 (xem change log). Tôi tin rằng đây là câu trả lời đúng, nhưng có không phải có thể xác nhận nó. Một số bằng chứng hạn chế có thể được tìm thấy bằng cách thay đổi dtype của mảng đầu vào và quan sát sự khác biệt tốc độ và thực tế là không phải ai cũng quan sát cùng một xu hướng về thời gian.

Answer 1

Giờ đây, phiên bản 1.8 được phát hành, theo các tài liệu, tất cả các ufuncs nên sử dụng SSE2, tôi muốn kiểm tra kỹ xem nhận xét của Seberg về SSE2 có hợp lệ không.

Để thực hiện thử nghiệm một trăn mới 2.7 cài đặt đã được tạo ra - numpy 1.7 và 1.8 được biên dịch với icc sử dụng các tùy chọn tiêu chuẩn trên lõi AMD Opteron chạy Ubuntu.

Đây là chạy thử cả trước và sau khi nâng cấp 1.8:

import numpy as np 
import timeit 

arr_1D=np.arange(5000,dtype=np.double) 
arr_2D=np.arange(500**2,dtype=np.double).reshape(500,500) 
arr_3D=np.arange(500**3,dtype=np.double).reshape(500,500,500) 

print 'Summation test:' 
print timeit.timeit('np.sum(arr_3D)', 
         'import numpy as np; from __main__ import arr_1D, arr_2D, arr_3D', 
         number=5)/5 
print timeit.timeit('np.einsum("ijk->", arr_3D)', 
         'import numpy as np; from __main__ import arr_1D, arr_2D, arr_3D', 
         number=5)/5 
print '----------------------\n' 


print 'Power test:' 
print timeit.timeit('arr_3D*arr_3D*arr_3D', 
         'import numpy as np; from __main__ import arr_1D, arr_2D, arr_3D', 
         number=5)/5 
print timeit.timeit('np.einsum("ijk,ijk,ijk->ijk", arr_3D, arr_3D, arr_3D)', 
         'import numpy as np; from __main__ import arr_1D, arr_2D, arr_3D', 
         number=5)/5 
print '----------------------\n' 


print 'Outer test:' 
print timeit.timeit('np.outer(arr_1D, arr_1D)', 
         'import numpy as np; from __main__ import arr_1D, arr_2D, arr_3D', 
         number=5)/5 
print timeit.timeit('np.einsum("i,k->ik", arr_1D, arr_1D)', 
         'import numpy as np; from __main__ import arr_1D, arr_2D, arr_3D', 
         number=5)/5 
print '----------------------\n' 


print 'Einsum test:' 
print timeit.timeit('np.sum(arr_2D*arr_3D)', 
         'import numpy as np; from __main__ import arr_1D, arr_2D, arr_3D', 
         number=5)/5 
print timeit.timeit('np.einsum("ij,oij->", arr_2D, arr_3D)', 
         'import numpy as np; from __main__ import arr_1D, arr_2D, arr_3D', 
         number=5)/5 
print '----------------------\n' 

NumPy 1.7.1:

Summation test: 
0.172988510132 
0.0934836149216 
---------------------- 

Power test: 
1.93524689674 
0.839519000053 
---------------------- 

Outer test: 
0.130380821228 
0.121401786804 
---------------------- 

Einsum test: 
0.979052495956 
0.126066613197 

NumPy 1.8:

Summation test: 
0.116551589966 
0.0920487880707 
---------------------- 

Power test: 
1.23683619499 
0.815982818604 
---------------------- 

Outer test: 
0.131808176041 
0.127472200394 
---------------------- 

Einsum test: 
0.781750011444 
0.129271841049 

Tôi nghĩ rằng đây là khá thuyết phục rằng SSE đóng một vai trò lớn trong sự khác biệt thời gian, cần lưu ý rằng lặp đi lặp lại các xét nghiệm này timings ve ry bởi chỉ ~ 0,003s. Sự khác biệt còn lại nên được đề cập trong các câu trả lời khác cho câu hỏi này.

Answer 2

Tôi nghĩ rằng những timings giải thích những gì đang xảy ra:

a = np.arange(1000, dtype=np.double) 
%timeit np.einsum('i->', a) 
100000 loops, best of 3: 3.32 us per loop 
%timeit np.sum(a) 
100000 loops, best of 3: 6.84 us per loop 

a = np.arange(10000, dtype=np.double) 
%timeit np.einsum('i->', a) 
100000 loops, best of 3: 12.6 us per loop 
%timeit np.sum(a) 
100000 loops, best of 3: 16.5 us per loop 

a = np.arange(100000, dtype=np.double) 
%timeit np.einsum('i->', a) 
10000 loops, best of 3: 103 us per loop 
%timeit np.sum(a) 
10000 loops, best of 3: 109 us per loop 

Vì vậy, về cơ bản bạn có một overhead 3us gần như liên tục khi gọi np.sum qua np.einsum, vì vậy họ về cơ bản chạy càng nhanh, nhưng người ta mất nhiều thời gian một chút để có được đang đi. Tại sao điều đó có thể? Tiền của tôi là về sau:

a = np.arange(1000, dtype=object) 
%timeit np.einsum('i->', a) 
Traceback (most recent call last): 
... 
TypeError: invalid data type for einsum 
%timeit np.sum(a) 
10000 loops, best of 3: 20.3 us per loop 

Không chắc những gì đang xảy ra một cách chính xác, nhưng có vẻ như np.einsum được bỏ qua một số kiểm tra để trích xuất loại chức năng cụ thể để thực hiện phép nhân và bổ sung, và đang xảy ra trực tiếp với * và + chỉ dành cho các loại C chuẩn.

---

Các trường hợp đa chiều không khác nhau:

n = 10; a = np.arange(n**3, dtype=np.double).reshape(n, n, n) 
%timeit np.einsum('ijk->', a) 
100000 loops, best of 3: 3.79 us per loop 
%timeit np.sum(a) 
100000 loops, best of 3: 7.33 us per loop 

n = 100; a = np.arange(n**3, dtype=np.double).reshape(n, n, n) 
%timeit np.einsum('ijk->', a) 
1000 loops, best of 3: 1.2 ms per loop 
%timeit np.sum(a) 
1000 loops, best of 3: 1.23 ms per loop 

Vì vậy, một overhead chủ yếu là liên tục, không phải là một nhanh hơn chạy khi họ bắt tay vào nó.

Answer 3

Trước hết, đã có rất nhiều cuộc thảo luận trong quá khứ về điều này trong danh sách khó xử. Ví dụ, xem: http://numpy-discussion.10968.n7.nabble.com/poor-performance-of-sum-with-sub-machine-word-integer-types-td41.html http://numpy-discussion.10968.n7.nabble.com/odd-performance-of-sum-td3332.html

Một số nắm thực tế là einsum là mới, và có lẽ đó là cố gắng để được tốt hơn về sự liên kết bộ nhớ cache và các vấn đề truy xuất bộ nhớ khác, trong khi rất nhiều các chức năng NumPy cũ tập trung vào một cách dễ dàng thực hiện di động qua một cái được tối ưu hóa rất nhiều. Tôi chỉ suy đoán, có, mặc dù.

---

Tuy nhiên, một số việc bạn đang làm không phải là một so sánh "táo" với "táo".

Ngoài những gì @Jamie đã nói, sum sử dụng một ắc thích hợp hơn cho mảng

Ví dụ, sum là cẩn thận hơn về việc kiểm tra các loại đầu vào và sử dụng một ắc thích hợp. Ví dụ, hãy xem xét những điều sau đây:

In [1]: x = 255 * np.ones(100, dtype=np.uint8) 

In [2]: x 
Out[2]: 
array([255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 
     255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 
     255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 
     255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 
     255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 
     255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 
     255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 255, 
     255, 255, 255, 255, 255, 255, 255, 255, 255], dtype=uint8) 

Lưu ý rằng sum là đúng:

In [3]: x.sum() 
Out[3]: 25500 

Trong khi einsum sẽ cho kết quả sai:

In [4]: np.einsum('i->', x) 
Out[4]: 156 

Nhưng nếu chúng ta sử dụng một ít hạn chế dtype, chúng tôi sẽ vẫn nhận được kết quả bạn mong đợi:

In [5]: y = 255 * np.ones(100) 

In [6]: np.einsum('i->', y) 
Out[6]: 25500.0