inv_erf_nsqrt(x::real)::real = begin
if is_nan(x) || is_inf(x) || x < 0 || x > 1
reject("inv_erf_nsqrt: x must be finite and between 0 and 1; found x = ", x)
end
t = log(fma(x, -x, 1.0))
p = 0.0
if abs(t) > 6.125
p = 3.03697567e-10
p = fma(p, t, 2.93243101e-8)
p = fma(p, t, 1.22150334e-6)
p = fma(p, t, 2.84108955e-5)
p = fma(p, t, 3.93552968e-4)
p = fma(p, t, 3.02698812e-3)
p = fma(p, t, 4.83185798e-3)
p = fma(p, t, -2.64646143e-1)
p = fma(p, t, 8.40016484e-1)
else
p = 5.43877832e-9
p = fma(p, t, 1.43285448e-7)
p = fma(p, t, 1.22774793e-6)
p = fma(p, t, 1.12963626e-7)
p = fma(p, t, -5.61530760e-5)
p = fma(p, t, -1.47697632e-4)
p = fma(p, t, 2.31468678e-3)
p = fma(p, t, 1.15392581e-2)
p = fma(p, t, -2.32015476e-1)
p = fma(p, t, 8.86226892e-1)
end
p * x
end