8

My simulation needs to implement

np.log(np.cosh(x))

This overflows for large x, i. e. I'm getting the RuntimeWarning: overflow encountered in cosh warning. In principle, as logarithm decreases the number in question, in some range of x, cosh should overflow while log(cosh()) should not.

Is there any solution for that in NumPy, for example similar in spirit to the np.log1p() function?

To provide more info: I am aware that a possible solution might be symbolic using SymPy https://github.com/sympy/sympy/issues/12671 however the simulation should be fast, and symbolic calculation AFAIK might slow it down significantly.

WojciechR
  • 323
  • 1
  • 6

1 Answers1

10

The following implementation of log(cosh(x)) should be numerically stable:

import numpy as np

def logcosh(x):
    # s always has real part >= 0
    s = np.sign(x) * x
    p = np.exp(-2 * s)
    return s + np.log1p(p) - np.log(2)

Explanation:

For real values you could use the following identity:

log(cosh(x)) = logaddexp(x, -x) - log(2)
             = abs(x) + log1p(exp(-2 * abs(x))) - log(2)

which is numerically stable because the argument to exp is always non-positive. For complex numbers we instead require that the argument to exp has non-positive real part, which we achieve by using -x when real(x) > 0 and x otherwise.

myrtlecat
  • 2,156
  • 12
  • 16