The problem asks you to carry out three tasks for each number in a range:
- Determine if it is square-free.
- Count the number of distinct prime divisors.
- Find the sum of its divisors.
When you have to do something involving the divisors of each number in a range, then the technique that you need is sieving. This is well-known when the problem is to finding prime numbers (you use the Sieve of Eratosthenes), but sieves can also be used to efficiently establish other facts about the divisibility of numbers in a range.
So let's look at task 1 (finding the square-free numbers in a range). If we had a list of prime numbers then we could take the square of each prime number in turn and then go through the multiples of that square in the range, marking them as non-square-free, like this:
def round_up(n, m):
"Return the smallest multiple of m that's greater than or equal to n."
n += m - 1
return n - n % m
def square_free(start, stop, primes):
"""Return a list of the square-free numbers in the range from start
(inclusive) to stop (exclusive). The argument primes must be an
iterable of prime numbers up to at least sqrt(stop - 1).
>>> square_free(1, 11, [2, 3]) # https://oeis.org/A005117
[1, 2, 3, 5, 6, 7, 10]
>>> square_free(100, 111, [2, 3, 5, 7])
[101, 102, 103, 105, 106, 107, 109, 110]
"""
square_free = [True] * (stop - start)
for p in primes:
p2 = p ** 2
if p2 >= stop:
break
for q in range(round_up(start, p2) - start, stop - start, p2):
square_free[q] = False
return [i for i, sf in enumerate(square_free, start) if sf]
Note:
I've written this code in the form of functions with documented arguments and results. This makes the code easier to understand (because you can consider each function in isolation, and you can read the documentation), easier to test (because you can call the function with whatever arguments you choose—and the examples in the docstrings can be automatically checked using the doctest module), and easier to reuse (we'll see below that we have more uses for the round_up function).
I tested the code against the On-Line Encyclopedia of Integer Sequences as shown by the link to sequence A005117 in the docstring. This encyclopedia is an invaluable resource when solving numerical problems: not just for sequences that can be used as test cases, but also for the mathematical discussion of each sequence.
Now, for task 2 (finding the sums of divisors of the square-free numbers in a range). If a number \$n\$ is square-free then it is the product of distinct primes, that is, \$n = pqr\dots\$ for distinct primes \$p, q, r, \ldots\$. Let's take an example, say \$n = 30\$, which is the product of the distinct primes \$2, 3, 5\$, and look at the sum of divisors: $$ σ(30) = 1 + 2 + 3 + 5 + 6 + 10 + 15 + 30. $$ Let's group these into odd and even divisors: $$ \eqalign{ σ(30) &= (1 + 3 + 5 + 15) + (2 + 6 + 10 + 30) \\ &= (1 + 2)(1 + 3 + 5 + 15)} $$ and then group the remainder \$ 1 + 3 + 5 + 15 \$ according to their divisibility by 3: $$ \eqalign{ σ(30) &= (1 + 2)((1 + 5) + (3 + 15)) \\ &= (1 + 2)(1 + 3)(1 + 5).} $$ In the general case, if \$n = pqr\dots\$, then $$ σ(n) = (1 + p)(1 + q)(1 + r)\dots. $$ (See this answer for a more detailed look at computing the sum of divisors.)
Now, if we have an iterable of prime numbers, then we can compute all the sums of square-free divisors of the numbers in a range like this:
def sum_square_free_divisors(start, stop, primes):
"""Return a list giving the sum of square-free divisors of every
number in the range from start (inclusive) to stop (exclusive).
The argument primes must be an iterable of prime numbers up to at
least stop - 1.
>>> sum_square_free_divisors(1, 11, [2, 3, 5, 7]) # https://oeis.org/A048250
[1, 3, 4, 3, 6, 12, 8, 3, 4, 18]
"""
result = [1] * (stop - start)
for p in primes:
if p >= stop:
break
for i in range(round_up(start, p) - start, stop - start, p):
result[i] *= (1 + p)
return result
There are two more things we need to be able to do before we can put together a solution to the whole problem. We need to be able to construct an iterable of primes. Again, we can do this using a sieve; see this answer for some suggestions about how to implement it. I'm going to use sieve3 here.
And finally, we need to be able to count the number of distinct primes dividing a number. Again, it's going to be most efficient if we write yet another sieve:
def distinct_prime_divisors(start, stop, primes):
"""Return a list of the number of distinct primes (from the iterable
primes) that divide the numbers in range from start (inclusive) to
stop (exclusive).
>>> distinct_prime_divisors(1, 11, [2, 3, 5, 7]) # https://oeis.org/A001221
[0, 1, 1, 1, 1, 2, 1, 1, 1, 2]
"""
result = [0] * (stop - start)
for p in primes:
if p >= stop:
break
for i in range(round_up(start, p) - start, stop - start, p):
result[i] += 1
return result
When putting all this together, we could just go through each of test cases one at a time and solve it separately. But there could be as many as 100,000 test cases, so it's important to avoid repetition of work. The set of prime numbers does not change from one test case to another, nor does the set of square-free numbers or their sums of divisors. So what we should do is to read all \$T\$ test cases, determine the smallest value for \$L\$ and the biggest value for \$R\$, and then run all our sieves just the once.
Even once we've run our sieves and have a list of all the good numbers, we still don't want to iterate over all the good numbers between \$L\$ and \$R\$, adding up the sums of their divisors (because there might be 100,000 test cases and there could be as many as 34,693 good numbers in the range, leading to more than 3 billion additions). Instead … well, this answer is getting quite long, so see if you can figure out how it works by reading the code:
from itertools import accumulate
from bisect import bisect_left, bisect_right
def good_numbers_solutions(cases):
"""Given a list of pairs (L, R) describing test cases for the Good
Numbers problem, return a list of solutions.
>>> good_numbers_solutions([(1, 5), (6, 10)])
[6, 30]
"""
# Minimum value of L for any test case.
L_min = min(L for L, _ in cases)
# Maximum value of R for any test case, plus 1.
R_limit = max(R for _, R in cases) + 1
# List and set of primes below R_limit.
primes = sieve3(R_limit)
prime_set = set(primes)
# Square-free numbers between L_min and R_limit.
sf = square_free(L_min, R_limit, primes)
# Sums of square-free divisors for n between L_min and R_limit.
ssfd = sum_square_free_divisors(L_min, R_limit, primes)
ssfd_limit = max(ssfd) + 1
# Counts of distinct prime divisors below ssfd_limit.
dpd = distinct_prime_divisors(0, ssfd_limit, sieve3(ssfd_limit))
# Good numbers between L_min and R_limit.
good = [n for n in sf if dpd[ssfd[n - L_min]] in prime_set]
# Running sum of sums-of-square-free divisors of good numbers.
running = [0] + list(accumulate(ssfd[n - L_min] for n in good))
# Solve the test cases.
return [running[bisect_right(good, R)] - running[bisect_left(good, L)]
for L, R in cases]
(You'll need the documentation for itertools.accumulate, bisect.bisect_left and bisect.bisect_right.)
6\n30running your program. \$\endgroup\$