Quadrat special primes: Difference between revisions

From Rosetta Code
Content added Content deleted
(Added Algol 68)
(julia example)
Line 108: Line 108:
10667 11243 11279 11423 12323 12647 12791 13367 13691 14591 14627 14771
10667 11243 11279 11423 12323 12647 12791 13367 13691 14591 14627 14771
15671
15671
</pre>


=={{header|Julia}}==
<lang julia>using Primes

function quadrat(N = 16000)
pmask = primesmask(1, 2N)
qprimes = [2]
while (n = qprimes[end]) < N
for i in 1:isqrt(N)
q = n + i * i
if q > N
return qprimes
elseif pmask[q] # got next qprime
push!(qprimes, q)
break
end
end
end
end

println("Quadrat special primes < 16000:")
foreach(p -> print(rpad(p[2], 6), p[1] % 10 == 0 ? "\n" : ""), enumerate(quadrat()))
</lang>{{out}}
<pre>
Quadrat special primes < 16000:
2 3 7 11 47 83 227 263 587 911
947 983 1019 1163 1307 1451 1487 1523 1559 2459
3359 4259 4583 5483 5519 5843 5879 6203 6779 7103
7247 7283 7607 7643 8219 8363 10667 11243 11279 11423
12323 12647 12791 13367 13691 14591 14627 14771 15671
</pre>
</pre>


Revision as of 00:04, 29 March 2021

Quadrat 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 squares of positive integers, where     n   <   16000.

ALGOL 68

Translation of: ALGOL W

<lang algol68>BEGIN # find some primes where the gap between the current prime and the next is a square # 

   # reurns a sieve of primes up to n #
   PROC eratosthenes = ( INT n )[]BOOL:
        BEGIN
           [ 1 : n ]BOOL p;
           p[ 1 ] := FALSE; p[ 2 ] := TRUE;
           FOR i FROM 3 BY 2 TO n DO p[ i ] := TRUE  OD;
           FOR i FROM 4 BY 2 TO n DO p[ i ] := FALSE OD;
           FOR i FROM 3 BY 2 TO ENTIER sqrt( n ) DO
               IF p[ i ] THEN FOR s FROM i * i BY i + i TO n DO p[ s ] := FALSE OD FI
           OD;
           p
        END # eratosthenes # ;
   # an array of squares #
   PROC get squares = ( INT n )[]INT:
        BEGIN
           [ 1 : n ]INT s;
           FOR i TO n DO s[ i ] := i * i OD;
           s
        END # get squares # ;
   INT max number = 16000;
   []BOOL prime = eratosthenes( max number );
   []INT square = get squares(  max number );
   INT   p count, this prime, next prime;
   # the first square gap is 1 (between 2 and 3) the gap between all other primes is even #
   # so we treat 2-3 as a special case                                                    #
   p count := 1; this prime := 2; next prime := 3;
   print( ( " ", whole( this prime, -5 ) ) );
   WHILE next prime < max number DO
       this prime := next prime;
       p count   +:= 1;
       print( ( " ", whole( this prime, -5 ) ) );
       IF p count MOD 12 = 0 THEN print( ( newline ) ) FI;
       INT sq pos := 2;
       WHILE next prime := this prime + square[ sq pos ];
             IF next prime < max number THEN NOT prime[ next prime ] ELSE FALSE FI
       DO sq pos +:= 2 OD
   OD

END</lang>

Output:
     2     3     7    11    47    83   227   263   587   911   947   983
  1019  1163  1307  1451  1487  1523  1559  2459  3359  4259  4583  5483
  5519  5843  5879  6203  6779  7103  7247  7283  7607  7643  8219  8363
 10667 11243 11279 11423 12323 12647 12791 13367 13691 14591 14627 14771
 15671

ALGOL W

