Consider the gamma density parameterized in terms of shape $(\alpha)$ and rate $(\beta)$: $$ f(x)=\frac{\beta^\alpha}{\Gamma(\alpha)}x^{\alpha-1}e^{-\beta x}\mathbf 1_{x>0}. $$ Direct computation of the density function from the above expression becomes difficult for large $\alpha$ so how do modern software packages compute it in this situation?
We do have ways of computing $\log\Gamma(\alpha)$ for large $\alpha$ so its seems reasonable to evaluate $$ \log f(x)=\alpha\log\beta-\log\Gamma(\alpha)+(\alpha-1)\log x-\beta x $$ and then exponentiate the result. But I'm skeptical this is how it's done. Can someone please elaborate?