Comment évaluer xe^x/(e^x-1) avec une stabilité numérique autour de zéro et quand x est très positif ou négatif? J'ai accès à toutes les fonctions mathématiques habituelles dans numpy
et scipy
.Comment évaluer xe^x/(e^x-1) avec la stabilité numérique en Python?
Répondre
def f(x):
if abs(x) > 0.1: return x*exp(x)/(exp(x)-1)
else: return 1/(1.-x/2.+x**2/6.-x**3/24.)
L'expansion dans la dernière ligne peut être prolongée de la façon évidente si une plus grande précision est nécessaire, et peut être plus rapidement en re-phrasé. En l'état, il se trompe autant que 1e-6.
Il y a une valeur pas si grande 'z' telle que pour' x> z' la valeur réelle de l'expression 'x * exp (x)/(exp (x) -1)' est très proche de 'x', alors que le calculer directement entraînera un débordement et nan. – Leon
Bon point, @Leon. Pour gérer cette possibilité, les cas de 'x> 0.1' et' x <-0.1' devraient être traités séparément. Si 'x> 0.1', retourner' x/(1-exp (-x)) '; si 'x <-0.1', retourner' x * exp (x)/(exp (x) -1) '; et sinon utiliser l'approximation polynomiale. –
Juste trouvé scipy.special.exprel, qui fonctionnerait. J'ai besoin que mon code fonctionne sur theano ou tensorflow aussi, donc une solution plus basique serait meilleure. –
Frustrant, 'np.expm1' semble donner' nan' au lieu de 'inf' pour les grandes entrées. Cela rend plus difficile d'éviter des choses comme «si» ou «où». – user2357112
Cela semble être un [bug dépendant de la plateforme] (https://github.com/numpy/numpy/issues/6818), toujours ouvert. – user2357112