Reciprocal-root inverse Gaussian sampler ($h(x)=\frac{\lambda (x-\mu {)}^2}{{\mu }^{2}x}$)

For a chi-squared-one variable $Y$, the positive inverse branches of $h(x)=\lambda (x-\mu {)}^2/({\mu }^{2}x)$ have product ${\mu }^2$. Weight a branch value $x$ by $\mu /(\mu +x)$; the two weights sum to one. Multiplying the chi-squared density by this weight and the branch Jacobian gives the inverse Gaussian density on both sides of $\mu$.

With $t=\mu Y/(2\lambda )$, the lower root is $\mu /(1+t+\sqrt{t(t+2)})$ and the upper root is ${\mu }^2/x_{-}$. Rationalizing the lower-root formula avoids subtracting nearly equal large numbers when $Y$ is large.