[Numpy-discussion] Matrix Expontial for differenr t.

[Robert Kern]

On Mon, Jan 28, 2013 at 5:31 PM, Till Stensitzki <mail.till at gmx.de> wrote:
> Hi group,
> is there a faster way to calculate the
> matrix exponential for different t's
> than this:
>
> def sol_matexp(A, tlist, y0):
>     w, v = np.linalg.eig(A)
>     out = np.zeros((tlist.size, y0.size))
>     for i, t in enumerate(tlist):
>         sol_t = np.dot(v,np.diag(np.exp(-w*t))).dot(np.linalg.inv(v)).dot(y0)
>         out[i, :] =  sol_t
>     return out
>
> This is the calculates exp(-Kt).dot(y0) for a list a ts.

You can precalculate the latter part of the expression and avoid the
inv() by using solve().

viy0 = np.linalg.solve(v, y0)
for i, t in enumerate(tlist):
    # And no need to dot() the first part. Broadcasting works just fine.
    sol_t = (v * np.exp(-w*t)).dot(viy0)
    ...

--
Robert Kern