Hi,
I was messing around with some electronics calculations on the train a
couple of days ago and (feeling very grown up) solved a first order ODE
in terms of a matrix exponential. Imagine my distress when I realised
that Maxima doesn't have a built-in facility for evaluating them in
floating point!
So I've spent quite some time over the last few evenings poring over
papers from the 1970's and I think I have some working code, which I'm
attaching as a patch: I'd love to get some comments.
In particular, if you think this is horribly broken and generally wrong,
I'd like to know about it! If not, I would really like to find a way to
hook into the simplifier so that
exp(matrix ([1,1],[1,1]));
is left alone (fair enough!), but
exp(matrix ([1,1],[1,1.0]));
or
exp(matrix ([1,1b0],[1,1]));
get computed without any more thought on the part of the user. But I
didn't want to spend time on that, or on writing proper documentation,
until I'd got some feedback.
There are some more modern papers about improving this calculation, but
they are basically all based on the same approach (scaled Pad?
approximants). In particular, MATLAB and Octave use this sort of thing,
so we're in good company. What I want to make sure we do better than
those two is getting "the right answer" with bigfloats.
Barton: I made the patch stick the results in the middle of your library
in share/linearalgebra, since it seems to belong well there. In
particular, it would be nice to integrate it with the existing
matrixexp.lisp, which deals with the rational version of this
problem. Also, I'm using a couple of share/linearalgebra's functions, so
it made sense from a dependency point of view too. Let me know if you
think there's a better place to put it?
Let the rotten fruit throwing begin :-)
Rupert
-------------- next part --------------
A non-text attachment was scrubbed...
Name: 0001-A-first-numerical-matrix-exponential-implementation.patch
Type: text/x-diff
Size: 30227 bytes
Desc: not available
URL: <http://www.math.utexas.edu/pipermail/maxima/attachments/20121214/d28458ce/attachment-0001.patch>;
-------------- next part --------------
A non-text attachment was scrubbed...
Name: not available
Type: application/pgp-signature
Size: 315 bytes
Desc: not available
URL: <http://www.math.utexas.edu/pipermail/maxima/attachments/20121214/d28458ce/attachment-0001.pgp>;