how to plot solution of differential system

asked 2012-09-21 23:24:40 +0200

5taras
1 ●1 ●1 ●4

updated 2012-09-22 08:49:18 +0200

calc314
4191 ●22 ●48 ●113

Hi, everyone! I got a question. I have a system of three Newton's equations

x''=-x/(x*x+y*y+z*z)^3/2-F(t)
y''=-x/(x*x+y*y+z*z)^3/2
z''=-x/(x*x+y*y+z*z)^3/2

I use ode_solver() to solve numerically the system by using y[0] for x, y[1] for y, y[2] for z like in example in the manual without using jacobian:

j_1= lambda t,y: [ y[3], y[4], y[5],  -y[0]/(y[0]*y[0]+y[1]*y[1]+y[2]*y[2])^(1.5)-F*cos(w*t), -y[1]/(y[0]*y[0]+y[1]*y[1]+y[2]*y[2])^(1.5), -y[2]/(y[0]*y[0]+y[1]*y[1]+y[2]*y[2])^(1.5)]

How may I plot the curve (x(t), y(t)) or (x(t), z(t))? T.plot_solution(i=0, filename=outfile) allows to plot only (x(t), t), (y(t), t) or (z(t), t) Thank you

edit retag flag offensive close merge delete

add a comment

2 Answers

Sort by » oldest newest most voted

answered 2012-09-22 08:48:45 +0200

calc314
4191 ●22 ●48 ●113

You can get access to the list of solution points using the .solution option.

a=T.solution
print a

The code above shows the solution points and you can see the data structure. To plot what you want, try these:

xysoln=[[soln[1][0],soln[1][1]] for soln in a]
list_plot(xysoln)

yzsoln=[[soln[1][1],soln[1][2]] for soln in a]
list_plot(yzsoln)

xzsoln=[[soln[1][0],soln[1][2]] for soln in a]
list_plot(xzsoln)

xyzsoln=[[soln[1][0],soln[1][1],soln[1][2]] for soln in a]
list_plot(xyzsoln)