Carmichael 3 strong pseudoprimes: Difference between revisions
Content added Content deleted
m (→numbers in sorted order: optimized displaying of sorted numbers. -- ~~~~) |
m (→numbers in order of calculation: optimized code. -- ~~~~) |
||
Line 282: | Line 282: | ||
<br>the '''modulus''' function was hard-coded below (when using a negative value). |
<br>the '''modulus''' function was hard-coded below (when using a negative value). |
||
===numbers in order of calculation=== |
===numbers in order of calculation=== |
||
<lang rexx>/*REXX program |
<lang rexx>/*REXX program calculates Carmichael 3-strong pseudoprimes (up to N).*/ |
||
numeric digits 30 /*in case user wants bigger nums.*/ |
numeric digits 30 /*in case user wants bigger nums.*/ |
||
parse arg N .; if N=='' then N=61 /*allow user to specify the limit*/ |
parse arg N .; if N=='' then N=61 /*allow user to specify the limit*/ |
||
Line 290: | Line 290: | ||
!.=0 /*a method of prime memoization. */ |
!.=0 /*a method of prime memoization. */ |
||
do p=3 to N by 2; if \isPrime(p) then iterate /*Not prime? Skip.*/ |
do p=3 to N by 2; if \isPrime(p) then iterate /*Not prime? Skip.*/ |
||
pm=p-1; |
pm=p-1; nps=-p*p; @.=0; min=1e9; max=0 /*some handy-dandy variables.*/ |
||
do h3=2 to pm; g=h3+p /*find Carmichael #s for this P. */ |
do h3=2 to pm; g=h3+p /*find Carmichael #s for this P. */ |
||
do d=1 to g-1 |
do d=1 to g-1 |
||
if g*pm//d\==0 |
if g*pm//d\==0 then iterate |
||
if ((nps//h3)+h3)//h3\==d//h3 then iterate |
|||
q=1+pm*g%d; if \isPrime(q) then iterate |
|||
r=1+p*q%h3; if q*r//pm\==1 then iterate |
|||
if \isPrime(r) then iterate |
|||
carms=carms+1 /*bump the Carmichael # counter. */ |
carms=carms+1 /*bump the Carmichael # counter. */ |
||
min=min(min,q); max=max(max,q); @.q=r /*build a list.*/ |
|||
end /*d*/ |
end /*d*/ |
||
end /*h3*/ |
end /*h3*/ |
||
/*display a list of some Carm #s.*/ |
|||
do j=min to max by 2; if @.j==0 then iterate /*one of the #s?*/ |
|||
say '──────── a Carmichael number: ' p times j times @.j |
|||
end /*j*/ |
|||
say /*show bueatification blank line.*/ |
say /*show bueatification blank line.*/ |
||
end |
end /*p*/ |
||
say; say carms ' Carmichael numbers found.' |
say; say carms ' Carmichael numbers found.' |
||
exit /*stick a fork in it, we're done.*/ |
exit /*stick a fork in it, we're done.*/ |
||
/*──────────────────────────────────ISPRIME subroutine──────────────────*/ |
/*──────────────────────────────────ISPRIME subroutine──────────────────*/ |
||
isPrime: procedure expose !.; parse arg x; if !.x then return 1 |
isPrime: procedure expose !.; parse arg x; if !.x then return 1 |
||
if wordpos(x,'2 3 5 7')\==0 then do; !.x=1; return 1; end |
if wordpos(x,'2 3 5 7 11 13')\==0 then do; !.x=1; return 1; end |
||
if x< |
if x<17 then return 0; if x//2==0 then return 0; if x//3==0 then return 0 |
||
if right(x,1)==5 then return 0; if x//7==0 then return 0 |
|||
do i=11 by 6 until i*i>x; if x// i ==0 then return 0 |
|||
if x//(i+2) ==0 then return 0 |
if x//(i+2) ==0 then return 0 |
||
end /*i*/ |
end /*i*/ |