Cubic special primes

From Rosetta Code
Revision as of 18:47, 20 June 2021 by Tigerofdarkness (talk | contribs) (Added Algol 68)
Cubic special primes is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.
Task

n   is smallest prime such that the difference of successive terms are the smallest cubics of positive integers, where     n   <   15000.

ALGOL 68

<lang algol68>BEGIN # find a sequence of primes where the members differ by a cube #

   INT max prime = 15 000;
   # sieve the primes to max prime #
   [ 1 : max prime ]BOOL prime;
   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( max prime ) DO
       IF prime[ i ] THEN FOR s FROM i * i BY i + i TO UPB prime DO prime[ s ] := FALSE OD FI
   OD;
   # 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</lang>

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

<lang AWK>

  1. syntax: GAWK -f CUBIC_SPECIAL_PRIMES.AWK
  2. 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)

} </lang>

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

F#

<lang fsharp> // 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 "" </lang>

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

Works with: Factor version 0.98

<lang 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</lang>

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

<lang freebasic>

  1. 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</lang>

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

Translation of: Wren

<lang 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.")

}</lang>

Output:
Same as Wren example.

Julia

<lang 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()))

</lang>

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

Nim

Translation of: Go

<lang 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."</lang>

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

Library: ntheory

<lang 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;</lang>

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

Raku

A two character difference from the Quadrat Special Primes entry. (And it could have been one.) <lang perl6>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)];</lang>
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

<lang rexx>/*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. */

                    @cbp= 'the smallest primes  < '    commas(hi)     " such that the" ,
                          'difference of successive terma are the smallest cubic numbers'

if cols>0 then say ' index │'center(@cbp , 1 + cols*(w+1) ) if cols>0 then say '───────┼'center("" , 1 + cols*(w+1), '─') cbp= 0; idx= 1 /*initialize number of cbp and index.*/ op= 1 $= /*a list of nice primes (so far). */

    do j=0  by 0
                  do k=1  until !.np;  np= op + k**3 /*find the next square + oldPrime.*/
                  if np>= hi  then leave j           /*Is newPrime too big?  Then quit.*/
                  end   /*k*/
    cbp= cbp + 1                                /*bump the number of   cbp's.          */
    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 nice prime ──► list, allow big#*/
    if cbp//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.*/ say say 'Found ' commas(cbp) " of " @cbp 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?             */
                                                /* [↑]  the above five lines saves time*/
              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</lang>
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

Ring

Also see Quadrat_Special_Primes#Ring <lang 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 </lang>

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...

Wren

Library: Wren-math
Library: Wren-fmt

<lang ecmascript>import "/math" for Int, Math import "/fmt" for Fmt

var isCube = Fn.new { |n|

   var c = Math.cbrt(n).round
   return c*c*c == n

}

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 (isCube.call(gap)) {
       Fmt.print("$,7d $,7d $,6d $4d", lastCubicSpecial, p, gap, Math.cbrt(gap).round)
       lastCubicSpecial = p
       count = count + 1
   }

} System.print("\n%(count+1) such primes found.")</lang>

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.