Ultra useful primes

An ultra-useful prime is a member of the sequence where each a(n) is the smallest positive integer k such that 2^(2ⁿ) - k is prime.

k must always be an odd number since 2 to any power is always even.

Task

Find and show here, on this page, the first 10 elements of the sequence.

Stretch

Find and show the next several elements. (The numbers get really big really fast. Only nineteen elements have been identified as of this writing.)

See also

OEIS:A058220 - Ultra-useful primes: smallest k such that 2^(2^n) - k is prime

ALGOL 68

Works with: ALGOL 68G version Any - tested with release 2.8.3.win32

Uses Algol 68G's LONG LONG INT which has programmer-specifiable precision. Uses Miller Rabin primality testing.

Library: ALGOL 68-primes

BEGIN # find members of the sequence a(n) = smallest k such that 2^(2^n) - k is prime #
    PR precision 650 PR # set number of digits for LONG LOMG INT       #
                        # 2^(2^10) has 308 digits but we need more for #
                        # Miller Rabin primality testing               #
    PR read "primes.incl.a68" PR # include the prime related utilities #
    FOR n TO 10 DO
        LONG LONG INT two up 2 up n = LONG LONG INT( 2 ) ^ ( 2 ^ n );
        FOR i BY 2
        WHILE IF is probably prime( two up 2 up n - i ) THEN
                  # found a sequence member #
                  print( ( " ", whole( i, 0 ) ) );
                  FALSE # stop looking #
              ELSE
                  TRUE # haven't found a sequence member yet #
              FI
        DO SKIP OD
    OD
END

Output:

 1 3 5 15 5 59 159 189 569 105

Arturo

ultraUseful: function [n][
    k: 1
    p: (2^2^n) - k
    while ø [
        if prime? p -> return k
        p: p-2
        k: k+2
    ]
]

print [pad "n" 3 "|" pad.right "k" 4]
print repeat "-" 10
loop 1..10 'x ->
    print [(pad to :string x 3) "|" (pad.right to :string ultraUseful x 4)]

Output:

C

Translation of: Wren

Library: GMP

#include <stdio.h>
#include <gmp.h>

int a(unsigned int n) {
    int k;
    mpz_t p;
    mpz_init_set_ui(p, 1);
    mpz_mul_2exp(p, p, 1 << n);
    mpz_sub_ui(p, p, 1);
    for (k = 1; ; k += 2) {
        if (mpz_probab_prime_p(p, 15) > 0) return k;
        mpz_sub_ui(p, p, 2);
    }
}

int main() {
    unsigned int n;
    printf(" n   k\n");
    printf("----------\n");
    for (n = 1; n < 15; ++n) printf("%2d   %d\n", n, a(n));
    return 0;
}

Output:

C#

using System.Numerics;
using System.Security.Cryptography;

internal class Program
{
    private static void Main(string[] args)
    {
        for (int i = 1; i <= 10; i++)
        {
            Console.WriteLine(GetUltraUsefulPrime(i));
        }

    }
    public static int GetUltraUsefulPrime(int n)
    {
        int k = 1;
        
        BigInteger num = BigInteger.Pow(2, (int)Math.Pow(2, n));

        while (!IsProbablePrime(num - k, 20))
        {
            k += 2;
        }

        return k;
    }

    // Miller-Rabin primality from this site
    public static bool IsProbablePrime(BigInteger source, int certainty)
    {
        if (source == 2 || source == 3)
            return true;
        if (source < 2 || source % 2 == 0)
            return false;

        BigInteger d = source - 1;
        int s = 0;

        while (d % 2 == 0)
        {
            d /= 2;
            s += 1;
        }

        // There is no built-in method for generating random BigInteger values.
        // Instead, random BigIntegers are constructed from randomly generated
        // byte arrays of the same length as the source.
        RandomNumberGenerator rng = RandomNumberGenerator.Create();
        byte[] bytes = new byte[source.ToByteArray().LongLength];
        BigInteger a;

        for (int i = 0; i < certainty; i++)
        {
            do
            {
                // This may raise an exception in Mono 2.10.8 and earlier.
                // http://bugzilla.xamarin.com/show_bug.cgi?id=2761
                rng.GetBytes(bytes);
                a = new BigInteger(bytes);
            }
            while (a < 2 || a >= source - 2);

            BigInteger x = BigInteger.ModPow(a, d, source);
            if (x == 1 || x == source - 1)
                continue;

            for (int r = 1; r < s; r++)
            {
                x = BigInteger.ModPow(x, 2, source);
                if (x == 1)
                    return false;
                if (x == source - 1)
                    break;
            }

            if (x != source - 1)
                return false;
        }

        return true;
    }
}

Output:

Craft Basic

for n = 1 to 10

	let k = -1

	do

		let k = k + 2
		wait

	loop prime(2 ^ (2 ^ n) - k) = 0

	print "n = ", n, " k = ", k

next n

Factor

Works with: Factor version 0.99 2021-06-02

USING: io kernel lists lists.lazy math math.primes prettyprint ;

