% Physics of the coupled oscillators (implementation of the ODEs)
function dydt = oscillators(y, t)

% Parameters m1, m2, k1, k2, k3, L1, L2, L3 are globals from driver program.
global m1 m2 k1 k2 k3 L1 L2 L3

% Just to make the notation clearer:
x1 = y(1);
x2 = y(2);
v1 = y(3);
v2 = y(4);

% Forces on mass 1 and mass 2:
F1 = -k1*(x1-L1) + k3*((x2-x1) - L3);
F2 = -k2*(x2-L2) - k3*((x2-x1) - L3);

% Resulting equations of motion:
dydt = [v1, v2, F1/m1, F2/m2];

