Dlaczego dwa algorytmy wyszukiwania liczb pierwszych różnią się tak bardzo, chociaż wydają się wykonywać tę samą liczbę iteracji?

Question

Mam dwa algorytmy wyszukiwania liczb pierwszych w Pythonie. Wewnętrzna pętla każdego z nich wydaje się być wykonywana tyle samo razy i jest równie prosta. Jednak jeden z nich zajmuje 10 razy więcej niż drugi. Moje pytanie brzmi:

Dlaczego? Czy to jakieś dziwactwo Pythona, które można zoptymalizować (jak?), Czy też brakuje mi czegoś innego?

Problem, który rozwiązuję, jest w istocie od http://www.spoj.pl/problems/PRIME1/. W moim przypadku mam N = 10 ** 9, delta = 10 ** 5, i chcę znaleźć wszystkie liczby pierwotne między N-trójkąt a deltą. Mam też smallprimes, pre-made lista wszystkich liczb pierwszych mniejsza lub równa pierwiastek kwadratowy z N. Pierwszy algorytm jest bardzo prosty:

def range_f1(lo, hi, smallprimes): 
    """Finds all primes p with lo <= p <= hi. 

    smallprimes is the sorted list of all primes up to (at least) square root of hi. 
    hi & lo might be large, but hi-lo+1 miust fit into a long.""" 

    primes =[] 
    for i in xrange(hi-lo+1): 
    n = lo + i 

    isprime = True 
    for p in smallprimes: 
     if n % p == 0: 
      isprime = False 
      break 

    if isprime: 
     primes.append(n) 

    return primes 

Wywołanie range_f1(N-delta,N,smallprimes) trwa długo (około 10 sekund). Wewnętrzna pętla nazywa się 195170 razy. Mam również wersję tego algorytmu, która zastępuje pętlę ze zrozumieniem listy (to jest ta, której faktycznie używam do profilowania, patrz koniec pytania), ale nie jest to szybsze.

Druga wersja (brzydki realizacja) sita Eratostenesa:

def range_sieve(lo, hi, smallprimes): 
    """Parameters as before""" 

    # So ugly! But SO FAST!! How?? 

    delta = hi-lo+1 
    iamprime = [True] * delta  # iamprime[i] says whether lo + i is prime 
    if lo<= 1: 
    iamprime[1 - lo] = False 

    def sillyfun():  # For profiling & speed comparison 
    pass 

    for p in smallprimes: 
    rem = lo % p 
    if rem == 0: 
     iamprime[0] = False 
    for i in xrange(p - rem, delta, p): 
     iamprime[i] = False 
     sillyfun() 

    if p >= lo and p <= hi: 
     iamprime[p - lo] = True 

    return [p + lo for (p, ami) in enumerate(iamprime) if ami] 

to jest około 10 razy szybciej, zajmuje mniej niż 2 sekundy. Jednak wewnętrzna pętla (sillyfun()) jest nazywana 259982 razy, więcej niż w ostatnim przypadku. Nie jestem w stanie wyjaśnić, dlaczego jest to szybkie.

Pomyślałem, że być może powodem jest to, że wewnętrzna pętla pierwszego algorytmu zawiera modularną arytmetykę, podczas gdy druga ma tylko przypisanie. Jednak dodaje wydaje się sugerować, że zadanie nie jest szybsze niż modułowej arytmetyki:

>>> from timeit import timeit 
>>> timeit("10**9 % 2341234") 
0.23445186469234613 
>>> timeit("a[5000]=False", "a = [True] * 10000") 
0.47924750212666822 

---

tu jest (mniej czytelny) wersja pierwszego algorytmu faktycznie używać:

def range_f2(lo, hi, smallprimes): 

    primes =[] 
    for i in xrange(hi-lo+1): 
    n = lo + i 

    try: 
     (1 for p in smallprimes if n % p ==0).next() 
    except StopIteration: 
     primes.append(n) 

    return primes 

Oto wynik wywołania profilera dla range_f2() (należy zwrócić uwagę na liczbę generujących czas wyrażenia):

>>> from cProfile import run as prof 
>>> prof("range_f2(N-delta,N,sp)") 
200005 function calls in 13.866 CPU seconds 

Ordered by: standard name 

ncalls tottime percall cumtime percall filename:lineno(function) 
     1 0.000 0.000 13.866 13.866 <string>:1(<module>) 
