 Stable root evaluation for inverse Gaussian sampling

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.