🚀 Feature

Implement the erfcx(x) special function, which computes exp(x^2) * erfc(x) in a numerically stable way. Also for convenience, add logerfc(x) = log(erfc(x)) and logerfcx(x) = log(erfcx(x)).

erfcx is available in many numerical packages, such as Matlab, Julia, SciPy R, and others.

From erfcx it is easy to obtain logerfc and logerfcx, but this involves a conditional which can be slow in pure Python code. So I recommend adding logerfc and logerfcx as well, which can be implemented as:

def logerfc(x): 
    if x > 0.0:
        return log(erfcx(x)) - x**2
    else:
        return log(erfc(x))

def logerfcx(x):
    if x < 0.0:
        return log(erfc(x)) + x^2
    else:
        return log(erfcx(x))

Motivation

These special functions are very useful whenever we have to work with truncated normal distributions.

Related: #2129, #32293

cc @mruberry @rgommers @heitorschueroff