Greatest prime dividing the n-th cubefree number: Difference between revisions
Content deleted Content added
Added FreeBASIC |
|||
Line 109: | Line 109: | ||
The 1000000th term of a[n] is 1202057 |
The 1000000th term of a[n] is 1202057 |
||
The 10000000th term of a[n] is 1202057</pre> |
The 10000000th term of a[n] is 1202057</pre> |
||
=={{header|Phix}}== |
|||
Quite possibly flawed, but it does finish in the blink of an eye. |
|||
<syntaxhighlight lang="phix"> |
|||
with javascript_semantics |
|||
-- Note this routine had a major hiccup at 27000 = 2^3*3^3*5^3 and another was predicted at 9261000 = 27000*7^3. |
|||
-- The "insufficient kludge" noted below makes it match Wren over than limit, but quite possibly only by chance. |
|||
function cubes_before(atom n) |
|||
-- nb: if n is /not/ cube-free it /is/ included in the result. |
|||
-- nth := n-cubes_before(n) means n is the nth cube-free integer, |
|||
-- but only if cubes_before(n-1)==cubes_before(n), |
|||
-- otherwise cubes_before(cubicate) isn't really very meaningful. |
|||
integer r = 0 |
|||
bool xpm = true -- extend prior multiples |
|||
sequence pm = {} |
|||
for i=1 to n do -- (or eg floor(cbrt(n))) |
|||
atom p3 = power(get_prime(i),3) |
|||
if p3>n then exit end if |
|||
integer k = floor(n/p3) |
|||
for mask=1 to power(2,length(pm))-1 do |
|||
integer m = mask, mi = 0 |
|||
atom kpm = p3 |
|||
while m do |
|||
mi += 1 |
|||
if odd(m) then |
|||
kpm *= pm[mi] |
|||
end if |
|||
m = floor(m/2) |
|||
end while |
|||
if kpm>n then |
|||
if count_bits(mask)=1 then |
|||
xpm = false |
|||
pm = pm[1..mi-1] |
|||
exit |
|||
end if |
|||
else |
|||
-- subtract? already counted multiples.. |
|||
integer l = floor(n/kpm) |
|||
-- DEV insufficient kludge... (probably) |
|||
--if count_bits(mask)=1 then |
|||
if odd(count_bits(mask)) then -- match Wren |
|||
k -= l |
|||
else |
|||
k += l |
|||
end if |
|||
end if |
|||
end for |
|||
r += k |
|||
if xpm then |
|||
pm &= p3 |
|||
end if |
|||
end for |
|||
return r |
|||
end function |
|||
function cube_free(atom nth) |
|||
-- get the nth cube-free integer |
|||
atom lo = nth, hi = lo*2, mid, cb, k |
|||
while hi-cubes_before(hi)<nth do |
|||
lo = hi |
|||
hi = lo*2 |
|||
end while |
|||
-- bisect until we have a possible... |
|||
while true do |
|||
mid = floor((lo+hi)/2) |
|||
cb = cubes_before(mid) |
|||
k = mid-cb |
|||
if k=nth then |
|||
-- skip omissions |
|||
while cubes_before(mid-1)!=cb do |
|||
mid -= 1 |
|||
cb -= 1 |
|||
end while |
|||
exit |
|||
elsif k<nth then |
|||
lo = mid |
|||
else |
|||
hi = mid |
|||
end if |
|||
end while |
|||
return mid |
|||
end function |
|||
function A370833(atom n) |
|||
if n=1 then return 1 end if |
|||
atom nth = cube_free(n) |
|||
sequence f = prime_powers(nth) |
|||
return f[$][1] |
|||
end function |
|||
atom t0 = time() |
|||
printf(1,"First 100 terms of a[n]:\n%s\n",join_by(apply(tagset(100),A370833),1,10," ",fmt:="%3d")) |
|||
for n in sq_power(10,tagset(7,3)) do |
|||
printf(1,"The %,dth term of a[n] is %,d.\n",{n,A370833(n)}) |
|||
end for |
|||
?elapsed(time()-t0) |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
First 100 terms of a[n]: |
|||
1 2 3 2 5 3 7 3 5 11 |
|||
3 13 7 5 17 3 19 5 7 11 |
|||
23 5 13 7 29 5 31 11 17 7 |
|||
3 37 19 13 41 7 43 11 5 23 |
|||
47 7 5 17 13 53 11 19 29 59 |
|||
5 61 31 7 13 11 67 17 23 7 |
|||
71 73 37 5 19 11 13 79 41 83 |
|||
7 17 43 29 89 5 13 23 31 47 |
|||
19 97 7 11 5 101 17 103 7 53 |
|||
107 109 11 37 113 19 23 29 13 59 |
|||
The 1,000th term of a[n] is 109. |
|||
The 10,000th term of a[n] is 101. |
|||
The 100,000th term of a[n] is 1,693. |
|||
The 1,000,000th term of a[n] is 1,202,057. |
|||
The 10,000,000th term of a[n] is 1,202,057. |
|||
"0.1s" |
|||
</pre> |
|||
=={{header|RPL}}== |
=={{header|RPL}}== |