Circular primes
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
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- 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
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.- Output:
First 19 circular primes: 2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933
AppleScript
on isPrime(n)
if (n < 4) then return (n > 1)
if ((n mod 2 is 0) or (n mod 3 is 0)) then return false
repeat with i from 5 to (n ^ 0.5) div 1 by 6
if ((n mod i is 0) or (n mod (i + 2) is 0)) then return false
end repeat
return true
end isPrime
on isCircularPrime(n)
if (not isPrime(n)) then return false
set temp to n
set c to 0
repeat while (temp > 9)
set temp to temp div 10
set c to c + 1
end repeat
set p to (10 ^ c) as integer
set temp to n
repeat c times
set temp to temp mod p * 10 + temp div p
if ((temp < n) or (not isPrime(temp))) then return false
end repeat
return true
end isCircularPrime
-- Return the first c circular primes.
-- Takes 2 as read and checks only odd numbers thereafter.
on circularPrimes(c)
if (c < 1) then return {}
set output to {2}
set n to 3
set counter to 1
repeat until (counter = c)
if (isCircularPrime(n)) then
set end of output to n
set counter to counter + 1
end if
set n to n + 2
end repeat
return output
end circularPrimes
return circularPrimes(19)
- Output:
{2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933}
Arturo
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
]
- Output:
2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933
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)
}
- Output:
first 19 circular primes: 2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933
BASIC
Rosetta Code problem: https://rosettacode.org/wiki/Circular_primes
by Jjuanhdez, 02/2023
BASIC256
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
end while
end
function isPrime(v)
if v < 2 then return False
if v mod 2 = 0 then return v = 2
if v mod 3 = 0 then return v = 3
d = 5
while d * d <= v
if v mod d = 0 then return False else d += 2
end while
return True
end function
function isCircularPrime(p)
n = floor(log(p)/log(10))
m = 10^n
q = p
for i = 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
- Output:
Same as FreeBASIC entry.
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
- Output:
Primeros 19 primos circulares: 2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933
PureBasic
Macro floor(x)
Round(x, #PB_Round_Down)
EndMacro
Procedure isPrime(v.i)
If v <= 1 : ProcedureReturn #False
ElseIf v < 4 : ProcedureReturn #True
ElseIf v % 2 = 0 : ProcedureReturn #False
ElseIf v < 9 : ProcedureReturn #True
ElseIf v % 3 = 0 : ProcedureReturn #False
Else
Protected r = Round(Sqr(v), #PB_Round_Down)
Protected f = 5
While f <= r
If v % f = 0 Or v % (f + 2) = 0
ProcedureReturn #False
EndIf
f + 6
Wend
EndIf
ProcedureReturn #True
EndProcedure
Procedure isCircularPrime(p.i)
n.i = floor(Log(p)/Log(10))
m.i = Pow(10, n)
q.i = p
For i.i = 0 To n
If q < p Or Not isPrime(q)
ProcedureReturn #False
EndIf
q = (q % m) * 10 + floor(q / m)
Next i
ProcedureReturn #True
EndProcedure
OpenConsole()
p.i = 2
dp.i = 1
cont.i = 0
PrintN("Primeros 19 primos circulares:")
While cont < 19
If isCircularPrime(p)
Print(Str(p) + " ")
cont + 1
EndIf
p + dp
dp = 2
Wend
PrintN(#CRLF$ + "--- terminado, pulsa RETURN---"): Input()
CloseConsole()
- Output:
Same as FreeBASIC entry.
Yabasic
p = 2
dp = 1
cont = 0
print("Primeros 19 primos circulares:")
while cont < 19
if isCircularPrime(p) then
print p," ";
cont = cont + 1
fi
p = p + dp
dp = 2
wend
end
sub isPrime(v)
if v < 2 return False
if mod(v, 2) = 0 return v = 2
if mod(v, 3) = 0 return v = 3
d = 5
while d * d <= v
if mod(v, d) = 0 then return False else d = d + 2 : fi
wend
return True
end sub
sub isCircularPrime(p)
n = floor(log(p)/log(10))
m = 10^n
q = p
for i = 0 to n
if (q < p or not isPrime(q)) return false
q = (mod(q, m)) * 10 + floor(q / m)
next i
return true
end sub
- Output:
Same as FreeBASIC entry.
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;
}
- 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++
#include <cstdint>
#include <algorithm>
#include <iostream>
#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);
}
bool is_circular_prime(const integer& p) {
if (!is_prime(p))
return false;
std::string str = p.get_str();
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 = p.get_str();
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;
}
- 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
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;
}
- Output:
First 19 circular primes: 2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933
Delphi
Again note how the code is broken up into easy to distinguish subroutines that can be reused, as opposed to mashing everything together into a mass of code that is difficult to understand, debug or enhance.
procedure ShowCircularPrimes(Memo: TMemo);
{Show list of the first 19, cicular primes}
var I,Cnt: integer;
var S: string;
procedure RotateStr(var S: string);
{Rotate characters in string}
var I: integer;
var C: char;
begin
C:=S[Length(S)];
for I:=Length(S)-1 downto 1 do S[I+1]:=S[I];
S[1]:=C;
end;
function IsCircularPrime(N: integer): boolean;
{Test if all rotations of number are prime and}
{A rotation of the number hasn't been used before}
var I,P: integer;
var NS: string;
begin
Result:=False;
NS:=IntToStr(N);
for I:=1 to Length(NS)-1 do
begin
{Rotate string and convert to integer}
RotateStr(NS);
P:=StrToInt(NS);
{Exit if number is not prime or}
{Is prime, is less than N i.e. we've seen it before}
if not IsPrime(P) or (P<N) then exit;
end;
Result:=True;
end;
begin
S:='';
Cnt:=0;
{Look for circular primes and display 1st 19}
for I:=0 to High(I) do
if IsPrime(I) then
if IsCircularPrime(I) then
begin
Inc(Cnt);
S:=S+Format('%7D',[I]);
if Cnt>=19 then break;
If (Cnt mod 5)=0 then S:=S+CRLF;
end;
Memo.Lines.Add(S);
Memo.Lines.Add('Count = '+IntToStr(Cnt));
end;
- Output:
2 3 5 7 11
13 17 37 79 113
197 199 337 1193 3779
11939 19937 193939 199933
Count = 19
Elapsed Time: 21.234 ms.
EasyLang
fastfunc prime n .
if n mod 2 = 0 and n > 2
return 0
.
i = 3
sq = sqrt n
while i <= sq
if n mod i = 0
return 0
.
i += 2
.
return 1
.
func cycle n .
m = n
p = 1
while m >= 10
p *= 10
m = m div 10
.
return m + n mod p * 10
.
func circprime p .
if prime p = 0
return 0
.
p2 = cycle p
while p2 <> p
if p2 < p or prime p2 = 0
return 0
.
p2 = cycle p2
.
return 1
.
p = 2
while count < 19
if circprime p = 1
print p
count += 1
.
p += 1
.F#
This task uses rUnitP and Extensible Prime Generator (F#)
// 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 ""
printf "The first 5 repunit primes are "; rUnitP(10)|>Seq.take 5|>Seq.iter(fun n->printf $"R(%d{n}) "); printfn ""
- 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).
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
- 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.
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
- 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
FutureBasic
begin globals
_last1 = 7
CFTimeInterval t
byte n(_last1 +1) //array of digits
uint64 list(20), p10 = 1
short count, digits, first1 = _last1
n(_last1) = 2 //First number to check
end globals
local fn convertToNumber as uint64
short i = first1
uint64 nr = n(i)
while i < _last1
i++ : nr = nr * 10 + n(i)
wend
end fn = nr
void local fn increment( dgt as short )
select n(dgt)
case 1, 7, 5 : n(dgt) += 2
case 3 : if digits then n(dgt) = 7 else n(dgt) = 5
case 9 : n(dgt) = 1
if dgt == first1 then first1-- : digits++ : p10 *= 10
fn increment( dgt -1 )
case else : n(dgt)++
end select
end fn
local fn isPrime( v as uint64 ) as bool
if v < 4 then exit fn = yes
if v mod 3 == 0 then exit fn = no
uint64 f = 5
while f*f <= v
if v mod f == 0 then exit fn = no
f += 2
if v mod f == 0 then exit fn = no
f += 4
wend
end fn = yes
void local fn circularPrimes
uint64 num, nr
short r
while ( count < 19 )
num = fn convertToNumber
if fn isPrime( num )
nr = num
for r = 1 to digits
nr = 10 * ( nr mod p10 ) + ( nr / p10 ) //Rotate
if nr < num then break //Rotation of lower number
if fn isPrime( nr ) == no then break //Not prime
next
if r > digits then list( count ) = num : count++
end if
fn increment( _last1 )
wend
end fn
window 1, @"Circular Primes in FutureBasic", (0, 0, 280, 320)
t = fn CACurrentMediaTime
fn circularPrimes
t = fn CACurrentMediaTime - t
for count = 0 to 18
printf @"%22u",list(count)
next
printf @"\n Compute time: %.3f ms", t * 1000
handleevents- Output:
File:Screenshot of Circular Primes in FutureBasic.png
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))
}
}
- 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
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
- 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
R=: 10x #. #&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
)
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=: 10x #. #&1 (;q:@R)&> 500 |limit error | (;q:@R)&>500 ) boxdraw_j_ 0 Filter=: (#~`)(`:6) R=: 10x #. #&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
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));
}
}
- 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.
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);The first task
limit(19; circular_primes)The second task (viewed independently):
last(limit(19; circular_primes)) as $max
| limit(4; repunits | select(. > $max and is_prime))
| "R(\(tostring|length))"- 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
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
- 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
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)
}
- 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
-- 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, ", "))
- Output:
2, 3, 5, 7, 11, 13, 17, 37, 79, 113, 197, 199, 337, 1193, 3779, 11939, 19937, 193939, 199933
Mathematica / Wolfram Language
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}]
- 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
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)
- 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 14 digit numbers.RepUnit test with gmp is to boring.
Using a base 4 downcounter to create the numbers with more than one digit.
Changed the manner of the counter, so that it counts that there is no digit smaller than the first.
199-> 333 and not 311 so a base4 counter 1_(1,3,7,9) changed to base3 3_( 3,7,9 )->base2 7_(7,9) -> base 1 = 99....99.
The main speed up is reached by testing with small primes within CycleNum.
program CircularPrimes;
//nearly the way it is done:
//http://www.worldofnumbers.com/circular.htm
//base 4 counter to create numbers with first digit is the samallest used.
//check if numbers are tested before and reduce gmp-calls by checking with prime 3,7
{$IFDEF FPC}
{$MODE DELPHI}{$OPTIMIZATION ON,ALL}
uses
Sysutils,gmp;
{$ENDIF}
{$IFDEF Delphi}
uses
System.Sysutils,?gmp?;
{$ENDIF}
{$IFDEF WINDOWS}
{$APPTYPE CONSOLE}
{$ENDIF}
const
MAXCNTOFDIGITS = 14;
MAXDGTVAL = 3;
conv : array[0..MAXDGTVAL+1] of byte = (9,7,3,1,0);
type
tDigits = array[0..23] of byte;
tUint64 = NativeUint;
var
mpz : mpz_t;
digits,
revDigits : tDigits;
CheckNum : array[0..19] of tUint64;
Found : array[0..23] of tUint64;
Pot_ten,Count,CountNumCyc,CountNumPrmTst : tUint64;
procedure CheckOne(MaxIdx:integer);
var
Num : Uint64;
i : integer;
begin
i:= MaxIdx;
repeat
inc(CountNumPrmTst);
num := CheckNum[i];
mpz_set_ui(mpz,Num);
If mpz_probab_prime_p(mpz,3)=0then
EXIT;
dec(i);
until i < 0;
Found[Count] := CheckNum[0];
inc(count);
end;
function CycleNum(MaxIdx:integer):Boolean;
//first create circular numbers to minimize prime checks
var
cycNum,First,P10 : tUint64;
i,j,cv : integer;
Begin
i:= MaxIdx;
j := 0;
First := 0;
repeat
cv := conv[digits[i]];
dec(i);
First := First*10+cv;
revDigits[j]:= cv;
inc(j);
until i < 0;
// if num is divisible by 3 then cycle numbers also divisible by 3 same sum of digits
IF First MOD 3 = 0 then
EXIT(false);
If First mod 7 = 0 then
EXIT(false);
//if one of the cycled number must have been tested before break
P10 := Pot_ten;
i := 0;
j := 0;
CheckNum[j] := First;
cycNum := First;
repeat
inc(CountNumCyc);
cv := revDigits[i];
inc(j);
cycNum := (cycNum - cv*P10)*10+cv;
//num was checked before
if cycNum < First then
EXIT(false);
if cycNum mod 7 = 0 then
EXIT(false);
CheckNum[j] := cycNum;
inc(i);
until i >= MaxIdx;
EXIT(true);
end;
var
T0: Int64;
idx,MaxIDx,dgt,MinDgt : NativeInt;
begin
T0 := GetTickCount64;
mpz_init(mpz);
fillchar(digits,Sizeof(digits),chr(MAXDGTVAL));
Count :=0;
For maxIdx := 2 to 10 do
if maxidx in[2,3,5,7] then
begin
Found[Count]:= maxIdx;
inc(count);
end;
Pot_ten := 10;
maxIdx := 1;
idx := 0;
MinDgt := MAXDGTVAL;
repeat
if CycleNum(MaxIdx) then
CheckOne(MaxIdx);
idx := 0;
repeat
dgt := digits[idx]-1;
if dgt >=0 then
break;
digits[idx] := MinDgt;
inc(idx);
until idx >MAXCNTOFDIGITS-1;
if idx > MAXCNTOFDIGITS-1 then
BREAK;
if idx<=MaxIDX then
begin
digits[idx] := dgt;
//change all to leading digit
if idx=MaxIDX then
Begin
For MinDgt := 0 to idx do
digits[MinDgt]:= dgt;
minDgt := dgt;
end;
end
else
begin
minDgt := MAXDGTVAL;
For maxidx := 0 to idx do
digits[MaxIdx] := MAXDGTVAL;
Maxidx := idx;
Pot_ten := Pot_ten*10;
writeln(idx:7,count:7,CountNumCyc:16,CountNumPrmTst:12,GetTickCount64-T0:8);
end;
until false;
writeln(idx:7,count:7,CountNumCyc:16,CountNumPrmTst:12,GetTickCount64-T0:8);
T0 := GetTickCount64-T0;
For idx := 0 to count-2 do
write(Found[idx],',');
writeln(Found[count-1]);
writeln('It took ',T0,' ms ','to check ',MAXCNTOFDIGITS,' decimals');
mpz_clear(mpz);
{$IFDEF WINDOWS}
readln;
{$ENDIF}
end.
- @ Tio.Run:
Digits found cycles created primetests time in ms
2 9 7 10 0
3 13 43 38 0
4 15 203 89 0
5 17 879 213 0
6 19 4209 816 1
7 19 16595 1794 2
8 19 67082 4666 6
9 19 270760 13108 19
10 19 1094177 39059 63
11 19 4421415 118787 220
12 19 23728859 1078484 1505
13 19 77952009 1869814 3562
14 19 296360934 4405393 11320
#### 14 19 754020918 28736408 26129 ##### before
2,3,5,7,11,13,17,37,79,113,197,199,337,1193,3779,11939,19937,193939,199933
It took 11320 ms to check 14 decimals
only creating numbers
14 4 0 0 184
creating and cycling numbers testing not ( MOD 3=0 OR MoD 7 = 0)
14 4 247209700 0 8842
that reduces the count of gmp-prime tests to 1/6
PascalABC.NET
##
function IsPrime(n: integer): boolean;
begin
result := false;
for var i := 2 to n.Sqrt.Floor do
if n mod i = 0 then
exit;
result := true;
end;
function rotations(n: integer): sequence of integer;
begin
var a := n.tostring;
for var i := 1 to a.Length do
yield (a[i:] + a[:i]).ToInteger;
end;
function isCircular(n: integer): boolean;
begin
var rot := rotations(n);
if rot.first = rot.min then
result := rot.All(x -> isPrime(x));
end;
(2..1000_000).Where(x -> isCircular(x)).Take(19).Println;
- Output:
2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933
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');
}
- 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
PicoLisp
For small primes, I only check numbers that are made up of the digits 1, 3, 7, and 9, and I also use a small prime checker to avoid the overhead of the Miller-Rabin Primality Test. For the large repunits, one can use the code from the Miller Rabin task.
(load "plcommon/primality.l") # see task: "Miller-Rabin Primality Test"
(de candidates (Limit)
(let Q (0)
(nth
(sort
(make
(while Q
(let A (pop 'Q)
(when (< A Limit)
(link A)
(setq Q
(cons
(+ (* 10 A) 1)
(cons
(+ (* 10 A) 3)
(cons
(+ (* 10 A) 7)
(cons (+ (* 10 A) 9) Q))))))))))
6)))
(de circular? (P0)
(and
(small-prime? P0)
(fully '((P) (and (>= P P0) (small-prime? P))) (rotations P0))))
(de rotate (L)
(let ((X . Xs) L)
(append Xs (list X))))
(de rotations (N)
(let L (chop N)
(mapcar
format
(make
(do (dec (length L))
(link (setq L (rotate L))))))))
(de small-prime? (N) # For small prime candidates only
(if (< N 2)
NIL
(let W (1 2 2 . (4 2 4 2 4 6 2 6 .))
(for (D 2 T (+ D (pop 'W)))
(T (> (* D D) N) T)
(T (=0 (% N D)) NIL)))))
(de repunit-primes (N)
(let (Test 111111 Remaining N K 6)
(make
(until (=0 Remaining)
(setq Test (inc (* 10 Test)))
(inc 'K)
(when (prime? Test)
(link K)
(dec 'Remaining))))))
(setq Circular
(conc
(2 3 5 7)
(filter circular? (candidates 1000000))
(mapcar '((X) (list 'R X)) (repunit-primes 4))))
(prinl "The first few circular primes:")
(println Circular)
(bye)- Output:
The first few 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))
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.
?- 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.
- 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
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.
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.
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'}"
}
- 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
Version 1
/*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
- 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
Version 2
Libraries: How to use
Library: Numbers
Library: Functions
Library: Settings
Library: Abend
Primes(n) retrieves all primes below 200000 using a sieve. Prime(n) is an optimized Miller-Rabin primality test. Part 2 of the task could be completed up to the 3rd Repunit. Since primality checking becomes very slow for numbers with over 1000 digits, part 3 of the task was out of reach. However, I verified R(1031) that was left over from part 2. Remember, all this is written in REXX, without relying on external libraries.
include Settings
say version; say 'Circular primes'; say
call First19
call Next3
call HigherReps
exit
First19:
procedure expose glob. prim.
call Time('r')
numeric digits 10
say 'First 19 circular primes:'
p = Primes(200000)
do i = 1 to p
a = prim.prime.i
if Verify(a,'024568','m') > 0 then
iterate i
b = a; l = Length(b)
do l-1
b = Right(b,l-1)||Left(b,1)
if b < a then
iterate i
if \ IsPrime(b) then
iterate i
end
call Charout ,a' '
end
say
say Format(Time('e'),,3) 'seconds'
say
return
Next3:
procedure expose glob.
call Time('r')
numeric digits 320
say 'Next 3 circular primes:'
do i = 7 to 320
r = Repunit(i)
if IsPrime(r) then
call charout ,'R('i') '
end
say
say Format(Time('e'),,3) 'seconds'
say
return
HigherReps:
procedure expose glob.
call Time('r')
numeric digits 1040
say 'Primality of R(1031):'
if IsPrime(Repunit(1031)) then
say 'R(1031) is probable prime'
else
say 'R(1031) is composite'
say Format(Time('e'),,3) 'seconds'
say
return
Repunit:
/* Repeat 1's function */
procedure expose glob.
arg x
/* Formula */
return Copies('1',x)
include Numbers
include Functions
include Abend
- Output:
REXX-ooRexx_5.0.0(MT)_64-bit 6.05 23 Dec 2022 Circular primes First 19 circular primes: 3 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933 0.126 seconds Next 3 circular primes: R(19) R(23) R(317) 62.924 seconds Primality of R(n): R(1031) is probable prime 110.665 seconds
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- 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...
RPL
To speed up execution, we use a generator of candidate numbers made only of 1, 3, 7 and 9.
≪ 1 SF DUP
DO
→STR TAIL LASTARG HEAD + STR→
CASE
DUP2 > OVER 3 MOD NOT OR THEN 1 CF END
DUP ISPRIME? NOT THEN 1 CF END
END
UNTIL DUP2 == 1 FC? OR END
DROP2 1 FS?
≫ 'CIRC?' STO
≪ →STR "9731" → prev digits
≪ 1 SF "" @ flag 1 set = carry
prev SIZE 1 FOR j
digits DUP
prev j DUP SUB POS 1 FS? -
IF DUP NOT THEN DROP digits SIZE ELSE 1 CF END
DUP SUB SWAP +
-1 STEP
IF 1 FS? THEN digits DUP SIZE DUP SUB SWAP + END
STR→
≫ ≫ 'NEXT1379' STO
≪ 2 CF { 2 3 5 7 } 7
DO
NEXT1379
IF DUP CIRC? THEN
SWAP OVER + SWAP
IF OVER SIZE 19 ≥ THEN 2 SF END
END
UNTIL 2 FS? END DROP
≫ 'TASK' STO
- Output:
1: { 2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933 }
Runs in 13 minutes 10 seconds on a HP-50g.
Ruby
It takes more then 25 minutes to verify that R49081 is probably prime - omitted here.
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)}" }
- 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
// [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);
}
- 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
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))
}
}
- 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
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:"
{|n| (10**n - 1)/9 -> is_prob_prime }.first(4, 4..Inf).each {|n|
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)")
}
- 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)
Uiua
The built-in factorisation °/× is fast enough that it's not worth doing anything clever. Just filter a large enough range of candidates for primality, and for each generate its rotations, check it's the smallest of these and if so test whether all are prime. Filter by outcome of this test.
# Find first 19 circular primes
IsP ← =1⧻°/×
Rots ← ≡(⋕↻)+1⊙¤⇡⧻.°⋕
Circ ← ⨬(⋅0|/↧≡IsP)⊸(=⊙/↧)⟜Rots
▽⊸≡Circ▽⊸≡IsP+1⇡1e6- Output:
[2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933]
Wren
Wren-cli
Second part is very slow - over 37 minutes to find all four.
import "./math" for Int
import "./big" for BigInt
import "./str" for Str
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(Str.repeat("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)
- 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. 11 minutes 19 seconds including the stretch goal.
/* Circular_primes_embedded.wren */
import "./gmp" for Mpz
import "./math" for Int
import "./fmt" for Fmt
import "./str" for Str
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.new()
for (p in primes[3..-1]) {
repunit.setStr(Str.repeat("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.setStr(Str.repeat("1", i), 10)
Fmt.print("R($-5d) : $s", i, repunit.probPrime(15) > 0)
}
- 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
XPL0
func IsPrime(N); \Return 'true' if N > 2 is a prime number
int N, I;
[if (N&1) = 0 \even number\ then return false;
for I:= 3 to sqrt(N) do
[if rem(N/I) = 0 then return false;
I:= I+1;
];
return true;
];
func CircPrime(N0); \Return 'true' if N0 is a circular prime
int N0, N, Digits, Rotation, I, R;
[N:= N0;
Digits:= 0; \count number of digits in N
repeat Digits:= Digits+1;
N:= N/10;
until N = 0;
N:= N0;
for Rotation:= 0 to Digits-1 do
[if not IsPrime(N) then return false;
N:= N/10; \rotate least sig digit into high end
R:= rem(0);
for I:= 0 to Digits-2 do
R:= R*10;
N:= N+R;
if N0 > N then \reject N0 if it has a smaller prime rotation
return false;
];
return true;
];
int Counter, N;
[IntOut(0, 2); ChOut(0, ^ ); \show first circular prime
Counter:= 1;
N:= 3; \remaining primes are odd
loop [if CircPrime(N) then
[IntOut(0, N); ChOut(0, ^ );
Counter:= Counter+1;
if Counter >= 19 then quit;
];
N:= N+2;
];
]- Output:
2 3 5 7 11 13 17 37 79 113 197 199 337 1193 3779 11939 19937 193939 199933
Zig
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.
const std = @import("std");
const math = std.math;
const heap = std.heap;
pub fn main() !void {
var arena = heap.ArenaAllocator.init(heap.page_allocator);
defer arena.deinit();
const allocator = arena.allocator();
var candidates = std.PriorityQueue(u32, void, u32cmp).init(allocator, {});
defer candidates.deinit();
const stdout = std.io.getStdOut().writer();
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(_: void, 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: u32 = @intFromFloat(@log10(@as(f32, @floatFromInt(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: u32 = @intFromFloat(@log10(@as(f32, @floatFromInt(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;
}
- Output:
The circular primes are:
2 3 5 7 11 13 17 37 79 113
197 199 337 1193 3779 11939 19937 193939 199933
Categories:
- Programming Tasks
- Prime Numbers
- ALGOL 68
- ALGOL W
- AppleScript
- Arturo
- AWK
- BASIC
- BASIC256
- FreeBASIC
- PureBasic
- Yabasic
- C
- GMP
- C++
- D
- Delphi
- SysUtils,StdCtrls
- EasyLang
- F Sharp
- Factor
- Forth
- FutureBasic
- Pages with broken file links
- Go
- GMP(Go wrapper)
- Haskell
- J
- Java
- Jq
- Julia
- Kotlin
- Lua
- Mathematica
- Wolfram Language
- Nim
- Bignum
- Pascal
- Free Pascal
- PascalABC.NET
- Perl
- Ntheory
- Phix
- PicoLisp
- Prolog
- Python
- Raku
- REXX
- Ring
- RPL
- Ruby
- Rust
- Scala
- Sidef
- Uiua
- Wren
- Wren-math
- Wren-big
- Wren-fmt
- Wren-str
- Wren-gmp
- XPL0
- Zig