with np.errstate(divide='ignore'):
        sla = _norm_logsf(a)
        slb = _norm_logsf(b)
        tab = np.log1p(-np.exp(slb - sla))
        slx = _norm_logsf(x)