Let \(A_i :=(\)the number of ways to get \(i\) rolls where all sides have shown up\()\). Additionally, let \(f(x) = \sum_{i=0}^{\infty} \frac{A_ix^i}{i!}\) and \(g(x)=\sum_{i=0}^{\infty} A_ix^i\).

We have \(f(x)= (e^x-1)^N\). Let us find \(g(x)\) from \(f(x)\). Here, note that \(e^{ax}\) in EGF corresponds to \(\frac{1}{1-ax}\) in OGF.

Thus, from \(f(x)=\sum_{i=0}^{N} \binom{N}{i} e^{ix} (-1)^{N-i}\), we have \(g(x) = \sum_{i=0}^{N} \frac{\binom{N}{i}(-1)^{N-i}}{1-ix}\).

Therefore, \(g(x)\) can be represented as \(g(x)=\frac{p(x)}{q(x)}\) using degree-\(n\) polynomials \(p(x)\) and \(q(x)\), which can be computed in \(\mathrm{O}(N\log^2 N)\) using the above formula of \(g(x)\). Here, note that \(q(x)\) has a factor \((1-Nx)\). Let \(q(x) = (1-Nx)r(x)\).

Here, let us define \(B_i:=(\)the probability that the die is rolled exactly \(i\) times before all sides have shown up\()\), and let \(h(x) = \sum_{i=0}^{\infty} B_ix^i\). Then, we have \(h(x) = g\left(\frac{x}{N}\right)(1-x)\), from which we get \(h(x) = \frac{p(\frac{x}{N})}{r(\frac{x}{N})}\).

Now, for the probability generating function \(F(x) = \sum_{i=0}^{\infty} a_ix^i\), let us consider \(F(e^x)\). We have \([x^k] F(e^x) = \sum_{i=0}^{\infty} [x^k] a_ie^{ix} = \sum_{i=0}^{\infty} \frac{a_i \times i^k}{k!}\), so the sought values are the coefficients of \(x^1, x^2, \ldots, x^M\) in \(h(e^x)\).

From \(h(x) = \frac{p(\frac{x}{N})}{r(\frac{x}{N})}\), the problem is reduced to finding the coefficients of degree up to \(M\) in \(C(e^x)\) for a degree-\(N\) polynomial \(C(x) = \sum_{i=0}^{N} c_ix^i\). This can be done by finding the coefficients of degree up to \(M\) in \(\sum_{i=0}^{N} \frac{c_i}{1-ix}\), in \(\mathrm{O}(N\log^2 N + M\log M)\) time.

Therefore, the problem can be solved in \(\mathrm{O}(N\log^2 N + M\log M)\) time, which is fast enough.