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}.

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";
}

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

PARI / GP, 4127 cyfr

(10 ⁴¹²⁷ -1) / 9 + 2 * ¹⁰⁵¹⁵

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.

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.

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))

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;
}

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.

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

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

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.)

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?)

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