As a follow up to this question I am also interested in a symbolic closed form for this integral
$$\int_0^\infty d r \,r^2\, j_{n_1}( k_1 r)\, j_{n_2}( k_2 r)\, j_{n_3}( k_3 r)\,, $$
where $j_n(r)$ is the $n^{\rm th}$ order spherical Bessel function, $k_1$,$k_2$ and $k_3$ are real positive numbers and $n_1,n_2$ and $n_3$ are positive integers.
The spherical Bessel function $j_n$ can be defined by
$$ j_n(x) = (-x)^n \left(\frac{1}{x}\frac{d}{dx}\right)^n\,\frac{\sin x}{x}.$$
Clue
As an answer to this question, @joriki provided a nice solution for $n_1=n_2=n_3=0$.
If I am to believe Mathematica again, for instance
$$\int_0^\infty d r \,r^2\, j_2( r)\, j_2( 2 r)\, j_2( 3 r)=-\frac{\pi}{48}$$
and
$$\int_0^\infty d r \,r^2\, j_2( r)\, j_4( 2 r)\, j_4( 3 r)=-\frac{\pi}{48}$$
so the integral seems possible. On the other hand, if some $n_i$ are odd the integral seems ill defined.
I would guess that for odd indices the answer is $\pi/(8k_1 k_2 k_3)$ times some function of the signs of $n_1$, $n_2$ and $n_3$.
Update
My guess seems to be wrong. Symbolic integration for the first $8\times8\times 8 $ values of $(n_1,n_2,n_3)$ yields (with $k_1=k_2=k_3=1$)
Best Answer
R. Mehrem gives this as equation 5.14 in the paper The Plane Wave Expansion, Infinite Integrals and Identities involving Spherical Bessel Functions (arXiv:0909.0494v4). It looks like this:
The symbols in curly braces are Wigner 6j-symbols, and the angle-bracketed letters are Clebsch-Gordan coefficients. The P_l is a Legendre polynomial, and Delta is k1^2+k2^2-k3^2/(2k1k2), coming from applying law of cosines to a triangle in k-space (it is just cos theta_12)). The beta function is zero unless the k1, k2, k3 combination is such that it forms a closed triangle (i.e. the side lengths satisfy the triangle inequality).