<lang algolw>begin % find some primes where the gap between the current prime and the next is a square % 

   % 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 ;
   % sets s( 1 :: n ) to the squares %
   procedure getSquares ( integer array s ( * ) ; integer value n ) ;
       for i := 1 until n do s( i ) := i * i;
   integer MAX_NUMBER;
   MAX_NUMBER := 16000;
   begin
       logical array prime(  1 :: MAX_NUMBER );
       integer array square( 1 :: MAX_NUMBER );
       integer       pCount, thisPrime, nextPrime;
       % sieve the primes to MAX_NUMBER %
       Eratosthenes( prime, MAX_NUMBER );
       % calculate the squares to MAX_NUMBER %
       getSquares( square, MAX_NUMBER );
       % the first gap is 1 (between 2 and 3) the gap between all other primes is even %
       % so we treat 2-3 as a special case                                             %
       pCount := 1; thisPrime := 2; nextPrime := 3;
       write( i_w := 6, s_w := 0, " ", thisPrime );
       while nextPrime < MAX_NUMBER do begin
           integer sqPos;
           thisPrime := nextPrime;
           pCount := pCount + 1;
           writeon( i_w := 6, s_w := 0, " ", thisPrime );
           if pCount rem 12 = 0 then write();
           sqPos := 2;
           while begin
               nextPrime := thisPrime + square( sqPos );
               nextPrime < MAX_NUMBER and not prime( nextPrime )
           end do sqPos := sqPos + 2;
       end while_thisPrime_lt_MAX_NUMBER
   end

end.</lang>

Output:
      2      3      7     11     47     83    227    263    587    911    947    983
   1019   1163   1307   1451   1487   1523   1559   2459   3359   4259   4583   5483
   5519   5843   5879   6203   6779   7103   7247   7283   7607   7643   8219   8363
  10667  11243  11279  11423  12323  12647  12791  13367  13691  14591  14627  14771
  15671


Julia

<lang julia>using Primes

function quadrat(N = 16000) 

   pmask = primesmask(1, 2N)
   qprimes = [2]
   while (n = qprimes[end]) < N
       for i in 1:isqrt(N)
           q = n + i * i
           if  q > N
               return qprimes
           elseif pmask[q]  # got next qprime
               push!(qprimes, q)
               break
           end
       end
   end

end

println("Quadrat special primes < 16000:") foreach(p -> print(rpad(p[2], 6), p[1] % 10 == 0 ? "\n" : ""), enumerate(quadrat()))

</lang>

Output:
Quadrat special primes < 16000:
2     3     7     11    47    83    227   263   587   911   
947   983   1019  1163  1307  1451  1487  1523  1559  2459
3359  4259  4583  5483  5519  5843  5879  6203  6779  7103
7247  7283  7607  7643  8219  8363  10667 11243 11279 11423
12323 12647 12791 13367 13691 14591 14627 14771 15671

Raku

<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: * > 16000, :k)];</lang>
Output:
49 matching numbers:
    2     3     7    11    47    83   227
  263   587   911   947   983  1019  1163
 1307  1451  1487  1523  1559  2459  3359
 4259  4583  5483  5519  5843  5879  6203
 6779  7103  7247  7283  7607  7643  8219
 8363 10667 11243 11279 11423 12323 12647
12791 13367 13691 14591 14627 14771 15671

REXX

<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>
output   when using the default inputs:
 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

Ring

<lang ring> load "stdlib.ring"

see "working..." + nl

Primes = [] limit1 = 50 oldPrime = 1 add(Primes,2)

for n = 1 to limit1 

   nextPrime = oldPrime + pow(n,2)
   if isprime(nextPrime)
      n = 1
      add(Primes,nextPrime)
      oldPrime = nextPrime
   else
      nextPrime = nextPrime - oldPrime
   ok

next

see "prime1 prime2 Gap" + nl for n = 1 to Len(Primes)-1 

   diff = Primes[n+1] - Primes[n]
   see ""+ Primes[n] + "      " + Primes[n+1] + "    " + diff + nl

next

see nl + "done..." + nl </lang>

Output:
working...
prime1 prime2 Gap
2      3    1
3      7    4
7      11    4
11      47    36
47      83    36
83      227    144
227      263    36
263      587    324
587      911    324
911      947    36
947      983    36
983      1019    36
1019      1163    144
1163      1307    144
1307      1451    144
1451      1487    36
1487      1523    36
1523      1559    36
1559      2459    900
2459      3359    900
3359      4259    900
4259      4583    324
4583      5483    900
5483      5519    36
5519      5843    324
5843      5879    36
5879      6203    324
6203      6779    576
6779      7103    324
7103      7247    144
7247      7283    36
7283      7607    324
7607      7643    36
7643      8219    576
8219      8363    144
8363      10667    2304
10667      11243    576
11243      11279    36
11279      11423    144
11423      12323    900
12323      12647    324
12647      12791    144
12791      13367    576
13367      13691    324
13691      14591    900
14591      14627    36
14627      14771    144
14771      15671    900
done...
Retrieved from "https://rosettacode.org/wiki/Quadrat_special_primes?oldid=310124"