: useful ( -- list )
    1 lfrom
    [ 2^ 2^ 1 lfrom [ - prime? ] with lfilter car ] lmap-lazy ;

10 useful ltake [ pprint bl ] leach nl

Output:

1 3 5 15 5 59 159 189 569 105

FreeBASIC

Translation of: Ring

#include "isprime.bas"

Dim As Longint n, k, limit = 10
Dim As ulongint num
For n = 1 To limit
    k = -1
    Do
        k += 2
        num = (2 ^ (2 ^ n)) - k
        If isPrime(num) Then
            Print "n = "; n; " k = "; k
            Exit Do
        End If
    Loop
Next

Sleep

Go

Library: GMP(Go wrapper)

package main

import (
    "fmt"
    big "github.com/ncw/gmp"
)

var two = big.NewInt(2)

func a(n uint) int {
    one := big.NewInt(1)
    p := new(big.Int).Lsh(one, 1 << n)
    p.Sub(p, one)
    for k := 1; ; k += 2 {
        if p.ProbablyPrime(15) {
            return k
        }
        p.Sub(p, two)
    }
}

func main() {
    fmt.Println(" n   k")
    fmt.Println("----------")
    for n := uint(1); n < 14; n++ {
        fmt.Printf("%2d   %d\n", n, a(n))
    }
}

Output:

J

Implementation:

uup=: {{
  ref=. 2^2x^y+1
  k=. 1
  while. -. 1 p: ref-k do. k=. k+2 end.
}}"0

I don't have the patience to get this little laptop to compute the first 10 such elements, so here I only show the first five:

   uup i.10
1 3 5 15 5 59 159 189 569 105

Java

import java.math.BigInteger;

public final class UltraUsefulPrimes {

	public static void main(String[] args) {
		for ( int n = 1; n <= 10; n++ ) {
			showUltraUsefulPrime(n);
		}		
	}
	
	private static void showUltraUsefulPrime(int n) {
		BigInteger prime = BigInteger.ONE.shiftLeft(1 << n);
		BigInteger k = BigInteger.ONE; 
		while ( ! prime.subtract(k).isProbablePrime(20) ) {
			k = k.add(BigInteger.TWO);
		}
		
		System.out.print(k + " ");		   
	}

}

Output:

1 3 5 15 5 59 159 189 569 105

Julia

using Primes

nearpow2pow2prime(n) = findfirst(k -> isprime(2^(big"2"^n) - k), 1:10000)

@time println([nearpow2pow2prime(n) for n in 1:12])

Output:

[1, 3, 5, 15, 5, 59, 159, 189, 569, 105, 1557, 2549]
  3.896011 seconds (266.08 k allocations: 19.988 MiB, 1.87% compilation time)

M2000 Interpreter

Change 10 to 6 if you don't have time to spend...

module Ultra_Useful {
	minusOne=biginteger("-1")	
	one=biginteger("1")
	two=biginteger("2")
	num=biginteger("0")
	k=num
	with num,"tostring" as numS
	with k,"tostring" as kS
	for n= 1 to 10
		k=minusOne
		kk=-1
		n1=biginteger(str$(n,""))
		method two,"intpower", n1 as n1
		do
			method k,"add", two as k
			kk+=2
			method two, "intpower", n1 as num
			method num, "subtract", k as num
			method num,"isProbablyPrime", 10 as ret
			if ret then
				Print n, kk : Refresh
				exit
			end if
		always
	next
}
Ultra_Useful

Output:

Mathematica /Wolfram Language

ClearAll[FindUltraUsefulPrimeK]
FindUltraUsefulPrimeK[n_] := Module[{num, tmp},
  num = 2^(2^n);
  Do[
   If[PrimeQ[num - k],
    tmp = k;
    Break[];
    ]
   ,
   {k, 1, \[Infinity], 2}
   ];
  tmp
  ]
res = FindUltraUsefulPrimeK /@ Range[13];
TableForm[res, TableHeadings -> Automatic]

Output:

Nim

Library: Nim-Integers

import std/strformat
import integers

let One = newInteger(1)

echo " n  k"
var count = 1
var n = 1
while count <= 13:
  var k = 1
  var p = One shl (1 shl n) - k
  while not p.isPrime:
    p -= 2
    k += 2
  echo &"{n:2}  {k}"
  inc n
  inc count

Output:

Perl

Library: ntheory

use strict;
use warnings;
use feature 'say';
use bigint;
use ntheory 'is_prime';

sub useful {
    my @n = @_;
    my @u;
    for my $n (@n) {
        my $p = 2**(2**$n);
        LOOP: for (my $k = 1; $k < $p; $k += 2) {
            is_prime($p-$k) and push @u, $k and last LOOP;
       }
    }
    @u
}

say join ' ', useful 1..13;

Output:

1 3 5 15 5 59 159 189 569 105 1557 2549 2439

Phix

with javascript_semantics
atom t0 = time()
include mpfr.e
mpz p = mpz_init()
 
function a(integer n)
    mpz_ui_pow_ui(p,2,power(2,n))
    mpz_sub_si(p,p,1)
    integer k = 1
    while not mpz_prime(p) do
        k += 2
        mpz_sub_si(p,p,2)
    end while
    return k
