62 lines
2.2 KiB
Text
62 lines
2.2 KiB
Text
scope # Legendre prime counting function - translation of EasyLang/Go non-memoised version
|
|
|
|
# generate and return a register containing the primes up to sieveSize
|
|
local proc primeSieve( sieveSize :: number ) :: register
|
|
local sieve; create register sieve( sieveSize ); # "vector" to be sieved
|
|
sieve[ 1 ] := false;
|
|
for sPos from 2 to sieveSize do sieve[ sPos ] := true od;
|
|
# sieve the primes
|
|
for sPos from 2 to entier( sqrt( sieveSize ) ) do
|
|
if sieve[ sPos ] then
|
|
for p from sPos * sPos to sieveSize by sPos do
|
|
sieve[ p ] := false
|
|
od
|
|
fi
|
|
od;
|
|
return sieve
|
|
end;
|
|
|
|
# construct a sieve of primes up to the maximum we will need
|
|
local constant psRoot1e9 := primeSieve( sqrt 1e9 );
|
|
|
|
# returns a sequence of primes extracted from a prime sieve
|
|
local proc primesUpTo( n :: number, sieve :: register ) :: sequence
|
|
local result := seq(); # sequence of primes - initially empty
|
|
for sPos from 1 to n do
|
|
if sieve[ sPos ] then insert sPos into result fi
|
|
od;
|
|
return result
|
|
end;
|
|
|
|
local proc pix( n :: integer ) :: integer
|
|
local result := 0;
|
|
if n = 2 then result := 1
|
|
elif n > 2 then
|
|
local constant rootN := entier sqrt n;
|
|
local constant primes := primesUpTo( rootN, psRoot1e9 );
|
|
local constant primeCount := size primes;
|
|
|
|
local proc phi( x :: integer, aIn :: integer ) :: integer
|
|
local a, sum, found1 := aIn, 0, false;
|
|
while a > 1 and not found1 do
|
|
local constant pa := primes[ a ];
|
|
if x <= pa then
|
|
found1 := true
|
|
else
|
|
a -:= 1;
|
|
sum +:= phi( x \ pa, a )
|
|
fi
|
|
od;
|
|
return if found1 then 1 else x - ( x \ 2 ) - sum fi
|
|
end;
|
|
|
|
result := phi( n, primeCount ) + primeCount - 1
|
|
fi;
|
|
return result
|
|
end;
|
|
|
|
for i from 0 to 9 do
|
|
local n := 10 ^ i;
|
|
printf( "10^%d %d\n", i, pix( n ) )
|
|
od
|
|
epocs
|