The product of two spherical harmonics can be written as the sum of spherical harmonics with coefficients related to Wigner 3j matrices (ref eqn 16) :
$$Y_{l_1m_1}(\theta,\phi)Y_{l_2m_2}(\theta,\phi) = \sum_{lm}W_{3j}(l_1,m_1,l_2,m_2,l,m)Y_{lm}(\theta,\phi)W_{3j}(l_1,0,l_2,0,l,0)$$
where $W_{3j}$ is the Wigner 3$j$ coefficient matrix.
I would like to generalize this to product of multiple spherical harmonics, an expression for :
$$\prod_{i=1}^N Y_{l_im_i}$$
Is there a nice expression for this? Or at least for $N=3$ and $N=4?$