[SciPy-dev] Matrix exponential

jason-sage@creativetra... jason-sage@creativetra...
Sat Feb 28 04:09:50 CST 2009

(It appears that this postscript from John (at the top of the message 
below) didn't make it to the scipy-dev list, probably because he's not 
subscribed, so I'm forwarding the message to the scipy-dev list.)


John Cremona wrote:
> PS An interesting quote from one of Higham's talks:
>                   The availability of expm(A) in
>                        early versions of MATLAB
>                     quite possibly contributed to
> the system’s technical and commercial success.”
>                          — Cleve Moler (2003)
> I get the impression that this is used a lot, though they only seem to
> want double precision (as opposed to multi) which is both fast and has
> predictably bounded error.  The method is a variant of a standard one
> (Pade approximations) with some nice tricks as some once-and-for all
> parameter tuning (for which he said they used Maple!)
> John Cremona
> 2009/2/27  <jason-sage@creativetrax.com>:
>> John Cremona posted the following message to the sage development list about
>> matrix exponentials.  I'm copying it to here since it asks about the scipy
>> matrix exponential method (we say numpy below, but we really mean scipy...)
>> John Cremona wrote:
>>>>> I have just been to a colloquium talk by numerical analyst Nick Higham
>>>>> (Manchester) called "How to compute and not to compute a matrix
>>>>> exponential".  He has new methods which are now in mathematica, matlab
>>>>> and NAG but (apparantly) nowhere else.  He only seemed interested in
>>>>> getting good speed & precision to 16 decimals but (when I asked)
>>>>> confirmed that the methods should apply to give arbitrary precision.
>>>>> I just checked and see that Sage's  matrix exp() uses something stupid
>>>>> except over RDF/CDF where it uses a pade approximation method via
>>>>> numpy.  The method of the talk was a variant of that, the main trick
>>>>> being to use exactly the right order of Pade approx. so maximise
>>>>> precision and speed.
>>>>> I would like to know how good the numpy method is, and  whether it can
>>>>> be improved to this "state of the art" version at least for RDF.  Then
>>>>> it could be another selling point for Sage.
>>> Could you CC the numpy devlist as well on this?  It sounds exciting!
>> I will if you give me the address (or you can perhaps?).  It might be
>> worth including Higham's URL:
>> http://www.maths.manchester.ac.uk/~higham/  as he has lots of his
>> talks up there including some which are similar to the one I heard.

More information about the Scipy-dev mailing list