function dYdt = kepler_ode(t,Y)
% assume 'Y' is a vector of form y = [ x ;  y ; vx ; vy ]

r = [Y(1);Y(2)];
v = [Y(3);Y(4)];

% use the following three variables from the main worksheet
% must declare them 'global' here and there
global mStar mPlanet G 

% assumes the sun is located at r=(0 0)
F = - mStar * mPlanet * G / norm(r)^3 .* r;
a = F ./ mPlanet;

dYdt = [ v(1) ; v(2) ; a(1) ; a(2)];  % must column vector 


        
