Given the densities $\varphi$ and $\psi$ of two independent random variables, the probability that their product is less than $z$ is
$$
\iint_{xy< z}\varphi(x)\psi(y)\,\mathrm{d}x\,\mathrm{d}y\tag{1}
$$
Letting $w=xy$ so that $x=w/y$ yields
$$
\iint_{w< z}\varphi\left(\frac{w}{y}\right)\psi(y)\,\mathrm{d}\frac{w}{y}\,\mathrm{d}y=\iint_{w< z}\varphi\left(\frac{w}{y}\right)\psi(y)\,\mathrm{d}w\,\frac{\mathrm{d}y}{y}\tag{2}
$$
Taking the derivative of $(2)$ with respect to $z$ gives the density of the product of the random variables to be
$$
\phi(z)=\int\varphi\left(\frac{z}{y}\right)\psi(y)\,\frac{\mathrm{d}y}{y}\tag{3}
$$
We can compute the expected value using this distribution as
$$
\begin{align}
\mathrm{E}(Z)
&=\int z\phi(z)\,\mathrm{d}z\\
&=\iint z\,\varphi\left(\frac{z}{y}\right)\psi(y)\,\frac{\mathrm{d}y}{y}\,\mathrm{d}z\\
&=\iint xy\,\varphi(x)\psi(y)\,\mathrm{d}y\,\mathrm{d}x\tag{4}
\end{align}
$$
which is exactly what one would expect when computing the expected value of the product directly.

In the same way, we can also compute
$$
\begin{align}
\mathrm{E}(Z^2)
&=\int z^2\phi(z)\,\mathrm{d}z\\
&=\iint z^2\,\varphi\left(\frac{z}{y}\right)\psi(y)\,\frac{\mathrm{d}y}{y}\,\mathrm{d}z\\
&=\iint x^2y^2\,\varphi(x)\psi(y)\,\mathrm{d}y\,\mathrm{d}x\tag{5}
\end{align}
$$
again getting the same result as when computing this directly.

The variance is then, as usual, $\mathrm{E}(Z^2)-\mathrm{E}(Z)^2$.