58 lines
1.6 KiB
Plaintext
58 lines
1.6 KiB
Plaintext
with javascript_semantics
|
|
constant LIMIT = 1e7, N = floor((floor(LIMIT/3)-1)/4)+1
|
|
|
|
function Sieve4n_3_Primes()
|
|
sequence sieve = repeat(0,N), p4n3 = {}
|
|
for idx=1 to N do
|
|
if sieve[idx]=0 then
|
|
integer n = idx*4-1
|
|
p4n3 &= n
|
|
if idx+n>N then
|
|
// collect the rest
|
|
for j=idx+1 to N do
|
|
if sieve[j]=0 then
|
|
p4n3 &= 4*j-1
|
|
end if
|
|
end for
|
|
exit
|
|
end if
|
|
for j=idx+n to N by n do
|
|
sieve[j] = 1
|
|
end for
|
|
end if
|
|
end for
|
|
return p4n3
|
|
end function
|
|
|
|
sequence p4n3 = Sieve4n_3_Primes(),
|
|
BlumField = repeat(false,LIMIT),
|
|
Blum50 = {}, counts = repeat(0,10)
|
|
|
|
for idx,n in p4n3 do
|
|
for bj in p4n3 from idx+1 do
|
|
atom k = n*bj
|
|
if k>LIMIT then exit end if
|
|
BlumField[k] = true
|
|
end for
|
|
end for
|
|
integer count = 0
|
|
for n,k in BlumField do
|
|
if k then
|
|
if count<50 then Blum50 &= n end if
|
|
counts[remainder(n,10)] += 1
|
|
count += 1
|
|
if count=50 then
|
|
printf(1,"First 50 Blum integers:\n%s\n",{join_by(Blum50,1,10," ",fmt:="%3d")})
|
|
elsif count=26828 or remainder(count,1e5)=0 then
|
|
printf(1,"The %,7d%s Blum integer is: %,9d\n", {count,ord(count),n})
|
|
if count=4e5 then exit end if
|
|
end if
|
|
end if
|
|
end for
|
|
printf(1,"\nPercentage distribution of the first 400,000 Blum integers:\n")
|
|
for i,n in counts do
|
|
if n then
|
|
printf(1," %6.3f%% end in %d\n", {n/4000, i})
|
|
end if
|
|
end for
|