　一般に $10^6$ を $2n$ とすれば, $M$ は多項式 $(x+y+z+w+v+u)^{2n}$ において各文字がすべて偶数べきであるような項の係数の和に等しいことに留意する. さらに
　　$$\displaystyle f(x,y,z,w,v,u)=\frac{1}{2^6}\sum_{\lbrace1,-1\rbrace^6} (ix+jy+kz+lw+mv+nu)^{2n}$$
とおけば, $M$ は $f$ の各項の係数の和に等しく, これはすなわち $f(1,1,1,1,1,1)$ で与えられる. したがって,
　　$$\begin{aligned} M&=\frac{1}{64}\sum_{\lbrace1,-1\rbrace^6} (i+j+k+l+m+n)^{2n}\\\\
　　&=\frac{1}{64}\left(\binom{6}{6}6^{2n}+\binom{6}{5}4^{2n}+\binom{6}{4}2^{2n}+\cdots+\binom{6}{0}(-6)^{2n}\right)\\\\
　　&=\frac{1}{32}\left(\binom{6}{6}6^{2n}+\binom{6}{5}4^{2n}+\binom{6}{4}2^{2n}\right)\\\\
　　&=\frac{1}{32}\left(6^{2n}+6\times 4^{2n}+15\times 2^{2n}\right)\end{aligned}$$
　特に $2n=10^6$ において $M$ を $1000$ で割った余りは, Eulerの定理などから $\textbf{696}$ と計算できる. 具体的には, 明らかに $M$ は $8$ の倍数であるから, あとは $M\equiv(1+6+15)\/32=11\/16\pmod{125}$ を計算すればよい.
公式解説
https://onlinemathcontest.com/contests/omc024/editorial/147
　$13$ 以下の**素数** $10^6$ 個からなる**順序付いた**組 $(a_{1},a_{2},\cdots,a_{10^6})$ であって, それらすべての積が平方数であるものは $M$ 個存在します. $M$ を $1000$ で割った余りを求めてください.