Efficiency Problem with project Euler 214 - Totient chains

Viewed 128

I am attempting project Euler question #214 (https://projecteuler.net/problem=214) and after a while playing with things I have the functionality working but not the efficiency.

The task is essentially to use Eulers totient function (phi), which determines the number of co-primes of a number, k, from 1 to k inclusive. Find phi(x), then x = phi(x) and so on until x = 1, forming a chain, and find the sum of all primes less than 40 million which form a chain of length 25.

I initially used a gcd calculator, but efficiency modifications I have made include removing this for a less computationally heavy method, finding just the primes at the start with a sieve to reduce runtime stress, using a dictionary so I don't have to recalculate phi values, and facts such as phi(p) where p is prime = p-1. However, the runtime of my system is still approximately 2 hours, where the program should run within a minute.

I'm also trying not to use library imports which solve the whole thing for me, I just have a timer to see how long the program takes, and the reduce function. Please tell me where my code is too slow and suggest a possible fix. Thanks

from timeit import default_timer as timer
from functools import reduce

def test_prime(n):
    if (n==1):
        return False
    elif (n==2):
        return True;
    else:
        for x in range(2,n):
            if(n % x==0):
                return False
        return True             

def phi(n, phi_values, primes):
    if n in phi_values:
        return phi_values[n]
    else:
        if test_prime(n):
            phi_values[n] = n - 1
            return n - 1
        else:
            factors = prime_factors(n)
            val = int(n * reduce((lambda x, y: x * y), factors))
            phi_values[n] = val
            return val

def prime_factors(n):
    i = 2
    factors = []
    while i * i <= n:
        if n % i:
            i += 1
        else:
            n //= i
            hold = 1 - (1 / i)
            if not hold in factors: factors.append(hold)
    if n > 1:
        hold = 1 - (1 / n)
        if not hold in factors: factors.append(hold)
    return factors

def primes_sieve(n):
    primes = []
    sieve = [True] * (n+1)
    for i in range(2, n+1):
        if (sieve[i]):
            primes.append(i)
            for j in range(i, n+1, i):
                sieve[j] = False
    return primes



def main(max, target):
    times = []
    phi_values = {}
    total = 0
    step = 1
    primes = primes_sieve(max)
    times.append(timer())
    print(times[-1])
    for i in primes:
        if i >= step*100000:
            print(i)
            times.append(timer())
            print(times[-1] - times [-2])
            step += 1
        count = 2
        x = i - 1
        while x > 1:
            x = phi(x, phi_values, primes)
            count += 1
            if count > target:
                break
        if count == target:
            total += i
    times.append(timer())
    print("Sum: " + str(total))
    print(times[-1] - times[0])
    print(times)

if __name__ == "__main__":
    main(40000000, 25)
2 Answers

interesting problem and nice structured code. You can improve the function test_prime by using sqrt(n). I think there is a huge potential.

def test_prime(n):
    if (n==1):
        return False
    elif (n==2):
        return True;
    else:
        for x in range(2,n**0.5):
            if(n % x==0):
                return False
        return True  

The obvious improvement is to note that phi is a multiplicative function, so you don't need to find ALL the prime factors of n, but just some prime factor p. If p**k divides n, (and k is maximal) then phi(n) = p**(k-1) (p-1) * phi(n/p**k). You can also precompute for prime powers.

Related