Circular primes

Circular primes is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

Definitions

A circular prime is a prime number with the property that the number generated at each intermediate step when cyclically permuting its (base 10) digits will also be prime.

For example: 1193 is a circular prime, since 1931, 9311 and 3119 are all also prime.

Note that a number which is a cyclic permutation of a smaller circular prime is not considered to be itself a circular prime. So 13 is a circular prime but 31 is not.

A repunit (denoted by R) is a number whose base 10 representation contains only the digit 1.

For example: R(2) = 11 and R(5) = 11111 are repunits.

Task

• Find the first 19 circular primes.

• If your language has access to arbitrary precision integer arithmetic, given that they are all repunits, find the next 4 circular primes.

• (Stretch) Determine which of the following repunits are probably circular primes: R(5003), R(9887), R(15073), R(25031), R(35317) and R(49081). The larger ones may take a long time to process so just do as many as you reasonably can.

See also

• Wikipedia article - Circular primes.
• Wikipedia article - Repunit.
• OEIS sequence A016114 - Circular primes.

Go

Library: GMP(Go wrapper)

<lang go>package main

import (

   "fmt"
   big "github.com/ncw/gmp"
   "strings"

)

// OK for 'small' numbers. func isPrime(n int) bool {

   switch {
   case n < 2:
       return false
   case n%2 == 0:
       return n == 2
   case n%3 == 0:
       return n == 3
   default:
       d := 5
       for d*d <= n {
           if n%d == 0 {
               return false
           }
           d += 2
           if n%d == 0 {
               return false
           }
           d += 4
       }
       return true
   }

}

func repunit(n int) *big.Int {

   ones := strings.Repeat("1", n)
   b, _ := new(big.Int).SetString(ones, 10)
   return b

}

var circs = []int{2, 3, 5, 7}

// binary search is overkill for a small number of elements func alreadyFound(n int) bool {

   for _, i := range circs {
       if i == n {
           return true
       }
   }
   return false

}

func isCircular(n int) bool {

   nn := n
   pow := 1 // will eventually contain 10 ^ d where d is number of digits in n
   for nn > 0 {
       pow *= 10
       nn /= 10
   }
   nn = n
   for {
       nn *= 10
       f := nn / pow // first digit
       nn += f * (1 - pow)
       if alreadyFound(nn) {
           return false
       }
       if nn == n {
           break
       }
       if !isPrime(nn) {
           return false
       }
   }
   return true

}

var digits = [5]int{0, 1, 3, 7, 9}

func to01379(n int) int {

   sum := 0
   pow := 1
   for n > 0 {
       d := n % 5
       sum += digits[d] * pow
       n /= 5
       pow *= 10
   }
   return sum

}

func main() {

   fmt.Println("The first 19 circular primes are:")
   count := 4
   for i := 6; count < 19; i++ {
       j := to01379(i)
       if j%10 == 0 || !isPrime(j) {
           continue
       }
       if !isPrime(j) {
           continue
       }
       if isCircular(j) {
           count++
           circs = append(circs, j)
       }
   }
   fmt.Println(circs)
   fmt.Println("\nThe next 4 circular primes, in repunit format, are:")
   count = 0
   var rus []string
   for i := 7; count < 4; i++ {
       if repunit(i).ProbablyPrime(10) {
           count++
           rus = append(rus, fmt.Sprintf("R(%d)", i))
       }
   }
   fmt.Println(rus)
   fmt.Println("\nThe following repunits are probably circular primes:")
   for _, i := range []int{5003, 9887, 15073, 25031, 35317, 49081} {
       fmt.Printf("R(%-5d) : %t\n", i, repunit(i).ProbablyPrime(10))
   }

}</lang>

Output:

The first 19 circular primes are:
[2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933]

The next 4 circular primes, in repunit format, are:
[R(19) R(23) R(317) R(1031)]

The following repunits are probably circular primes:
R(5003 ) : false
R(9887 ) : false
R(15073) : false
R(25031) : false
R(35317) : false
R(49081) : true

Phix

<lang Phix>function circular(integer p)

   integer len = length(sprintf("%d",p)),
           pow = power(10,len-1),
           p0 = p
   for i=1 to len-1 do
       p = pow*remainder(p,10)+floor(p/10)
       if p<p0 or not is_prime(p) then return false end if
   end for
   return true

end function

sequence c = {} integer n = 1 while length(c)<19 do

   integer p = get_prime(n)
   if circular(p) then c &= p end if
   n += 1

end while printf(1,"The first 19 circular primes are:\n%v\n\n",{c})

include mpfr.e procedure repunit(mpz z, integer n)

   mpz_set_str(z,repeat('1',n))

end procedure

c = {} n = 7 mpz z = mpz_init() randstate state = gmp_randinit_mt() while length(c)<4 do

   repunit(z,n)
   if mpz_probable_prime_p(z,state) then
       c = append(c,sprintf("R(%d)",n))
   end if
   n += 1

end while printf(1,"The next 4 circular primes, in repunit format, are:\n%s\n\n",{join(c)})</lang>

Output:

The first 19 circular primes are:
{2,3,5,7,11,13,17,37,79,113,197,199,337,1193,3779,11939,19937,193939,199933}

The next 4 circular primes, in repunit format, are:
R(19) R(23) R(317) R(1031)

stretch

<lang Phix>constant t = {5003, 9887, 15073, 25031, 35317, 49081} printf(1,"The following repunits are probably circular primes:\n") for i=1 to length(t) do

   integer ti = t[i]
   atom t0 = time()
   repunit(z,ti)
   bool bPrime = mpz_probable_prime_p(z,state,1)
   printf(1,"R(%d) : %t (%s)\n", {ti, bPrime, elapsed(time()-t0)})

