Zoptymalizowany produkt dotowy w Pythonie

Question

Produkt punktowy dwóch n-wymiarowych wektorów u=[u1,u2,...un] i v=[v1,v2,...,vn] jest podany przez u1*v1 + u2*v2 + ... + un*vn.

Pytanie posted yesterday zachęciło mnie do znalezienia najszybszego sposobu obliczania produktów dot w Pythonie przy użyciu tylko standardowej biblioteki, bez modułów firm trzecich lub wywołań C/Fortran/C++.

Wyliczyłem cztery różne podejścia; jak dotąd najszybszym wydaje się być sum(starmap(mul,izip(v1,v2))) (gdzie starmap i izip pochodzą z modułu itertools).

Dla kodu przedstawionego poniżej, są minęły czasy (w sekundach, na milion przebiegów):

d0: 12.01215 
d1: 11.76151 
d2: 12.54092 
d3: 09.58523 

można myśleć szybszy sposób to zrobić?

import timeit # module with timing subroutines                
import random # module to generate random numnbers               
from itertools import imap,starmap,izip 
from operator import mul 

def v(N=50,min=-10,max=10): 
    """Generates a random vector (in an array) of dimension N; the           
    values are integers in the range [min,max].""" 
    out = [] 
    for k in range(N): 
     out.append(random.randint(min,max)) 
    return out 

def check(v1,v2): 
    if len(v1)!=len(v2): 
     raise ValueError,"the lenght of both arrays must be the same" 
    pass 

def d0(v1,v2): 
    """                          
    d0 is Nominal approach:                     
    multiply/add in a loop                     
    """ 
    check(v1,v2) 
    out = 0 
    for k in range(len(v1)): 
     out += v1[k] * v2[k] 
    return out 

def d1(v1,v2): 
    """                          
    d1 uses an imap (from itertools)                   
    """ 
    check(v1,v2) 
    return sum(imap(mul,v1,v2)) 

def d2(v1,v2): 
    """                          
    d2 uses a conventional map                    
    """ 
    check(v1,v2) 
    return sum(map(mul,v1,v2)) 

def d3(v1,v2): 
    """                          
    d3 uses a starmap (itertools) to apply the mul operator on an izipped (v1,v2)       
    """ 
    check(v1,v2) 
    return sum(starmap(mul,izip(v1,v2))) 

# generate the test vectors                     
v1 = v() 
v2 = v() 

if __name__ == '__main__': 

    # Generate two test vectors of dimension N                

    t0 = timeit.Timer("d0(v1,v2)","from dot_product import d0,v1,v2") 
    t1 = timeit.Timer("d1(v1,v2)","from dot_product import d1,v1,v2") 
    t2 = timeit.Timer("d2(v1,v2)","from dot_product import d2,v1,v2") 
    t3 = timeit.Timer("d3(v1,v2)","from dot_product import d3,v1,v2") 

    print "d0 elapsed: ", t0.timeit() 
    print "d1 elapsed: ", t1.timeit() 
    print "d2 elapsed: ", t2.timeit() 
    print "d3 elapsed: ", t3.timeit() 

Należy zauważyć, że nazwa pliku musi być dot_product.py dla skryptu do uruchomienia; Użyłem Pythona 2.5.1 na Mac OS X wersja 10.5.8.

EDIT:

wpadłem skryptu dla n = 1000 i są to wyniki (w sekundach, dla jednego miliona przejazdów):

d0: 205.35457 
d1: 208.13006 
d2: 230.07463 
d3: 155.29670 

myślę, że można bezpiecznie przyjąć, że rzeczywiście , opcja trzecia jest najszybsza, a opcja druga najwolniejsza (z czterech prezentowanych).

Answer 1

Tak dla zabawy napisałem "D4", który wykorzystuje numpy:

from numpy import dot 
def d4(v1, v2): 
    check(v1, v2) 
    return dot(v1, v2) 

Moje wyniki (Python 2.5.1, XP Pro SP3, 2 GHz Core2 Duo T7200):