end function
 
for i=1 to 10 do
    printf(1,"%d ",a(i))
end for
if machine_bits()=64 then
    ?elapsed(time()-t0)
    for i=11 to 13 do
        printf(1,"%d ",a(i))
    end for
end if
?elapsed(time()-t0)

Output:

1 3 5 15 5 59 159 189 569 105 "0.0s"
1557 2549 2439 "1 minute and 1s"

Raku

The first 10 take less than a quarter second. 11 through 13, a little under 30 seconds. Drops off a cliff after that.

sub useful ($n) {
    (|$n).map: {
        my $p = 1 +< ( 1 +< $_ );
        ^$p .first: ($p - *).is-prime
    }
}

put useful 1..10;

put useful 11..13;

Output:

1 3 5 15 5 59 159 189 569 105
1557 2549 2439

Python

# useful.py by xing216
from gmpy2 import is_prime
def useful(n):
    k = 1
    is_useful = False
    while is_useful == False:
        if is_prime(2**(2**n) - k):
            is_useful = True
            break
        k += 2
    return k
if __name__ == "__main__":
    print("n | k")
    for i in range(1,14):
        print(f"{i:<4}{useful(i)}")

Output:

Quackery

prime is defined at Miller–Rabin primality test#Quackery.

  12 times 
    [ i^ 1+ bit bit
      1 
      [ 2dup - 
        prime not while 
        2 + again ]
      echo sp drop ]

Output:

1 3 5 15 5 59 159 189 569 105 1557 2549

Ring

see "works..." + nl
limit = 10

for n = 1 to limit
    k = -1
    while true
          k = k + 2
          num = pow(2,pow(2,n)) - k
          if isPrime(num)
             ? "n = " + n + " k = " + k
             exit
          ok
    end
next
see "done.." + nl

func isPrime num
     if (num <= 1) return 0 ok
     if (num % 2 = 0 and num != 2) return 0 ok
     for i = 3 to floor(num / 2) -1 step 2
         if (num % i = 0) return 0 ok
     next
     return 1

Output:

works...
n = 1 k = 1
n = 2 k = 3
n = 3 k = 5
n = 4 k = 15
n = 5 k = 5
n = 6 k = 59
n = 7 k = 159
n = 8 k = 189
n = 9 k = 569
n = 10 k = 105
done...

Ruby

The 'prime?' method in Ruby's OpenSSL (standard) lib is much faster than the 'prime' lib. 0.05 sec. for this, about a minute for the next three.

require 'openssl'

(1..10).each do |n|
   pow = 2 ** (2 ** n)
   print "#{n}:\t"
   puts (1..).step(2).detect{|k| OpenSSL::BN.new(pow-k).prime?}
end

Output:

1:      1
2:      3
3:      5
4:      15
5:      5
6:      59
7:      159
8:      189
9:      569
10:     105

Sidef

say(" n   k")
say("----------")

for n in (1..13) {
    var t = 2**(2**n)
    printf("%2d   %d\n", n, {|k| t - k -> is_prob_prime }.first)
}

Output:

(takes ~20 seconds)

V (Vlang)

import math

fn main() {
	limit := 10 // depending on computer, higher numbers = longer times
	mut num, mut k := i64(0), i64(0)
	println("n   k\n-------")
	for n in 1..limit {
		k = -1
		for n < limit {
			  k = k + 2
			  num = math.powi(2, math.powi(2 , n)) - k
			  if is_prime(num) == true {
				 println("${n}   ${k}")
				 break
			  }
		}
	}
}

fn is_prime(num i64) bool {
     if num <= 1 {return false}
     if num % 2 == 0 && num != 2 {return false}
	 for idx := 3; idx <= math.floor(num / 2) - 1; idx += 2 {
         if num % idx == 0 {return false}
     }
     return true
}

Output:

Wren

CLI

Library: Wren-big

Library: Wren-fmt

Manages to find the first ten but takes 84 seconds to do so.

import "./big" for BigInt
import "./fmt" for Fmt

var one = BigInt.one
var two = BigInt.two

var a = Fn.new { |n|
    var p = (BigInt.one << (1 << n)) - one
    var k = 1
    while (true) {
        if (p.isProbablePrime(5)) return k
        p = p - two
        k = k + 2
    }
}

System.print(" n   k")
System.print("----------")
for (n in 1..10) Fmt.print("$2d   $d", n, a.call(n))

Output:

Embedded

Library: Wren-gmp

The following takes about 17 seconds to get to n = 13 but 7 minutes 10 seconds to reach n = 14. I didn't bother after that.

import "./gmp" for Mpz
import "./fmt" for Fmt

var one = Mpz.one
var two = Mpz.two

var a = Fn.new { |n|
    var p = Mpz.one.lsh(1 << n).sub(one)
    var k = 1
    while (true) {
        if (p.probPrime(15) > 0) return k
        p.sub(two)
        k = k + 2
    }
}

System.print(" n   k")
System.print("----------")
for (n in 1..14) Fmt.print("$2d   $d", n, a.call(n))

Output: