[SciPy-dev] Re: [SciPy-user] Matrix square root

Nils Wagner nwagner at mecha.uni-stuttgart.de
Tue Oct 5 05:47:54 CDT 2004


Travis Oliphant wrote:

> Pearu Peterson wrote:
>
>>
>>
>> On Mon, 4 Oct 2004, Nils Wagner wrote:
>>
>>> AFAIK there is no direct support for sqrtm but you may use funm.
>>
>>
>>
>> I just discovered that there might be a bug in funm. I have added
>> a test to scipy.linalg.test_matfuncs that demonstrates the bug.
>> So, be aware to check the results of funm(..,sqrt).
>
>
> As far as I undertand, the idea of a general purpose matrix function 
> has precision problems.   The algorithm in scipy is based on a 
> published algorithm using the Schur decomposition.   That is why it is 
> usually a good idea to write separate algorithms for special matrix 
> functions if possible.
>
> Is this really a bug, or just a problem with using the generic funm 
> algorithm when a more accurate algorithm exists.
>
> The fact that you can return a reasonable square root for non-special 
> matrices makes me believe that it is not a bug in funm, but a 
> demonstration that the generic funm algorithm does not always produce 
> the best results for specific functions of matrices. 
>
>
> -Travis O.
>
>
> _______________________________________________
> Scipy-dev mailing list
> Scipy-dev at scipy.net
> http://www.scipy.net/mailman/listinfo/scipy-dev


Travis,

Please find attached a quick and dirty implementation of sqrtm using
Denman-Beavers and Schulz iteration as suggested by Higham.
This might be a starting point for sqrtm in scipy...
Any comments ?

Nils



 


-------------- next part --------------
A non-text attachment was scrubbed...
Name: narep.py
Type: text/x-python
Size: 1100 bytes
Desc: not available
Url : http://projects.scipy.org/pipermail/scipy-dev/attachments/20041005/5d208142/attachment.py 


More information about the Scipy-dev mailing list