OP points to a particular implementation of the mish activation function for accuracy specifications, so I had to characterize this first. That implementation uses single precision (float), and is stable and accurate in the positive half-plane. In the negative half-plane, because it uses logf instead of log1pf, relative error quickly grows a $x\to-\infty$. Loss of accuracy starts around $-1$ and already at $-16.6355324$ the implementation falsely returns $0$, because $\exp(-16.6355324) = 2^{-24}$.
The same accuracy and behavior can be achieved by using a simple mathematical transformation that eliminates $\mathrm{tahn}$, and considering that GPUs typically offer a fused multiply-add (FMA) as well as a fast reciprocal, which one would want to utilize. Exemplary CUDA code looks as follows:
__device__ float my_mishf (float x)
{
float r;
float e = expf (x);
r = 1.0f / fmaf (fmaf (-0.5f, e, -1.0f), e, -1.0f);
r = fmaf (r, x, x);
return r;
}
As with the reference implementation pointed to by OP, this has excellent accuracy in the positive half-plane, and in the negative half-plane error increases rapidly so that at $-16.6355324$ the implementation falsely returns $0$.
If there is a desire to address these accuracy issues, we apply the following observations. For sufficiently small $x$, $f(x) = x \exp(x)$ to within floating-point accuracy. for float computation this holds for $x < -15$. For the interval $[-15,-1]$, we can use a rational approximation $R(x)$ to compute $f(x) := R(x)\exp(x)$. Exemplary CUDA code looks as follows:
__device__ float my_mishf (float x)
{
float r;
float e = expf (x);
if (x >= -1.0f) {
r = 1.0f / fmaf (fmaf (-0.5f, e, -1.0f), e, -1.0f);
r = fmaf (r, x, x);
} else {
float eh = expf (0.5f * x);
float p = 1.03628484e-3f;
p = fmaf (p, x, -7.28869578e-3f);
p = fmaf (p, x, 3.47027853e-2f);
p = fmaf (p, x, -3.54762197e-1f);
p = fmaf (p, x, 8.58785510e-1f);
p = fmaf (p, x, -1.38066018e+0f);
p = fmaf (p, x, 5.97694337e-1f);
float q = 1.03527203e-3f;
q = fmaf (q, x, -7.35638570e-3f);
q = fmaf (q, x, 3.28683965e-2f);
q = fmaf (q, x, -3.79927367e-1f);
q = fmaf (q, x, 6.86127186e-1f);
q = fmaf (q, x, -1.81509328e+0f);
q = fmaf (q, x, 1.00000000e+0f);
r = (1.0f / q) * p;
if (x < -15.0f) r = 1.0f;
r = r * x * eh * eh;
}
return r;
}