Efficient way to reduce a binomial coefficient as a fraction Here is the full problem.

You need to calculate Euler's totient function of a binomial coefficient $C_n^k$.
Input
The first line contains two integers: $n$ and $k$ $(0 \le k \le n \le 500000)$.
Output
Print one number $\varphi (C_n^k)$ modulo $10^9+7$.

My thoughts:
It is known that $$\varphi(a)=a \prod_{p|a}(1-\frac{1}{p}) $$
where $p$ are prime numbers divide $a$.
Hence, if we can obtain somehow vector&lt;int&gt; multipliers that contains divisors of $C_n^k$ then we can easily do the following steps in order to calculate $\varphi(C_n^k)$:

Multiply all elements of that vector modulo $10^9+7$. Let's call the result by result
Then we can iterate through all prime numbers that divide any element of multipliers(these prime numbers can be obtained by a minor modification of sieve of Eratosthenes). Since $1-\frac{1}{p}=\frac{p-1}{p}$ we can update the result by:

result = divideMod(multiplyMod(result, p-1), p)

where divideMod and multiplyMod are functions doing corresponding operations modulo $10^9+7$.
And yes, we can do modulus division since $10^9+7$ is prime.
By doing all that stuff we get what we needed: $\varphi(C_n^k)$ modulo $10^9+7$. This all idea now requires just a vector multipliers. Here is my attempt to get it:
I need to write a function calculates the combinations number $C_n^k$. The function shouldn't return the total result of the operation(because it can be too large since $(0 \le k \le n \le 500000)$). It should return the vector&lt;int&gt; which contains divisors of that number. Let's do some math:
$$
C_n^k = \frac{n!}{(n-k)! k!} \\
=\frac{n(n-1)(n-2)...(n-k+1)}{k(k-1)(k-2)...1}
$$
So now I need to reduce this fraction. And the question is: what is the most efficient way to do this(in terms of time)?
I've tried the following. Consider the numerator and denominator are represented by vector&lt;int&gt; numerator={n, n-1, ..., n-k+1} and vector&lt;int&gt; denominator={k, k-1, ..., 1} respectively.
    vector&lt;long&gt; numerator(k);
    vector&lt;long&gt; denominator(k);
    for (int i = 0; i&lt;k; i++) {
        numerator[i] = n-i;
        denominator[i] = k-i;
    }

    vector&lt;long&gt; multipliers;
    for (int i = 0; i &lt; k; i++) {
        for (int j = 0; j &lt; k; j++) {
            if (numerator[i] == 1) 
                break;
            
            long greatest_common_divisor = gcd(numerator[i], denominator[j]);
            numerator[i] /= greatest_common_divisor;
            denominator[j] /= greatest_common_divisor;
        }
        if (numerator[i] != 1) 
            multipliers.push_back(numerator[i]);
    }

As you can see I just go through all numbers in numerator and denominator and divide them by their greatest common divisor.
Time complexity of this algorithm is $O( k^2 log(nk) )$
It's too big and for this solution contest system returns time limit exceeded.($0 \le k \le n \le 500000$)
Does there exist more efficient way?
You store the numbers n-k+1 to n in an array.

Then for each prime number p ≤ k: Find which power of p is a factor of k! (Thats k/p + k/p^2 + k/p^3 ... ) Then remove that power from the array: Find the first number divisible by p (that would be the number at index 0 if n-k+1 is divisible by p, otherwise at index p - ((n-k+1) modulo p)). That number is divisible by p, possible by p^2 etc., and the next number divisible by p is at the index p higher). 

You are now left with an array of k numbers that you need to multiply.