123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778 |
- #!/usr/bin/ruby
- # Author: Trizen
- # Date: 28 February 2021
- # Edit: 23 February 2024
- # https://github.com/trizen
- # Fast recursive algorithm for generating all the odd 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
- # See also:
- # https://oeis.org/A062739
- func odd_powerful_numbers(A, B, k=2) {
- var odd_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) {
- odd_powerful << (m * ipow(j, r)) if j.is_odd
- }
- return nil
- }
- for j in (from .. upto) {
- j.is_even && next
- j.is_coprime(m) || next
- j.is_squarefree || next
- __FUNC__(m * ipow(j, r), r-1)
- }
- }(1, 2*k - 1)
- odd_powerful.sort
- }
- var a = 1e5.irand
- var b = 1e7.irand
- for k in (2..5) {
- say "Testing: k = #{k}"
- assert_eq(
- odd_powerful_numbers(a, b, k),
- k.powerful(a, b).grep{.is_odd}
- )
- }
- say odd_powerful_numbers(1e6 - 2e4, 1e6, 2) #=> [982081, 984987, 985527, 986049, 990025, 990125, 994009, 998001]
|