为何直接传NumPy数组调用函数遇exp溢出警告，间接调用却正常？

背景与实现代码 #

我正在用Python结合SciPy和NumPy实现生存函数（对应doi: 10.1063/1.447954的公式49），实现的函数如下：

from scipy.special import erfc
import numpy as np

def fokker_plank_solution_1_sur(t, k, c, x_0):
    if k == 0 & c == 0:
        pass
    
    else:
        print("Case 2")
        t_1 = (k - c)/(k - 2*c)
        t_2 = np.exp((k * x_0 + k*(k - 2*c)*t))
        t_3 = erfc((x_0 + 2 * (k - c) * t)/(np.sqrt(4*t)))
        t_4 = 0.5
        t_5 = erfc((- x_0 + 2 * c * t)/np.sqrt(4 * t))
        t_6 = k/(k - (2 * c))
        t_7 = np.exp(2 * c * x_0)
        t_8 = erfc((x_0 + 2*c*t)/(np.sqrt(4*t)))
        
        y = t_1 * t_2 * t_3 + t_4 * (t_5 - t_6 * t_7 * t_8)
        
        return y

两种调用情况 #

情况1：间接传入数组（无溢出） #

使用以下代码调用时，能得到正常数值（t_2=1.94e+05）并正常绘图：

import matplotlib.pyplot as plt

D = .25
x_val = np.linspace(1,500,250)
c = 0
x_0 = 5

for k in [.8]:
    y_val = fokker_plank_solution_1_sur(x_val * D, k, c, x_0)
    plt.plot(x_val, y_val)

情况2：直接传入数组（触发溢出警告） #

直接在控制台执行以下调用时：

fokker_plank_solution_1_sur(np.array((50.0)) * .25, 10, 0, 0.8)

触发警告：

RuntimeWarning: overflow encountered in exp

排查发现t_2的计算值为e^1258.0，超出了浮点数的处理范围。

疑问 #

为何间接使用NumPy数组调用正常，直接传入却触发溢出警告？

解答 #

根本原因不是间接/直接传数组的差异，而是两次调用的参数取值完全不同：

• 情况1的参数：k=0.8，t的范围是0.25~125（x_val*D），x_0=5，c=0。计算t_2的指数为：

  k*x_0 + k*(k-2c)*t = 0.8*5 + 0.8*(0.8)*t = 4 + 0.64*t
  

  即使t取最大值125，指数也仅为4+0.64*125=84，e^84约为1.4e+36，远小于64位浮点数的上限（约8e+307），不会溢出。
• 情况2的参数：k=10，t=50*0.25=12.5，x_0=0.8，c=0。计算t_2的指数为：

  k*x_0 + k*(k-2c)*t = 10*0.8 + 10*(10)*12.5 = 8 + 1250 = 1258
  

  e^1258远超64位浮点数的上限，直接触发溢出警告。

简单说，两次调用的参数组合导致指数项的大小天差地别，和数组的传入方式无关。

内容的提问来源于stack exchange，提问作者eB_Bo