I am trying to compute
$$ \begin{align*} \mathrm{E}[X^2] &= \lim_{t\to\infty} \int_{0}^{t} x^2 \frac{\lambda^rx^{r-1}\exp(-\lambda x)}{\Gamma (r)}dx \\[2em] &= \frac{\lambda^r}{\Gamma (r)} \lim_{t\to\infty} \int_{0}^{t} x^{r+1}\exp(-\lambda x)dx \end{align*} $$
I have tried integration by parts for the improper integral, but it is rather messy and does not match my numerical result when I fix $r$ and $\lambda$ to specific values. I am messing up somewhere and any help would be appreciated.
Unless I'm misinterpreting the question, we have $$ \begin{align} \int_0^\infty x^{r+1}\exp(-\lambda x)\text{d} x & = \lambda^{-r-2}\int_0^\infty u^{r+1}\exp(-u)\text{d} u \\ & = \lambda^{-r-2}\Gamma(r+2) \end{align} $$ By substitution $x\mapsto \frac x\lambda$. Thus $$ \begin{align} \Bbb{E}[X^2] & = \frac{\lambda^r}{\Gamma(r)}\frac{\Gamma(r+2)}{\lambda^{r+2}} \\ & = \lambda^{-2}\frac{\Gamma(r+2)}{\Gamma(r+1)}\frac{\Gamma(r+1)}{\Gamma(r)} \\ & = \frac{r(r+1)}{\lambda^2} \end{align} $$