195170 12.632 0.000 12.632 0.000 prime1.py:101(<genexpr>) 
     1 1.224 1.224 13.865 13.865 prime1.py:90(range_f2) 
    4832 0.009 0.000 0.009 0.000 {method 'append' of 'list' objects} 
     1 0.000 0.000 0.000 0.000 {method 'disable' of '_lsprof.Profiler' objects} 

Oto wynik profiler dla range_sieve():

>>> prof("range_sieve(N-delta,N,sp)") 
259985 function calls in 1.370 CPU seconds 

Ordered by: standard name 
ncalls tottime percall cumtime percall filename:lineno(function) 
    1 0.003 0.003 1.370 1.370 <string>:1(<module>) 
    1 0.877 0.877 1.367 1.367 prime1.py:39(range_sieve) 
259982 0.490 0.000 0.490 0.000 prime1.py:51(sillyfun) 
    1 0.000 0.000 0.000 0.000 {method 'disable' of '_lsprof.Profiler' objects} 

Wreszcie, tutaj jest on kompletny kod, który generuje listę małych liczb pierwszych (w bardzo głupi sposób), dzięki czemu można sprawdzić, jakie wyniki można uzyskać: http://pastebin.com/7sfN4mG4

UPDATE Według popularnego zapotrzebowania, dane profilowania dla pierwszego fragmentu kodu. Brak danych o tym, ile razy wykonywana jest pętla wewnętrzna, ale wydaje się całkiem jasne, że jest taka sama jak trzecia.

>>> prof("range_f1(N-delta,N,sp)") 
     4835 function calls in 14.184 CPU seconds 
Ordered by: standard name 
ncalls tottime percall cumtime percall filename:lineno(function) 
    1 0.000 0.000 14.184 14.184 <string>:1(<module>) 
    1 14.162 14.162 14.183 14.183 prime1.py:69(range_f1) 
    4832 0.021 0.000 0.021 0.000 {method 'append' of 'list' objects} 
    1 0.000 0.000 0.000 0.000 {method 'disable' of '_lsprof.Profiler' objects} 

Answer 1

Różnica polega na algorytmie.W pierwszej wersji, w wersji próbnej, testujesz każdego kandydata na wszystkie małe liczby pierwsze - nie zatrzymujesz się, gdy mała prime przekracza candidate ** 0.5, nie ma znaczenia dla range(10**9 - 10**5 ,10**9), jeśli smallprimes ma dobrą górną granicę, ale byłaby, gdyby długość zakresu były znacznie większe w stosunku do górnej granicy. W przypadku kompozytów nie wiąże się to z dużymi kosztami, ponieważ większość z nich ma co najmniej jeden mały dzielnik główny. Ale w liczbach pierwszych musisz przejść całą drogę do N**0.5. W tym przedziale jest około 10**5/log(10**9) liczb pierwszych, a każda z nich podzielona jest na próbę przez około 10**4.5/log(10**4.5) liczb pierwszych, co oznacza około 1.47*10**7 podziałów próbnych na liczby pierwsze.

Z drugiej strony za pomocą sita wykreśla się tylko kompozyty, każdy kompozyt jest przekreślany tyle razy, ile ma głównych dzielników, podczas gdy liczby pierwsze nie są w ogóle przekreślane (więc liczby pierwsze są wolne). Liczba głównych dzielników n jest ograniczona (wielokrotnością) log(n) (to surowe górne ograniczenie, zwykle znacznie przeszacowane), więc daje górną granicę 10**5*log(10**9) (razy mała stała) przekroczenia, około 2*10**6. Oprócz tego, przekroczenie może być mniej pracy niż podział (nie wiem o Pythonie, to jest dla tablic C). Więc robisz mniej pracy z sitkiem.

Edytuj: zebrał rzeczywiste liczby dla 10**9-10**5 na 10**9.

Ticks: 259987 
4832 
Divisions: 20353799 
4832 

Sito ma tylko 259.987 Przejścia-off, widać, że surowy górna granica powyżej jest przecenianie przez dużą czynnik. Podział próbny wymaga więcej niż 20 milionów dywizji, 16433632 tych dla liczb pierwszych (x/log x niedoszacowuje liczby liczb pierwszych, dla x = 10**4.5 około 10%), 3435183 są używane dla 3297 liczb w tym przedziale, którego najmniejszy czynnik główny jest większy niż n**(1/3).