Why expm(2*A) != expm(A) @ expm(A)

1049 views python
7

According to Matrix exponential, if XY = YX, then exp(X)exp(Y) = exp(X+Y). However when I run the following in Python:

import numpy as np
from scipy.linalg import expm

A = np.arange(1,17).reshape(4,4)

print(expm(2*A))
[[ 306.63168024  344.81465009  380.01335176  432.47730444]
 [ 172.59336774  195.36562731  214.19453937  243.76985501]
 [ -35.40485583  -39.87705598  -42.94545895  -50.01324379]
 [-168.44316833 -190.32607875 -209.76427134 -237.72069322]]

print(expm(A) @ expm(A))
[[1.87271814e+30 2.12068332e+30 2.36864850e+30 2.61661368e+30]
 [4.32685652e+30 4.89977229e+30 5.47268806e+30 6.04560383e+30]
 [6.78099490e+30 7.67886126e+30 8.57672762e+30 9.47459398e+30]
 [9.23513328e+30 1.04579502e+31 1.16807672e+31 1.29035841e+31]]

I get two very different results. Note that @ is just the dot product.

I also tried it in Matlab and the two results are the same as expected. What am I missing here?

answered question

Maybe this is related to the rounding errors?

What version of numpy? on '1.15.0', (expm(2*A) == expm(A)@expm(A)).all() returns TRUE

I have Numpy '1.15.3'

1 Answer

11

To me the root cause of this issue should be similar to the one described in the following NumPy bug report.

posted this

Have an answer?

JD

Please login first before posting an answer.