Sequence of primorial primes
The sequence of primorial primes is given as the increasing values of n where primorial(n) ± 1 is prime.
Noting that the n'th primorial is defined as the multiplication of the smallest n primes, the sequence is of the number of primes, in order that when multiplied together is one-off being a prime number itself.
The task is to generate and show here the first ten values of the sequence. An optional extended task is to show the first twenty members of the series.
Note:
- This task asks for the primorial indices that create the final primorial prime numbers, so there should be no ten-or-more digit numbers in the program output (although extended precision integers will be needed for intermediate results).
- There is some confusion in the references, but for the purposes of this task the sequence begins with n = 1.
- Cf.
- Primorial prime Wikipedia.
- Primorial prime from The Prime Glossary.
- Sequence A088411 from The On-Line Encyclopedia of Integer Sequences
Related tasks:
EchoLisp
<lang lisp> (lib 'timer) ;; for (every (proc t) interval) (lib 'bigint)
- memoize primorial
(define p1000 (cons 1 (primes 1000))) ; remember first 1000 primes (define (primorial n)
(if(zero? n) 1 (* (list-ref p1000 n) (primorial (1- n)))))
(remember 'primorial)
(define N 0)
- search one at a time
(define (search t) ;; time parameter, set by (every), not used (set! N (1+ N)) (if (or (prime? (1+ (primorial N))) (prime? (1- (primorial N)))) (writeln 'HIT N )
(writeln N (date->time-string (current-date )))))
</lang>
- Output:
<lang lisp>
- run in batch mode to make the browser happy
- search next N every 2000 msec
(every 2000 search)
→ HIT 1 HIT 2 HIT 3 HIT 4 HIT 5 HIT 6 7 "13:47:03" HIT 11 HIT 13 HIT 24 HIT 66 HIT 68 HIT 75 HIT 167 HIT 171 HIT 172 HIT 287 HIT 310 HIT 352 HIT 384 455 "14:07:08" 456 "14:07:12" HIT 457 458 "14:07:25" 459 "14:07:29"
</lang>
Haskell
<lang Haskell> import Data.List (scanl1, findIndices, nub)
primes :: [Integer] primes = 2 : filter isPrime [3,5..]
isPrime :: Integer -> Bool isPrime = isPrime_ primes
where isPrime_ :: [Integer] -> Integer -> Bool
isPrime_ (p:ps) n
| p * p > n = True
| n `mod` p == 0 = False
| otherwise = isPrime_ ps n
primorials :: [Integer] primorials = 1 : scanl1 (*) primes
primorialsPlusMinusOne :: [Integer] primorialsPlusMinusOne = concatMap (\n -> [n-1, n+1]) primorials
sequenceOfPrimorialPrimes :: [Int] sequenceOfPrimorialPrimes = tail $ nub $ map (`div` 2) $ findIndices (==True) bools
where bools = map isPrime primorialsPlusMinusOne
main :: IO () main = mapM_ print $ take 10 sequenceOfPrimorialPrimes </lang>
- Output:
1 2 3 4 5 6 11 13 24 66
PARI/GP
This program, in principle, generates arbitrarily many primirial prime indices. Practically speaking it depends on your patience; it generates them up to 643 in about 5 seconds. <lang parigp>n=0; P=1; forprime(p=2,, P*=p; n++; if(ispseudoprime(P+1) || ispseudoprime(P-1), print1(n", ")))</lang>
- Output:
1, 2, 3, 4, 5, 6, 11, 13, 24, 66, 68, 75, 167, 171, 172, 287, 310, 352, 384, 457, 564, 590, 616, 620, 643, 849,
Perl
<lang perl>use ntheory ":all"; my $i = 0; for (1..1e6) {
my $n = pn_primorial($_);
if (is_prime($n-1) || is_prime($n+1)) {
say;
last if ++$i >= 20;
}
}</lang>
- Output:
1 2 3 4 5 6 11 13 24 66 68 75 167 171 172 287 310 352 384 457
Perl 6
<lang perl6>constant @primes = grep *.is-prime, 2..*; constant @primorials = [\*] 1, @primes;
my @pp_indexes := @primorials.pairs.map: {
.key if ( .value + any(1, -1) ).is-prime
};
say ~ @pp_indexes[ 0 ^.. 20 ]; # Skipping bogus first element.</lang>
- Output:
1 2 3 4 5 6 11 13 24 66 68 75 167 171 172 287 310 352 384 457
Python
Uses the pure python library pyprimes . Note that pyprimes.isprime switches to a fast and good probabilistic algorithm for larger integers.
<lang python>import pyprimes
def primorial_prime(_pmax=500):
isprime = pyprimes.isprime
n, primo = 0, 1
for prime in pyprimes.nprimes(_pmax):
n, primo = n+1, primo * prime
if isprime(primo-1) or isprime(primo+1):
yield n
if __name__ == '__main__':
# Turn off warning on use of probabilistic formula for prime test
pyprimes.warn_probably = False
for i, n in zip(range(20), primorial_prime()):
print('Primorial prime %2i at primorial index: %3i' % (i+1, n))</lang>
- Output:
Primorial prime 1 at primorial index: 1 Primorial prime 2 at primorial index: 2 Primorial prime 3 at primorial index: 3 Primorial prime 4 at primorial index: 4 Primorial prime 5 at primorial index: 5 Primorial prime 6 at primorial index: 6 Primorial prime 7 at primorial index: 11 Primorial prime 8 at primorial index: 13 Primorial prime 9 at primorial index: 24 Primorial prime 10 at primorial index: 66 Primorial prime 11 at primorial index: 68 Primorial prime 12 at primorial index: 75 Primorial prime 13 at primorial index: 167 Primorial prime 14 at primorial index: 171 Primorial prime 15 at primorial index: 172 Primorial prime 16 at primorial index: 287 Primorial prime 17 at primorial index: 310 Primorial prime 18 at primorial index: 352 Primorial prime 19 at primorial index: 384 Primorial prime 20 at primorial index: 457
Racket
We use a memorized version of primordial to make it faster.
<lang Racket>#lang racket
(require math/number-theory
racket/generator)
(define-syntax-rule (define/cache (name arg) body ...)
(begin
(define cache (make-hash))
(define (name arg)
(hash-ref! cache arg (lambda () body ...)))))
(define/cache (primorial n)
(if (zero? n)
1
(* (nth-prime (sub1 n))
(primorial (sub1 n)))))
(for ([i (in-range 20)]
[n (in-generator (for ([i (in-naturals 1)])
(define pr (primorial i))
(when (or (prime? (add1 pr)) (prime? (sub1 pr)))
(yield i))))])
(displayln n))</lang>
- Output:
1 2 3 4 5 6 11 13 24 66 68 75 167 171 172 287 310 352 384 457