[SciPy-user] Reliability of linalg.eig for complex matrices A,B

Nils Wagner nwagner at mecha.uni-stuttgart.de
Fri Aug 30 02:43:55 CDT 2002


Pearu Peterson schrieb:
> 
> On Thu, 29 Aug 2002, My VDI Freemail wrote:
> 
> > Hi,
> >
> > I tried to solve a transcendental eigenvalue problem with scipy.
> > However It works more or less for r e a l arguments x, but
> > for complex arguments I cannot observe convergency to any
> > eigenvalue. I guess its a problem of the eigenvalue solver.
> > Please find attached my small example.
> >
> > Any ideas or suggestions ?
> 
> The problem is with your algorithm (that was also broken for real
> arguments). The following code tends to converge with arbitrary (but
> reasonable) initial data:
> 
Sorry I made a mistake. It should be x = x + e instead of x = x + 0.9*e
This is the actual output

(1+1j)
1 (1.13701208592+0.19792389747j)
2 (1.29658201804+0.000731255467871j)
3 (1.30037778364+2.78727309871e-07j)
4 (1.30037848976+2.04704281954e-14j)
5 (1.30037848976+4.10207670715e-29j)

Nils

> x = 3+4j
> print x
> e = 1.0
> iter = 0
> while abs(e) > 1.e-10:
>   iter = iter + 1
>   A = array(([x,cos(x)],[x*x+3.0,sin(x)]))
>   B = array(([-1.,sin(x)],[-2.0*x,-cos(x)]))
>   s = linalg.eig(A,B,left=0,right=0)
>   as = abs(s)
>   e = s[list(as).index(min(as))]
>   x = x + 0.9*e
>   print iter, x
> 
> Here is the output:
> (3+4j)
> 1 (3.00008114791+3.09923713594j)
> 2 (3.00077696838+2.1949847382j)
> 3 (3.00510351351+1.27097555405j)
> 4 (3.03297969798+0.218683158764j)
> 5 (0.949729120741+0.269605657507j)
> 6 (1.25450225927+0.0142850278504j)
> 7 (1.29572155364+0.001465000977j)
> 8 (1.29991190209+0.000147112117971j)
> 9 (1.30033182184+1.47176102089e-05j)
> 10 (1.30037382288+1.47182528726e-06j)
> 11 (1.30037802307+1.47183171673e-07j)
> 12 (1.30037844309+1.4718323597e-08j)
> 13 (1.3003784851+1.471832424e-09j)
> 14 (1.3003784893+1.47183243043e-10j)
> 15 (1.30037848972+1.47183243107e-11j)
> 16 (1.30037848976+1.47183243114e-12j)
> 
> Regards,
>         Pearu
> 
> PS: No need to send messages to both scipy-user and scipy-dev mailing
> lists. Just pick one and stick to it.
> 
> _______________________________________________
> SciPy-user mailing list
> SciPy-user at scipy.net
> http://www.scipy.net/mailman/listinfo/scipy-user



More information about the SciPy-user mailing list