Relation between the Hilbert-Hankel operator and Laplace transform on $L^2(\mathbb R_{>0})$

Question

Let $\mathbb R_{>0} = (0, \infty)$. The Hilbert-Hankel operator $H$ is the integral kernel operator on $L^2(\mathbb R_{>0})$ defined as $$(Hf)(x) = \int_0^\infty \frac{f(y)}{x+y}\, dy$$ Prove that $H = \mathcal L^2$, where $\mathcal L$ is (the Laplace transform) the integral kernel operator on $L^2(\mathbb R_{>0})$ defined as $$(\mathcal L f)(x) = \int_0^\infty e^{-xy} f(y)\, dy$$ Thus, conclude that $\|H\| = \pi$.

---

Questions:

1. The task is to show that $\mathcal L$ is a square root of the Hilbert-Hankel operator, i.e. $H = \mathcal L^2$. Is this square root unique, i.e. does it make sense to define $H^{1/2}:= \mathcal L$?

2. Is my proof (in particular the applicability of Fubini's theorem) correct? I worry because $\frac{1}{\sqrt x} \to \infty$ as $x\to 0$, so I am not sure if Fubini's theorem is really applicable? Although, we only have to fix $x \in (0, \infty)$ and then evaluate the integrals - in which case, $\frac{1}{\sqrt x}$ is certainly just a finite number.

3. So I wonder, are there any interesting results known about $\mathcal L^n$, for $n\in \mathbb N$? Also, can we find the $n$th root of the Hilbert-Hankel operator, i.e. what is $H^{1/n}$ for $n\in \mathbb N$?

---

My work:

If $H = \mathcal L^2$, it is easy to see that $\|H\| = \|\mathcal L^2\| = \|\mathcal L\|^2 = \pi$. This follows since $\|\mathcal L\| = \sqrt\pi$ and $\mathcal L$ is self-adjoint, as shown in this post. For self-adjoint operators $A$ (and more generally, for normal operators) it is true that $\|A\|^2 = \|A^2\|$. The main effort lies in showing $H = \mathcal L^2$. Directly from definitions, we have $$\begin{align} \mathcal L^2f(x) &= \int_0^\infty \int_0^\infty e^{-xy-yz} f(z)\, dz\, dy\\ &= \int_0^\infty\int_0^\infty e^{-xy-yz} f(z)\, dy\, dz\\ &= \int_0^\infty f(z) \int_0^\infty e^{-xy-yz} \, dy \, dz\\ &= \int_0^\infty \frac{f(z)}{x+z}\, dz\\ &= Hf(x) \end{align}$$ where we have used Fubini's theorem. Fubini's theorem is applicable because $$\begin{align} \int_0^\infty\int_0^\infty e^{-xy-yz} |f(z)| dz dy &\le \|f\|_2 \int_0^\infty \left(\int_0^\infty e^{-2xy-2yz} \, dz\right)^{1/2}\, dy\\ &= \|f\|_2 \int_0^\infty \frac{e^{-xy}}{\sqrt{2y}}\, dy \\ &= \sqrt{\frac{\pi}{2x}} \|f\|_2 < \infty \end{align}$$ for $0 < x < \infty$.

Thanks a lot!