d0 elapsed: 12.1977242918 
d1 elapsed: 13.885232341 
d2 elapsed: 13.7929552499 
d3 elapsed: 11.0952246724 

d4 upłynął: 56.3278584289 # idź numpy!

I jeszcze więcej zabawy, włączyłem psyco:

d0 elapsed: 0.965477735299 
d1 elapsed: 12.5354792299 
d2 elapsed: 12.9748163524 
d3 elapsed: 9.78255448667 

d4 upłynął: +54,4599059378

na tej podstawie, oświadczam d0 zwycięzcę :)

---

Aktualizacja

@kizer.SE: I prawdopodobnie należy wspomnieć, że zrobiłem wszystko, aby przekształcić numpy tablice pierwszy:

from numpy import array 
v3 = [array(vec) for vec in v1] 
v4 = [array(vec) for vec in v2] 

# then 
t4 = timeit.Timer("d4(v3,v4)","from dot_product import d4,v3,v4") 

I zawarte check(v1, v2) ponieważ jest to ujęte w innych testach. Pozostawienie go dawałoby numpy nieuczciwą przewagę (choć wygląda na to, że może z niej skorzystać). Konwersja macierzy ogoliła się około sekundy (znacznie mniej, niż myślałem).

Wszystkie moje testy zostały przeprowadzone z N = 50.

@nikow: Używam numpy 1.0.4, która jest niewątpliwie stara, z pewnością możliwe, że poprawili wydajność przez ostatnie półtora roku od czasu, gdy go zainstalowałem.

---

Aktualizacja # 2

@ kaiser.se Wow, jesteś całkowicie w prawo. Pewnie myślałem, że są to listy list lub coś takiego (naprawdę nie mam pojęcia, o czym myślałem ... +1 dla programowania parami).

Jak to wygląda:

v3 = array(v1) 
v4 = array(v2) 

Nowe wyniki:

d4 elapsed: 3.22535741274 

Z psyco:

d4 elapsed: 2.09182619579 

d0 wciąż wygrywa z psyco, ale numpy jest prawdopodobnie lepszy ogólny, zwłaszcza z większymi zbiorami danych.

Wczoraj byłam trochę zaniepokojona moim powolnym, mizernym rezultatem, ponieważ prawdopodobnie numpy jest używany do partii obliczeń i miał dużo optymalizacji. Oczywiście, nie przejął się wystarczająco, aby sprawdzić mój wynik :)

Answer 2

nie wiem o szybsze, ale sugeruję:

sum(i*j for i, j in zip(v1, v2)) 

to o wiele łatwiejsze do odczytania i nie wymaga nawet moduły standardowej bibliotece.

Answer 3

Oto porównanie z numpy. Porównamy szybki Starmap podejście numpy.dot

pierwsze, iteracji nad normalnych list Python:

$ python -mtimeit "import numpy as np; r = range(100)" "np.dot(r,r)" 
10 loops, best of 3: 316 usec per loop 

$ python -mtimeit "import operator; r = range(100); from itertools import izip, starmap" "sum(starmap(operator.mul, izip(r,r)))" 
10000 loops, best of 3: 81.5 usec per loop 

Następnie numpy ndarray:

$ python -mtimeit "import numpy as np; r = np.arange(100)" "np.dot(r,r)" 
10 loops, best of 3: 20.2 usec per loop 

$ python -mtimeit "import operator; import numpy as np; r = np.arange(100); from itertools import izip, starmap;" "sum(starmap(operator.mul, izip(r,r)))" 
10 loops, best of 3: 405 usec per loop 

Widząc to, wydaje numpy na NumPy tablic jest najszybszy, a następnie konstrukcje funkcjonalne Pythona pracujące z listami.

Answer 4

Proszę porównać tę funkcję "d2a" i porównaj ją z funkcją "d3".

from itertools import imap, starmap 
from operator import mul 

def d2a(v1,v2): 
    """ 
    d2a uses itertools.imap 
    """ 
    check(v1,v2) 
    return sum(imap(mul, v1, v2)) 

map (na Pythona 2.x, która jest co zakładam używasz) niepotrzebnie tworzy fikcyjny wykaz przed obliczeń.