Cubic special primes
- Task
n is smallest prime such that the difference of successive terms are the smallest cubics of positive integers,
where n < 15000.
11l
F is_prime(a)
I a == 2
R 1B
I a < 2 | a % 2 == 0
R 0B
L(i) (3 .. Int(sqrt(a))).step(2)
I a % i == 0
R 0B
R 1B
V p = 2
V n = 1
print(2, end' ‘ ’)
L p + n ^ 3 < 15000
I is_prime(p + n ^ 3)
p += n ^ 3
n = 1
print(p, end' ‘ ’)
E
n++
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
Action!
INCLUDE "H6:SIEVE.ACT"
PROC Main()
DEFINE MAX="14999"
BYTE ARRAY primes(MAX+1)
INT i=[2],next,count=[1],n=[1]
Put(125) PutE() ;clear the screen
Sieve(primes,MAX+1)
PrintI(i)
DO
next=i+n*n*n
IF next>=MAX THEN
EXIT
FI
IF primes(next) THEN
n=1 i=next count==+1
Put(32) PrintI(i)
ELSE
n==+1
FI
OD
PrintF("%E%EThere are %I cubic special primes",count)
RETURN
- Output:
Screenshot from Atari 8-bit computer
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419 There are 23 cubic special primes
ALGOL 68
BEGIN # find a sequence of primes where the members differ by a cube #
INT max prime = 15 000;
# sieve the primes to max prime #
PR read "primes.incl.a68" PR
[]BOOL prime = PRIMESIEVE max prime;
# construct a table of cubes, we will need at most the cube root of max prime #
[ 1 : ENTIER exp( ln( max prime ) / 3 ) ]INT cube;
FOR i TO UPB cube DO cube[ i ] := i * i * i OD;
# find the prime sequence #
print( ( "2" ) );
INT p := 2;
WHILE p < max prime DO
# find a prime that is p + a cube #
INT q := 0;
FOR i WHILE q := p + cube[ i ];
IF q > max prime THEN FALSE ELSE NOT prime[ q ] FI
DO SKIP OD;
IF q <= max prime THEN print( ( " ", whole( q, 0 ) ) ) FI;
p := q
OD;
print( ( newline ) )
END
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
ALGOL W
begin % find a sequence of primes where the members differ by a cube %
integer MAX_PRIME, MAX_CUBE;
MAX_PRIME := 15000;
MAX_CUBE := entier( exp( ln( MAX_PRIME ) / 3 ) );
begin
logical array prime ( 1 :: MAX_PRIME );
integer array cube ( 1 :: MAX_CUBE );
integer p, q, c;
% sieve the primes to MAX_PRIME %
prime( 1 ) := false; prime( 2 ) := true;
for i := 3 step 2 until MAX_PRIME do prime( i ) := true;
for i := 4 step 2 until MAX_PRIME do prime( i ) := false;
for i := 3 step 2 until entier( sqrt( MAX_PRIME ) ) do begin
integer ii; ii := i + i;
if prime( i ) then for s := i * i step ii until MAX_PRIME do prime( s ) := false
end for_i;
% construct a table of cubes, we will need at most the cube root of max prime %
for i := 1 until MAX_CUBE do cube( i ) := i * i * i;
% find the prime sequence %
writeon( "2" );
p := 2;
while p < MAX_PRIME do begin
% find a prime that is p + a cube %
c := 1;
q := p + cube( c );
while q <= MAX_PRIME and not prime( q ) do begin
c := c + 1;
q := p + cube( c )
end while_not_prime_q ;
if q <= MAX_PRIME then writeon( i_w := 1, s_w := 0, " ", q );
p := q
end while_p_lt_MAX_PRIME
end;
write()
end.
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
Arturo
[p n]: [2 1]
prints -> 2
until -> 15000 =< p + n^3 [
if? prime? p + n^3 [
'p + n^3
n: 1
prints -> p
]
else -> 'n+1
]
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
AWK
# syntax: GAWK -f CUBIC_SPECIAL_PRIMES.AWK
# converted from FreeBASIC
BEGIN {
start = p = 2
stop = 15000
n = 1
printf("%5d ",p)
count = 1
do {
if (is_prime(p + n^3)) {
p += n^3
n = 1
printf("%5d%1s",p,++count%10?"":"\n")
}
else {
n++
}
} while (p + n^3 <= stop)
printf("\nCubic special primes %d-%d: %d\n",start,stop,count)
exit(0)
}
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:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419 Cubic special primes 2-15000: 23
C
#include <stdio.h>
#include <stdlib.h>
#include <stdbool.h>
#include <math.h>
#include <locale.h>
bool *sieve(int limit) {
int i, p;
limit++;
// True denotes composite, false denotes prime.
bool *c = calloc(limit, sizeof(bool)); // all false by default
c[0] = true;
c[1] = true;
for (i = 4; i < limit; i += 2) c[i] = true;
p = 3; // Start from 3.
while (true) {
int p2 = p * p;
if (p2 >= limit) break;
for (i = p2; i < limit; i += 2 * p) c[i] = true;
while (true) {
p += 2;
if (!c[p]) break;
}
}
return c;
}
bool isCube(int n) {
int s = (int)cbrt((double)n);
return s*s*s == n;
}
int main() {
const int limit = 14999;
int i, j, p, pc = 0, gap = 1, count = 1, lastCubicSpecial = 3;
const char *fmt = "%'7d %'7d %'6d %'4d\n";
bool *c = sieve(limit);
for (i = 0; i < limit; ++i) {
if (!c[i]) ++pc;
}
int *primes = (int *)malloc(pc * sizeof(int));
for (i = 0, j = 0; i < limit; ++i) {
if (!c[i]) primes[j++] = i;
}
free(c);
printf("Cubic special primes under 15,000:\n");
printf(" Prime1 Prime2 Gap Cbrt\n");
setlocale(LC_NUMERIC, "");
printf(fmt, 2, 3, 1, 1);
for (i = 2; i < pc; ++i) {
p = primes[i];
gap = p - lastCubicSpecial;
if (isCube(gap)) {
printf(fmt, lastCubicSpecial, p, gap, (int)cbrt((double)gap));
lastCubicSpecial = p;
++count;
}
}
printf("\n%d such primes found.\n", count+1);
free(primes);
return 0;
}
- Output:
Same as Wren example.
C++
#include <cstdint>
#include <iostream>
bool is_prime(uint32_t number) {
if ( number % 2 == 0 ) {
return number == 2;
}
int k = 3;
while ( k * k <= number ) {
if ( number % k == 0 ) {
return false;
}
k += 2;
}
return true;
}
int main() {
uint32_t p = 2, n = 1;
std::cout << 2 << " ";
while ( p + n * n * n < 15'000 ) {
if ( is_prime(p + n * n * n) ) {
p += n * n * n;
n = 1;
std::cout << p << " ";
} else {
n += 1;
}
}
}
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
Delphi
This program makes extensive use of the Delphi Prime-Generator Object
procedure CubicSpecialPrimes(Memo: TMemo);
const Limit = 15000;
var Sieve: TPrimeSieve;
var N,I,PP, Count: integer;
var S: string;
begin
Sieve:=TPrimeSieve.Create;
try
{Get all primes up to limit}
Sieve.Intialize(Limit);
N:=1;
Count:=0;
I:=1;
S:='';
while true do
begin
{Find first cube increment that is prime}
PP:=N +I*I*I;
if PP>Limit then break;
if Sieve.Flags[PP] then
begin
Inc(Count);
S:=S+Format('%6D',[PP]);
if (Count mod 5) = 0 then S:=S+CRLF;
{Step to next cube position}
I:=1;
N:=PP;
end
else Inc(I);
end;
Memo.Lines.Add(Format('There are %d cubic special primes',[count]));
Memo.Lines.Add(S);
finally Sieve.Free; end;
end;
- Output:
There are 23 cubic special primes 2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419 Elapsed Time: 2.178 ms.
EasyLang
fastfunc isprim num .
i = 2
while i <= sqrt num
if num mod i = 0
return 0
.
i += 1
.
return 1
.
p = 2
n = 1
write 2 & " "
repeat
h = p + n * n * n
until h >= 15000
if isprim h = 1
p = h
n = 1
write p & " "
else
n += 1
.
.
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
F#
// Cubic Special Primes: Nigel Galloway. March 30th., 2021
let fN=let N=[for n in [0..25]->n*n*n] in let mutable n=2 in (fun g->match List.contains(g-n)N with true->n<-g; true |_->false)
primes32()|>Seq.takeWhile((>)16000)|>Seq.filter fN|>Seq.iter(printf "%d "); printfn ""
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
Factor
USING: fry io kernel lists lists.lazy math math.functions
math.primes prettyprint ;
2 [ 1 lfrom swap '[ 3 ^ _ + ] lmap-lazy [ prime? ] lfilter car ]
lfrom-by [ 15000 < ] lwhile [ pprint bl ] leach nl
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
FreeBASIC
#include once "isprime.bas"
dim as integer p = 2, n = 1
print 2;" ";
do
if isprime(p + n^3) then
p += n^3
n = 1
print p;" ";
else
n += 1
end if
loop until p + n^3 >= 15000
print
end
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
FutureBasic
local fn IsPrime( n as NSUInteger ) as BOOL
—
BOOL isPrime = YES
NSUInteger i
if n < 2 then exit fn = NO
if n = 2 then exit fn = YES
if n mod 2 == 0 then exit fn = NO
for i = 3 to int(n^.5) step 2
if n mod i == 0 then exit fn = NO
next
end fn = isPrime
int p, n, counter
p = 2
n = 1
counter = 2
printf @"%5d \b", 2
do
if fn IsPrime ( p + n^3 )
p += n^3
n = 1
printf @"%5d \b", p
if counter mod 8 == 0 then print
counter++
else
n++
end if
until ( p + n^3 >= 15000 )
HandleEvents
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
Go
package main
import (
"fmt"
"math"
)
func sieve(limit int) []bool {
limit++
// True denotes composite, false denotes prime.
c := make([]bool, limit) // all false by default
c[0] = true
c[1] = true
// no need to bother with even numbers over 2 for this task
p := 3 // Start from 3.
for {
p2 := p * p
if p2 >= limit {
break
}
for i := p2; i < limit; i += 2 * p {
c[i] = true
}
for {
p += 2
if !c[p] {
break
}
}
}
return c
}
func isCube(n int) bool {
s := int(math.Cbrt(float64(n)))
return s*s*s == n
}
func commas(n int) string {
s := fmt.Sprintf("%d", n)
if n < 0 {
s = s[1:]
}
le := len(s)
for i := le - 3; i >= 1; i -= 3 {
s = s[0:i] + "," + s[i:]
}
if n >= 0 {
return s
}
return "-" + s
}
func main() {
c := sieve(14999)
fmt.Println("Cubic special primes under 15,000:")
fmt.Println(" Prime1 Prime2 Gap Cbrt")
lastCubicSpecial := 3
gap := 1
count := 1
fmt.Printf("%7d %7d %6d %4d\n", 2, 3, 1, 1)
for i := 5; i < 15000; i += 2 {
if c[i] {
continue
}
gap = i - lastCubicSpecial
if isCube(gap) {
cbrt := int(math.Cbrt(float64(gap)))
fmt.Printf("%7s %7s %6s %4d\n", commas(lastCubicSpecial), commas(i), commas(gap), cbrt)
lastCubicSpecial = i
count++
}
}
fmt.Println("\n", count+1, "such primes found.")
}
- Output:
Same as Wren example.
J
cubes=. 3 ^~ 1 , 10 + i: 8j8
nextp=. ({.@#~ 1&p:)@(+&cubes)
nextp^:(14e3&>)^:a: 2
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
Java
public final class CubicSpecialPrimes {
public static void main(String[] args) {
int p = 2;
int n = 1;
System.out.print(2 + " ");
while ( p + n * n * n < 15_000 ) {
if ( isPrime(p + n * n * n) ) {
p += n * n * n;
n = 1;
System.out.print(p + " ");
} else {
n += 1;
}
}
}
private static boolean isPrime(int number) {
if ( number % 2 == 0 ) {
return number == 2;
}
int k = 3;
while ( k * k <= number ) {
if ( number % k == 0 ) {
return false;
}
k += 2;
}
return true;
}
}
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
jq
Adapted from Julia
Works with gojq, the Go implementation of jq
For the definition of `is_prime` used here, see https://rosettacode.org/wiki/Additive_primes
# Output: an array
def cubic_special_primes:
. as $N
| sqrt as $maxidx
| {cprimes: [2], q: 0}
| until( .cprimes[-1] >= $N or .q >= $N;
label $out
| foreach range(1; $maxidx + 1) as $i (.;
.q = (.cprimes[-1] + ($i * $i * $i))
| if .q >= $N
then .emit = true
elif .q | is_prime
then .cprimes = .cprimes + [.q]
| .emit = true
else .
end;
select(.emit)) | {cprimes, q}, break $out )
| .cprimes ;
15000
| ("Cubic special primes < \(.):",
cubic_special_primes)
- Output:
[2,3,11,19,83,1811,2027,2243,2251,2467,2531,2539,3539,3547,4547,5059,10891,12619,13619,13627,13691,13907,14419]
Julia
using Primes
function cubicspecialprimes(N = 15000)
pmask = primesmask(1, N)
cprimes, maxidx = [2], isqrt(N)
while (n = cprimes[end]) < N
for i in 1:maxidx
q = n + i * i * i
if q > N
return cprimes
elseif pmask[q] # got next qprime
push!(cprimes, q)
break
end
end
end
end
println("Cubic special primes < 15000:")
foreach(p -> print(rpad(p[2], 6), p[1] % 10 == 0 ? "\n" : ""), enumerate(cubicspecialprimes()))
- Output:
Cubic special primes < 16000: 2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
Lua
do -- find a sequence of primes where the members differ by a cube
local maxPrime = 15000
-- sieve the primes to maxPrime
local isPrime = {}
for i = 1, maxPrime do isPrime[ i ] = i % 2 ~= 0 end
isPrime[ 1 ] = false
isPrime[ 2 ] = true
for s = 3, math.floor( math.sqrt( maxPrime ) ), 2 do
if isPrime[ s ] then
for p = s * s, maxPrime, s do isPrime[ p ] = false end
end
end
-- construct a table of cubes, we will need at most the cube root of maxPrime
local cube = {}
for i = 1, math.floor( ( maxPrime ^ ( 1 / 3 ) ) ) do cube[ i ] = i * i * i end
-- find the prime sequence
io.write( "2" )
local p = 2
while p < maxPrime do
-- find a prime that is p + a cube
local q, i = 0, 1
repeat
q, i = p + cube[ i ], i + 1
until q > maxPrime or isPrime[ q ]
if q <= maxPrime then io.write( " ", q ) end
p = q
end
io.write( "\n" )
end
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
Mathematica /Wolfram Language
start = {2};
Do[
If[PrimeQ[i],
last = Last[start];
If[IntegerQ[Surd[i - last, 3]],
AppendTo[start, i]
]
]
,
{i, 3, 15000}
]
start
- Output:
{2,3,11,19,83,1811,2027,2243,2251,2467,2531,2539,3539,3547,4547,5059,10891,12619,13619,13627,13691,13907,14419}
Nim
import math, strformat
func sieve(limit: Positive): seq[bool] =
# True denotes composite, false denotes prime.
result = newSeq[bool](limit + 1) # All false by default.
result[0] = true
result[1] = true
# No need to bother with even numbers over 2 for this task.
var p = 3
while true:
let p2 = p * p
if p2 > limit: break
for i in countup(p2, limit, 2 * p):
result[i] = true
while true:
inc p, 2
if not result[p]: break
func isCube(n: int): bool =
let s = cbrt(n.toFloat).int
result = s * s * s == n
let c = sieve(14999)
echo "Cubic special primes under 15_000:"
echo " Prime1 Prime2 Gap Cbrt"
var lastCubicSpecial = 3
var count = 1
echo &"{2:7} {3:7} {1:6} {1:4}"
for n in countup(5, 14999, 2):
if c[n]: continue
let gap = n - lastCubicSpecial
if gap.isCube:
let gapCbrt = cbrt(gap.toFloat).int
echo &"{lastCubicSpecial:7} {n:7} {gap:6} {gapCbrt:4}"
lastCubicSpecial = n
inc count
echo &"\n{count + 1} such primes found."
- Output:
Cubic special primes under 15_000: Prime1 Prime2 Gap Cbrt 2 3 1 1 3 11 8 2 11 19 8 2 19 83 64 4 83 1811 1728 12 1811 2027 216 6 2027 2243 216 6 2243 2251 8 2 2251 2467 216 6 2467 2531 64 4 2531 2539 8 2 2539 3539 1000 10 3539 3547 8 2 3547 4547 1000 10 4547 5059 512 8 5059 10891 5832 18 10891 12619 1728 12 12619 13619 1000 10 13619 13627 8 2 13627 13691 64 4 13691 13907 216 6 13907 14419 512 8 23 such primes found.
Perl
use strict;
use warnings;
use feature 'say';
use ntheory 'is_prime';
my @sp = my $previous = 2;
do {
my($next,$n);
while () { last if is_prime( $next = $previous + ++$n**3 ) }
push @sp, $next;
$previous = $next;
} until $sp[-1] >= 15000;
pop @sp and say join ' ', @sp;
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
Phix
See Quadrat_Special_Primes#Phix
Python
#!/usr/bin/python
def isPrime(n):
for i in range(2, int(n**0.5) + 1):
if n % i == 0:
return False
return True
if __name__ == '__main__':
p = 2
n = 1
print("2",end = " ")
while True:
if isPrime(p + n**3):
p += n**3
n = 1
print(p,end = " ")
else:
n += 1
if p + n**3 >= 15000:
break
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
Raku
A two character difference from the Quadrat Special Primes entry. (And it could have been one.)
my @sqp = 2, -> $previous {
my $next;
for (1..∞).map: *³ {
$next = $previous + $_;
last if $next.is-prime;
}
$next
} … *;
say "{+$_} matching numbers:\n", $_».fmt('%5d').batch(7).join: "\n" given
@sqp[^(@sqp.first: * > 15000, :k)];
- Output:
23 matching numbers: 2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
REXX
Version 1 Inline code
/*REXX program finds the smallest primes such that the difference of successive terms */
/*───────────────────────────────────────────────────── are the smallest cubic numbers. */
parse arg hi cols . /*obtain optional argument from the CL.*/
if hi=='' | hi=="," then hi= 15000 /* " " " " " " */
if cols=='' | cols=="," then cols= 10 /* " " " " " " */
call genP /*build array of semaphores for primes.*/
w= 10 /*width of a number in any column. */
title= 'the smallest primes < ' commas(hi) " such that the" ,
'difference of successive terma are the smallest cubic numbers'
if cols>0 then say ' index │'center(title, 1 + cols*(w+1) ) /*display the title.*/
if cols>0 then say '───────┼'center("" , 1 + cols*(w+1), '─') /* " " sep.*/
found= 0; idx= 1 /*initialize number of cbp and index.*/
op= 1
$= /*a list of cubic primes (so far). */
do j=0 by 0
do k=1 until !.np
np= op + k**3 /*find the next cube plus the oldPrime.*/
if np>= hi then leave j /*Is newPrime too big? Then quit. */
end /*k*/
found= found + 1 /*bump number of primes of this type. */
op= np /*assign the newPrime to the oldPrime*/
if cols<0 then iterate /*Build the list (to be shown later)? */
c= commas(np) /*maybe add commas to the number. */
$= $ right(c, max(w, length(c) ) ) /*add a cubic prime──►list, allow big#.*/
if found//cols\==0 then iterate /*have we populated a line of output? */
say center(idx, 7)'│' substr($, 2); $= /*display what we have so far (cols). */
idx= idx + cols /*bump the index count for the output*/
end /*j*/
if $\=='' then say center(idx, 7)"│" substr($, 2) /*possible display residual output.*/
if cols>0 then say '───────┴'center("" , 1 + cols*(w+1), '─') /*display foot sep. */
say
say 'Found ' commas(found) " of " title
exit 0 /*stick a fork in it, we're all done. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
commas: parse arg ?; do jc=length(?)-3 to 1 by -3; ?=insert(',', ?, jc); end; return ?
/*──────────────────────────────────────────────────────────────────────────────────────*/
genP: !.= 0 /*placeholders for primes (semaphores).*/
@.1=2; @.2=3; @.3=5; @.4=7; @.5=11 /*define some low primes. */
!.2=1; !.3=1; !.5=1; !.7=1; !.11=1 /* " " " " flags. */
#=5; s.#= @.# **2 /*number of primes so far; prime². */
/* [↓] generate more primes ≤ high.*/
do j=@.#+2 by 2 to hi /*find odd primes from here on. */
parse var j '' -1 _; if _==5 then iterate /*J divisible by 5? (right dig)*/
if j// 3==0 then iterate /*" " " 3? */
if j// 7==0 then iterate /*" " " 7? */
do k=5 while s.k<=j /* [↓] divide by the known odd primes.*/
if j // @.k == 0 then iterate j /*Is J ÷ X? Then not prime. ___ */
end /*k*/ /* [↑] only process numbers ≤ √ J */
#= #+1; @.#= j; s.#= j*j; !.j= 1 /*bump # of Ps; assign next P; P²; P# */
end /*j*/; return
- output when using the default inputs:
index │ the smallest primes < 15,000 such that the difference of successive terma are the smallest cubic numbers ───────┼─────────────────────────────────────────────────────────────────────────────────────────────────────────────── 1 │ 2 3 11 19 83 1,811 2,027 2,243 2,251 2,467 11 │ 2,531 2,539 3,539 3,547 4,547 5,059 10,891 12,619 13,619 13,627 21 │ 13,691 13,907 14,419 ───────┴─────────────────────────────────────────────────────────────────────────────────────────────────────────────── Found 23 of the smallest primes < 15,000 such that the difference of successive terma are the smallest cubic numbers 0.040 seconds
Version 2 Using REXX libraries and simplified
Libraries: How to use
Library: Sequences
Library: Functions
Library: Settings
Library: Abend
Below program shows how simple this task can be. Compare it with Version 1, which essentially uses the same technique! Procedure Xpon() (in Functions) returns the exponent of a number, used to control digits and layout. Procedure Primes() (in Sequences, using a sieve) returns an array with primes and a sparse array with numbers flagged as prime; the latter one is needed. So no functions like IsPrime? or IsCube? used! Some other entries employ the same fast solution. Version 2 works up to t = 100 million within acceptable time.
call Time('r')
say 'Cubic special primes - Using REXX libraries'
parse version version; say version; say
t = 150000; d = Xpon(t)+2; numeric digits d
call Primes t
call ShowCubics t,d
say
say Format(Time('e'),,3) 'seconds'
exit
ShowCubics:
procedure expose prim.
arg x,y
call Charout ,Right(2,y)
a = 2; b = 3; k = 1; n = 1
do while b <= x
b = a + k*k*k
if prim.flag.b = 1 then do
n = n+1
call Charout ,Right(b,y)
if n//10 = 0 then
say
a = b; k = 1
end
k = k+1
end
say; say
say n 'cubic special primes found below' x
return
include Sequences
include Functions
- Output:
REXX-Regina_3.9.6(MT) 5.00 29 Apr 2024 Cubic special primes 2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419 23 cubic special primes found below 15000 0.006 seconds
Ring
Also see Quadrat_Special_Primes#Ring
load "stdlib.ring"
see "working..." + nl
Primes = []
limit1 = 50
oldPrime = 2
add(Primes,2)
for n = 1 to limit1
nextPrime = oldPrime + pow(n,3)
if isprime(nextPrime)
n = 1
add(Primes,nextPrime)
oldPrime = nextPrime
else
nextPrime = nextPrime - oldPrime
ok
next
see "prime1 prime2 Gap Cbrt" + nl
for n = 1 to Len(Primes)-1
diff = Primes[n+1] - Primes[n]
for m = 1 to diff
if pow(m,3) = diff
cbrt = m
exit
ok
next
see ""+ Primes[n] + " " + Primes[n+1] + " " + diff + " " + cbrt + nl
next
see "Found " + Len(Primes) + " of the smallest primes < 15,000 such that the difference of successive terma are the smallest cubic numbers" + nl
see "done..." + nl
- Output:
working... prime1 prime2 Gap Cbrt 2 3 1 1 3 11 8 2 11 19 8 2 19 83 64 4 83 1811 1728 12 1811 2027 216 6 2027 2243 216 6 2243 2251 8 2 2251 2467 216 6 2467 2531 64 4 2531 2539 8 2 2539 3539 1000 10 3539 3547 8 2 3547 4547 1000 10 4547 5059 512 8 5059 10891 5832 18 10891 12619 1728 12 12619 13619 1000 10 13619 13627 8 2 13627 13691 64 4 13691 13907 216 6 13907 14419 512 8 Found 23 of the smallest primes < 15,000 such that the difference of successive terma are the smallest cubic numbers done...
RPL
« { } 2 1
DO DUP2 3 ^ +
IF DUP ISPRIME?
THEN 4 ROLL 4 ROLL + ROT DROP SWAP 1
ELSE DROP 1 + END
UNTIL OVER 15000 > END
DROP2
» 'TASK' STO
- Output:
1: { 2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419 }
Ruby
require 'prime'
res = [2]
until res.last > 15000 do
res << (1..).detect{|n| (res.last + n**3).prime? } ** 3 + res.last
end
puts res[..-2].join(" ")
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
Rust
use num_prime::{nt_funcs, Primality, PrimalityTestConfig};
const LIMIT: usize = 15_000;
fn main() {
let mut p = 2_usize;
let mut n = 1;
print!("2 ");
loop {
let new_p = p + n * n * n;
if nt_funcs::is_prime(&new_p, Some(PrimalityTestConfig::strict())) == Primality::Yes {
p = new_p;
n = 1;
print!("{} ", p);
} else {
n += 1;
if new_p >= LIMIT {
break;
}
}
}
}
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
Sidef
func cubic_primes(callback) {
var prev = 2
callback(prev)
loop {
var curr = (1..Inf -> lazy.map { prev + _**3 }.first { .is_prime })
callback(curr)
prev = curr
}
}
say gather {
cubic_primes({|k|
break if (k >= 15000)
take(k)
})
}
- Output:
[2, 3, 11, 19, 83, 1811, 2027, 2243, 2251, 2467, 2531, 2539, 3539, 3547, 4547, 5059, 10891, 12619, 13619, 13627, 13691, 13907, 14419]
Wren
import "./math" for Int
import "./fmt" for Fmt
var primes = Int.primeSieve(14999)
System.print("Cubic special primes under 15,000:")
System.print(" Prime1 Prime2 Gap Cbrt")
var lastCubicSpecial = 3
var gap = 1
var count = 1
Fmt.print("$,7d $,7d $,6d $4d", 2, 3, 1, 1)
for (p in primes.skip(2)) {
gap = p - lastCubicSpecial
if (Int.isCube(gap)) {
Fmt.print("$,7d $,7d $,6d $4d", lastCubicSpecial, p, gap, gap.cbrt.truncate)
lastCubicSpecial = p
count = count + 1
}
}
System.print("\n%(count+1) such primes found.")
- Output:
Cubic special primes under 15,000: Prime1 Prime2 Gap Cbrt 2 3 1 1 3 11 8 2 11 19 8 2 19 83 64 4 83 1,811 1,728 12 1,811 2,027 216 6 2,027 2,243 216 6 2,243 2,251 8 2 2,251 2,467 216 6 2,467 2,531 64 4 2,531 2,539 8 2 2,539 3,539 1,000 10 3,539 3,547 8 2 3,547 4,547 1,000 10 4,547 5,059 512 8 5,059 10,891 5,832 18 10,891 12,619 1,728 12 12,619 13,619 1,000 10 13,619 13,627 8 2 13,627 13,691 64 4 13,691 13,907 216 6 13,907 14,419 512 8 23 such primes found.
XPL0
func IsPrime(N); \Return 'true' if N is prime
int N, I;
[if N <= 2 then return N = 2;
if (N&1) = 0 then \even >2\ return false;
for I:= 3 to sqrt(N) do
[if rem(N/I) = 0 then return false;
I:= I+1;
];
return true;
];
int N, C, T;
[N:= 2;
repeat C:= 1;
loop [T:= N + C*C*C;
if IsPrime(T) then
[IntOut(0, N);
ChOut(0, ^ );
N:= T;
quit;
]
else C:= C+1;
];
until N > 15000;
]
- Output:
2 3 11 19 83 1811 2027 2243 2251 2467 2531 2539 3539 3547 4547 5059 10891 12619 13619 13627 13691 13907 14419
- Draft Programming Tasks
- Prime Numbers
- 11l
- Action!
- Action! Sieve of Eratosthenes
- ALGOL 68
- ALGOL 68-primes
- ALGOL W
- Arturo
- AWK
- C
- C++
- Delphi
- SysUtils,StdCtrls
- EasyLang
- F Sharp
- Factor
- FreeBASIC
- FutureBasic
- Go
- J
- Java
- Jq
- Julia
- Lua
- Mathematica
- Wolfram Language
- Nim
- Perl
- Ntheory
- Phix
- Python
- Raku
- REXX
- Ring
- RPL
- Ruby
- Rust
- Sidef
- Wren
- Wren-math
- Wren-fmt
- XPL0