-
Notifications
You must be signed in to change notification settings - Fork 0
/
prog.sf
50 lines (41 loc) · 1.09 KB
/
prog.sf
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
#!/usr/bin/ruby
# The number of n-almost primes less than or equal to 9^n, starting with a(0)=1.
# https://oeis.org/A116429
# Known terms:
# 1, 4, 26, 181, 1095, 6416, 35285, 187929, 973404, 4934952, 24628655, 121375817, 592337729, 2868086641, 13798982719, 66043675287, 314715355788
# Corrected terms:
# a(16) = 314715355786
# New terms:
# a(17) = 1494166794434
# a(18) = 7071357084444
# a(19) = 33374079939405
#`{
# PARI/GP program:
almost_prime_count(N, k) = if(k==1, return(primepi(N))); (f(m, p, k, j=0) = my(c=0, s=sqrtnint(N\m, k)); if(k==2, forprime(q=p, s, c += primepi(N\(m*q))-j; j += 1), forprime(q=p, s, c += f(m*q, q, k-1, j); j += 1)); c); f(1, 2, k);
a(n) = if(n == 0, 1, almost_prime_count(9^n, n)); \\ ~~~~
}
#for n in (0..100) {
for n in (0..100) {
say [n, n.almost_prime_count(9**n)]
}
__END__
[0, 1]
[1, 4]
[2, 26]
[3, 181]
[4, 1095]
[5, 6416]
[6, 35285]
[7, 187929]
[8, 973404]
[9, 4934952]
[10, 24628655]
[11, 121375817]
[12, 592337729]
[13, 2868086641]
[14, 13798982719]
[15, 66043675287]
[16, 314715355786]
[17, 1494166794434]
[18, 7071357084444]
[19, 33374079939405]