包除原理を用いて解きます。

\(i\) のうち、決められた \(k_{j}\) 個が \(\left\lfloor{\frac{i}{M}}\right\rfloor=\left\lfloor{\frac{P_{i}}{M}}\right\rfloor=j\) を満たす状態を考えます。

このとき、 \(k_{j}\) 個の決め方は \(\dbinom{M}{k_{j}}\) 通りで、これらの index の \(P_{i}\) の決め方は \(\frac{M!}{(M-k_{j})!}\) 通りです。

そして、制約の課されていない残りの \((NM-\sum_{j=0}^{N-1}k_{j})\) 個の index の決め方は \((NM-\sum_{j=0}^{N-1}k_{j})!\) 通りです。

以上より、 \((k_{0},k_{1},...,k_{N-1})\) が決まっているときの場合の数は以下のようになります。

\[(NM-\sum_{j=0}^{N-1}k_{j})!\prod_{j=0}^{N-1}\frac{M!}{((M-k_{j})!)^2k_{j}!}\]

よって答えは

\[ans=\sum_{(k_{0},k_{1},...,k_{N-1})}(-1)^{\sum_{j=0}^{N-1}k_{j}}(NM-\sum_{j=0}^{N-1}k_{j})!\prod_{j=0}^{N-1}\frac{M!}{((M-k_{j})!)^2k_{j}!}\]

となります。

この式は \(k\) の総和に依存した式と、 \(k_{j}\) に依存した式の積の総和なので、まとめて計算することができます。

具体的には

\[f(x)=\sum_{k=0}^{M}\frac{M!}{((M-k)!)^2k!}(-x)^k,\ g(x)=(f(x))^{N}\]

としたとき、

\[ans=\sum_{s=0}^{NM}(NM-s)!([x^s]g(x))\]

となります。

\(g(x)\) は、FFTと繰り返し二乗法を用いることで \(\mathrm{O}(NM \log(NM))\) で求めることができるので、全体の計算量も \(\mathrm{O}(NM\log(NM))\) です。

実装例(C++)

なお、この問題では \(NM\lt 2^{20}\) であって \(g(x)\) のすべての項の係数を計算できることから、長さ \(2^{20}\) の係数列の FFT / 逆 FFT を \(1\) 回ずつ行うことで \(g(x)\) の係数列を計算することができます。詳しくは実装例をご覧ください。

実装例(C++)