Dlaczego pow (a, d, n) jest dużo szybsze niż a ** d% n?


110

Próbowałem zaimplementować test pierwszości Millera-Rabina i byłem zdziwiony, dlaczego trwa to tak długo (> 20 sekund) dla średnich liczb (~ 7 cyfr). Ostatecznie znalazłem następujący wiersz kodu jako źródło problemu:

x = a**d % n

(gdzie a,d i nsą podobne, ale nierówne, średnie liczby, **to operator potęgowanie i %jest operatorem modulo)

Następnie próbowałem zastąpić go następującym:

x = pow(a, d, n)

i to w porównaniu jest prawie natychmiastowe.

W kontekście, oto oryginalna funkcja:

from random import randint

def primalityTest(n, k):
    if n < 2:
        return False
    if n % 2 == 0:
        return False
    s = 0
    d = n - 1
    while d % 2 == 0:
        s += 1
        d >>= 1
    for i in range(k):
        rand = randint(2, n - 2)
        x = rand**d % n         # offending line
        if x == 1 or x == n - 1:
            continue
        for r in range(s):
            toReturn = True
            x = pow(x, 2, n)
            if x == 1:
                return False
            if x == n - 1:
                toReturn = False
                break
        if toReturn:
            return False
    return True

print(primalityTest(2700643,1))

Przykładowe obliczenia czasowe:

from timeit import timeit

a = 2505626
d = 1520321
n = 2700643

def testA():
    print(a**d % n)

def testB():
    print(pow(a, d, n))

print("time: %(time)fs" % {"time":timeit("testA()", setup="from __main__ import testA", number=1)})
print("time: %(time)fs" % {"time":timeit("testB()", setup="from __main__ import testB", number=1)})

Wyjście (uruchom z PyPy 1.9.0):

2642565
time: 23.785543s
2642565
time: 0.000030s

Dane wyjściowe (uruchomione z Pythonem 3.3.0, 2.7.2 zwraca bardzo podobne czasy):

2642565
time: 14.426975s
2642565
time: 0.000021s

I powiązane pytanie, dlaczego te obliczenia są prawie dwa razy szybsze, gdy są uruchamiane z Pythonem 2 lub 3 niż z PyPy, podczas gdy zwykle PyPy jest znacznie szybsze ?

Odpowiedzi:


164

Zobacz artykuł w Wikipedii dotyczący potęgowania modularnego . Zasadniczo, kiedy to zrobisz a**d % n, musisz obliczyć a**d, co może być dość duże. Ale są sposoby obliczania a**d % nbez konieczności obliczania a**dsamego siebie i tak właśnie powjest. **Operator nie może tego zrobić, ponieważ nie może „widzieć w przyszłości”, aby wiedzieć, że idziesz do natychmiastowego odbioru modułu.


14
+1 to właściwie to, co sugeruje >>> print pow.__doc__ pow(x, y[, z]) -> number With two arguments, equivalent to x**y. With three arguments, equivalent to (x**y) % z, but may be more efficient (e.g. for longs).
— ciąg dokumentów

6
W zależności od wersji Pythona może to być prawdą tylko w określonych warunkach. IIRC, w wersjach 3.x i 2.7, możesz używać tylko trzyargumentowej formy z typami całkowitymi (i nieujemną potęgą) i zawsze otrzymasz modułowe potęgowanie z inttypem natywnym , ale niekoniecznie z innymi typami całkowitymi. Ale w starszych wersjach istniały zasady dotyczące dopasowania do C long, dozwolona była forma z trzema argumentami floatitp. (Mam nadzieję, że nie używasz wersji 2.1 lub starszej i nie używasz żadnych niestandardowych typów całkowitych z modułów C, więc żadne ma to dla ciebie znaczenie.)
— abarnert

13
Z Twojej odpowiedzi wynika, że ​​kompilator nie może zobaczyć wyrażenia i zoptymalizować go, co nie jest prawdą. Tak się składa, że żaden z obecnych kompilatorów Pythona tego nie robi.
— danielkza

5
@danielkza: To prawda, nie chciałem sugerować, że jest to teoretycznie niemożliwe. Może „nie patrzy w przyszłość” byłoby bardziej trafne niż „nie patrzy w przyszłość”. Należy jednak pamiętać, że optymalizacja może być niezwykle trudna lub w ogóle niemożliwa. Dla stałych argumentów może być zoptymalizowane, ale x ** y % n, xmoże być obiektem narzędzia __pow__i na podstawie liczby losowej, powraca się z kilku różnych obiektów wykonawcze __mod__w sposób, który również zależeć od liczby losowe, itp
— BrenBarn

