|
=={{header|Phix}}==
{{libheader|mpfr}}
Horribly slow...
<lang Phix>include primes.e
Uses primes[] and add_block() from [[Extensible_prime_generator#Phix|Extensible_prime_generator]]
include mpfr.e
and Miller_Rabin() from [[Miller–Rabin_primality_test#Phix|Miller–Rabin_primality_test]]
<lang Phix>while length(primes)<100 do add_block() end while
constant limit = 9999
integer primorial = 1
mpz p = mpz_init(1),
p1 = mpz_init()
randstate state = gmp_randinit_mt()
atom t0 = time()
constant integer base = iff(machine_bits()=32?1000:1000000000)
constant integer digits_perfound = length(sprint(base-1))0, i
constant string dpfmt = sprintf("%%0%dd",digits_per)
mpz_mul_si(p, p, get_prime(n))
for i= -1 to length(bigint)+1 do ▼-- start at the back, number grows to the right (like little endian)
mpz_add_si(p1, p, i)
sequence bigint = {primorial}
if mpz_probable_prime_p(p1,state,25) then
atom result
integer l = mpz_sizeinbase(p,10)
string bisps = iff(l>20?sprintf(" [big (%d digits )]", length(bis)l) ▼
procedure bi_mul(integer pn)
:mpz_get_str(p))
integer carry = 0
printf(1,"%d (%s) is a primorial prime\n",{n, bisps}) ▼▲ for i=1 to length(bigint) do
result = pn*bigint[i]+carry
carry = floor(result/base)
bigint[i] = result-carry*base
while carry <> 0 do
result = carry
carry = floor(carry/base)
bigint &= result-carry*base
end while
end procedure
procedure bi_add(integer carry)
for i=1 to length(bigint) do
result = bigint[i]+carry
carry = floor(result/base)
bigint[i] = result-carry*base
if carry=0 then return end if ▼ end for
if carry!=0 then
if carry<0
or carry!=floor(carry/base) then
?9/0 -- sanity check
end if
bigint &= carry
end if
end procedure
integer found = 0
bi_mul(primes[n])
sequence save_bigint = bigint
for pm1=-1 to 2 by 3 do -- (ie n-1 then n+1)
bi_add(pm1)
string bis = sprint(bigint[$])
for i=length(bigint)-1 to 1 by -1 do
bis &= sprintf(dpfmt,bigint[i])
end for
printf(1,"working [%d]...\r",{n})
if Miller_Rabin(bis,1)=PROBABLY_PRIME then
if length(bis)>20 then
▲ bis = sprintf("[big (%d digits)]",length(bis))
end if
▲ printf(1,"%d (%s) is a primorial prime\n",{n,bis})
found += 1
exit
end if
end for
▲ if carryfound>= 020 then returnexit end if
bigint = save_bigint
if found>=12 then exit end if
{p,p1} = mpz_clear({p,p1})
end for</lang>
state = gmp_randclear(state)</lang>
{{out}}
<pre>
1 (32) is a primorial prime
2 (56) is a primorial prime
3 (2930) is a primorial prime
4 (211210) is a primorial prime
5 (23092310) is a primorial prime
6 (3002930030) is a primorial prime
11 (200560490131200560490130) is a primorial prime
13 (304250263527209304250263527210) is a primorial prime
24 ([big (35 digits)]) is a primorial prime
66 ([big (131 digits)]) is a primorial prime
68 ([big (136 digits)]) is a primorial prime
75 ([big (154 digits)]) is a primorial prime
167 (413 digits) is a primorial prime
171 (425 digits) is a primorial prime
172 (428 digits) is a primorial prime
287 (790 digits) is a primorial prime
310 (866 digits) is a primorial prime
352 (1007 digits) is a primorial prime
384 (1116 digits) is a primorial prime
457 (1368 digits) is a primorial prime
</pre>
|