Как работать с scipy и чрезвычайно большими числами


Я хотел бы использовать специальные функции scipy для работы с очень большими числами

alpha = 9999
def y(t):
    return 1 / (special.lambertw(alpha * math.exp(alpha-t)) + 1)

math.exp выдает ошибку переполнения, что неудивительно. Поэтому я попытался использовать десятичный модуль вместо

alpha = 9999
def y(t):
    exp = decimal.Decimal(math.exp(1))
    exp = exp ** alpha
    exp = exp * decimal.Decimal(math.exp(-t))
    return 1 / (special.lambertw(alpha * math.exp(alpha-t)) + 1)

Но получаем следующую ошибку:

TypeError: ufunc '_lambertw' not supported for the input types, and the 
inputs could not be safely coerced to any supported types according to 
the casting rule ''safe''

special.lambertw происходит от Сципиона

Каков был бы правильный способ справиться с этим?

1 4

1 ответ:

Один из вариантов-использовать mpmath. Она включает в себя реализацию lambertw.

Например,

In [20]: import mpmath

In [21]: mpmath.mp.dps = 30

In [22]: alpha = 9999

In [23]: def y(t):
    ...:     return 1 / (mpmath.lambertw(alpha * mpmath.exp(alpha-t)) + 1)
    ...: 

In [24]: y(1.5)
Out[24]: mpf('0.000100015000749774938119797735952206')

В общем случае вы не сможете использовать специальные функции scipy с очень большими значениями. Большая часть кода scipy реализована на C, C++ или Fortran и ограничена 64-битными значениями с плавающей запятой, с максимальным значением около 1. 8e308:

In [11]: np.finfo(np.float64).max
Out[11]: 1.7976931348623157e+308