Python中支持向量化大数值sinh的mpmath更快替代方案

你遇到的核心矛盾是：NumPy的向量化sinh速度快但会溢出，mpmath的高精度sinh无溢出但循环遍历速度极慢。以下是几种兼顾速度与精度的解决方案：

方案1：使用NumPy扩展精度类型 #

NumPy的longdouble（或float128，依平台支持）提供比float64更大的数值范围，能避免大部分溢出问题，同时保持向量化计算的速度。

import numpy as np

# 将数组转换为扩展精度类型
A1_ext = A1.astype(np.longdouble)
# 执行向量化sinh计算
sinhZ1 = np.sinh(A1_ext)
# 若后续计算允许，可转回float64（视数值范围而定）
sinhZ1 = sinhZ1.astype(np.float64)

注意：不同硬件平台对longdouble的支持不同（如x86为80位，部分ARM平台为64位），需提前测试兼容性。

方案2：用Numba加速mpmath循环 #

Numba的JIT编译能大幅优化Python循环的执行效率，无需修改原有逻辑即可提升mpmath计算速度。

首先安装Numba：

pip install numba

然后修改代码：

import numpy as np
import mpmath as mp
from numba import jit

@jit(nopython=False)  # 因调用mpmath无法使用nopython模式，仍有显著加速效果
def compute_mp_sinh(arr):
    result = np.empty(arr.shape, dtype=object)
    for i in range(arr.shape[0]):
        for j in range(arr.shape[1]):
            result[i,j] = mp.sinh(arr[i,j])
    return result

# 调用加速后的函数
sinhZ1 = compute_mp_sinh(A1)

该方案的速度通常能接近NumPy向量化的水平，同时保留mpmath的高精度特性。

方案3：SymPy向量化符号计算 #

SymPy支持对矩阵进行逐元素的符号函数操作，既能处理大数值无溢出，也能通过向量化减少循环开销。

import numpy as np
import sympy as sp

# 将NumPy数组转为SymPy矩阵
A1_sym = sp.Matrix(A1)
# 逐元素应用sinh
sinhZ1_sym = A1_sym.applyfunc(sp.sinh)
# 转换回NumPy数值数组
sinhZ1 = np.array(sinhZ1_sym).astype(np.float64)

SymPy的速度介于NumPy与纯mpmath循环之间，适合需要符号计算与数值结果结合的场景。

方案4：大数值sinh近似计算 #

当x极大时，sinh(x) ≈ exp(x)/2（因exp(-x)趋近于0可忽略），直接用NumPy的exp计算，完全保持向量化速度，且能处理更大的数值范围。

import numpy as np

# 设定阈值：float64中exp(709)为最大有限值，对应sinh(x)的阈值约700
mask = np.abs(A1) > 700
sinhZ1 = np.empty_like(A1)

# 普通数值用标准sinh
sinhZ1[~mask] = np.sinh(A1[~mask])
# 大数值用近似公式
sinhZ1[mask] = np.exp(A1[mask]) / 2

此方案速度最快，若你的级数计算中极大值项的近似误差可接受，这是最优选择。

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