I try to find an upper estimate of the number of integer tuples $(i_1,\ldots,i_M)$ such that $i_1!\cdots i_M!\leq s$ for a given real number $s$. I'm especially interested in asymptotics of this number for $M\rightarrow\infty$. Are there any well known results?

First I will assume that you don't count $0!$ and $1!$ as different. If $s$ is a fixed number, and $M\to\infty$, the asymptotic number of solutions is $$\binom{M}{\lfloor \log_2 s\rfloor}.$$ Proof: Note that at most $\lfloor \log_2 s\rfloor$ of the $i_j$ values can be 2 or more. The rest, which is most of them, must be 1. Let $m(\ell)$ be the number of distinct products $i_1!\cdots i_\ell!\le s$ with each factor at least 2. Then the total number of solutions is $$\sum_{\ell=0}^{\lfloor \log_2 s\rfloor} \binom{M}{\ell} m(\ell).$$ Since $s$ is bounded, so $m(\ell)$ is uniformly bounded over all $\ell$, so the sum is asymptotically determined by its last term $\ell=\lfloor \log_2 s\rfloor$. Since $m(\lfloor \log_2 s\rfloor)=1$ (the only case is a lot of $2!$s), the claim follows. If you want $0!$ and $1!$ to be counted as different, there are 2 possibilities for each value not at least 2. The same argument gives the asymptotic value as $$2^{M\lfloor \log_2 s\rfloor}\binom{M}{\lfloor \log_2 s\rfloor}.$$ In both cases the relative error is $O(1/M)$ for fixed $s$. The question becomes more interesting if $s$ is not fixed but increases with $M$. 


If you just care about asymptotics, then you can use central limit theorem or LDP techniques to get an exponential estimate. Here is how: let $x_1, \ldots, x_M$ be iid uniform random variables in $[0,\log s]$, and consider the random sum $$ \sum_{i=1}^M x_i (1 + \log x_i) + \frac{1}{2} \log 2\pi x_i$$, which should converge to some Gaussian with nonzero mean and some variance. You basically want to ask what is the probability that this is less than $\log s$. Since then you can just multiply the probability by $(\log s)^M$ to get an estimate of the form $\exp(f(M,s))$. I am pretty sure if you use Cramer's theorem in Large deviations principle you can compute the tail probability of the sum above. In fact there is a version of Stirling I believe that is a strict upper bound, i.e., something like $n! \le c\exp (n (\log n  1) + \frac{1}{2} \log n)$ for some universal $c$ , which shouldn't matter much. If you want I can try to compute for you. 

