almost_primes.sf 2.0 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778
  1. #!/usr/bin/ruby
  2. # a(n) is the smallest n-gonal number with exactly n prime factors (counted with multiplicity).
  3. # https://oeis.org/A358863
  4. # Known terms:
  5. # 28, 16, 176, 4950, 8910, 1408, 346500, 277992, 7542080, 326656, 544320, 120400000, 145213440, 48549888, 4733575168, 536813568, 2149576704, 3057500160, 938539560960, 1358951178240
  6. # PARI/GP program:
  7. #`(
  8. bigomega_polygonals(A, B, n, k) = A=max(A, 2^n); (f(m, p, n) = my(list=List()); if(n==1, forprime(q=max(p,ceil(A/m)), B\m, my(t=m*q); if(ispolygonal(t,k), listput(list, t))), forprime(q = p, sqrtnint(B\m, n), my(t=m*q); if(ceil(A/t) <= B\t, list=concat(list, f(t, q, n-1))))); list); vecsort(Vec(f(1, 2, n)));
  9. a(n, k=n) = if(k < 3, return()); my(x=2^n, y=2*x); while(1, my(v=bigomega_polygonals(x, y, n, k)); if(#v >= 1, return(v[1])); x=y+1; y=2*x); \\ ~~~~
  10. )
  11. func upper_bound(n, from = 2, upto = 2*from) {
  12. say "\n:: Searching an upper-bound for a(#{n})\n"
  13. loop {
  14. var count = n.almost_prime_count(from, upto)
  15. if (count > 0) {
  16. say "Sieving range: [#{from}, #{upto}]"
  17. say "This range contains: #{count.commify} elements\n"
  18. n.almost_primes_each(from, upto, {|v|
  19. if (v.is_polygonal(n)) {
  20. say "a(#{n}) = #{v}"
  21. return v
  22. }
  23. })
  24. }
  25. from = upto+1
  26. upto *= 2
  27. }
  28. }
  29. upper_bound(29)
  30. __END__
  31. a(3) = 28
  32. a(4) = 16
  33. a(5) = 176
  34. a(6) = 4950
  35. a(7) = 8910
  36. a(8) = 1408
  37. a(9) = 346500
  38. a(10) = 277992
  39. a(11) = 7542080
  40. a(12) = 326656
  41. a(13) = 544320
  42. a(14) = 120400000
  43. a(15) = 145213440
  44. a(16) = 48549888
  45. a(17) = 4733575168
  46. a(18) = 536813568
  47. a(19) = 2149576704
  48. a(20) = 3057500160
  49. a(21) = 938539560960
  50. a(22) = 1358951178240
  51. a(23) = 36324805836800
  52. a(24) = 99956555776
  53. a(25) = 49212503949312
  54. a(26) = 118747221196800
  55. a(27) = 59461613912064
  56. a(28) = 13749193801728
  57. a(29) = 7526849672380416
  58. a(30) = 98516240758210560
  59. a(31) = 4969489493917696
  60. a(32) = 78673429816934400
  61. a(33) = 4467570822566903808
  62. a(34) = 1013309912383488000