Mersenne primes: Difference between revisions

Content added Content deleted
m (→‎{{header|REXX}}: changed a glyph in a comment.)
m (→‎{{header|REXX}}: changed some comments, simplified the first DO loop.)
Line 615: Line 615:
=={{header|REXX}}==
=={{header|REXX}}==
This REXX version   (using a 32-bit Regina REXX interpreter)   will find those Mersenne primes which are less than
This REXX version   (using a 32-bit Regina REXX interpreter)   will find those Mersenne primes which are less than
<br>8 million decimal digits &nbsp; (which would be '''M38''').
<br>8 million decimal digits &nbsp; (which would be '''M43''').
<lang rexx>/*REXX program uses exponent─and─mod operator to test possible Mersenne numbers. */
<lang rexx>/*REXX program uses exponent─and─mod operator to test possible Mersenne numbers. */
do j=1; z=j /*process a range (or run out of time).*/
do j=1; /*process a range, or run out of time.*/
if \isPrime(z) then iterate /*if Z isn't a prime, keep plugging.*/
if \isPrime(j) then iterate /*if J isn't a prime, keep plugging.*/
r= testMer(z) /*If Z is prime, give Z the 3rd degree.*/
r= testMer(j) /*If J is prime, give J the 3rd degree.*/
if r==0 then say right('M'z, 10) "──────── is a Mersenne prime."
if r==0 then say right('M'j, 10) "──────── is a Mersenne prime."
else say right('M'z, 50) "is composite, a factor:" r
else say right('M'j, 50) "is composite, a factor:" r
end /*j*/
end /*j*/
exit /*stick a fork in it, we're all done. */
exit /*stick a fork in it, we're all done. */
Line 631: Line 631:
end /*j*/ /* └─◄ Is j>√ x ? Then return 1*/
end /*j*/ /* └─◄ Is j>√ x ? Then return 1*/
/*──────────────────────────────────────────────────────────────────────────────────────*/
/*──────────────────────────────────────────────────────────────────────────────────────*/
iSqrt: procedure; parse arg x; #= 1; r= 0; do while # <= x; #= # * 4; end
iSqrt: procedure; parse arg x; #= 1; r= 0; do while #<=x; #=#*4; end
do while #>1; #= #%4; _= x-r-#; r= r%2; if _>=0 then do; x=_; r=r+#; end
do while #>1; #=#%4; _= x-r-#; r= r%2; if _>=0 then do; x=_; r=r+#; end
end /*while*/ /*iSqrt ≡ integer square root.*/
end /*while*/ /*iSqrt ≡ integer square root.*/
return r /*───── ─ ── ─ ─ */
return r /*───── ─ ── ─ ─ */
/*──────────────────────────────────────────────────────────────────────────────────────*/
/*──────────────────────────────────────────────────────────────────────────────────────*/
testMer: procedure; parse arg x; p= 2**x /* [↓] do we have enough digits?*/
testMer: procedure; parse arg x; p =2**x /* [↓] do we have enough digits?*/
$$= x2b( d2x(x) ) + 0
$$=x2b( d2x(x) ) + 0
if pos('E',p)\==0 then do; parse var p "E" _; numeric digits _+2; p=2**x; end
if pos('E',p)\==0 then do; parse var p "E" _; numeric digits _+2; p=2**x; end
!.=1; !.1=0; !.7=0 /*array used for a quicker test. */
!.=1; !.1=0; !.7=0 /*array used for a quicker test. */
R=iSqrt(p) /*obtain integer square root of P*/
R=iSqrt(p) /*obtain integer square root of P*/
do k=2 by 2; q= k*x + 1 /*(shortcut) compute value of Q. */
do k=2 by 2; q=k*x + 1 /*(shortcut) compute value of Q. */
m=q // 8 /*obtain the remainder when ÷ 8. */
m=q // 8 /*obtain the remainder when ÷ 8. */
if !.m then iterate /*M must be either one or seven.*/
if !.m then iterate /*M must be either one or seven.*/
parse var q '' -1 _; if _==5 then iterate /*last digit a five ?*/
parse var q '' -1 _; if _==5 then iterate /*last digit a five?*/
if q// 3==0 then iterate /* ÷ by three ?*/
if q// 3==0 then iterate /* ÷ by three? */
if q// 7==0 then iterate /* " " seven ?*/
if q// 7==0 then iterate /* " " seven? */
if q//11==0 then iterate /* " " eleven ?*/
if q//11==0 then iterate /* " " eleven?*/
/* ____ */
/* ____ */
if q>R then return 0 /*Is q>√2**x ? A Mersenne prime*/
if q>R then return 0 /*Is q>√2**x ? A Mersenne prime*/
sq= 1; $= $$ /*obtain binary version from $. */
sq=1; $=$$ /*obtain binary version from $. */
do until $==''; sq= sq * sq
do until $==''; sq=sq*sq
parse var $ _ 2 $ /*obtain 1st digit and the rest. */
parse var $ _ 2 $ /*obtain 1st digit and the rest. */
if _ then sq= (sq+sq) // q
if _ then sq=(sq+sq) // q
end /*until*/
end /*until*/
if sq==1 then return q /*Not a prime? Return a factor.*/
if sq==1 then return q /*Not a prime? Return a factor.*/