擅长:python、mysql、java
<p>我认为,@askewchan的解决方案和<a href="http://docs.scipy.org/doc/scipy/reference/generated/scipy.special.log_ndtr.html" rel="nofollow">^{<cd1>}</a>的组合可以做到:</p>
<pre><code>from scipy.special import log_ndtr
_log2 = np.log(2)
_sqrt2 = np.sqrt(2)
def my_func(x):
return np.exp(x ** 2) * 2 * ndtr(x * np.sqrt(2))
def my_func2(x):
return np.exp(x * x + _log2 + log_ndtr(x * _sqrt2))
print(my_func(-150))
# nan
print(my_func2(-150)
# 0.0037611803122451198
</code></pre>
<p>对于<code>x <= -20</code>,<code>log_ndtr(x)</code><a href="https://github.com/scipy/scipy/blob/07c263874ca93312d9a54c1249ff3a7b9a88226e/scipy/special/cephes/ndtr.c#L296-L338" rel="nofollow">uses a Taylor series expansion of the error function to iteratively compute the log CDF directly</a>,这比简单地取{<cd4>}要稳定得多。在</p>
<hr/>
<h2>更新</h2>
<p>正如您在评论中提到的,如果<code>x</code>足够大,<code>exp</code>也可能溢出。虽然您可以使用<code>mpmath.exp</code>来解决这个问题,但是一个更简单、更快的方法是转换成一个<code>np.longdouble</code>,在我的机器上,它可以表示高达1.189731495357231765e+4932的值:</p>
^{pr2}$