Ask Your Question
1

Plotting question

asked 2011-01-05 10:20:55 -0500

mhfrey gravatar image

updated 2011-05-12 09:16:46 -0500

Kelvin Li gravatar image

I am trying to plot Sage calculation listed below. I get the following error:

verbose 0 (3989: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 200 points. verbose 0 (3989: plot.py, generate_plot_points) Last error message: 'unable to simplify to float approximation'

What am I doing wrong?

Thank you for your help.

Gather all for Full H off resonance

Ix=1/2*matrix(SR,2,2,[[0,1],[1,0]])
Iy=1/(2*I)*matrix(SR,2,2,[[0,1],[-1,0]])
Iz=1/2*matrix(SR,2,2,[[1,0],[0,-1]])
Al=vector(SR,[1,0])
Be=vector(SR,[0,1])

g = 26.7522128e7 # 1H Gamma rad/s/T
rf_field = 100 # B1 rf field in kHz
Brf = rf_field*1e3*2*pi/g*2 # B1 in T for Wnut/2pi= 100 kHz
Thetarf = pi/2 # angle between BO and Brf

Wnut = (1/2*g*Brf*sin(Thetarf))
t360 = 1/((1/2*g*Brf*sin(Thetarf))/(2*pi))# 360 pulse
tip = 180 # Tip angle
tp = tip/360*t360
Bo = 500e6*2*pi/g# T Spectrometer Field Strength
wref = g*Bo # spectrometer reference frequency for Protons, rad/S
Bo.n(),(wref/(2*pi)).simplify()/1e6;(Wnut/(2*pi)).simplify()/1e3,tp.simplify()*1e6

Ratio = var ('Ratio')# Wnut/O0, Ratio creates offset, 0=on res
Boo=(Wnut*Ratio+wref)/g
woo = g*Boo
O0=woo-wref

Phip = var('Phip') # Phip = phase of rf pulse
Hrfrf = (O0)*Iz + Wnut*(Ix*cos(Phip)+Iy*sin(Phip))
Pabrf=((abs(Be*exp(-i*Hrfrf(Phip=0)*tp)*Al))^2)



for Ratio in srange (-1,1,0.25,include_endpoint=True):
    Pabrf(Ratio=Ratio).n(digits=5) #(pi)x on Alpha off resonance Full H

plot (Pabrf,(Ratio, -1,1))
edit retag flag offensive close merge delete

Comments

looks like nmr to me

Evgeny gravatar imageEvgeny ( 2011-01-05 10:41:52 -0500 )edit

two quick tips that don't address your actual problem - you can use "var('Ratio')" instead of "Ratio = var('Ratio')", and the loop in the second and third last lines is unnecessary if you want to use plot()

Felix Lawrence gravatar imageFelix Lawrence ( 2011-01-05 11:51:30 -0500 )edit

5 answers

Sort by ยป oldest newest most voted
1

answered 2011-05-24 12:03:33 -0500

nate_iverson gravatar image

Try

def Pabrf:
    return N((abs(Be*exp(-i*Hrfrf(Phip=0)*tp)*Al))^2)
plot(Pabrf)
edit flag offensive delete link more

Comments

you could also use a lambda function

niles gravatar imageniles ( 2011-05-24 19:46:42 -0500 )edit
1

answered 2011-01-05 11:59:35 -0500

Felix Lawrence gravatar image

One solution is to use list_plot instead of plot. plot takes a function, whereas list_plot accepts a list of (x,y) coordinates to plot. list_plot therefore gives you a bit more control and flexibility. If you replace your last three lines by:

Pabrfs= []
for r in srange (-1,1,0.25,include_endpoint=True):
    #(pi)x on Alpha off resonance Full H
    Pabrfs.append((r,Pabrf(Ratio=r).n(digits=5)))

list_plot (Pabrfs)

then it should work. Alternatively you can make this a bit shorter using a list comprehension:

#(pi)x on Alpha off resonance Full H
Pabrfs = [(r,Pabrf(Ratio=r).n()) for r in srange(-1,1,0.25,include_endpoint=True)]

list_plot(Pabrfs)
edit flag offensive delete link more
1

answered 2011-01-05 11:58:06 -0500

jk gravatar image

'unable to simplify to float approximation'

Came here because I have a similar problem

Apparently the plotting command provided by Sage evaluate expressions by coercing to float, i.e. plot calls

float(Pabrf)

This fails on your expression because it's too complicated.

Can someone explain why Sage doesn't evaluate expressions numerically when plotting? What would be a straightforward why to do this?

edit flag offensive delete link more
1

answered 2011-01-06 03:58:03 -0500

mhfrey gravatar image

updated 2011-01-06 05:04:40 -0500

I have been able to get most things working with:

list_plot([(r,Pabrf(Ratio=r).n()) for r in srange(-1, 1, 0.1, include_endpoint=True)], plotjoined=True, frame=true).show(ymin=0,axes_labels=('$\Omega _{0}/\omega _{nut}$',"$ P _{a\Rightarrow b}$"))

I have an additional questions about controlling plot features:

  1. How do I turn off or turn on the lines at the origin(o)?

  2. How do I turn off the minor ticks?

  3. How do I get the \rightarrow to display instead of the \Rightarrow? Only the \Rightarrow displays properly.

Thank you again for your help.

edit flag offensive delete link more

Comments

Not really sure about any of these. Have a look through the plotting documentation http://sagemath.org/doc/reference/sage/plot/plot.html . If you want lots of control over your graphs then it might be worth checking out matplotlib/pylab, which is included with Sage but has separate documentation. Personally I've found sage's plotting functions convenient for quick and dirty figures, but use matplotlib when I want finesse.

Felix Lawrence gravatar imageFelix Lawrence ( 2011-01-06 18:53:47 -0500 )edit
0

answered 2011-01-06 02:27:00 -0500

mhfrey gravatar image

updated 2011-01-06 02:27:27 -0500

Thank you this works for me. I was able to simplify it further to:

list_plot([(r,Pabrf(Ratio=r).n()) for r in srange(-1, 1, 0.1, include_endpoint=True)], plotjoined=True)

edit flag offensive delete link more

Your Answer

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

Add Answer

Question Tools

1 follower

Stats

Asked: 2011-01-05 10:20:55 -0500

Seen: 2,292 times

Last updated: May 24 '11