Ask Your Question

# How to plot field lines in sage

I want to plot the field lines of an electric field, take for example that of an electric dipole. How can I do this in sage?

Edit: Perhaps my question was formulated too much in physicist slang... Mathematically I just want to plot the integral curves of a vector field for example: let $\vec{r}(x,y) = (x,y)$, $\vec{E_1} = \frac{\vec{r}}{|r|^{3}}$ and $\vec{E_2}(x,y) = \vec{E_1}(x-1,y)$

Then I want to plot the field lines of $\vec{E} = \vec{E_1} - \vec{E_2}$

The field line plot should look similar to this one:

I think in mathematica this is called a stream plot: http://reference.wolfram.com/mathemat...

edit retag close merge delete

## 5 answers

Sort by ยป oldest newest most voted

I don't think there's a way to automate this totally and avoid putting in all of the initial conditions. Here is another version.

var('x y')
E1=vector([x,3*y])/(sqrt(x^2+y^2))^(3/2)
E2=E1.subs(x=x-1)
F=(E1-E2)/norm(E1-E2)
int_conds=[[0.05,0.001],[0.05,0.005],[0.05,-0.001],[0.05,-0.005]]
solns=[desolve_odeint(F,init,srange(0,1.5,0.05),[x,y]) for init in int_conds]
nhalf=len(solns[0])/3
a=[arrow(solns[i][nhalf],solns[i][nhalf+1]) for i in range(0,len(solns))]
b=map(lambda x: line(x,thickness=2),solns)
show(sum(a)+sum(b)+point((0,0),size=50)+point((1,0),size=50)+plot_vector_field(F,(x,-0.1,1.1),(y,-.5,.5)),axes=false)


more

Thanks for the edit. See if something like this does what you want.

var('x y')
E1=vector([x,y])/(sqrt(x^2+y^2))^3
E2=E1.subs(x=x-1)
F=(E1-E2)

tstop=0.4
step=0.005
int_conds=[[0.5,0.1*k] for k in [1..3,-3..-1]]
solns=sum([line(desolve_odeint(F,init,srange(0,tstop,step),[x,y]),thickness=2) for init in int_conds])
solns2=sum([line(desolve_odeint(F,init,srange(0,-1*tstop,-1*step),[x,y]),thickness=2) for init in int_conds])

F=(E1-E2)/norm(E1-E2)
solns+solns2+plot_vector_field(F,(x,-0.2,1.2),(y,-1,1))

more

## Comments

Thanks, that's almost what I want. The missing point is that I want not only one or three lines but something like I added in the original post above.

( 2012-09-29 03:28:25 -0600 )edit

Sorry I had two typos in my Edit, which I fixed just now, see above. You should not normalize F! Perhaps you could include those changes in your examples. Thanks.

( 2012-09-29 05:21:04 -0600 )edit

I've updated my code a good bit to include your edits and some solving backwards and forwards from the initial condition.

( 2012-09-30 04:10:01 -0600 )edit

In the attached file one can see the maxima version (code in my previous answer)

sdfplot.pdf

more

I'm not sure whether you already have equations to use or not. How you plot these all depends, of course, on how the equations are given. But, here is something. Using information in section D.73.3 of this webpage on quantum mechanics, we need to plot ${(x^2+z^2)}^{(3/2)}=c x^2$ for various values of $c$. The following code produces a plot.

var('x z')
sum([implicit_plot((x^2+z^2)^(3/2)==c*x^2,(x,-5,5),(z,-5,5)) for c in [1,2,3,4]])


more

## Comments

That are the equipotential surfaces and not the field lines. You can assume that I know an analytic expression for the electric field strength at every point.

( 2012-09-28 21:40:28 -0600 )edit

Can you post the expression?

( 2012-09-29 01:37:47 -0600 )edit

See my edit

( 2012-09-29 02:29:04 -0600 )edit

People interested in this question will want to perhaps help out in making our ticket for streamline plots a reality.

more

## Your Answer

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

Add Answer

## Stats

Asked: 2012-09-28 11:23:45 -0600

Seen: 1,936 times

Last updated: Sep 30 '12