end for</lang>

Output:

64-bit can only cope with the first five (it terminates abruptly on the sixth)

The following repunits are probably circular primes:
R(5003) : false (2.0s)
R(9887) : false (13.5s)
R(15073) : false (45.9s)
R(25031) : false (1 minute and 19s)
R(35317) : false (3 minutes and 04s)

32-bit is much slower and can only cope with the first four

The following repunits are probably circular primes:
R(5003) : false (10.2s)
R(9887) : false (54.9s)
R(15073) : false (2 minutes and 22s)
R(25031) : false (7 minutes and 45s)
diag looping, error code is 1, era is #00644651

Raku

Most of the repunit testing is relatively speedy using the ntheory library. The really slow ones are R(25031), at ~42 seconds and R(49081) at 922(!!) seconds.

<lang perl6>#!/usr/bin/env raku

1. 20200406 Raku programming solution

sub isCircular(\n) {

  return False unless n.is-prime;
  my @circular = n.comb;
  return False if @circular.min < @circular[0];
  for (1..^@circular.elems) {
     return False if n > my $rotated = @circular.rotate($_).join;
     return False unless $rotated.is-prime;
  }
  True

}

say "The first 19 circular primes are:"; say ((1..∞).hyper.grep: { isCircular($_) })[^19];

say "\nThe next 4 circular primes, in repunit format, are:"; loop ( my $i = 7, my $count = 0; $count < 4; $i++ ) {

  my $target = 1 x $i;
  if $target.is-prime  {
     say "R(",$i,")";
     $count++
  }

}

use ntheory:from<Perl5> qw[is_prime];

say "\nRepunit testing:";

(5003, 9887, 15073, 25031, 35317, 49081).map: {

   my $now = now;
   say "R($_): Prime? ", ?is_prime("{1 x $_}"), "  {(now - $now).fmt: '%.2f'}"

}</lang>

Output:

The first 19 circular primes are:
(2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933)

The next 4 circular primes, in repunit format, are:
R(19)
R(23)
R(317)
R(1031)

Repunit testing:
R(5003): Prime? False  0.00
R(9887): Prime? False  0.01
R(15073): Prime? False  0.02
R(25031): Prime? False  41.40
R(35317): Prime? False  0.32
R(49081): Prime? True  921.73

REXX

<lang rexx>/*REXX program finds & displays circular primes (with a title & in a horizontal format).*/ parse arg N hp . /*obtain optional arguments from the CL*/ if N== | N=="," then N= 19 /* " " " " " " */ if hp== | hp=="," then hp= 1000000 /* " " " " " " */ call genP /*gen primes up to hp (200,000). */ q= 024568 /*digs that most circular P can't have.*/

1. =0; $= /*#: circular P count; $: is a list.*/

     do j=2  until #==N                         /* [↓]   traipse through the number(s).*/
     if \!.j                       then iterate /*Is  J  not prime?   Then skip number.*/
     if j>9 & verify(j, q, 'M')>0  then iterate /*Does J contain forbidden digs?  Skip.*/
     if \circP(j)                  then iterate /*Not circular?         "    "    "    */
     #= # + 1                                   /*bump the count of circular primes.   */
     $= $ j                                     /*add this prime number  ──►  $  list. */
     end   /*j*/

say center(' first' # "circular primes ", 79, '─') say strip($) exit /*stick a fork in it, we're all done. */ /*──────────────────────────────────────────────────────────────────────────────────────*/ circP: procedure expose @. !.; parse arg x 1 ox /*obtain a prime number to be examined.*/

               do  length(x)-1                  /*parse  X  number, rotating the digits*/
               parse var  x    2  y  1  f  2    /*get rightmost digs, append first dig.*/
               x= y || f                        /*construct a new possible circular P. */
               if x<ox  then return 0           /*is the number less than the original?*/
               if \!.x  then return 0           /*this version of rotated P isn't prime*/
               end   /*len···*/
      return 1                                  /*passed all tests,  X is a circular P.*/

/*──────────────────────────────────────────────────────────────────────────────────────*/ genP: @.1=2; @.2=3; @.3=5; @.4=7; @.5=11; @.6= 13; nP=6 /*assign low primes; # primes. */

     !.= 0; !.2=1; !.3=1; !.5=1; !.7=1; !.11=1; !.13=1  /*assign primality to numbers. */
      do j=@.nP+4  by 2  to hp                  /*only find odd primes from here on.   */
      if j// 3==0  then iterate                 /*is J divisible by  #3  Then not prime*/
      parse var j  -1 _;if _==5  then iterate /*Is last digit a "5"?     "   "    "  */
      if j// 7==0  then iterate                 /*is J divisible by  7?    "   "    "  */
      if j//11==0  then iterate                 /* " "     "      " 11?    "   "    "  */
      if j//13==0  then iterate                 /*is "     "      " 13?    "   "    "  */
          do k=7  while k*k<=j                  /*divide by some generated odd primes. */
          if j // @.k==0  then iterate j        /*Is J divisible by  P?  Then not prime*/
          end   /*k*/                           /* [↓]  a prime  (J)  has been found.  */
      nP= nP+1;           !.j=1;       @.nP= j  /*bump P cnt;  assign P to @.  and  !. */
      end       /*j*/;    return</lang>

output   when using the default inputs:

────────────────────────── first 19 circular primes ───────────────────────────
2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933