Totient maximum¶
Euler's Totient function, φ(n) [sometimes called the phi function], is used to determine the number of numbers less than n which are relatively prime to n. For example, as 1, 2, 4, 5, 7, and 8, are all less than nine and relatively prime to nine, φ(9)=6.
n Relatively Prime φ(n) n/φ(n) 2 1 1 2 3 1,2 2 1.5 4 1,3 2 2 5 1,2,3,4 4 1.25 6 1,5 2 3 7 1,2,3,4,5,6 6 1.1666... 8 1,3,5,7 4 2 9 1,2,4,5,7,8 6 1.5 10 1,3,7,9 4 2.5 It can be seen that n=6 produces a maximum n/φ(n) for n ≤ 10.
Find the value of n ≤ 1,000,000 for which n/φ(n) is a maximum.
As we know that to find the value of Euler Totient Function $\phi(n)$ we can use the following formula: $$\phi(n) = n\prod_{i = 1}^{k}\left(1 - \frac{1}{p^k}\right)$$
Where $p_1, p_2, p_3.....p_k$ are the prime factors of a given number $n$.
But if we had to generate prime factors for each and every number using brute force method would increase time and space complexity. So I wrote an algorithm taking inspiration from Sieve of Eratosthenes.
I will explain my code in bits so that it will be very easy for you to understand.
First thing is to generate prime factors for all the numbers upto 100. We will assume that each and every number until 100 are prime. This assumptions tells us that the prime factors of the numbers are themseleves. For example, if we assume that 20 is a prime number, then the prime factors of 20 will be [20]
.
# we are using 101, so that we can also
# include 0 for the sake of simplicity
# original_code: prime_factors = [[i] for i in xrange(n+1)]
prime_factors = [[i] for i in xrange(101)]
print(prime_factors)
If you have observed, our assumption is completely flawed. We will correct our mistake by finding prime factors for each and every composite number. We will also make a Boolean list which will tell us if the given number is a prime or composite. The list(composite
) will have the following properties.
- For a given number, if the value is
False
then the number is prime. - If the value is
True
for a given number, then the number is a composite number(not a prime number).
# Boolean list for composite numbers
# 101 is to make indexing easy.
composite = [False]*(101)
We will further break down the problem of finding the prime factors of all the numbers. We will first find the prime factors of all the odd composite numbers and later on we will find the prime factors of all the even composite numbers. Remember, the largest prime factor for a given number will be less than or equal to the square root of the number if the number is composite. So we can limit our loop till $\sqrt{100} = 10$.
for i in range(3, int(100**0.5)+1, 2):
if not composite[i]:
for j in range(i, 100/i+1, 2):
mul = i*j
prime_factors[mul] = [i] + prime_factors[j]
composite[mul] = True
print prime_factors
print composite
To make it simple, we will take an example,
i = 3Now
composite[3] = False
, so we pass the if condition. Moving on for first iteration in the nested for loop
j = 3 mul = 9 prime_factors[9] = [3, 3] # we are reassigning the value composite[9] = TrueMoving on to second iteration in the nested for loop,
j = 5 mul = 15 prime_factors[15] = [3, 5] # value of prime_factors[5] = [5] composite[15] = TrueFor third iteration,
j = 7 mul = 21 prime_factors[21] = [3, 7] # value of prime_factors[7] = [7] composite[21] = TrueFor fourth iteration,
j = 9 mul = 27 prime_factors[27] = [3, 3, 3] # value of prime_factors[9] = [3, 3] composite[15] = TrueAnd these iterations will continue...
As you can see, we have only a part of the flawed assumption. Only the odd composite number values are correct, but all the even numbers still seem to follow the assumption. We will correct that mistake.
for j in range(2, 100/2+1):
mul = 2*j
prime_factors[mul] = [2] + prime_factors[j]
composite[mul] = True
print prime_factors
print composite
Now everything is correct, Yay! we have found the prime factors for all the numbers less than 100.
The reason why we have found the prime factors of odd numbers first and later on even numbers will be better understood with the following example:
Lets go back in time, our prime_factors
list had the following values:
[[0], [1], [2], [3], [4], [5], [6], [7], [8], [9], [10], [11], [12], [13], [14], [15], [16], [17], [18], [19],....]First finding the prime factors for all the even composite numbers will result in the
prime_factors
as follows:
[[0], [1], [2], [3], [2, 2], [5], [2, 3], [7], [2, 2, 2], [9], [2, 5], [11], [2, 2, 3], [13], [2, 7], [15], [2, 2, 2, 2], [17], [2, 9], [19],...]That will be a blunder! So we found the prime factors for all the odd composite numbers first and then even composite numbers.
With the help of the prime factors for each number we can find the value of $\phi(n)$. But our requirement is $n/\phi(n)$ $$\phi(n) = n\prod_{i = 1}^{k}\left(1 - \frac{1}{p^k}\right)$$ $$\implies \frac{n}{\phi(n)} = \frac{1}{\prod_{i = 1}^{k}\left(1 - \frac{1}{p^k}\right)}$$
To make it simple for us, we will write a function to do this work. Remember that in the above formula, we can only use a given prime number once. If we had to find the value for 8, we will only use 2 once.
# function to calculate the value of
# n/phi(n)
def totient_function(prime_factors):
multiplication = 1
for factor in set(prime_factors):
multiplication *= (1.0 - 1.0/factor)
return 1.0/multiplication
print totient_function(prime_factors[6]) # lets find the value for 6
And turns out that our function is also working correctly.
Time to see the code for the problem
Program¶
# http://radiusofcircle.blogspot.com
# import time
import time
# start time
start = time.time()
# use xrange function if available
try:
zrange = xrange
except:
zrange = range
# function to calculate the value of
# n/phi(n)
def totient_function(prime_factors):
multiplication = 1
for factor in set(prime_factors):
multiplication *= (1.0 - 1.0/factor)
return 1.0/multiplication
# function to generate prime factors of numbers
# upto given number n
def prime_factors_generator(n):
prime_factors = [[i] for i in zrange(n+1)]
composite = [False]*(n+1)
for i in zrange(3, int(n**0.5)+1, 2):
if not composite[i]:
for j in zrange(i, n/i+1, 2):
mul = i*j
prime_factors[mul] = [i] + prime_factors[j]
composite[mul] = True
for j in zrange(2, n/2+1):
mul = 2*j
prime_factors[mul] = [2] + prime_factors[j]
composite[mul] = True
return prime_factors
# generate prime factors upto 1000000
p_f = prime_factors_generator(1000000)
# remove 0, 1 from the list and find the
# value of n/phi(n)
p_f_nums = map(totient_function, p_f[2:])
# find the index of the maximum of p_h_nums
# add 2 to compensate the removal of 0, 1
print p_f_nums.index(max(p_f_nums))+2
# time at the end of the program
end = time.time()
# total time taken to execute the program
print end - start
The code is very simple and its just an extention of the example which I have explained in detail.
You can download the program from Github Gist pep69.py
Summary¶
Eventhough this problem was less difficult when compared to few other problems, it took a lot of time for me to write a good algorithm to solve this problem. The main problem was a question on whether to use prime numbers or not. After thinking for some time, I mixed up the concept of generating the prime numbers and finding the prime factors and have written a solution. There is a lot to improve in this code. I think we can even write a different algorithm. But I am satisfied with this solution and have sticked to this solution.
If you have any doubt or didn't understand anything then comment in the comment box below and I will be glad to help you. You can also Conact me for any suggestions or for any other help.