2
@danielkza: Ponadto funkcje nie mają tej samej domeny: .3 ** .4 % .5jest całkowicie legalne, ale jeśli kompilator przekształciłby to w pow(.3, .4, .5)to, wygenerowałoby rozszerzenie TypeError. Kompilator musiałby być w stanie wiedzieć, że a, di nsą gwarancją wartości integralną typu (a może właśnie konkretnie typu int, ponieważ transformacja nie pomaga w inny sposób), a dna pewno będzie nieujemna. To jest coś, co JIT mógłby zrobić, ale statyczny kompilator dla języka z typami dynamicznymi i bez wnioskowania po prostu nie może.
— abarnert

37

BrenBarn odpowiedział na twoje główne pytanie. Poza tym:

dlaczego jest prawie dwa razy szybszy, gdy działa z Pythonem 2 lub 3 niż PyPy, podczas gdy zwykle PyPy jest znacznie szybszy?

Jeśli przeczytasz stronę wydajności PyPy , to jest dokładnie w tym rodzaju rzeczy, w których PyPy nie jest dobry - w rzeczywistości jest to pierwszy przykład, który podają:

Złe przykłady obejmują wykonywanie obliczeń z dużymi długimi wartościami - które są wykonywane przez nieoptymalizowany kod pomocniczy.

Teoretycznie, przekształcenie ogromnego potęgowania, po którym następuje mod, w modułowe potęgowanie (przynajmniej po pierwszym przejściu) jest transformacją, którą JIT może wykonać… ale nie JIT PyPy.

Na marginesie, jeśli musisz wykonywać obliczenia z dużymi liczbami całkowitymi, możesz spojrzeć na moduły innych firm, takie jak gmpy, które czasami mogą być znacznie szybsze niż natywna implementacja CPythona w niektórych przypadkach poza głównymi zastosowaniami, a także ma wiele dodatkowej funkcjonalności, którą w innym przypadku musiałbyś napisać samodzielnie, kosztem mniejszej wygody.


2
długie zostały naprawione. wypróbuj pypy 2.0 beta 1 (nie będzie szybszy niż CPython, ale nie powinien też być wolniejszy). gmpy nie ma sposobu na obsłużenie MemoryError :(
— fijal

@fijal: Tak, gmpyw kilku przypadkach jest wolniejszy zamiast szybszego i sprawia, że ​​wiele prostych rzeczy jest mniej wygodnych. Nie zawsze taka jest odpowiedź - ale czasami tak jest. Warto więc przyjrzeć się temu, jeśli masz do czynienia z dużymi liczbami całkowitymi, a typ natywny Pythona nie wydaje się wystarczająco szybki.
— abarnert

1
a jeśli nie obchodzi cię, czy twoje liczby są duże, twój program nie
— działa

1
Jest to czynnik, który sprawił, że PyPy nie używał biblioteki GMP przez długi czas. To może być w porządku dla ciebie, nie jest w porządku dla programistów Python VM. Malloc może zawieść bez użycia dużej ilości pamięci RAM, po prostu umieść tam bardzo dużą liczbę. Zachowanie GMP od tego momentu jest nieokreślone i Python nie może na to pozwolić.
— fijal

1
@fijal: Całkowicie zgadzam się, że nie powinno się go używać do implementowania wbudowanego typu Pythona. To nie znaczy, że nigdy nie powinno być używane do niczego.
— abarnert

11

Istnieją skróty do wykonywania modularnego potęgowania: na przykład możesz znaleźć a**(2i) mod ndla każdego iod 1do log(d)i pomnożyć razem (mod n) potrzebne wyniki pośrednie. Dedykowana funkcja potęgowania modularnego, taka jak 3-argumentowa, pow()może wykorzystać takie sztuczki, ponieważ wie, że wykonujesz arytmetykę modularną. Parser Pythona nie może tego rozpoznać na podstawie samego wyrażenia a**d % n, więc wykona pełne obliczenie (co zajmie znacznie więcej czasu).


3

Sposób x = a**d % nobliczania to podniesienie ado dpotęgi, a następnie modulo to z n. Po pierwsze, jeśli ajest duża, tworzy to ogromną liczbę, która jest następnie obcięta. Jednak x = pow(a, d, n)najprawdopodobniej jest zoptymalizowany tak, że nśledzone są tylko ostatnie cyfry, które są wszystkim, co jest wymagane do obliczenia mnożenia modulo liczba.


6
„wymaga mnożenia d do obliczenia x ** d” - niepoprawne. Możesz to zrobić w mnożeniach O (log d) (bardzo szerokich). Potęgowanie przez podniesienie do kwadratu może być używane bez modułu. Sam rozmiar mnożników jest tym, co tu prowadzi.
— John Dvorak

@JanDvorak Prawda, nie jestem pewien, dlaczego myślałem, że Python nie użyłby tego samego algorytmu potęgowania **co dla pow.
— Yuushi,

5
Nie ostatnie cyfry „n” ... po prostu przechowuje obliczenia w Z / nZ.
— Thomas
Korzystając z naszej strony potwierdzasz, że przeczytałeś(-aś) i rozumiesz nasze zasady używania plików cookie i zasady ochrony prywatności.
Licensed under cc by-sa 3.0 with attribution required.