This post will give the analysis to Project Euler #242 from hackerrank, which is an extended version from the original.
Given 5 integers $m$, $r$, $n$, $k$, and $M$, count the number of k-subsets of $\{1, 2, \ldots, n\}$ such that the sum of the subset is $r \pmod m$. Let the answer be $S$, output $S \times m$ modulo $M$.
This is actually an extended version of IMO 1995 P6 (having $m = p$, $n = 2p$, $k = p$, $r = 0$). A discussion to this problem can be found in AOPS.
Back to the original problem, consider the generating function $(1 + xy)(1 + xy^2) \ldots (1 + xy^n)$. Our answer will be the sum of coefficients of $x^k y^r, x^k y^{r+m}, x^k y^{r+2m}, \ldots$.
To get this sum, we use the Series Multisection formula on the y coefficients, then our generating function become:
$$\frac{1}{m} \sum_{l=0}^{m-1} \left( \omega^{-lr} \prod_{i=1}^n (1 + x\omega^{li}) \right)$$
where $\omega$ is the primitive $m$-th root of unity. Since the final answer will be multiplied by $m$, the $\frac{1}{m}$ will be ignored for the rest of this post. Now our goal is to find the coefficient of $x^k$ from that equation.
Using the Cauchy Binomial Theorem, we convert our equation above to:
$$\sum_{l=0}^{m-1} \left( \omega^{-lr} \sum_{i=0}^n x^i \omega^{li(i+1)/2} {n \brack i}_{\omega^l} \right)$$
Since we only care with the $x^k$ coefficient, then the answer to our original problem is the evaluation of:
$$ \sum_{i=0}^{m-1} \left( (\omega^{i})^{k(k+1)/2 - r} {n \brack k}_{\omega^i} \right) $$
We will compute the sum for each primitive of $\omega^i$ separately. Denote $\omega_d$ as the primitive $d$-th root of unity. Also note the q-Lucas theorem: ${n \brack k}_{\omega_d} = {\lfloor n/d \rfloor \choose \lfloor k/d \rfloor } {n \mod d \brack k \mod d}_{\omega_d}$.
Since $n \mod d \leq k \mod d$, the value of ${n \mod d \brack k \mod d} = 1$ if $n \equiv k \pmod d$, $0$ otherwise.
With these in mind, we rewrite our equation to:
$$ \sum_{d|m} \left( \sum_{(d, i) = 1}^d (\omega_d^i)^{k(k+1)/2 - r} {n \choose k}[n \equiv k \pmod d] \right) $$
where $[]$ is the Iverson bracket. Recall that the sum of $k$-th power of the $d$-th root of unity is $\frac{\mu(d / \gcd(k, d)) \varphi(d)}{\varphi(d / \gcd(k, d))}$, where $\mu$ is the Möbius function and $\varphi$ is the Euler’s Totient function.
Therefore, our final answer is:
$$ \sum_{d|m} \left( \frac{\mu\left(\frac{d}{\gcd(k(k+1)/2 - r, d)}\right) \varphi(d)}{\varphi\left(\frac{d}{\gcd(k(k+1)/2 - r, d)}\right)} {n \choose k}[n \equiv k \pmod d] \right) $$
And now, for the computational part. Compute the final sum for every prime factor of $M$, then combine the answers using Chinese Remainder Theorem. The computation of binomial theorem modulo $p^e$ can be computed in $O(e(p + e + \log n))$, which is described in this post.
The computation of the divisors of $m$ and all its mobius and euler totient values can be computed in $O(\sqrt m)$
Therefore, the overall complexity is $O(\sqrt m + \log M \log n)$.