1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677787980 |
- #!/usr/bin/ruby
- # Daniel "Trizen" Șuteu
- # Date: 28 February 2021
- # https://github.com/trizen
- # Fast recursive algorithm for generating all the k-powerful numbers in a given range [a,b].
- # A positive integer n is considered k-powerful, if for every prime p that divides n, so does p^k.
- # Example:
- # 2-powerful = a^2 * b^3, for a,b >= 1
- # 3-powerful = a^3 * b^4 * c^5, for a,b,c >= 1
- # 4-powerful = a^4 * b^5 * c^6 * d^7, for a,b,c,d >= 1
- # OEIS:
- # https://oeis.org/A001694 -- 2-powerful numbers
- # https://oeis.org/A036966 -- 3-powerful numbers
- # https://oeis.org/A036967 -- 4-powerful numbers
- # https://oeis.org/A069492 -- 5-powerful numbers
- # https://oeis.org/A069493 -- 6-powerful numbers
- func powerful_numbers(A, B, k=2) {
- var powerful = []
- func (m,r) {
- var from = 1
- var upto = iroot(idiv(B, m), r)
- if (r <= k) {
- if (A > m) {
- # Optimization by Dana Jacobsen (from Math::Prime::Util::PP)
- with (idiv_ceil(A,m)) {|l|
- if ((l >> r) == 0) {
- from = 2
- }
- else {
- from = l.iroot(r)
- from++ if (ipow(from, r) != l)
- }
- }
- }
- return nil if (from > upto)
- for j in (from .. upto) {
- powerful << (m * ipow(j, r))
- }
- return nil
- }
- for j in (from .. upto) {
- j.is_coprime(m) || next
- j.is_squarefree || next
- __FUNC__(m * ipow(j, r), r-1)
- }
- }(1, 2*k - 1)
- powerful.sort
- }
- var a = 1e5.irand
- var b = 1e7.irand
- for k in (2..5) {
- say "Testing: k = #{k}"
- assert_eq(
- powerful_numbers(a, b, k),
- k.powerful(a, b)
- )
- }
- say powerful_numbers(1e6 - 1e4, 1e6, 2) #=> [990025, 990125, 990584, 991232, 992016, 994009, 995328, 996004, 996872, 998001, 998784, 1000000]
|