Ask Your Question
1

Plotting eigenvalues as a function of a variable

asked 2013-01-01 15:30:29 +0200

wrighteap gravatar image

updated 2013-01-01 16:04:23 +0200

I am attempting to plot the the eigenvalues for the matrix:

H=matrix([[1+x, 0, x],[0, 2, 2*x],[x, 2*x, 1+x]])

which were calculated using the following method:

a=H.eigenvalues()

in the following way:

h=fast_callable(real_part(a[2]),vars=[x],domain=CC)

plot(h,(0,0.5))

but unfortunately this does not work for a[0] and a[1],

f=fast_callable(real_part(a[0]),vars=[x],domain=CC)
g=fast_callable(real_part(a[1]),vars=[x],domain=CC)

plot(f,(0,0.5))+plot(g,(0,0.5))

Any suggestions would be greatly appreciated.

edit retag flag offensive close merge delete

2 Answers

Sort by ยป oldest newest most voted
2

answered 2013-01-02 09:35:50 +0200

achrzesz gravatar image

In this approach you can check that the eigenvalues are real:

reset()
f=lambda x:matrix([[1+x, 0, x],[0, 2, 2*x],[x, 2*x, 1+x]]).eigenvalues()
#the first ***real*** eigenvalue versus x 
l0=[[x,f(x)[0]] for x in srange(0,0.5,0.01)]
p0=list_plot(l0,plotjoined=true)
#the second eigenvalue
l1=[[x,f(x)[1]] for x in srange(0,0.5,0.01)]
p1=list_plot(l1,plotjoined=true)
#the third eigenvalue
l2=[[x,f(x)[2]] for x in srange(0,0.5,0.01)]
p2=list_plot(l2,plotjoined=true)
# all
show(p0+p1+p2)
edit flag offensive delete link more

Comments

Very elegant thanks.

wrighteap gravatar imagewrighteap ( 2013-01-02 12:50:36 +0200 )edit
1

answered 2013-01-01 22:25:51 +0200

calc314 gravatar image

updated 2013-01-02 09:37:08 +0200

Something odd is definitely happening here, and I am not the one to help on this situation. But, here is some code that will get the plot you want. It's a bit brute force and not as elegant as your approach, but it does work.

H=matrix([[1+x, 0, x],[0, 2, 2*x],[x, 2*x, 1+x]])
pts0=[]
pts1=[]
pts2=[]
xlist=srange(0,0.5,0.01)
ans=[map(lambda x: n(real_part(x)),H.subs(x=x0).eigenvalues()) for x0 in xlist]
ans=[sorted(a) for a in ans]
for i in range(0,len(xlist)):
    pts0.append([xlist[i],ans[i][0]])
    pts1.append([xlist[i],ans[i][1]])
    pts2.append([xlist[i],ans[i][2]])  
list_plot(pts0,plotjoined=true,color='red')+list_plot(pts1,plotjoined=true,color='blue')+list_plot(pts2,plotjoined=true,color='green')

image description

Using the characteristic polynomial and the implicit_plot command also gives a really nice plot.

var('y')
H=matrix([[1+x, 0, x],[0, 2, 2*x],[x, 2*x, 1+x]])
p(x,y)=H.charpoly('y')
implicit_plot(p(x,y),(x,-3,5),(y,0,5))

image description

You can also use contour_plot to color the regions in the plane by sign.

contour_plot(p(x,y),(x,-3,5),(y,0,5),contours=0,cmap=['red','blue'])
edit flag offensive delete link more

Comments

Thank you very much, at least now I have a working method. The real problem here is my knowledge of programming is somewhat limited so I don't understand what the eigenvalues() method is doing under the hood. One thing that struck me as odd is that both a[0] and a[1] take on complex values for different values of x when H is a hermitian matrix. Granted this could be a computational thing because the complex parts are of order 10^-15 and -19 for the values of x I tried, i.e. effectively 0.

wrighteap gravatar imagewrighteap ( 2013-01-02 07:23:02 +0200 )edit

If you are getting such small complex parts, then you are correct that this is a computational issue. @achrzesz gives some very nice code in his answer.

calc314 gravatar imagecalc314 ( 2013-01-02 09:38:35 +0200 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

Stats

Asked: 2013-01-01 15:30:29 +0200

Seen: 705 times

Last updated: Jan 02 '13