# Writing a system of ODEs in vector form

I have the following system:

$$\frac{dx}{dt} = p_x \\ \frac{dy}{dt} = p_y \\ \frac{dp_x}{dt}=-\frac{\partial V}{\partial x} \\ \frac{dp_y}{dt}=-\frac{\partial V}{\partial y}$$

,where $$V(x,y) = \frac{1}{2} (\omega_x x^2 + \omega_y y^2)$$, and $$\omega_x = 1, \omega_y = 2$$. I would like to write this system in vector form $$\frac{d \vec{u}}{dt} = \vec{f} (t, \vec{u})$$, where $$\vec{u} = \{x,y, p_x, p_y\}$$, in order to solve it with a RK2 method, but how do I define the function $$\vec{f} (t, \vec{u})$$ and how do I write the system as a vector?

The general solution can be obtained by e.g.:

u[t_] = {x[t], y[t], px[t], py[t]};
V = 1/2 (wx^2 x[t]^2 + wy y[t]^2) /. {wx -> 1, wy -> 2};
eq = {D[u[t], t] == {u[t][[3]], u[t][[4]], -D[V, x[t]], -D[V, y[t]]}};
DSolve[eq, u[t], t]


However if you want a specific solution you must define initial conditions u[0]. E.g.:

u[t_] = {x[t], y[t], px[t], py[t]};
V = 1/2 (wx^2 x[t]^2 + wy y[t]^2) /. {wx -> 1, wy -> 2};
eq = {D[u[t], t] == {u[t][[3]], u[t][[4]], -D[V, x[t]], -D[V, y[t]]},
u[0] == {1, 1, 0, 0}};
sol = DSolve[eq, u[t], t][[1]]

ParametricPlot[Evaluate[{x[t], y[t]} /. sol], {t, 0, 10}]


• I see, thank you! Can the vector eq be used to solve the system with a numerical method, like Runge-Kutta 2 or the trapezium rule? Jan 12 at 16:14
• Of course. Instead of "DSolve" you would then need "NDSolve" like e.g.: sol = NDSolve[eq, u[t], {t, 0, 10}] and the option "Method" of "NDSolve" Jan 12 at 16:25