Optimizing code for sequence A064604

Question

I am implementing A064604 - OEIS for positive integers up to 10 billion.

I am finding the divisors in \$O(\sqrt N)\$. So, the overall time complexity of running the formula right now is \$O(N\sqrt N)\$. How do I improve on this?

Pastebin

import math

def factors(n):    
    return set(reduce(list.__add__, 
                ([i, n//i] for i in range(1, int(n**0.5) + 1) if n % i == 0)))

def sigma_4(n):
    l = factors(n)
    ans=0
    for factor in l:
        ans += (pow(factor,4))
    return ans

n=int(raw_input())
ans=0
for i in xrange(1,n+1):
    ans+=sigma_4(i)
print ans

Answer 1

That's a very naive bruteforce algorithm. To optimize this sort of calculation you can usually "turn it around" by doing something like: for each possible divisor, computing the number of times it'll feature as a contributor and adding them up. And if you can get it to work by just iterating over primes, it'll be even better/faster.

Answer 2

Use List Comprehensions for sigma_4

def sigma_4(n):
    return sum([i**4 for i in factors(n)])