Numerical Integration of System + x/y Plot

asked 1 year ago

Jack Zuffante gravatar image

updated 1 year ago

How do I numerically integrate such a system for a certain range of t and then make an x/y plot of the curve?

x''(t)=A(t,x,y,x',y')

y''(t)=B(t,x,y,x',y')

A and B are functions t, x, y, dx/dt, and dy/dt.

Preview: (hide)

Comments

I'll post my specific code tomorrow.

Jack Zuffante gravatar imageJack Zuffante ( 1 year ago )

Well first off all theres the fact that just typing the double quotes there doesn't work. But my main problem is that in the example achrzesz posted, the dependent variable (which it looks like it is also t in his case) does not appear in the equations, so I guess I need some way to define t as a variable and x and y as functions of t in a way that would work with the example. My instinct is to do (using simpler equations for now):

t=var('t')

x=function('x')(t)

y=function('y')(t)

u=diff(x,t)

v=diff(y,t)

diff(u,t)=y*u

diff(v,t)=x*v

T=ode_solver()

f = lambda t,g:[u,v,diff(u,t),diff(v,t)]

T.function=f

T.ode_solve(g_0=[t0,x0,y0,dx0,dy0],t_span=[0,20],num_points=1000)

I would just ...(more)

Jack Zuffante gravatar imageJack Zuffante ( 1 year ago )

Yeah, I'm not sure how to make it work

Jack Zuffante gravatar imageJack Zuffante ( 1 year ago )

The first part of that comment is showing the procedure that you need to do on the paper (reducing the second-order system to a system of the first-order).

The following part might be frustrating. After having your first-order equations as

y1'=y3 y2'=y4 y3'=-y2-1 y4'=-y1

you write in the code:

f = lambda t,y:[y[2],y[3],-y[1]-1,-y[0]]

This is because the indices start with zero in SageMath. So y1 is y[0], y2 is y[1], etc. (You think of y1, y2, etc as element of the array y)

The rest is the code. You will run this:

T=ode_solver()
f = lambda t,y:[y[2],y[3],-y[1]-1,-y[0]]
T.function=f
T.ode_solve(y_0=[1,1,0,0],t_span=[0,20],num_points=1000)
f = T.interpolate_solution()
plot(f,0,5).show()

I hope that helps.

tolga gravatar imagetolga ( 1 year ago )