Circular primes: Difference between revisions
(→{{header|Perl}}: perpend =={{header|Pascal}}== test up to 10 digits) |
m (→{{header|Free Pascal}}: doesn't take 2.6 s for the first 19 but 0.073s til 999999) |
||
| Line 1,656: | Line 1,656: | ||
const |
const |
||
MAXCNTOFDIGITS = 10; |
MAXCNTOFDIGITS = 6;//10;//11 |
||
MAXDGTVAL = 3; |
MAXDGTVAL = 3; |
||
conv : array[0..MAXDGTVAL+1] of byte = (9,7,3,1,0); |
conv : array[0..MAXDGTVAL+1] of byte = (9,7,3,1,0); |
||
| Line 1,814: | Line 1,814: | ||
19 199933 |
19 199933 |
||
Real time: |
Real time: 0.073 s digit Count MAXCNTOFDIGITS = 6; |
||
Real time: 2.691 s digit Count MAXCNTOFDIGITS = 10; |
|||
Real time: 24.545 s |
Real time: 24.545 s digit Count MAXCNTOFDIGITS = 11; |
||
</pre> |
|||
=={{header|Perl}}== |
=={{header|Perl}}== |
||
{{trans|Raku}} |
{{trans|Raku}} |
||
Revision as of 13:09, 18 November 2021
You are encouraged to solve this task according to the task description, using any language you may know.
- 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.
ALGOL 68
<lang algol68>BEGIN # find circular primes - primes where all cyclic permutations #
# of the digits are also prime #
# genertes a sieve of circular primes, only the first #
# permutation of each prime is flagged as TRUE #
OP CIRCULARPRIMESIEVE = ( INT n )[]BOOL:
BEGIN
[ 0 : n ]BOOL prime;
prime[ 0 ] := prime[ 1 ] := FALSE;
prime[ 2 ] := TRUE;
FOR i FROM 3 BY 2 TO UPB prime DO prime[ i ] := TRUE OD;
FOR i FROM 4 BY 2 TO UPB prime DO prime[ i ] := FALSE OD;
FOR i FROM 3 BY 2 TO ENTIER sqrt( UPB prime ) DO
IF prime[ i ] THEN
FOR s FROM i * i BY i + i TO UPB prime DO prime[ s ] := FALSE OD
FI
OD;
INT first digit multiplier := 10;
INT max with multiplier := 99;
# the 1 digit primes are non-curcular, so start at 10 #
FOR i FROM 10 TO UPB prime DO
IF i > max with multiplier THEN
# starting a new power of ten #
first digit multiplier *:= 10;
max with multiplier *:= 10 +:= 9
FI;
IF prime[ i ] THEN
# have a prime #
# cycically permute the number until we get back #
# to the original - flag all the permutations #
# except the original as non-prime #
INT permutation := i;
WHILE permutation := ( permutation OVER 10 )
+ ( ( permutation MOD 10 ) * first digit multiplier )
;
permutation /= i
DO
IF NOT prime[ permutation ] THEN
# the permutation is not prime #
prime[ i ] := FALSE
ELIF permutation > i THEN
# haven't permutated e.g. 101 to 11 #
IF NOT prime[ permutation ] THEN
# i is not a circular prime #
prime[ i ] := FALSE
FI;
prime[ permutation ] := FALSE
FI
OD
FI
OD;
prime
END # CIRCULARPRIMESIEVE # ;
# construct a sieve of circular primes up to 999 999 #
# only the first permutation is included #
[]BOOL prime = CIRCULARPRIMESIEVE 999 999;
# print the first 19 circular primes #
INT c count := 0;
print( ( "First 19 circular primes: " ) );
FOR i WHILE c count < 19 DO
IF prime[ i ] THEN
print( ( " ", whole( i, 0 ) ) );
c count +:= 1
FI
OD;
print( ( newline ) )
END</lang>
- Output:
First 19 circular primes: 2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933
ALGOL W
<lang algolw>begin % find circular primes - primes where all cyclic permutations %
% of the digits are also prime %
% sets p( 1 :: n ) to a sieve of primes up to n %
procedure Eratosthenes ( logical array p( * ) ; integer value n ) ;
begin
p( 1 ) := false; p( 2 ) := true;
for i := 3 step 2 until n do p( i ) := true;
for i := 4 step 2 until n do p( i ) := false;
for i := 3 step 2 until truncate( sqrt( n ) ) do begin
integer ii; ii := i + i;
if p( i ) then for pr := i * i step ii until n do p( pr ) := false
end for_i ;
end Eratosthenes ;
% find circular primes in p in the range lo to hi, if they are circular, flag the %
% permutations as non-prime so we do not consider them again %
% non-circular primes are also flageed as non-prime %
% lo must be a power of ten and hi must be at most ( lo * 10 ) - 1 %
procedure keepCircular ( logical array p ( * ); integer value lo, hi ) ;
for n := lo until hi do begin
if p( n ) then begin
% have a prime %
integer c, pCount;
logical isCircular;
integer array permutations ( 1 :: 10 );
c := n;
isCircular := true;
pCount := 0;
% cyclically permute c until we get back to p or find a non-prime value for c %
while begin
integer first, rest;
first := c div lo;
rest := c rem lo;
c := ( rest * 10 ) + first;
isCircular := p( c );
c not = n and isCircular
end do begin
pCount := pCount + 1;
permutations( pCount ) := c
end while_have_another_prime_permutation ;
if not isCircular
then p( n ) := false
else begin
% have a circular prime - flag the permutations as non-prime %
for i := 1 until pCount do p( permutations( i ) ) := false
end if_not_isCircular__
end if_p_n
end keepCircular ;
integer cCount;
% sieve the primes up to 999999 %
logical array p ( 1 :: 999999 );
Eratosthenes( p, 999999 );
% remove non-circular primes from the sieve %
% the single digit primes are all circular so we start at 10 %
keepCircular( p, 10, 99 );
keepCircular( p, 100, 999 );
keepCircular( p, 1000, 9999 );
keepCircular( p, 10000, 99999 );
keepCircular( p, 100000, 200000 );
% print the first 19 circular primes %
cCount := 0;
write( "First 19 circular primes: " );
for i := 1 until 200000 do begin
if p( i ) then begin
writeon( i_w := 1, s_w := 1, i );
cCount := cCount + 1;
if cCount = 19 then goto end_circular
end if_p_i
end for_i ;
end_circular: end.</lang>
- Output:
First 19 circular primes: 2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933
Arturo
<lang rebol>perms: function [n][
str: repeat to :string n 2
result: new []
lim: dec size digits n
loop 0..lim 'd ->
'result ++ slice str d lim+d
return to [:integer] result
]
circulars: new []
circular?: function [x][
if not? prime? x -> return false
loop perms x 'y [
if not? prime? y -> return false
if contains? circulars y -> return false
]
'circulars ++ x
return true
]
i: 2 found: 0 while [found < 19][
if circular? i [
print i
found: found + 1
]
i: i + 1
]</lang>
- Output:
2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933
AWK
<lang AWK>
- syntax: GAWK -f CIRCULAR_PRIMES.AWK
BEGIN {
p = 2
printf("first 19 circular primes:")
for (count=0; count<19; p++) {
if (is_circular_prime(p)) {
printf(" %d",p)
count++
}
}
printf("\n")
exit(0)
} function cycle(n, m,p) { # E.G. if n = 1234 returns 2341
m = n
p = 1
while (m >= 10) {
p *= 10
m /= 10
}
return int(m+10*(n%p))
} function is_circular_prime(p, p2) {
if (!is_prime(p)) {
return(0)
}
p2 = cycle(p)
while (p2 != p) {
if (p2 < p || !is_prime(p2)) {
return(0)
}
p2 = cycle(p2)
}
return(1)
} function is_prime(x, i) {
if (x <= 1) {
return(0)
}
for (i=2; i<=int(sqrt(x)); i++) {
if (x % i == 0) {
return(0)
}
}
return(1)
} </lang>
- Output:
first 19 circular primes: 2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933
C
<lang c>#include <stdbool.h>
- include <stdint.h>
- include <stdio.h>
- include <stdlib.h>
- include <string.h>
- include <gmp.h>
bool is_prime(uint32_t n) {
if (n == 2)
return true;
if (n < 2 || n % 2 == 0)
return false;
for (uint32_t p = 3; p * p <= n; p += 2) {
if (n % p == 0)
return false;
}
return true;
}
// e.g. returns 2341 if n = 1234 uint32_t cycle(uint32_t n) {
uint32_t m = n, p = 1;
while (m >= 10) {
p *= 10;
m /= 10;
}
return m + 10 * (n % p);
}
bool is_circular_prime(uint32_t p) {
if (!is_prime(p))
return false;
uint32_t p2 = cycle(p);
while (p2 != p) {
if (p2 < p || !is_prime(p2))
return false;
p2 = cycle(p2);
}
return true;
}
void test_repunit(uint32_t digits) {
char* str = malloc(digits + 1);
if (str == 0) {
fprintf(stderr, "Out of memory\n");
exit(1);
}
memset(str, '1', digits);
str[digits] = 0;
mpz_t bignum;
mpz_init_set_str(bignum, str, 10);
free(str);
if (mpz_probab_prime_p(bignum, 10))
printf("R(%u) is probably prime.\n", digits);
else
printf("R(%u) is not prime.\n", digits);
mpz_clear(bignum);
}
int main() {
uint32_t p = 2;
printf("First 19 circular primes:\n");
for (int count = 0; count < 19; ++p) {
if (is_circular_prime(p)) {
if (count > 0)
printf(", ");
printf("%u", p);
++count;
}
}
printf("\n");
printf("Next 4 circular primes:\n");
uint32_t repunit = 1, digits = 1;
for (; repunit < p; ++digits)
repunit = 10 * repunit + 1;
mpz_t bignum;
mpz_init_set_ui(bignum, repunit);
for (int count = 0; count < 4; ) {
if (mpz_probab_prime_p(bignum, 15)) {
if (count > 0)
printf(", ");
printf("R(%u)", digits);
++count;
}
++digits;
mpz_mul_ui(bignum, bignum, 10);
mpz_add_ui(bignum, bignum, 1);
}
mpz_clear(bignum);
printf("\n");
test_repunit(5003);
test_repunit(9887);
test_repunit(15073);
test_repunit(25031);
test_repunit(35317);
test_repunit(49081);
return 0;
}</lang>
- Output:
With GMP 6.2.0, execution time on my system is about 13 minutes (3.2 GHz Quad-Core Intel Core i5, macOS 10.15.4).
First 19 circular primes: 2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933 Next 4 circular primes: R(19), R(23), R(317), R(1031) R(5003) is not prime. R(9887) is not prime. R(15073) is not prime. R(25031) is not prime. R(35317) is not prime. R(49081) is probably prime.
C++
<lang cpp>#include <cstdint>
- include <algorithm>
- include <iostream>
- include <sstream>
- include <gmpxx.h>
typedef mpz_class integer;
bool is_prime(const integer& n, int reps = 50) {
return mpz_probab_prime_p(n.get_mpz_t(), reps);
}
std::string to_string(const integer& n) {
std::ostringstream out; out << n; return out.str();
}
bool is_circular_prime(const integer& p) {
if (!is_prime(p))
return false;
std::string str(to_string(p));
for (size_t i = 0, n = str.size(); i + 1 < n; ++i) {
std::rotate(str.begin(), str.begin() + 1, str.end());
integer p2(str, 10);
if (p2 < p || !is_prime(p2))
return false;
}
return true;
}
integer next_repunit(const integer& n) {
integer p = 1;
while (p < n)
p = 10 * p + 1;
return p;
}
integer repunit(int digits) {
std::string str(digits, '1'); integer p(str); return p;
}
void test_repunit(int digits) {
if (is_prime(repunit(digits), 10))
std::cout << "R(" << digits << ") is probably prime\n";
else
std::cout << "R(" << digits << ") is not prime\n";
}
int main() {
integer p = 2;
std::cout << "First 19 circular primes:\n";
for (int count = 0; count < 19; ++p) {
if (is_circular_prime(p)) {
if (count > 0)
std::cout << ", ";
std::cout << p;
++count;
}
}
std::cout << '\n';
std::cout << "Next 4 circular primes:\n";
p = next_repunit(p);
std::string str(to_string(p));
int digits = str.size();
for (int count = 0; count < 4; ) {
if (is_prime(p, 15)) {
if (count > 0)
std::cout << ", ";
std::cout << "R(" << digits << ")";
++count;
}
p = repunit(++digits);
}
std::cout << '\n';
test_repunit(5003);
test_repunit(9887);
test_repunit(15073);
test_repunit(25031);
test_repunit(35317);
test_repunit(49081);
return 0;
}</lang>
- Output:
First 19 circular primes: 2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933 Next 4 circular primes: R(19), R(23), R(317), R(1031) R(5003) is not prime R(9887) is not prime R(15073) is not prime R(25031) is not prime R(35317) is not prime R(49081) is probably prime
D
<lang D>import std.bigint; import std.stdio;
immutable PRIMES = [
2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97, 101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197, 199, 211, 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 271, 277, 281, 283, 293, 307, 311, 313, 317, 331, 337, 347, 349, 353, 359, 367, 373, 379, 383, 389, 397, 401, 409, 419, 421, 431, 433, 439, 443, 449, 457, 461, 463, 467, 479, 487, 491, 499, 503, 509, 521, 523, 541, 547, 557, 563, 569, 571, 577, 587, 593, 599, 601, 607, 613, 617, 619, 631, 641, 643, 647, 653, 659, 661, 673, 677, 683, 691, 701, 709, 719, 727, 733, 739, 743, 751, 757, 761, 769, 773, 787, 797, 809, 811, 821, 823, 827, 829, 839, 853, 857, 859, 863, 877, 881, 883, 887, 907, 911, 919, 929, 937, 941, 947, 953, 967, 971, 977, 983, 991, 997
];
bool isPrime(BigInt n) {
if (n < 2) {
return false;
}
foreach (p; PRIMES) {
if (n == p) {
return true;
}
if (n % p == 0) {
return false;
}
if (p * p > n) {
return true;
}
}
for (auto m = BigInt(PRIMES[$ - 1]); m * m <= n ; m += 2) {
if (n % m == 0) {
return false;
}
}
return true;
}
// e.g. returns 2341 if n = 1234 BigInt cycle(BigInt n) {
BigInt m = n;
BigInt p = 1;
while (m >= 10) {
p *= 10;
m /= 10;
}
return m + 10 * (n % p);
}
bool isCircularPrime(BigInt p) {
if (!isPrime(p)) {
return false;
}
for (auto p2 = cycle(p); p2 != p; p2 = cycle(p2)) {
if (p2 < p || !isPrime(p2)) {
return false;
}
}
return true;
}
BigInt repUnit(int len) {
BigInt n = 0;
while (len > 0) {
n = 10 * n + 1;
len--;
}
return n;
}
void main() {
writeln("First 19 circular primes:");
int count = 0;
foreach (p; PRIMES) {
if (isCircularPrime(BigInt(p))) {
if (count > 0) {
write(", ");
}
write(p);
count++;
}
}
for (auto p = BigInt(PRIMES[$ - 1]) + 2; count < 19; p += 2) {
if (isCircularPrime(BigInt(p))) {
if (count > 0) {
write(", ");
}
write(p);
count++;
}
}
writeln;
}</lang>
- Output:
First 19 circular primes: 2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933
F#
<lang fsharp> // Circular primes - Nigel Galloway: September 13th., 2021 let fG n g=let rec fG y=if y=g then true else if y>g && isPrime y then fG(10*(y%n)+y/n) else false in fG(10*(g%n)+g/n) let rec fN g l=seq{let g=[for n in g do for g in [1;3;7;9] do let g=n*10+g in yield g] in yield! g|>List.filter(fun n->isPrime n && fG l n); yield! fN g (l*10)} let circP()=seq{yield! [2;3;5;7]; yield! fN [1;3;7;9] 10} circP()|> Seq.take 19 |>Seq.iter(printf "%d "); printfn "" let isPrimeI g=Open.Numeric.Primes.MillerRabin.IsProbablePrime(&g) printf "The first 5 repunit primes are "; Seq.initInfinite((+)1)|>Seq.filter(fun n->isPrimeI((10I**n-1I)/9I))|>Seq.take 5|>Seq.iter(fun n->printf $"R(%d{n}) "); printfn "" </lang>
- Output:
2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933 The first 5 repunit primes are R(2) R(19) R(23) R(317) R(1031)
Factor
Unfortunately Factor's miller-rabin test or bignums aren't quite up to the task of finding the next four circular prime repunits in a reasonable time. It takes ~90 seconds to check R(7)-R(1031).
<lang factor>USING: combinators.short-circuit formatting io kernel lists lists.lazy math math.combinatorics math.functions math.parser math.primes sequences sequences.extras ;
! Create an ordered infinite lazy list of circular prime ! "candidates" -- the numbers 2, 3, 5 followed by numbers ! composed of only the digits 1, 3, 7, and 9.
- candidates ( -- list )
L{ "2" "3" "5" "7" } 2 lfrom
[ "1379" swap selections >list ] lmap-lazy lconcat lappend ;
- circular-prime? ( str -- ? )
all-rotations {
[ [ infimum ] [ first = ] bi ]
[ [ string>number prime? ] all? ]
} 1&& ;
- circular-primes ( -- list )
candidates [ circular-prime? ] lfilter ;
- prime-repunits ( -- list )
7 lfrom [ 10^ 1 - 9 / prime? ] lfilter ;
"The first 19 circular primes are:" print 19 circular-primes ltake [ write bl ] leach nl nl
"The next 4 circular primes, in repunit format, are:" print 4 prime-repunits ltake [ "R(%d) " printf ] leach nl</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)
Forth
Forth only supports native sized integers, so we only implement the first part of the task. <lang Forth> create 235-wheel 6 c, 4 c, 2 c, 4 c, 2 c, 4 c, 6 c, 2 c,
does> swap 7 and + c@ ;
0 1 2constant init-235 \ roll 235 wheel at position 1
- next-235 over 235-wheel + swap 1+ swap ;
\ check that n is prime excepting multiples of 2, 3, 5.
- sq dup * ;
- wheel-prime? ( n -- f )
>r init-235 begin
next-235
dup sq r@ > if rdrop 2drop true exit then
r@ over mod 0= if rdrop 2drop false exit then
again ;
- prime? ( n -- f )
dup 2 < if drop false exit then dup 2 mod 0= if 2 = exit then dup 3 mod 0= if 3 = exit then dup 5 mod 0= if 5 = exit then wheel-prime? ;
- log10^ ( n -- 10^[log n], log n )
dup 0<= abort" log10^: argument error."
1 0 rot
begin dup 9 > while
>r swap 10 * swap 1+ r> 10 /
repeat drop ;
- log10 ( n -- n ) log10^ nip ;
- rotate ( n -- n )
dup log10^ drop /mod swap 10 * + ;
- prime-rotation? ( p0 p -- f )
tuck <= swap prime? and ;
- circular? ( n -- f ) \ assume n is not a multiple of 2, 3, 5
dup wheel-prime? invert
if drop false exit
then dup >r true
over log10 0 ?do
swap rotate j over prime-rotation? rot and
loop nip rdrop ;
- .primes
2 . 3 . 5 .
16 init-235 \ -- count, [n1 n2] as 2,3,5 wheel
begin
next-235 dup circular?
if dup . rot 1- -rot
then
third 0= until 2drop drop ;
." The first 19 circular primes are:" cr .primes cr bye </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
FreeBASIC
<lang freebasic>#define floor(x) ((x*2.0-0.5)Shr 1)
Function isPrime(Byval p As Integer) As Boolean
If p < 2 Then Return False
If p Mod 2 = 0 Then Return p = 2
If p Mod 3 = 0 Then Return p = 3
Dim As Integer d = 5
While d * d <= p
If p Mod d = 0 Then Return False Else d += 2
If p Mod d = 0 Then Return False Else d += 4
Wend
Return True
End Function
Function isCircularPrime(Byval p As Integer) As Boolean
Dim As Integer n = floor(Log(p)/Log(10))
Dim As Integer m = 10^n, q = p
For i As Integer = 0 To n
If (q < p Or Not isPrime(q)) Then Return false
q = (q Mod m) * 10 + floor(q / m)
Next i
Return true
End Function
Dim As Integer p = 2, dp = 1, cont = 0 Print("Primeros 19 primos circulares:") While cont < 19
If isCircularPrime(p) Then Print p;" "; : cont += 1 p += dp: dp = 2
Wend Sleep</lang>
- Output:
Primeros 19 primos circulares: 2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933
Go
<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{}
// 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
}
func main() {
fmt.Println("The first 19 circular primes are:")
digits := [4]int{1, 3, 7, 9}
q := []int{1, 2, 3, 5, 7, 9} // queue the numbers to be examined
fq := []int{1, 2, 3, 5, 7, 9} // also queue the corresponding first digits
count := 0
for {
f := q[0] // peek first element
fd := fq[0] // peek first digit
if isPrime(f) && isCircular(f) {
circs = append(circs, f)
count++
if count == 19 {
break
}
}
copy(q, q[1:]) // pop first element
q = q[:len(q)-1] // reduce length by 1
copy(fq, fq[1:]) // ditto for first digit queue
fq = fq[:len(fq)-1]
if f == 2 || f == 5 { // if digits > 1 can't contain a 2 or 5
continue
}
// add numbers with one more digit to queue
// only numbers whose last digit >= first digit need be added
for _, d := range digits {
if d >= fd {
q = append(q, f*10+d)
fq = append(fq, fd)
}
}
}
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
Haskell
Uses arithmoi Library: http://hackage.haskell.org/package/arithmoi-0.10.0.0 <lang haskell>import Math.NumberTheory.Primes (Prime, unPrime, nextPrime) import Math.NumberTheory.Primes.Testing (isPrime, millerRabinV) import Text.Printf (printf)
rotated :: [Integer] -> [Integer] rotated xs
| any (< head xs) xs = [] | otherwise = map asNum $ take (pred $ length xs) $ rotate xs where rotate [] = [] rotate (d:ds) = ds <> [d] : rotate (ds <> [d])
asNum :: [Integer] -> Integer asNum [] = 0 asNum n@(d:ds)
| all (==1) n = read $ concatMap show n | otherwise = (d * (10 ^ length ds)) + asNum ds
digits :: Integer -> [Integer] digits 0 = [] digits n = digits d <> [r]
where (d, r) = n `quotRem` 10
isCircular :: Bool -> Integer -> Bool isCircular repunit n
| repunit = millerRabinV 0 n | n < 10 = True | even n = False | null rotations = False | any (<n) rotations = False | otherwise = all isPrime rotations where rotations = rotated $ digits n
repunits :: [Integer] repunits = go 2
where go n = asNum (replicate n 1) : go (succ n)
asRepunit :: Int -> Integer asRepunit n = asNum $ replicate n 1
main :: IO () main = do
printf "The first 19 circular primes are:\n%s\n\n" $ circular primes printf "The next 4 circular primes, in repunit format are:\n" mapM_ (printf "R(%d) ") $ reps repunits printf "\n\nThe following repunits are probably circular primes:\n" mapM_ (uncurry (printf "R(%d) : %s\n") . checkReps) [5003, 9887, 15073, 25031, 35317, 49081] where primes = map unPrime [nextPrime 1..] circular = show . take 19 . filter (isCircular False) reps = map (sum . digits). tail . take 5 . filter (isCircular True) checkReps = (,) <$> id <*> show . isCircular True . asRepunit</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 ./circular_primes 277.56s user 1.82s system 97% cpu 4:47.91 total
J
<lang J>
R=: [: ". 'x' ,~ #&'1' assert 11111111111111111111111111111111x -: R 32
Filter=: (#~`)(`:6)
rotations=: (|."0 1~ i.@#)&.(10&#.inv) assert 123 231 312 -: rotations 123
primes_less_than=: i.&.:(p:inv) assert 2 3 5 7 11 -: primes_less_than 12
NB. circular y --> y is the order of magnitude.
circular=: monad define
P25=: ([: -. (0 e. 1 3 7 9 e.~ 10 #.inv ])&>)Filter primes_less_than 10^y NB. Q25 are primes with 1 3 7 9 digits
P=: 2 5 , P25
en=: # P
group=: en # 0
next=: 1
for_i. i. # group do.
if. 0 = i { group do. NB. if untested
j =: P i. rotations i { P NB. j are the indexes of the rotated numbers in the list of primes
if. en e. j do. NB. if any are unfound
j=: j -. en NB. prepare to mark them all as searched, and failed.
g=: _1
else.
g=: next NB. mark the set as found in a new group. Because we can.
next=: >: next
end.
group=: g j} group NB. apply the tested mark
end.
end.
group </. P
) </lang> J lends itself to investigation. Demonstrate and play with the definitions.
circular 3 NB. the values in the long list have a composite under rotation
┌─┬─┬─┬─┬──┬─────┬─────┬──────────────────────────────────────────────────────────────────────────┬─────┬─────┬───────────┬───────────┬───────────┬───────────┐
│2│5│3│7│11│13 31│17 71│19 137 139 173 179 191 193 313 317 331 379 397 739 773 797 911 937 977 997│37 73│79 97│113 131 311│197 719 971│199 919 991│337 373 733│
└─┴─┴─┴─┴──┴─────┴─────┴──────────────────────────────────────────────────────────────────────────┴─────┴─────┴───────────┴───────────┴───────────┴───────────┘
circular 2 NB. VV
┌─┬─┬─┬─┬──┬─────┬─────┬──┬─────┬─────┐
│2│5│3│7│11│13 31│17 71│19│37 73│79 97│
└─┴─┴─┴─┴──┴─────┴─────┴──┴─────┴─────┘
q: 91 NB. factor the lone entry
7 13
RC=: circular 8
{: RC NB. tail
┌─────────────────────────────────────────┐
│199933 319993 331999 933199 993319 999331│
└─────────────────────────────────────────┘
(< >./)@:(#&>) Filter circular 3 NB. remove the box containing most items
┌─┬─┬─┬─┬──┬─────┬─────┬─────┬─────┬───────────┬───────────┬───────────┬───────────┐
│2│5│3│7│11│13 31│17 71│37 73│79 97│113 131 311│197 719 971│199 919 991│337 373 733│
└─┴─┴─┴─┴──┴─────┴─────┴─────┴─────┴───────────┴───────────┴───────────┴───────────┘
] CPS=: {.&> (< >./)@:(#&>) Filter RC NB. first 19 circular primes
2 5 3 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933
# CPS NB. yes, 19 found.
19
Brief investigation into repunits.
Note 'The current Miller-Rabin test implemented in c is insufficient for this task' R=: ([: ". 'x' ,~ #&'1')&> (;q:@R)&> 500 |limit error | (;q:@R)&>500 ) boxdraw_j_ 0 Filter=: (#~`)(`:6) R=: ([: ". 'x' ,~ #&'1')&> (; q:@R)&> (0 ~: 3&|)Filter >: i. 26 NB. factor some repunits ┌──┬─────────────────────────────────┐ │1 │ │ ├──┼─────────────────────────────────┤ │2 │11 │ ├──┼─────────────────────────────────┤ │4 │11 101 │ ├──┼─────────────────────────────────┤ │5 │41 271 │ ├──┼─────────────────────────────────┤ │7 │239 4649 │ ├──┼─────────────────────────────────┤ │8 │11 73 101 137 │ ├──┼─────────────────────────────────┤ │10│11 41 271 9091 │ ├──┼─────────────────────────────────┤ │11│21649 513239 │ ├──┼─────────────────────────────────┤ │13│53 79 265371653 │ ├──┼─────────────────────────────────┤ │14│11 239 4649 909091 │ ├──┼─────────────────────────────────┤ │16│11 17 73 101 137 5882353 │ ├──┼─────────────────────────────────┤ │17│2071723 5363222357 │ ├──┼─────────────────────────────────┤ │19│1111111111111111111 │ ├──┼─────────────────────────────────┤ │20│11 41 101 271 3541 9091 27961 │ ├──┼─────────────────────────────────┤ │22│11 11 23 4093 8779 21649 513239 │ ├──┼─────────────────────────────────┤ │23│11111111111111111111111 │ ├──┼─────────────────────────────────┤ │25│41 271 21401 25601 182521213001 │ ├──┼─────────────────────────────────┤ │26│11 53 79 859 265371653 1058313049│ └──┴─────────────────────────────────┘ NB. R(2) R(19), R(23) are probably prime.
Java
<lang java>import java.math.BigInteger; import java.util.Arrays;
public class CircularPrimes {
public static void main(String[] args) {
System.out.println("First 19 circular primes:");
int p = 2;
for (int count = 0; count < 19; ++p) {
if (isCircularPrime(p)) {
if (count > 0)
System.out.print(", ");
System.out.print(p);
++count;
}
}
System.out.println();
System.out.println("Next 4 circular primes:");
int repunit = 1, digits = 1;
for (; repunit < p; ++digits)
repunit = 10 * repunit + 1;
BigInteger bignum = BigInteger.valueOf(repunit);
for (int count = 0; count < 4; ) {
if (bignum.isProbablePrime(15)) {
if (count > 0)
System.out.print(", ");
System.out.printf("R(%d)", digits);
++count;
}
++digits;
bignum = bignum.multiply(BigInteger.TEN);
bignum = bignum.add(BigInteger.ONE);
}
System.out.println();
testRepunit(5003);
testRepunit(9887);
testRepunit(15073);
testRepunit(25031);
}
private static boolean isPrime(int n) {
if (n < 2)
return false;
if (n % 2 == 0)
return n == 2;
if (n % 3 == 0)
return n == 3;
for (int p = 5; p * p <= n; p += 4) {
if (n % p == 0)
return false;
p += 2;
if (n % p == 0)
return false;
}
return true;
}
private static int cycle(int n) {
int m = n, p = 1;
while (m >= 10) {
p *= 10;
m /= 10;
}
return m + 10 * (n % p);
}
private static boolean isCircularPrime(int p) {
if (!isPrime(p))
return false;
int p2 = cycle(p);
while (p2 != p) {
if (p2 < p || !isPrime(p2))
return false;
p2 = cycle(p2);
}
return true;
}
private static void testRepunit(int digits) {
BigInteger repunit = repunit(digits);
if (repunit.isProbablePrime(15))
System.out.printf("R(%d) is probably prime.\n", digits);
else
System.out.printf("R(%d) is not prime.\n", digits);
}
private static BigInteger repunit(int digits) {
char[] ch = new char[digits];
Arrays.fill(ch, '1');
return new BigInteger(new String(ch));
}
}</lang>
- Output:
First 19 circular primes: 2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933 Next 4 circular primes: R(19), R(23), R(317), R(1031) R(5003) is not prime. R(9887) is not prime. R(15073) is not prime. R(25031) is not prime.
jq
Works with gojq, the Go implementation of jq
jq's integer-arithmetic accuracy is only sufficient for the first task; gojq has the accuracy for the second task but is not well-suited for it. Nevertheless, the code for solving both tasks is shown.
For an implementation of `is_prime` suitable for the first task, see for example Erdős-primes#jq. <lang jq> def is_circular_prime:
def circle: range(0;length) as $i | .[$i:] + .[:$i]; tostring as $s | [$s|circle|tonumber] as $c | . == ($c|min) and all($c|unique[]; is_prime);
def circular_primes:
2, (range(3; infinite; 2) | select(is_circular_prime));
- Probably only useful with unbounded-precision integer arithmetic:
def repunits:
1 | recurse(10*. + 1);
</lang> The first task <lang jq>limit(19; circular_primes)</lang> The second task (viewed independently): <lang jq> last(limit(19; circular_primes)) as $max | limit(4; repunits | select(. > $max and is_prime)) | "R(\(tostring|length))"</lang>
- Output:
For the first task:
2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933
Julia
Note that the evalpoly function used in this program was added in Julia 1.4 <lang julia>using Lazy, Primes
function iscircularprime(n)
!isprime(n) && return false dig = digits(n) return all(i -> (m = evalpoly(10, circshift(dig, i))) >= n && isprime(m), 1:length(dig)-1)
end
filtcircular(n, rang) = Int.(collect(take(n, filter(iscircularprime, rang)))) isprimerepunit(n) = isprime(evalpoly(BigInt(10), ones(Int, n))) filtrep(n, rang) = collect(take(n, filter(isprimerepunit, rang)))
println("The first 19 circular primes are:\n", filtcircular(19, Lazy.range(2))) print("\nThe next 4 circular primes, in repunit format, are: ",
mapreduce(n -> "R($n) ", *, filtrep(4, Lazy.range(6))))
println("\n\nChecking larger repunits:") for i in [5003, 9887, 15073, 25031, 35317, 49081]
println("R($i) is ", isprimerepunit(i) ? "prime." : "not prime.")
end
</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) Checking larger repunits: R(5003) is not prime. R(9887) is not prime. R(15073) is not prime. R(25031) is not prime. R(35317) is not prime. R(49081) is prime.
Kotlin
<lang scala>import java.math.BigInteger
val SMALL_PRIMES = listOf(
2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97, 101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197, 199, 211, 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 271, 277, 281, 283, 293, 307, 311, 313, 317, 331, 337, 347, 349, 353, 359, 367, 373, 379, 383, 389, 397, 401, 409, 419, 421, 431, 433, 439, 443, 449, 457, 461, 463, 467, 479, 487, 491, 499, 503, 509, 521, 523, 541, 547, 557, 563, 569, 571, 577, 587, 593, 599, 601, 607, 613, 617, 619, 631, 641, 643, 647, 653, 659, 661, 673, 677, 683, 691, 701, 709, 719, 727, 733, 739, 743, 751, 757, 761, 769, 773, 787, 797, 809, 811, 821, 823, 827, 829, 839, 853, 857, 859, 863, 877, 881, 883, 887, 907, 911, 919, 929, 937, 941, 947, 953, 967, 971, 977, 983, 991, 997
)
fun isPrime(n: BigInteger): Boolean {
if (n < 2.toBigInteger()) {
return false
}
for (sp in SMALL_PRIMES) {
val spb = sp.toBigInteger()
if (n == spb) {
return true
}
if (n % spb == BigInteger.ZERO) {
return false
}
if (n < spb * spb) {
//if (n > SMALL_PRIMES.last().toBigInteger()) {
// println("Next: $n")
//}
return true
}
}
return n.isProbablePrime(10)
}
fun cycle(n: BigInteger): BigInteger {
var m = n
var p = 1
while (m >= BigInteger.TEN) {
p *= 10
m /= BigInteger.TEN
}
return m + BigInteger.TEN * (n % p.toBigInteger())
}
fun isCircularPrime(p: BigInteger): Boolean {
if (!isPrime(p)) {
return false
}
var p2 = cycle(p)
while (p2 != p) {
if (p2 < p || !isPrime(p2)) {
return false
}
p2 = cycle(p2)
}
return true
}
fun testRepUnit(digits: Int) {
var repUnit = BigInteger.ONE
var count = digits - 1
while (count > 0) {
repUnit = BigInteger.TEN * repUnit + BigInteger.ONE
count--
}
if (isPrime(repUnit)) {
println("R($digits) is probably prime.")
} else {
println("R($digits) is not prime.")
}
}
fun main() {
println("First 19 circular primes:")
var p = 2
var count = 0
while (count < 19) {
if (isCircularPrime(p.toBigInteger())) {
if (count > 0) {
print(", ")
}
print(p)
count++
}
p++
}
println()
println("Next 4 circular primes:")
var repUnit = BigInteger.ONE
var digits = 1
count = 0
while (repUnit < p.toBigInteger()) {
repUnit = BigInteger.TEN * repUnit + BigInteger.ONE
digits++
}
while (count < 4) {
if (isPrime(repUnit)) {
print("R($digits) ")
count++
}
repUnit = BigInteger.TEN * repUnit + BigInteger.ONE
digits++
}
println()
testRepUnit(5003) testRepUnit(9887) testRepUnit(15073) testRepUnit(25031) testRepUnit(35317) testRepUnit(49081)
}</lang>
- Output:
First 19 circular primes: 2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933 Next 4 circular primes: R(19) R(23) R(317) R(1031) R(5003) is not prime. R(9887) is not prime. R(15073) is not prime. R(25031) is not prime. R(35317) is not prime.
Lua
<lang lua>-- Circular primes, in Lua, 6/22/2020 db local function isprime(n)
if n < 2 then return false end if n % 2 == 0 then return n==2 end if n % 3 == 0 then return n==3 end for f = 5, math.sqrt(n), 6 do if n % f == 0 or n % (f+2) == 0 then return false end end return true
end
local function iscircularprime(p)
local n = math.floor(math.log10(p)) local m, q = 10^n, p for i = 0, n do if (q < p or not isprime(q)) then return false end q = (q % m) * 10 + math.floor(q / m) end return true
end
local p, dp, list, N = 2, 1, {}, 19 while #list < N do
if iscircularprime(p) then list[#list+1] = p end p, dp = p + dp, 2
end print(table.concat(list, ", "))</lang>
- Output:
2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933
Mathematica / Wolfram Language
<lang Mathematica>ClearAll[RepUnit, CircularPrimeQ] RepUnit[n_] := (10^n - 1)/9 CircularPrimeQ[n_Integer] := Module[{id = IntegerDigits[n], nums, t},
AllTrue[
Range[Length[id]]
,
Function[{z},
t = FromDigits[RotateLeft[id, z]];
If[t < n,
False
,
PrimeQ[t]
]
]
]
]
Select[Range[200000], CircularPrimeQ]
res = {}; Dynamic[res] Do[
If[CircularPrimeQ[RepUnit[n]], AppendTo[res, n]]
,
{n, 1000}
]
Scan[Print@*PrimeQ@*RepUnit, {5003, 9887, 15073, 25031, 35317, 49081}]</lang>
- Output:
{2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933}
{2, 19, 23, 317}
False
False
False
False
False
True
Nim
<lang Nim>import bignum import strformat
const SmallPrimes = [
2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97, 101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197, 199, 211, 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 271, 277, 281, 283, 293, 307, 311, 313, 317, 331, 337, 347, 349, 353, 359, 367, 373, 379, 383, 389, 397, 401, 409, 419, 421, 431, 433, 439, 443, 449, 457, 461, 463, 467, 479, 487, 491, 499, 503, 509, 521, 523, 541, 547, 557, 563, 569, 571, 577, 587, 593, 599, 601, 607, 613, 617, 619, 631, 641, 643, 647, 653, 659, 661, 673, 677, 683, 691, 701, 709, 719, 727, 733, 739, 743, 751, 757, 761, 769, 773, 787, 797, 809, 811, 821, 823, 827, 829, 839, 853, 857, 859, 863, 877, 881, 883, 887, 907, 911, 919, 929, 937, 941, 947, 953, 967, 971, 977, 983, 991, 997]
let
One = newInt(1) Two = newInt(2) Ten = newInt(10)
- ---------------------------------------------------------------------------------------------------
proc isPrime(n: Int): bool =
if n < Two: return false
for sp in SmallPrimes: # let spb = newInt(sp) if n == sp: return true if (n mod sp).isZero: return false if n < sp * sp: return true
result = probablyPrime(n, 25) != 0
- ---------------------------------------------------------------------------------------------------
proc cycle(n: Int): Int =
var m = n var p = 1 while m >= Ten: p *= 10 m = m div 10 result = m + Ten * (n mod p)
- ---------------------------------------------------------------------------------------------------
proc isCircularPrime(p: Int): bool =
if not p.isPrime(): return false
var p2 = cycle(p)
while p2 != p:
if p2 < p or not p2.isPrime():
return false
p2 = cycle(p2)
result = true
- ---------------------------------------------------------------------------------------------------
proc testRepunit(digits: int) =
var repunit = One
var count = digits - 1
while count > 0:
repunit = Ten * repunit + One
dec count
if repunit.isPrime():
echo fmt"R({digits}) is probably prime."
else:
echo fmt"R({digits}) is not prime."
- ---------------------------------------------------------------------------------------------------
echo "First 19 circular primes:" var p = 2 var line = "" var count = 0 while count < 19:
if newInt(p).isCircularPrime():
if count > 0: line.add(", ")
line.add($p)
inc count
inc p
echo line
echo "" echo "Next 4 circular primes:" var repunit = One var digits = 1 while repunit < p:
repunit = Ten * repunit + One inc digits
line = "" count = 0 while count < 4:
if repunit.isPrime():
if count > 0: line.add(' ')
line.add(fmt"R({digits})")
inc count
repunit = Ten * repunit + One
inc digits
echo line
echo "" testRepUnit(5003) testRepUnit(9887) testRepUnit(15073) testRepUnit(25031) testRepUnit(35317) testRepUnit(49081)</lang>
- Output:
First 19 circular primes: 2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933 Next 4 circular primes: R(19) R(23) R(317) R(1031) R(5003) is not prime. R(9887) is not prime. R(15073) is not prime. R(25031) is not prime. R(35317) is not prime. R(49081) is probably prime.
Pascal
Free Pascal
Only test up to 10 Digit numbers.RepUnit test with gmp to boring.
Using a base 4 downcounter to create the numbers with more than one digit.
<lang pascal>{$IFDEF FPC}
{$MODE DELPHI}{$OPTIMIZATION ON,ALL}
{$ENDIF} {$IFDEF WINDOWS}
{$APPTYPE CONSOLE}
{$ENDIF}
const
MAXCNTOFDIGITS = 6;//10;//11 MAXDGTVAL = 3; conv : array[0..MAXDGTVAL+1] of byte = (9,7,3,1,0);
type
tDigits = array[0..47] of byte; tUint64 = NativeUint;
var
digits, rotdigits : tDigits; Pot10 : array[0..19] of tUint64; CheckNum : array[0..19] of tUint64; Found : array[0..24] of tUint64; Count : tUint64;
function isPrimeSpecial(n:tUint64):boolean; // no multiples of 2,5 const
wheeldiff : array[0..7] of Uint32 = (+6,+4,+2,+4,+2,+4,+6,+2);
var
p: tUint64; flipflop : Int32;
begin
if n< 32 then
EXIT(n in [2,3,5,7,11,13,17,19,23,29,31])
else
begin
IF (n mod 3 = 0 )then
EXIT(false);
result := true;
p := 1;
flipflop := 6;
while result do
Begin
p += wheeldiff[flipflop];
if p*p>n then
BREAK;
result := n mod p <> 0;
flipflop -= 1;
if flipflop<0 then
flipflop :=7;
end
end
end;
function CheckInList(num:tUint64;MaxIdx:integer):Boolean; Begin
CheckInList := true;
end; procedure CheckOne(MaxIdx:integer); var
pRot : pbyte; i,j,cv : integer; Num,First : tUint64;
begin
pRot := @rotDigits[0]; i:= MaxIdx; inc(maxIdx); j := 0; Num := 0; repeat cv := conv[digits[i]]; dec(i); pRot[j]:= cv; num := num*10+cv; pRot[j+MaxIdx]:= cv; inc(j); until i < 0;
//first create circular numbers to minimize prime checks
i := MaxIdx;
First := num;
j := 0;
repeat
Num := (Num - rotDigits[i-MaxIdx]*Pot10[MaxIdx-1])*10+rotDigits[i];
//num was checked before
if num < First then
EXIT;
CheckNum[j] := num;
inc(j);
inc(i);
until i >= 2*MaxIdx-1;
If isPrimeSpecial(First) then
Begin
repeat
dec(j);
IF Not(isPrimeSpecial(CheckNum[j])) then
EXIT;
until j = 0;
cv := Count;
Found[cv] := First;
inc(count);
end;
end;
var
idx,MaxIDx,dgt : NativeInt; Pot_ten : tUint64;
begin
fillchar(digits,Sizeof(digits),chr(MAXDGTVAL));
Pot_ten := 1;
For idx := Low(Pot10) to High(Pot10) do
begin
Pot10[idx] := Pot_ten;
Pot_ten := Pot_ten*10;
end;
For idx := 2 to 10 do
if isPrimeSpecial(idx) then
begin
Found[Count]:= idx;
inc(count);
end;
maxIdx := 1;
repeat
CheckOne(MaxIdx);
idx := 0;
repeat
dgt := digits[idx]-1;
if dgt >=0 then
break;
digits[idx] := MAXDGTVAL;
inc(idx);
until idx >MAXCNTOFDIGITS-1;
if idx > MAXCNTOFDIGITS-1 then
BREAK;
if idx<=MaxIDX then
digits[idx] := dgt
else
Maxidx := idx;
until false;
CheckOne(MaxIdx);
For idx := 0 to count-1 do
writeln(idx+1:5,Found[idx]:12);
end.</lang>
- @ Tio.Run:
1 2
2 3
3 5
4 7
5 11
6 13
7 17
8 37
9 79
10 113
11 197
12 199
13 337
14 1193
15 3779
16 11939
17 19937
18 193939
19 199933
Real time: 0.073 s digit Count MAXCNTOFDIGITS = 6;
Real time: 2.691 s digit Count MAXCNTOFDIGITS = 10;
Real time: 24.545 s digit Count MAXCNTOFDIGITS = 11;
Perl
<lang perl>use strict; use warnings; use feature 'say'; use List::Util 'min'; use ntheory 'is_prime';
sub rotate { my($i,@a) = @_; join , @a[$i .. @a-1, 0 .. $i-1] }
sub isCircular {
my ($n) = @_;
return 0 unless is_prime($n);
my @circular = split //, $n;
return 0 if min(@circular) < $circular[0];
for (1 .. scalar @circular) {
my $r = join , rotate($_,@circular);
return 0 unless is_prime($r) and $r >= $n;
}
1
}
say "The first 19 circular primes are:"; for ( my $i = 1, my $count = 0; $count < 19; $i++ ) {
++$count and print "$i " if isCircular($i);
}
say "\n\nThe next 4 circular primes, in repunit format, are:"; for ( my $i = 7, my $count = 0; $count < 4; $i++ ) {
++$count and say "R($i)" if is_prime 1 x $i
}
say "\nRepunit testing:";
for (5003, 9887, 15073, 25031, 35317, 49081) {
say "R($_): Prime? " . (is_prime 1 x $_ ? 'True' : 'False');
}</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 R(9887): Prime? False R(15073): Prime? False R(25031): Prime? False R(35317): Prime? False R(49081): Prime? True
Phix
with javascript_semantics 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() while length(c)<4 do repunit(z,n) if mpz_prime(z) 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)})
- 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
(It is probably quite unwise to throw this in the direction of pwa/p2js, I am not even going to bother trying.)
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_prime(z,1) printf(1,"R(%d) : %t (%s)\n", {ti, bPrime, elapsed(time()-t0)}) end for
- Output:
64-bit can only cope with the first five (it terminates abruptly on the sixth)
For comparison, the above Julia code (8/4/20, 64 bit) manages 1s, 5.6s, 15s, 50s, 1 min 50s (and 1 hour 45 min 40s) on the same box.
And Perl (somehow) manages 0/0/0/55s/0/21 mins 35 secs...
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
Prolog
Uses a miller-rabin primality tester that I wrote which includes a trial division pass for small prime factors. One could substitute with e.g. the Miller Rabin Primaliy Test task.
The circular(P) predicate generates all circular primes; for those > 1e6, it returns it as a term r(K) for repunit K.
Also the code is smart in that it only checks primes > 9 that are composed of the digits 1, 3, 7, and 9. <lang Prolog> ?- use_module(library(primality)).
circular(N) :- member(N, [2, 3, 5, 7]). circular(N) :-
limit(15, (
candidate(N),
N > 9,
circular_prime(N))).
circular(r(K)) :-
between(6, inf, K), N is (10**K - 1) div 9, prime(N).
candidate(0). candidate(N) :-
candidate(M), member(D, [1, 3, 7, 9]), N is 10*M + D.
circular_prime(N) :-
K is floor(log10(N)) + 1, circular_prime(N, N, K).
circular_prime(_, _, 0) :- !. circular_prime(P0, P, K) :-
P >= P0, prime(P), rotate(P, Q), succ(DecK, K), circular_prime(P0, Q, DecK).
rotate(N, M) :-
D is floor(log10(N)), divmod(N, 10, Q, R), M is R*10**D + Q.
main :-
findall(P, limit(23, circular(P)), S),
format("The first 23 circular primes:~n~w~n", [S]),
halt.
?- main. </lang>
- Output:
The first 23 circular primes: [2,3,5,7,11,13,17,37,79,113,197,199,337,1193,3779,11939,19937,193939,199933,r(19),r(23),r(317),r(1031)]
Python
<lang Python> import random
def is_Prime(n):
""" Miller-Rabin primality test.
A return value of False means n is certainly not prime. A return value of
True means n is very likely a prime.
"""
if n!=int(n):
return False
n=int(n)
#Miller-Rabin test for prime
if n==0 or n==1 or n==4 or n==6 or n==8 or n==9:
return False
if n==2 or n==3 or n==5 or n==7:
return True
s = 0
d = n-1
while d%2==0:
d>>=1
s+=1
assert(2**s * d == n-1)
def trial_composite(a):
if pow(a, d, n) == 1:
return False
for i in range(s):
if pow(a, 2**i * d, n) == n-1:
return False
return True
for i in range(8):#number of trials
a = random.randrange(2, n)
if trial_composite(a):
return False
return True
def isPrime(n: int) -> bool:
https://www.geeksforgeeks.org/python-program-to-check-whether-a-number-is-prime-or-not/
# Corner cases
if (n <= 1) :
return False
if (n <= 3) :
return True
# This is checked so that we can skip
# middle five numbers in below loop
if (n % 2 == 0 or n % 3 == 0) :
return False
i = 5
while(i * i <= n) :
if (n % i == 0 or n % (i + 2) == 0) :
return False
i = i + 6
return True
def rotations(n: int)-> set((int,)):
>>> {123, 231, 312} == rotations(123)
True
a = str(n)
return set(int(a[i:] + a[:i]) for i in range(len(a)))
def isCircular(n: int) -> bool:
>>> [isCircular(n) for n in (11, 31, 47,)]
[True, True, False]
return all(isPrime(int(o)) for o in rotations(n))
from itertools import product
def main():
result = [2, 3, 5, 7]
first = '137'
latter = '1379'
for i in range(1, 6):
s = set(int(.join(a)) for a in product(first, *((latter,) * i)))
while s:
a = s.pop()
b = rotations(a)
if isCircular(a):
result.append(min(b))
s -= b
result.sort()
return result
assert [2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933] == main()
repunit = lambda n: int('1' * n)
def repmain(n: int) -> list:
returns the first n repunit primes, probably.
result = []
i = 2
while len(result) < n:
if is_Prime(repunit(i)):
result.append(i)
i += 1
return result
assert [2, 19, 23, 317, 1031] == repmain(5)
- because this Miller-Rabin test is already on rosettacode there's no good reason to test the longer repunits.
</lang>
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
- 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 -> $i {
return False unless .is-prime and $_ >= n given @circular.rotate($i).join;
}
True
}
say "The first 19 circular primes are:"; say ((2..*).hyper.grep: { isCircular $_ })[^19];
say "\nThe next 4 circular primes, in repunit format, are:"; loop ( my $i = 7, my $count = 0; $count < 4; $i++ ) {
++$count, say "R($i)" if (1 x $i).is-prime
}
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 hip= 1000000 /* " " " " " " */ call genP /*gen primes up to hp (200,000). */ q= 024568 /*digs that most circular P can't have.*/ found= 0; $= /*found: circular P count; $: a list.*/
do j=1 until found==N; p= @.j /* [↓] traipse through all the primes.*/
if p>9 & verify(p, q, 'M')>0 then iterate /*Does J contain forbidden digs? Skip.*/
if \circP(p) then iterate /*Not circular? Then skip this number.*/
found= found + 1 /*bump the count of circular primes. */
$= $ p /*add this prime number ──► $ list. */
end /*j*/ /*at this point, $ has a leading blank.*/
say center(' first ' found " circular primes ", 79, '─') say strip($) exit 0 /*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 var x f 2 y /*parse X number, rotating the digits*/
x= y || f /*construct a new possible circular P. */
if x<ox then return 0 /*is number < the original? ¬ circular*/
if \!.x then return 0 /* " " not prime? ¬ circular*/
end /*length(x)···*/
return 1 /*passed all tests, X is a circular P.*/
/*──────────────────────────────────────────────────────────────────────────────────────*/ genP: @.1=2; @.2=3; @.3=5; @.4=7; @.5=11; @.6=13; @.7=17; @.8=19 /*assign Ps; #Ps*/
!.= 0; !.2=1; !.3=1; !.5=1; !.7=1; !.11=1; !.13=1; !.17=1; !.19=1 /* " primality*/
#= 8; sq.#= @.# **2 /*number of primes so far; prime square*/
do j=@.#+4 by 2 to hip; parse var j -1 _ /*get last decimal digit of J. */
if _==5 then iterate; if j// 3==0 then iterate; if j// 7==0 then iterate
if j//11==0 then iterate; if j//13==0 then iterate; if j//17==0 then iterate
do k=8 while sq.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. */
#= #+1; !.j= 1; sq.#= j*j; @.#= 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
Ring
<lang ring> see "working..." + nl see "First 19 circular numbers are:" + nl n = 0 row = 0 Primes = []
while row < 19
n++
flag = 1
nStr = string(n)
lenStr = len(nStr)
for m = 1 to lenStr
leftStr = left(nStr,m)
rightStr = right(nStr,lenStr-m)
strOk = rightStr + leftStr
nOk = number(strOk)
ind = find(Primes,nOk)
if ind < 1 and strOk != nStr
add(Primes,nOk)
ok
if not isprimeNumber(nOk) or ind > 0
flag = 0
exit
ok
next
if flag = 1
row++
see "" + n + " "
if row%5 = 0
see nl
ok
ok
end
see nl + "done..." + nl
func isPrimeNumber(num)
if (num <= 1) return 0 ok
if (num % 2 = 0) and (num != 2) return 0 ok
for i = 2 to sqrt(num)
if (num % i = 0) return 0 ok
next
return 1
</lang>
- Output:
working... First 19 circular numbers are: 2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933 done...
Ruby
It takes more then 25 minutes to verify that R49081 is probably prime - omitted here. <lang ruby>require 'gmp' require 'prime' candidate_primes = Enumerator.new do |y|
DIGS = [1,3,7,9]
[2,3,5,7].each{|n| y << n.to_s}
(2..).each do |size|
DIGS.repeated_permutation(size) do |perm|
y << perm.join if (perm == min_rotation(perm)) && GMP::Z(perm.join).probab_prime? > 0
end
end
end
def min_rotation(ar) = Array.new(ar.size){|n| ar.rotate(n)}.min
def circular?(num_str)
chars = num_str.chars
return GMP::Z(num_str).probab_prime? > 0 if chars.all?("1")
chars.size.times.all? do
GMP::Z(chars.rotate!.join).probab_prime? > 0
# chars.rotate!.join.to_i.prime?
end
end
puts "First 19 circular primes:" puts candidate_primes.lazy.select{|cand| circular?(cand)}.take(19).to_a.join(", "),"" puts "First 5 prime repunits:" reps = Prime.each.lazy.select{|pr| circular?("1"*pr)}.take(5).to_a puts reps.map{|r| "R" + r.to_s}.join(", "), "" [5003, 9887, 15073, 25031].each {|rep| puts "R#{rep} circular_prime ? #{circular?("1"*rep)}" } </lang>
- Output:
First 19 circular primes: 2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933 First 5 prime repunits: R2, R19, R23, R317, R1031 R5003 circular_prime ? false R9887 circular_prime ? false R15073 circular_prime ? false R25031 circular_prime ? false
Rust
<lang rust>// [dependencies] // rug = "1.8"
fn is_prime(n: u32) -> bool {
if n < 2 {
return false;
}
if n % 2 == 0 {
return n == 2;
}
if n % 3 == 0 {
return n == 3;
}
let mut p = 5;
while p * p <= n {
if n % p == 0 {
return false;
}
p += 2;
if n % p == 0 {
return false;
}
p += 4;
}
true
}
fn cycle(n: u32) -> u32 {
let mut m: u32 = n;
let mut p: u32 = 1;
while m >= 10 {
p *= 10;
m /= 10;
}
m + 10 * (n % p)
}
fn is_circular_prime(p: u32) -> bool {
if !is_prime(p) {
return false;
}
let mut p2: u32 = cycle(p);
while p2 != p {
if p2 < p || !is_prime(p2) {
return false;
}
p2 = cycle(p2);
}
true
}
fn test_repunit(digits: usize) {
use rug::{integer::IsPrime, Integer};
let repunit = "1".repeat(digits);
let bignum = Integer::from_str_radix(&repunit, 10).unwrap();
if bignum.is_probably_prime(10) != IsPrime::No {
println!("R({}) is probably prime.", digits);
} else {
println!("R({}) is not prime.", digits);
}
}
fn main() {
use rug::{integer::IsPrime, Integer};
println!("First 19 circular primes:");
let mut count = 0;
let mut p: u32 = 2;
while count < 19 {
if is_circular_prime(p) {
if count > 0 {
print!(", ");
}
print!("{}", p);
count += 1;
}
p += 1;
}
println!();
println!("Next 4 circular primes:");
let mut repunit: u32 = 1;
let mut digits: usize = 1;
while repunit < p {
repunit = 10 * repunit + 1;
digits += 1;
}
let mut bignum = Integer::from(repunit);
count = 0;
while count < 4 {
if bignum.is_probably_prime(15) != IsPrime::No {
if count > 0 {
print!(", ");
}
print!("R({})", digits);
count += 1;
}
digits += 1;
bignum = bignum * 10 + 1;
}
println!();
test_repunit(5003);
test_repunit(9887);
test_repunit(15073);
test_repunit(25031);
test_repunit(35317);
test_repunit(49081);
}</lang>
- Output:
Execution time is about 805 seconds on my system (3.2 GHz Quad-Core Intel Core i5, macOS 10.15.4).
First 19 circular primes: 2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933 Next 4 circular primes: R(19), R(23), R(317), R(1031) R(5003) is not prime. R(9887) is not prime. R(15073) is not prime. R(25031) is not prime. R(35317) is not prime. R(49081) is probably prime.
Scala
<lang scala>object CircularPrimes {
def main(args: Array[String]): Unit = {
println("First 19 circular primes:")
var p = 2
var count = 0
while (count < 19) {
if (isCircularPrime(p)) {
if (count > 0) {
print(", ")
}
print(p)
count += 1
}
p += 1
}
println()
println("Next 4 circular primes:")
var repunit = 1
var digits = 1
while (repunit < p) {
repunit = 10 * repunit + 1
digits += 1
}
var bignum = BigInt.apply(repunit)
count = 0
while (count < 4) {
if (bignum.isProbablePrime(15)) {
if (count > 0) {
print(", ")
}
print(s"R($digits)")
count += 1
}
digits += 1
bignum = bignum * 10
bignum = bignum + 1
}
println()
testRepunit(5003) testRepunit(9887) testRepunit(15073) testRepunit(25031) }
def isPrime(n: Int): Boolean = {
if (n < 2) {
return false
}
if (n % 2 == 0) {
return n == 2
}
if (n % 3 == 0) {
return n == 3
}
var p = 5
while (p * p <= n) {
if (n % p == 0) {
return false
}
p += 2
if (n % p == 0) {
return false
}
p += 4
}
true
}
def cycle(n: Int): Int = {
var m = n
var p = 1
while (m >= 10) {
p *= 10
m /= 10
}
m + 10 * (n % p)
}
def isCircularPrime(p: Int): Boolean = {
if (!isPrime(p)) {
return false
}
var p2 = cycle(p)
while (p2 != p) {
if (p2 < p || !isPrime(p2)) {
return false
}
p2 = cycle(p2)
}
true
}
def testRepunit(digits: Int): Unit = {
val ru = repunit(digits)
if (ru.isProbablePrime(15)) {
println(s"R($digits) is probably prime.")
} else {
println(s"R($digits) is not prime.")
}
}
def repunit(digits: Int): BigInt = {
val ch = Array.fill(digits)('1')
BigInt.apply(new String(ch))
}
}</lang>
- Output:
First 19 circular primes: 2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933 Next 4 circular primes: R(19), R(23), R(317), R(1031) R(5003) is not prime. R(9887) is not prime. R(15073) is not prime. R(25031) is not prime.
Sidef
<lang ruby>func is_circular_prime(n) {
n.is_prime || return false
var circular = n.digits circular.min < circular.tail && return false
for k in (1 ..^ circular.len) {
with (circular.rotate(k).digits2num) {|p|
(p.is_prime && (p >= n)) || return false
}
}
return true
}
say "The first 19 circular primes are:" say 19.by(is_circular_prime)
say "\nThe next 4 circular primes, in repunit format, are:"
say "R(#{n})" } say "\nRepunit testing:" [5003, 9887, 15073, 25031, 35317, 49081].each {|n| var now = Time.micro say ("R(#{n}) -> ", is_prob_prime((10**n - 1)/9) ? 'probably prime' : 'composite', " (took: #{'%.3f' % Time.micro-now} sec)") }</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) -> composite (took: 0.024 sec) R(9887) -> composite (took: 0.006 sec) R(15073) -> composite (took: 0.389 sec) R(25031) -> composite (took: 54.452 sec) R(35317) -> composite (took: 0.875 sec)
Wren
Wren-cli
Second part is very slow - over 37 minutes to find all four. <lang ecmascript>import "/math" for Int import "/big" for BigInt
var circs = []
var isCircular = Fn.new { |n|
var nn = n
var pow = 1 // will eventually contain 10 ^ d where d is number of digits in n
while (nn > 0) {
pow = pow * 10
nn = (nn/10).floor
}
nn = n
while (true) {
nn = nn * 10
var f = (nn/pow).floor // first digit
nn = nn + f * (1 - pow)
if (circs.contains(nn)) return false
if (nn == n) break
if (!Int.isPrime(nn)) return false
}
return true
}
var repunit = Fn.new { |n| BigInt.new("1" * n) }
System.print("The first 19 circular primes are:") var digits = [1, 3, 7, 9] var q = [1, 2, 3, 5, 7, 9] // queue the numbers to be examined var fq = [1, 2, 3, 5, 7, 9] // also queue the corresponding first digits var count = 0 while (true) {
var f = q[0] // peek first element
var fd = fq[0] // peek first digit
if (Int.isPrime(f) && isCircular.call(f)) {
circs.add(f)
count = count + 1
if (count == 19) break
}
q.removeAt(0) // pop first element
fq.removeAt(0) // ditto for first digit queue
if (f != 2 && f != 5) { // if digits > 1 can't contain a 2 or 5
// add numbers with one more digit to queue
// only numbers whose last digit >= first digit need be added
for (d in digits) {
if (d >= fd) {
q.add(f*10+d)
fq.add(fd)
}
}
}
} System.print(circs)
System.print("\nThe next 4 circular primes, in repunit format, are:") count = 0 var rus = [] var primes = Int.primeSieve(10000) for (p in primes[3..-1]) {
if (repunit.call(p).isProbablePrime(1)) {
rus.add("R(%(p))")
count = count + 1
if (count == 4) break
}
} System.print(rus)</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)]
Embedded
A massive speed-up, of course, when GMP is plugged in for the 'probably prime' calculations. Around 11.5 minutes including the stretch goal. <lang ecmascript>/* circular_primes_embedded.wren */
import "./math" for Int import "./fmt" for Fmt
foreign class Mpz {
construct init() {}
foreign setString(str, base)
foreign probPrime(reps)
}
var circs = []
var isCircular = Fn.new { |n|
var nn = n
var pow = 1 // will eventually contain 10 ^ d where d is number of digits in n
while (nn > 0) {
pow = pow * 10
nn = (nn/10).floor
}
nn = n
while (true) {
nn = nn * 10
var f = (nn/pow).floor // first digit
nn = nn + f * (1 - pow)
if (circs.contains(nn)) return false
if (nn == n) break
if (!Int.isPrime(nn)) return false
}
return true
}
System.print("The first 19 circular primes are:") var digits = [1, 3, 7, 9] var q = [1, 2, 3, 5, 7, 9] // queue the numbers to be examined var fq = [1, 2, 3, 5, 7, 9] // also queue the corresponding first digits var count = 0 while (true) {
var f = q[0] // peek first element
var fd = fq[0] // peek first digit
if (Int.isPrime(f) && isCircular.call(f)) {
circs.add(f)
count = count + 1
if (count == 19) break
}
q.removeAt(0) // pop first element
fq.removeAt(0) // ditto for first digit queue
if (f != 2 && f != 5) { // if digits > 1 can't contain a 2 or 5
// add numbers with one more digit to queue
// only numbers whose last digit >= first digit need be added
for (d in digits) {
if (d >= fd) {
q.add(f*10+d)
fq.add(fd)
}
}
}
} System.print(circs)
System.print("\nThe next 4 circular primes, in repunit format, are:") count = 0 var rus = [] var primes = Int.primeSieve(10000) var repunit = Mpz.init() for (p in primes[3..-1]) {
repunit.setString("1" * p, 10)
if (repunit.probPrime(10) > 0) {
rus.add("R(%(p))")
count = count + 1
if (count == 4) break
}
} System.print(rus) System.print("\nThe following repunits are probably circular primes:") for (i in [5003, 9887, 15073, 25031, 35317, 49081]) {
repunit.setString("1" * i, 10)
Fmt.print("R($-5d) : $s", i, repunit.probPrime(15) > 0)
}</lang>
and the C program to embed the above script in:
<lang c>#include <stdio.h>
- include <stdlib.h>
- include <string.h>
- include <gmp.h>
- include "wren.h"
/* C <=> Wren interface functions */
void C_mpzAllocate(WrenVM* vm) {
mpz_t *pz = (mpz_t*)wrenSetSlotNewForeign(vm, 0, 0, sizeof(mpz_t)); mpz_init(*pz);
}
void C_setString(WrenVM* vm) {
mpz_t *pz = (mpz_t*)wrenGetSlotForeign(vm, 0); const char *str = wrenGetSlotString(vm, 1); int base = (int)wrenGetSlotDouble(vm, 2); mpz_set_str(*pz, str, base);
}
void C_probPrime(WrenVM* vm) {
mpz_t *pz = (mpz_t*)wrenGetSlotForeign(vm, 0); int reps = (int)wrenGetSlotDouble(vm, 1); int ret = mpz_probab_prime_p(*pz, reps); wrenSetSlotDouble(vm, 0, (double)ret);
}
WrenForeignClassMethods bindForeignClass(WrenVM* vm, const char* module, const char* className) {
WrenForeignClassMethods methods;
methods.allocate = NULL;
methods.finalize = NULL;
if (strcmp(module, "main") == 0) {
if (strcmp(className, "Mpz") == 0) {
methods.allocate = C_mpzAllocate;
}
}
return methods;
}
WrenForeignMethodFn bindForeignMethod(
WrenVM* vm,
const char* module,
const char* className,
bool isStatic,
const char* signature) {
if (strcmp(module, "main") == 0) {
if (strcmp(className, "Mpz") == 0) {
if (!isStatic && strcmp(signature, "setString(_,_)") == 0) return C_setString;
if (!isStatic && strcmp(signature, "probPrime(_)") == 0) return C_probPrime;
}
}
return NULL;
}
static void writeFn(WrenVM* vm, const char* text) {
printf("%s", text);
}
void errorFn(WrenVM* vm, WrenErrorType errorType, const char* module, const int line, const char* msg) {
switch (errorType) {
case WREN_ERROR_COMPILE:
printf("[%s line %d] [Error] %s\n", module, line, msg);
break;
case WREN_ERROR_STACK_TRACE:
printf("[%s line %d] in %s\n", module, line, msg);
break;
case WREN_ERROR_RUNTIME:
printf("[Runtime Error] %s\n", msg);
break;
}
}
char *readFile(const char *fileName) {
FILE *f = fopen(fileName, "r"); fseek(f, 0, SEEK_END); long fsize = ftell(f); rewind(f); char *script = malloc(fsize + 1); fread(script, 1, fsize, f); fclose(f); script[fsize] = 0; return script;
}
static void loadModuleComplete(WrenVM* vm, const char* module, WrenLoadModuleResult result) {
if( result.source) free((void*)result.source);
}
WrenLoadModuleResult loadModule(WrenVM* vm, const char* name) {
WrenLoadModuleResult result = {0};
if (strcmp(name, "random") != 0 && strcmp(name, "meta") != 0) {
result.onComplete = loadModuleComplete;
char fullName[strlen(name) + 6];
strcpy(fullName, name);
strcat(fullName, ".wren");
result.source = readFile(fullName);
}
return result;
}
int main(int argc, char **argv) {
WrenConfiguration config;
wrenInitConfiguration(&config);
config.writeFn = &writeFn;
config.errorFn = &errorFn;
config.bindForeignClassFn = &bindForeignClass;
config.bindForeignMethodFn = &bindForeignMethod;
config.loadModuleFn = &loadModule;
WrenVM* vm = wrenNewVM(&config);
const char* module = "main";
const char* fileName = "circular_primes_embedded.wren";
char *script = readFile(fileName);
WrenInterpretResult result = wrenInterpret(vm, module, script);
switch (result) {
case WREN_RESULT_COMPILE_ERROR:
printf("Compile Error!\n");
break;
case WREN_RESULT_RUNTIME_ERROR:
printf("Runtime Error!\n");
break;
case WREN_RESULT_SUCCESS:
break;
}
wrenFreeVM(vm);
free(script);
return 0;
}</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
Zig
As of now (Zig release 0.8.1), Zig has large integer support, but there is not yet a prime test in the standard library. Therefore, we only find the circular primes < 1e6. As with the Prolog version, we only check numbers composed of 1, 3, 7, or 9. <lang Zig> const std = @import("std"); const math = std.math; const heap = std.heap; const stdout = std.io.getStdOut().writer();
pub fn main() !void {
var arena = heap.ArenaAllocator.init(heap.page_allocator); defer arena.deinit();
var candidates = std.PriorityQueue(u32).init(&arena.allocator, u32cmp); defer candidates.deinit();
try stdout.print("The circular primes are:\n", .{});
try stdout.print("{:10}" ** 4, .{ 2, 3, 5, 7 });
var c: u32 = 4;
try candidates.add(0);
while (true) {
var n = candidates.remove();
if (n > 1_000_000)
break;
if (n > 10 and circular(n)) {
try stdout.print("{:10}", .{n});
c += 1;
if (c % 10 == 0)
try stdout.print("\n", .{});
}
try candidates.add(10 * n + 1);
try candidates.add(10 * n + 3);
try candidates.add(10 * n + 7);
try candidates.add(10 * n + 9);
}
try stdout.print("\n", .{});
}
fn u32cmp(a: u32, b: u32) math.Order {
return math.order(a, b);
}
fn circular(n0: u32) bool {
if (!isprime(n0))
return false
else {
var n = n0;
var d = @floatToInt(u32, @log10(@intToFloat(f32, n)));
return while (d > 0) : (d -= 1) {
n = rotate(n);
if (n < n0 or !isprime(n))
break false;
} else true;
}
}
fn rotate(n: u32) u32 {
if (n == 0)
return 0
else {
const d = @floatToInt(u32, @log10(@intToFloat(f32, n))); // digit count - 1
const m = math.pow(u32, 10, d);
return (n % m) * 10 + n / m;
}
}
fn isprime(n: u32) bool {
if (n < 2)
return false;
inline for ([3]u3{ 2, 3, 5 }) |p| {
if (n % p == 0)
return n == p;
}
const wheel235 = [_]u3{
6, 4, 2, 4, 2, 4, 6, 2,
};
var i: u32 = 1;
var f: u32 = 7;
return while (f * f <= n) {
if (n % f == 0)
break false;
f += wheel235[i];
i = (i + 1) & 0x07;
} else true;
} </lang>
- Output:
The circular primes are:
2 3 5 7 11 13 17 37 79 113
197 199 337 1193 3779 11939 19937 193939 199933