Quadrat special primes: Difference between revisions

<lang rexx>/*REXX program finds the smallest primes such that the difference of successive terms */

/*─────────────────────────────────────────────────── are the smallest quadrat numbers. */

parse arg hi cols . /*obtain optional argument from the CL.*/

if hi=='' | hi=="," then hi= 16000 /* " " " " " " */

if cols=='' | cols=="," then cols= 10 /* " " " " " " */

call genP /*build array of semaphores for primes.*/

w= 10 /*width of a number in any column. */

@sqp= 'the smallest primes < ' commas(hi) " such that the" ,

'difference of successive terma are the smallest quadrat numbers'

if cols>0 then say ' index │'center(@sqp , 1 + cols*(w+1) )

if cols>0 then say '───────┼'center("" , 1 + cols*(w+1), '─')

sqp= 0; idx= 1 /*initialize number of sqp and index.*/

op= 1

$= /*a list of nice primes (so far). */

do j=0 by 0

do k=1 until !.np; np= op + k*k /*find the next square + oldPrime.*/

if np>= hi then leave j /*Is newPrime too big? Then quit.*/

end /*k*/

sqp= sqp + 1 /*bump the number of sqp'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 sqp//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(sqp) " of " @sqp

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>

{{out|output|text=&nbsp; when using the default inputs:}}

<pre>

index │ the smallest primes < 16,000 such that the difference of successive terma are the smallest quadrat numbers

───────┼───────────────────────────────────────────────────────────────────────────────────────────────────────────────

1 │ 2 3 7 11 47 83 227 263 587 911

11 │ 947 983 1,019 1,163 1,307 1,451 1,487 1,523 1,559 2,459

21 │ 3,359 4,259 4,583 5,483 5,519 5,843 5,879 6,203 6,779 7,103

31 │ 7,247 7,283 7,607 7,643 8,219 8,363 10,667 11,243 11,279 11,423

41 │ 12,323 12,647 12,791 13,367 13,691 14,591 14,627 14,771 15,671

Found 49 of the smallest primes < 16,000 such that the difference of successive terma are the smallest quadrat numbers

</pre>