Znajdź największą liczbę pierwszą, której długość, suma i iloczyn jest liczbą pierwszą

37

Liczba 113jest pierwszą liczbą pierwszą, której długość 3jest liczbą pierwszą, suma cyfrowa 5 = 1 + 1 + 3jest liczbą pierwszą, a produkt cyfrowy 3 = 1 * 1 * 3jest liczbą pierwszą.

Liczba pierwsza, która ma te 3 właściwości, będzie nazywana najwyższą liczbą pierwszą . Liczby pierwsze 11117i 1111151inne przykłady.

Cel

Napisz program, który może znaleźć największą możliwą liczbę pierwszą w niespełna godzinę na przyzwoitym nowoczesnym komputerze osobistym (takim jak preferowana specyfikacja tutaj ).

Nie powinieneś po prostu dać nam dużej najwyższej liczby. Musisz pokazać nam swój proces wyszukiwania za pomocą kodu, który faktycznie działa. Możesz opierać się na rozwiązaniach swoich lub innych osób, ale pamiętaj, aby dać im uznanie. Wspólnie staramy się znaleźć największą najwyższą liczbę pierwszą na zwykłym komputerze w ciągu godziny.

Punktacja

Wygrana, która znajdzie największą najwyższą liczbę pierwszą. Jeśli okaże się, że istnieje ostatecznie wiele najwyższych liczb pierwszych, wygrywa pierwsze podporządkowanie, które generuje najwyższą liczbę pierwszą.

(Jeśli możesz matematycznie udowodnić, że istnieje albo nie ma nieskończenie wiele najwyższych liczb pierwszych, dam ci 200 nagród za to tylko dlatego :)))

Detale

  • Możesz użyć dowolnego źródła do wygenerowania liczb pierwszych (np. Internetu).
  • Możesz użyć probabilistycznych metod testowych.
  • Wszystko jest w bazie 10.
  • Zero i jeden NIE są uważane za pierwsze.
  • Liczby pierwsze zawierające 0mają tak cyfrowy produkt, 0że oczywiście nie mogą być najwyższe.
  • Aby strona była mniej zaśmiecona, umieść duże (ponad 100 cyfr) najwyższe liczby pierwsze w postaci:

    {[number of 1's before the prime digit]}[prime digit]{[number of 1's after the prime digit]}
    

    Więc 1111151może być wyrażona jako {5}5{1}.

Hobby Calvina
źródło
Czy możemy zacząć od listy liczb pierwszych lub pobrać listę z Internetu i spędzić godzinę na sprawdzaniu supremacji?
Sparr
2
Jeśli możesz zacząć od najwyższej znanej najwyższej liczby pierwszej, staje się to wyzwaniem dla tego, kto może napisać program, który spędza dokładnie godzinę, obejmując największą możliwą lukę między dwoma najwyższymi liczbami pierwszymi. :(
Sparr
8
Poza niezawieraniem zera, każda potencjalna najwyższa liczba pierwsza musi oczywiście mieć postać 1 ^ n [3 | 5 | 7] 1 ^ m, tj. Niektóre 1, dowolna liczba pierwsza poniżej 10 i niektóre 1. Istnieje więcej ograniczeń, które możesz nałożyć od razu.
Ingo Bürk
3
Ryan rozpoczął powiązane pytanie dotyczące MSE dotyczące istnienia nieskończenie wielu najwyższych liczb pierwszych. Jeśli masz jakieś spostrzeżenia na temat tego pytania, proszę się zastanowić!
Semiclassical
1
Mogę łatwo wykazać, że obecnie nie ma dowodu na nieskończoną liczbę najwyższych liczb pierwszych (i że poświęcono temu znaczną część pracy). Według michaelnielsen.org/polymath1/… wiemy, że liczby pierwsze występują nieskończenie z przerwami tak małymi jak 246, ale dla dowodu nieskończonych najwyższych liczb pierwszych potrzebujemy przerwy 2, 4 lub 6 (odpowiadającej liczbom pierwszym z gdzieś w nich 3, 5 lub 7).
Tim S.

Odpowiedzi:

9

Perl, 15101 cyfr, {83} 7 {15017}, 8 minut. Max znaleziono: 72227 cyfr

Korzystanie z mojego modułu Math :: Prime :: Util i jego zaplecza GMP . Posiada szereg testów złożoności, w tym is_prob_prime () z testem ES BPSW (nieco bardziej rygorystycznym niż ispseudoprime Pari), is_prime (), który dodaje jeden losowy MR, i is_provable_prime (), który uruchomi BLS75 T5 lub ECPP. Przy tych rozmiarach i typach wykonanie dowodu zajmie dużo czasu. Wrzuciłem kolejny test MR w sub weryfikatorze pedantycznym. Czasy na Core2 E7500, który zdecydowanie nie jest moim najszybszym komputerem (zajmuje to 2,5 minuty na moim i7-4770K).

Jak zauważa Tim S., moglibyśmy szukać większych wartości, aż do momentu, gdy pojedynczy test zajmie godzinę. Przy ~ 15000 cyfr na tym E7500 potrzeba około 26 sekund na test MR i 2 minuty na pełny is_prime (podział próbny plus podstawowa 2 MR plus ES Lucas plus jedna losowa podstawowa MR). Mój i7-4770K jest ponad 3 razy szybszy. Próbowałem kilku rozmiarów, głównie widząc, jak to zrobiło na wynikach innych ludzi. Próbowałem 8k, 20k i 16k, zabijając każde po ~ 5 minutach. Następnie próbowałem 15k w progresji na ~ 10m każdy i miałem szczęście na czwartym.

Testy PRP OpenPFGW są z pewnością szybsze po około 4000 cyfr i znacznie szybsze w zakresie 50k +. Jego test ma jednak sporo fałszywych wyników pozytywnych, co czyni go doskonałym testem wstępnym, ale nadal chcielibyśmy zweryfikować wyniki za pomocą czegoś innego.

Można to zrównoleglić za pomocą wątków perla lub przy użyciu MCE podobnego do równoległych przykładów Fibonacciego w module.

Czas i wyniki na bezczynnym komputerze i7-4770K z pojedynczym rdzeniem:

  • wprowadź 3000, 16 sekund, 3019 cyfr, {318} 5 {2700}
  • wprowadź 4000, 47 sekund, 4001 cyfr, {393} 7 {3607}
  • wprowadź 4100, 5 sekund, 4127 cyfr, {29} 7 {4097}
  • wprowadź 6217, 5 sekund, 6217 cyfr, {23} 5 {6193}
  • wprowadź 6500, 5 minut, 6547 cyfr, {598} 5 {5948}
  • wprowadź 7000, 15 minut, 7013 cyfr, {2411} 7 {4601}
  • wprowadź 9000, 11 minut, 9001 cyfr, {952} 7 {8048}
  • wprowadź 12000, 10 minut, 12007 cyfr, {652} 5 {11354}
  • wprowadź 15100, 2,5 minuty, 15101 cyfr, {83} 7 {15017}
  • wprowadź 24600, 47 minut, 24671 cyfr, {621} 7 {24049}
  • wprowadź 32060, 18 minut, 32063 cyfr, {83} 7 {31979}
  • wprowadź 57000, 39 minut, 57037 cyfr, {112} 5 {56924}
  • wprowadź 72225, 42 minuty, 72227 cyfr, {16} 3 {72210}

Aby uzyskać wynik 32 000 cyfr, uruchomiłem 6 skryptów działających jednocześnie z każdym kolejnym argumentem zaczynającym się od 32 000. Po 26,5 minutach jeden skończył się z wyświetlonym wynikiem 32063 cyfr. Za 57k pozwalam kolejnym skryptom uruchamiać 6 jednocześnie na godzinę z przyrostem wejściowym 500, aż wynik 57k powróci za 57 minut. 72-cyfrowy wynik został znaleziony przez wykonanie kolejnych liczb pierwszych od 70k w górę, więc na pewno nie został znaleziony w ciągu godziny (chociaż kiedy już wiesz, od czego zacząć, to jest).

Scenariusz:

#!/usr/bin/env perl
use warnings;
use strict;
use Math::Prime::Util qw/:all/;
use Math::Prime::Util::GMP;  # Just to ensure it is used.

my $l = shift || 1000;  $l--;

while (1) {
  $l = next_prime($l);
  my @D = grep { is_prime($l-1 + $_) } (3,5,7);
  next unless scalar @D > 0;
  for my $s (0 .. $l-1) {
    my $e = $l-$s-1;
    warn "   checking $l $s\n" unless $s % 100;
    for my $d (@D) {
      my $n = "1"x$s . $d . "1"x$e;
      die unless length($n) == $l;
      verify_supreme($n,$s,$d,$e) if is_prime($n);  # ES BPSW + 1 rand-base M-R
    }
  }
}
sub verify_supreme {  # Be pedantic and verify the result
  my($n,$s,$d,$e) = @_;
  die "Incorrect length" unless is_prime(length($n));
  die "Incorrect sum" unless is_prime(vecsum(split(//,$n)));
  my $prod = 1; $prod *= $_ for split(//,$n);
  die "Incorrect product" unless is_prime($prod);
  die "n is not a prime!" unless miller_rabin_random($n,1);  # One more M-R test
  die "{$s} $d {$e}\n";
}
DanaJ
źródło
+1 za wprowadzenie mnie do tej biblioteki! Czasy na moim komputerze do iteracji liczb pierwszych mniejszych niż 10 ^ 7 (w porównaniu do CPython z gmpy2i PyPy z my_math): codepad.org/aSzc0esT
primo
Cieszę się ze to lubisz! Istnieją inne sposoby, w tym forprimes { ...do stuff... } 1e7;10-krotnie szybsze (podziękowania dla Pari / GP za wiele świetnych pomysłów). Zawsze doceniam opinie, więc daj mi znać, jeśli coś nie działa tak, jak chcesz.
DanaJ
21

Python 2.7 na PyPy, {2404} 3 {1596} (~ 10 ^ 4000)

11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111113111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111

Znalazłem go około 50 minut po uruchomieniu od 4000. Dlatego oszacowałbym, że jest to górna granica tego podejścia do kodu.

Zmiana: zauważyłem, że niektóre długości wydają się bardziej owocne dla generowania tego rodzaju liczby pierwszej niż inne, dlatego postanowiłem sprawdzić tylko 50 losowych lokalizacji cyfry, która nie jest 1 zamiast wszystkich możliwych lokalizacji, przed przeniesieniem na. Nie jestem do końca pewien, czy poprawi to wydajność lub czy 50 jest poprawne, ale zobaczymy.

Lista możliwości jest generowana na podstawie faktu, że aby spełniony był wymóg dotyczący produktu, liczba musi być równa wszystkim oprócz liczby pierwszej. Ponadto liczba pierwsza nie może być równa 2 ze względu na stosunek sumy i długości, a suma cyfrowa nie może być podzielna przez trzy, co daje wymaganie% 3.

is_prime pochodzi z http://codepad.org/KtXsydxK , napisane przez @primo

Uwaga: ta funkcja is_prime jest tak naprawdę pseudopierwszym testem Baillie-PSW, ale nie ma znanych przeciw-przykładów, więc nie będę się martwił o to rozróżnienie.

#http://codepad.org/KtXsydxK
from my_math import is_prime
import time,random
LLIMIT=2748
time.clock()
start_time=time.time()
checked=0
while time.time()-start_time<3600:
    small_primes = [a for a in range(LLIMIT,2*LLIMIT) if is_prime(a)]
    leng,dig=(0,0)
    for a in small_primes:
        if a+2 in small_primes:
            leng,dig=(a,3)
            break
        if a+4 in small_primes:
            leng,dig=(a,5)
            break
        if a+6 in small_primes:
            leng,dig=(a,7)
            break
    start=time.clock()
    print leng,dig,time.clock(),checked
    for loc in random.sample(range(leng),50):
        checked+=1
        if is_prime(int('1'*loc+str(dig)+'1'*(leng-loc-1))):
            print leng-1,loc,dig,time.clock(),time.clock()-start, \
                  int('1'*loc+str(dig)+'1'*(leng-loc-1))
            break
    LLIMIT=leng+1
isaacg
źródło
Niestety nie znam nic poza linkiem. Znalazłem link tutaj: codegolf.stackexchange.com/questions/10739/… Pierwsza odpowiedź
isaacg
No więc. Uznam cię.
isaacg
10
To jak ASCII, gdzie jest wally ...
trichoplax
5
Może powinieneś zmienić nazwę funkcji is_very_very_very_very_very_very_very_probably_prime()...
trichoplax
2
Mathmatica i Klon używają tej samej metody, więc nie może być aż tak źle.
primo
13

PARI / GP, 4127 cyfr

(10 4127 -1) / 9 + 2 * 10515

Jest to dość proste wyszukiwanie: sprawdź tylko długości liczb pierwszych, a następnie oblicz możliwe liczby pierwsze do użycia, a następnie powtórz wszystkie możliwości. Specjalnie opisałem typowe przypadki, w których do użycia jest 0 lub 1 odpowiednia liczba pierwsza.

supreme(lim,startAt=3)={
    forprime(d=startAt,lim,
        my(N=10^d\9, P=select(p->isprime(d+p),[1,2,4,6]), D, n=1);
        if(#P==0, next);
        if(#P==1,
            for(i=0,d-1,
                if (ispseudoprime(D=N+n*P[1]), print(D));
                n*=10
            );
            next
        );
        D=vector(#P-1,i,P[i+1]-P[i]);
        for(i=0,d-1,
            forstep(k=N+n*P[1],N+n*P[#P],n*D,
                if (ispseudoprime(k), print(k))
            );
            n*=10
        )
    )
};
supreme(4200, 4100)

Obliczenie zajęło 36 minut na jednym rdzeniu dość starej maszyny. Jestem pewien, że nie byłoby problemu ze znalezieniem takiej liczby ponad 5000 cyfr w ciągu godziny, ale jestem też niecierpliwy.

Lepszym rozwiązaniem byłoby użycie dowolnego rozsądnego języka do zrobienia wszystkiego poza najbardziej wewnętrzną pętlą, a następnie zbudowanie pliku abc dla pierwszej formy, który jest zoptymalizowany dla tego rodzaju obliczeń. Powinno to umożliwić zwiększenie obliczeń do co najmniej 10 000 cyfr.

Edycja : wdrożyłem opisane powyżej rozwiązanie hybrydowe, ale na mojej starej maszynie nie mogę znaleźć pierwszego terminu z> = 10 000 cyfr w mniej niż godzinę. O ile nie uruchomię tego na czymś szybszym, będę musiał zmienić cel na mniej wzniosły.

Charles
źródło
Skąd wiedziałeś, że możesz zacząć od 4100?
isaacg
@isaacg: Właśnie próbowałem być większy niż (niepoprawne) rozwiązanie Mathematica, które było nieco ponad 4000. Po prostu poszedłem do następnej wielokrotności 100 jako „nic-up-my-sleeve”. Właściwie wydaje się, że było to niefortunne miejsce startu, ponieważ musiałem iść dłużej niż się spodziewałem (i dłużej niż Mathematica!), Aby znaleźć liczbę pierwszą.
Charles,
Nie, miałeś niewiarygodne szczęście. (4127,3) jest pierwszą parą po 4100 i przypadkiem ma pierwszą. Wiele par w ogóle nie ma liczb pierwszych.
isaacg
@isaacg: Może tak. Moja heurystyka jest wyraźnie wyłączona, ponieważ oczekiwałbym prawdopodobieństwa ~ 80% znalezienia liczby pierwszej w danej parze: 1 - exp (-15 / (4 * log 10)), ale wydaje się, że są rzadsze, więc nie zachowuj się jak losowe {2, 3, 5} - gładkie liczby ich wielkości (chyba że wygłuszyłem obliczenia).
Charles
@isaacg: W każdym razie pracuję nad „lepszym rozwiązaniem”, o którym wspominałem: przenosząc ciężką pracę na pfgw. Przeszukano pierwsze 20 par powyżej 10 ^ 10000, nie znajdując niczego, ale zajęło to tylko ~ 15 minut.
Charles,
7

Mathematica 3181 cyfr

Aktualizacja: W moim pierwszym zgłoszeniu wystąpiły poważne błędy. Byłem w stanie poświęcić trochę czasu na sprawdzenie wyników dla tego. Dane wyjściowe są sformatowane jako lista cyfr. Dzięki temu można łatwo sprawdzić każdy z warunków.

f[primeDigitLength_]:=
Module[{id=ConstantArray[1,primeDigitLength-1]},
permutations=Reverse@Sort@Flatten[Table[Insert[id,#,pos],{pos,primeDigitLength}]&/@{3,5,7},1];
Flatten[Select[permutations,PrimeQ[FromDigits[#]]\[And]PrimeQ[Plus@@#]&,1],1]]

Przykład

To był mój pierwszy test, poszukiwanie rozwiązania z 3181 cyframi. Pierwszy przypadek znaleziono w 26 sekund.

Przejdźmy do rozumowania. Następnie przejdziemy do samego programu.

Zacznijmy, tak jak ja, „jaka jest czwarta liczba pierwsza?” Czy możemy znaleźć rozwiązanie z tyloma cyframi (3181)?

primeDigits = Prime[450]

3181


Numer można znaleźć, łącząc cyfry.

number = FromDigits[digits];

Zamiast wyświetlać, zamiast tego możemy zapytać, jakie są cyfry i gdzie się znajdują.

DigitCount[number]

{3180, 0, 0, 0, 0, 0, 1, 0, 0, 0}

Oznacza to, że było 3180 wystąpień cyfry 1 i jedno wystąpienie cyfry 7.

Na której pozycji jest cyfra 7?

Position[digits, 7][[1, 1]]

142

Cyfra 7 to 142. cyfra. Wszystkie pozostałe to 1.


Oczywiście iloczyn cyfr musi być liczbą pierwszą, a mianowicie 7.

digitProduct = Times @@ digits

7


A suma cyfr jest również liczbą pierwszą.

digitSum = Plus @@ digits
PrimeQ[digitSum]

3187
Prawda


I wiemy, że liczba cyfr jest liczbą pierwszą. Pamiętaj, że wybraliśmy 450. pierwszą, a mianowicie 3118.

Tak więc wszystkie warunki zostały spełnione.

DavidC
źródło
3
Jeśli się nie mylę, jego suma wynosi 4009, co nie jest liczbą pierwszą.
gerw
Jedno: czy nie powinna to być suma liczb pierwszych, a nie liczba cyfr? W twoim przypadku musisz przetestować, 4002 * 1 + 7 = 4009a nie 4003 zgodnie ze specyfikacją.
Johnride
2
@Johnride: Oba powinny być pierwsze.
gerw
@ gerw ma rację. Liczba cyfr ORAZ suma cyfr ORAZ iloczyn cyfr muszą być liczbą pierwszą.
Calvin's Hobbies
Wszyscy mieliście rację. W moim wcześniejszym zgłoszeniu zapomniałem sprawdzić sumę cyfr pod kątem pierwszeństwa. Odbywa się to teraz poprzez sprawdzenie, czy jedna z poniższych (nie ma znaczenia, która) jest liczbą pierwszą: długość cyfry + 2, długość cyfry _4 lub długość cyfry +6.
DavidC
7

Python 2.7, 6217 cyfr: {23} 5 {6193} 6 minut 51 sekund

Pracowałem nad własną wersją i byłem rozczarowany, widząc, że @issacg pobił mnie do sedna bardzo podobnym podejściem, aczkolwiek z is_ (bardzo_prawdopodobnie) _prime (). Widzę jednak, że mam pewne znaczące różnice, które skutkują lepszą odpowiedzią w krótszym czasie (gdy używam również is_prime). Aby to wyjaśnić, kiedy zaczynam od 4000, uzyskuję lepszą 4001-cyfrową odpowiedź ({393} 7 {3607}) w zaledwie 26 minut i 37 sekund przy użyciu standardowego interpretera Pythona (również w wersji 2.7), a nie PyPy wersja. Poza tym nie sprawdzam liczb na miejscu. Wszyscy kandydaci są sprawdzani.

Oto główne ulepszenia:

  1. Użyj generatora liczb pierwszych ( https://stackoverflow.com/questions/567222/simple-prime-generator-in-python/568618#568618 ), aby utworzyć listę liczb pierwszych do sprawdzenia i (jego wersja „małych liczb pierwszych”) i do generowania dopuszczalnych długości liczb.

  2. Chcemy spędzać czas na szukaniu największej liczby pierwszej o danej długości, a nie najmniejszej, dlatego najpierw sprawdzam możliwie największe liczby, a nie najmniejsze. Następnie, gdy zostanie znaleziony, możemy natychmiast przejść do następnej długości.

EDYCJA: Teraz z wieloprocesowym przetwarzaniem

Jest to znacząca zmiana w stosunku do poprzednich wersji. Wcześniej zauważyłem, że moja 8-rdzeniowa maszyna prawie nie działa, więc postanowiłem spróbować swoich sił w wieloprocesowym przetwarzaniu w Pythonie (pierwszy raz). Wyniki są bardzo miłe!

W tej wersji odradza się 7 procesów potomnych, które chwytają „zadanie” z kolejki potencjalnych możliwości (num_length + odpowiednie cyfry). Przebijają się, próbując różnych [7,5,3] pozycji, aż znajdzie jedną. Jeśli tak, informuje proces główny o nowej najdłuższej znalezionej długości. Jeśli dzieci pracują na długości num_l, która jest krótsza, po prostu za kaucją idą do następnej długości.

Rozpocząłem ten bieg z 6000 i nadal działa, ale jak dotąd jestem bardzo zadowolony z wyników.

Program nie zatrzymuje się jeszcze poprawnie, ale dla mnie nie jest to wielka sprawa.

Teraz kod:

#!/usr/bin/env python
from __future__ import print_function

import sys
from multiprocessing import Pool, cpu_count, Value
from datetime import datetime, timedelta

# is_prime() et al from: http://codepad.org/KtXsydxK - omitted for brevity
# gen_primes() from: https://stackoverflow.com/questions/567222/simple-prime-generator-in-python/568618#568618 - ommitted for brevity
from external_sources import is_prime, gen_primes


def gen_tasks(start_length):
    """
    A generator that produces a stream of eligible number lengths and digits
    """
    for num_length in gen_primes():
        if num_length < start_length:
            continue

        ns = [ n for n in [7,5,3] if num_length + n - 1 in prime_list ]
        if ns:
            yield (num_length, ns)


def hunt(num_length, ns):
    """
    Given the num_length and list of eligible digits to try, build combinations
    to try, and try them.
    """

    if datetime.now() > end_time or num_length <= largest.value:
        return

    print('Tasked with: {0}, {1}'.format(num_length, ns))
    sys.stdout.flush()
    template = list('1' * num_length)
    for offset in range(num_length):
        for n in ns:
            if datetime.now() > end_time or num_length <= largest.value:
                return

            num_list = template[:]
            num_list[offset] = str(n)
            num = int(''.join(num_list))

            if is_prime(num):
                elapsed = datetime.now() - start_time
                largest.value = num_length
                print('\n{0} - "{1}"\a'.format(elapsed, num))


if __name__ == '__main__':
    start_time = datetime.now()
    end_time = start_time + timedelta(seconds=3600)

    print('Starting @ {0}, will stop @ {1}'.format(start_time, end_time))

    start_length = int(sys.argv[1])

    #
    # Just create a list of primes for checking. Up to 20006 will cover the first
    # 20,000 digits of solutions
    #
    prime_list = []
    for prime in gen_primes():
        prime_list.append(prime)
        if prime > 20006:
            break;
    print('prime_list is primed.')

    largest = Value('d', 0)

    task_generator = gen_tasks(start_length)

    cores = cpu_count()
    print('Number of cores: {0}'.format(cores))


    #
    # We reduce the number of cores by 1 because __main__ is another process
    #
    pool = Pool(processes=cores - 1)

    while datetime.now() < end_time:
        pool.apply_async(hunt, next(task_generator))
mkoistinen
źródło
czytałby bardziej czytelnie, jeśli reprezentujesz link do codepad jako [uszkodzony, jeśli to konieczne] import
Sparr
Myślę, że byłoby to mylące, ponieważ kod na drugim końcu tak naprawdę nie jest tak importowalny.
mkoistinen
użyj składni isaacga. skomentuj adres URL, a następnie zaimportuj z nieistniejącej paczki (w moim przypadku moja_math)
Sparr
1
Właściwie generuję też liczby od największej do najmniejszej. Nie sądzę, aby nasze różnice w kodzie były bardzo znaczące. Spodziewałbym się, że różnica tkwi bardziej w naszych komputerach niż w kodzie. Niemniej jednak dobrze zrobione i witamy na stronie.
isaacg
my_mathmoże również służyć do generowania listy liczb pierwszych, à la while prime < 20006: prime = next_prime(prime). Wydaje się być około 3 razy szybszy gen_primesi znacznie wydajniejszy.
primo
6

C, GMP - {7224} 5 {564} = 7789

Uznanie dla @issacg i was wszystkich za inspiracje i triki.
A także mistrzowskie pytanie zadające @ Calvin's Hobbies dla tego pytania.

Skompilować: gcc -I/usr/local/include -o p_out p.c -pthread -L/usr/local/lib -lgmp

Jeśli masz ochotę przekazać swoją moc obliczeniową lub jesteś ciekawy wydajności, skopiuj kod i skompiluj. ;) Będziesz musiał zainstalować GMP.

#include<stdio.h>
#include<stdlib.h>
#include<sys/time.h>
#include<gmp.h>
#include<pthread.h>

#define THREAD_COUNT 1
#define MAX_DIGITS   7800
#define MIN_DIGITS   1000

static void huntSupremePrime(int startIndex) {

    char digits[MAX_DIGITS + 1];

    for (int i = 0; i < MAX_DIGITS; digits[i++] = '1');

    digits[MAX_DIGITS] = '\0';
    mpz_t testPrime, digitSum, digitCount, increment;

    for (int j = 0; j < MAX_DIGITS - startIndex; digits[j++] = '0');

    int step = THREAD_COUNT * 2;

    for (int i = startIndex, l = MAX_DIGITS - startIndex; i > MIN_DIGITS - 1; 
        i -= step, l += step) {
        fprintf(stderr, "Testing for %d digits.\n", i);
        mpz_init_set_ui(digitCount, i);
        if (mpz_millerrabin(digitCount, 10)) {
            for (int j = 3; j < 8; j += 2) {
                mpz_init_set_ui(digitSum, i - 1 + j);
                if (mpz_millerrabin(digitSum, 10)) {
                    gmp_printf("%Zd \n", digitSum);
                    digits[MAX_DIGITS - 1] = j + 48;
                    mpz_init_set_str(testPrime, digits, 10);
                    mpz_init_set_ui(increment, (j - 1) * 99);
                    for (int k = 0; k < i/20; ++k) {
                        if (mpz_millerrabin(testPrime, 25)) {
                            i = 0;
                            j = 9;
                            k = l;
                            gmp_printf("%Zd\n", testPrime);
                            break;
                        }
                        mpz_add(testPrime, testPrime, increment);
                        mpz_mul_ui(increment, increment, 100);
                        fprintf(stderr, "TICK %d\n", k);
                    }

                }
            }
        }
        for (int j = 0; j < step; digits[l + j++] = '0');

    }
}

static void *huntSupremePrimeThread(void *p) {
    int* startIndex = (int*) p;
    huntSupremePrime(*startIndex);
    pthread_exit(NULL);
}

int main(int argc, char *argv[]) {

    int  startIndexes[THREAD_COUNT];
    pthread_t threads[THREAD_COUNT];

    int startIndex = MAX_DIGITS;
    for (int i = 0; i < THREAD_COUNT; ++i) {
        for (;startIndex % 2 == 0; --startIndex);
        startIndexes[i] = startIndex;
        int rc = pthread_create(&threads[i], NULL, huntSupremePrimeThread, (void*)&startIndexes[i]); 
        if (rc) { 
            printf("ERROR; return code from pthread_create() is %d\n", rc);
            exit(-1);
        }
        --startIndex;
    }

    for (int i = 0; i < THREAD_COUNT; ++i) {
        void * status;
        int rc = pthread_join(threads[i], &status);
        if (rc) {
            printf("ERROR: return code from pthread_join() is %d\n", rc);
            exit(-1);
        }
    }

    pthread_exit(NULL);
    return 0;
}
Vectorized
źródło
5

PFGW, 6067 cyfr, {5956} 7 {110}

Uruchom program PFGW z następującym plikiem wejściowym i -f100wprowadź liczby wstępne. Po około 2-3 minutach procesora na moim komputerze (i5 Haswell) znajduje PRP (10 ^ (6073-6) -1) / 9 + 6 * 10 ^ 110 lub {5956} 7 {110} . Jako punkt początkowy wybrałem 6000 cyfr jako liczbę „nic w mojej kieszeni”, która jest nieco wyższa niż we wszystkich poprzednich zgłoszeniach.

ABC2 $a-$b & (10^($a-$b)-1)/9+$b*10^$c
a: primes from 6000 to 6200
b: in { 2 4 6 }
c: from 0 to 5990

Na podstawie tego, jak szybko udało mi się go znaleźć, mogłem z łatwością zwiększyć liczbę cyfr i nadal znaleźć PRP w ciągu godziny. Dzięki temu, jak napisane są reguły, mogę nawet znaleźć rozmiar, w którym mój procesor, działający na wszystkich 4 rdzeniach, może ukończyć jeden test PRP w ciągu godziny, poświęcić dużo czasu na znalezienie PRP i mieć moje „wyszukiwanie” wyłącznie jednego PRP.

PS W pewnym sensie nie jest to rozwiązanie „kodowe”, ponieważ nie napisałem niczego poza plikiem wejściowym ... ale wiele rozwiązań Mathematica dla problemów matematycznych w jednym wierszu można opisać w ten sam sposób, jak można używając biblioteki, która wykonuje dla Ciebie ciężką pracę. W rzeczywistości myślę, że trudno jest narysować dobrą linię między nimi. Jeśli chcesz, mogę napisać skrypt, który tworzy plik wejściowy PFGW i wywołuje PFGW. Skrypt może nawet wyszukiwać równolegle, aby użyć wszystkich 4 rdzeni i przyspieszyć wyszukiwanie ~ 4 razy (na moim CPU).

PPS Myślę, że LLR może wykonać testy PRP dla tych liczb i spodziewałbym się, że będzie to znacznie szybsze niż PFGW . Dedykowany program przesiewania mógłby również lepiej uwzględniać te liczby niż jednorazowy PFGW. Jeśli je połączysz, jestem pewien, że możesz przekroczyć granice znacznie wyżej niż obecne rozwiązania.

Tim S.
źródło
4

Python 2.7, 17-19 cyfr

11111111171111111

Znaleziono 5111111111111 (13 cyfr) w 3 sekundy, a ta 17-cyfrowa najwyższa liczba pierwsza w 3 minuty. Domyślam się, że maszyna docelowa może to uruchomić i uzyskać 19-cyfrową najwyższą liczbę pierwszą w niecałą godzinę. Podejście to nie skaluje się dobrze, ponieważ utrzymuje liczby pierwsze w pamięci do połowy liczby cyfr docelowych. Wyszukiwanie 17-cyfrowe wymaga przechowywania tablicy 100 milionów boolanów. 19 cyfr wymagałoby tablicy elementów 1B, a pamięć byłaby wyczerpana przed przejściem do 23 cyfr. Prawdopodobnie też będzie to środowisko uruchomieniowe.

Podejścia oparte na testach pierwszeństwa, które nie obejmują absurdalnie dużej liczby liczb pierwszych dzielników, wypadną znacznie lepiej.

#!/usr/bin/env python
import math
import numpy as np
import sys

max_digits = int(sys.argv[1])
max_num = 10**max_digits

print "largest supreme prime of " + str(max_digits) + " or fewer digits"

def sum_product_digits(num):
    add = 0
    mul = 1
    while num:
         add, mul, num = add + num % 10, mul * (num % 10), num / 10
    return add, mul

def primesfrom2to(n):
    # http://stackoverflow.com/questions/2068372/fastest-way-to-list-all-primes-below-n-in-python/3035188#3035188
    """ Input n>=6, Returns a array of primes, 2 <= p < n """
    sieve = np.ones(n/3 + (n%6==2), dtype=np.bool)
    sieve[0] = False
    for i in xrange(int(n**0.5)/3+1):
        if sieve[i]:
            k=3*i+1|1
            sieve[      ((k*k)/3)      ::2*k] = False
            sieve[(k*k+4*k-2*k*(i&1))/3::2*k] = False
    return np.r_[2,3,((3*np.nonzero(sieve)[0]+1)|1)]

def checkprime(n):
    for divisor in primes:
        if (divisor>math.sqrt(n)):
            break
        if n%divisor==0:
            return False
    return True

# make an array of all primes we need to check as divisors of our max_num
primes = primesfrom2to(math.sqrt(max_num))
# only consider digit counts that are prime
for num_digits in primes:
    if num_digits > max_digits:
        break
    for ones_on_right in range(0,num_digits):
        for mid_prime in ['3','5','7']:
            # assemble a number of the form /1*[357]1*/
            candidate = int('1'*(num_digits-ones_on_right-1)+mid_prime+'1'*ones_on_right)
            # check for primeness of digit sum first digit product first
            add, mul = sum_product_digits(candidate)
            if add in primes and mul in primes:
                # check for primality next
                if checkprime(candidate):
                    # supreme prime!
                    print candidate
Sparr
źródło
3

Mathematica 4211 4259 cyfr

Z numerem: {1042} 7 {3168} {388} 3 {3870}

Który został wygenerowany przez następujący kod:

TimeConstrained[
 Do[
  p = Prime[n];
  curlargest = Catch[
    If[PrimeQ[p + 6],
     list = ConstantArray[1, p];
     Do[
      temp = FromDigits[ReplacePart[list, i -> 7]];
      If[PrimeQ[temp],
       Throw[temp]
       ], {i, p}]
     ];

    If[PrimeQ[p + 4],
     list = ConstantArray[1, p];
     Do[
      temp = FromDigits[ReplacePart[list, i -> 5]];
      If[PrimeQ[temp],
       Throw[temp]
       ], {i, p}]
     ];
    If[PrimeQ[p + 2],
     list = ConstantArray[1, p];
     Do[
      temp = FromDigits[ReplacePart[list, i -> 3]];
      If[PrimeQ[temp],
       Throw[temp]
       ], {i, p}]
     ];
    Throw[curlargest];
    ]

  , {n, 565, 10000}]
 , 60*60]

Rzuty powodują, że przestaje on testować inne liczby z tymi samymi cyframi, co aktualnie znalezione. ponieważ rozpoczyna testowanie od najbardziej znaczącej cyfry, będzie to oznaczać, że zawsze zwraca największą liczbę, chyba że liczba cyfr jest członkiem pierwszej trójki.

Po prostu zacząłem testować tuż poniżej wartości jednej z poprzednich odpowiedzi :)

Po zakończeniu liczba jest przechowywana w zmiennej najbardziej aktualnej

Zestawienie
źródło
2

JavaScript, 3019 cyfr, {2 273} 5 {745}

To wykorzystuje test MillerRabin zawarty w BigInteger.js przez Toma Wu.

Począwszy od 0 => 2,046 cyfr = {1799} 7 {263} w ciągu godziny .

Począwszy od 3000 => 3019 cyfr = {2273} 5 {745} w ciągu godziny, mniej niż 3 sekundy .

Gdy zaczął się od zera, program przeskoczył do przodu i zaczął ponownie szukać na długości 1,5 x długości ostatniej znalezionej s-prime. Potem, gdy zobaczyłem, jak szybko działa, pomyślałem, że w ciągu godziny znajdzie taki, który zaczyna się od 3000 - co wystarczyło na zaledwie 3 sekundy.

Możesz go wypróbować tutaj: http://goo.gl/t3TmTk
(Ustaw, aby obliczyć wszystkie s-liczby pierwsze lub przejść do przodu.)

wprowadź opis zdjęcia tutaj wprowadź opis zdjęcia tutaj
wprowadź opis zdjęcia tutaj

Program działa poprzez konstruowanie ciągów wszystkich „1”, ale z jednym „3”, „5” lub „7”. Dodałem szybkie sprawdzenie w funkcji IsStrPrime, aby odrzucić liczby kończące się na „5”.

if (IsIntPrime(length)) {

    var do3s = IsIntPrime(length + 2);
    var do5s = IsIntPrime(length + 4);
    var do7s = IsIntPrime(length + 6);

    if (do3s || do5s || do7s) {

        // loop through length of number
        var before, digit, after;

        for (var after = 0; after <= length - 1; after++) {

            before = length - after - 1;
            beforeStr = Ones(before);
            afterStr = Ones(after);

            if (do3s && IsStrPrime(beforeStr + (digit = "3") + afterStr)) { RecordAnswer(); if (brk) break; }
            if (AreDone()) break;

            if (do5s && IsStrPrime(beforeStr + (digit = "5") + afterStr)) { RecordAnswer(); if (brk) break; }
            if (AreDone()) break;

            if (do7s && IsStrPrime(beforeStr + (digit = "7") + afterStr)) { RecordAnswer(); if (brk) break; }
            if (AreDone()) break;

            if (after % 10 == 0) document.title = "b=" + bestLength + ", testing=" + length + "-" + after;
        }
    }
}

To była niezła zabawa. Przypomina mi zagadkę, którą zrobiłem wiele lat temu, aby obliczyć, co nazywa się cyfrą pierwszą usuniętą . Jest to liczba pierwsza, która po usunięciu dowolnej cyfry pozostanie liczbą pierwszą. Na przykład 1037 jest liczbą pierwszą usuniętą, ponieważ 1037, 037, 137, 107 i 103 są liczbą pierwszą. Znalazłem jedną 84 cyfry, a najdłuższa, jaką znam, ma 332 cyfry. Jestem pewien, że moglibyśmy znaleźć go znacznie dłużej dzięki technikom stosowanym w tej układance. (Ale wybór numerów próbnych jest nieco trudniejszy - może?)

JeffSB
źródło
RE: cyfra usunięta liczba pierwsza, mieliśmy ją tutaj . Wygrałyby również 332 cyfry.
primo
0

Axiom, 3019 cyfr {318} 5 {2700}

)time on

-- Return true if n is pseudo prime else return false
-- **Can Fail**
prime1(n:PI):Boolean==
     n=1=>false
     n<4=>true
     i:=5;sq:=sqrt(1.*n)
     repeat
       if i>sq or i>50000 then break
       if n rem i=0       then return false
       i:=i+2
     if i~=50001        then return true
     --output("i")
     if powmod(3,n,n)=3 then return true
     --output("e")
     false

-- input  'n': must be n>1 prime
-- output [0] if not find any number, else return 
-- [n,a,b,c,z] 'n' digits of solution, 
-- 'a' number of '1', 'b' central digit, 'b' number of final digit '1'
-- 'z' the number found
g(n:PI):List NNI==
    x:=b:=z:=1
    for i in 1..n-1 repeat
        z:=z*10+1
        b:=b*10
    repeat
       --output b
       k:=0    -- 3 5 7 <-> 2 4 6
       for i in [2,4,6] repeat
           ~prime?(n+i)=>0 --somma
           k:=k+1
           t:=z+b*i
           if prime1(t) then return [n,x-1,i+1,n-x,t]
       --if x=1 then output ["k=", k] 
       if k=0  then break
       x:=x+1
       b:=b quo 10
       if b<=0 then break
    [0]

-- start from number of digits n
-- and return g(i) with i prime i>=n 
find(n:PI):List NNI==
    i:=n
    if i rem 2=0 then i:=i+1 
    repeat
        if prime?(i) then --solo le lunghezze prime sono accettate
             output i 
             a:=g(i)
             if a~=[0] then return a
        i:=i+2

wynik z wartości początkowej 3000 w 529 sek

(4) -> find(3000)
   3001
   3011
   3019

   (4)
   [3019, 318, 5, 2700, Omissis]
                                            Type: List NonNegativeInteger
       Time: 0.02 (IN) + 525.50 (EV) + 0.02 (OT) + 3.53 (GC) = 529.07 sec
RosLuP
źródło