% getF.m
%
% Evaluate state derivatives for Central Force Motion Example
%
% X_dot = F(X)
% X = [x-position; y-position; x-velocity; y-velocity]
%
% Assume a gravitational force given by F/m = -mu/r^2
% x component: Fx/m = -mu/r^2*cos(theta) = -mu*x/r^3
% y component: Fy/m = -mu/r^2*sin(theta) = -mu*y/r^3 
% t: current time given by the integrator
% X: current state vector given by the integrator
% dX: returned vector of state derivatives
%

function dX = getF(t, X)

mu = 3.;

% unpack the state vector (follow the definition)
x     = X(1);
y     = X(2);
x_dot = X(3);
y_dot = X(4);
    
r     = sqrt(x^2 + y^2);

x_dot_dot = -(mu/r^3)*x; % x component of attractive force per unit mass
y_dot_dot = -(mu/r^3)*y; % y component of attractive force per unit mass

% pack up the vector of state derivatives (follow the definition)
dX = [ x_dot;
       y_dot;
       x_dot_dot;
       y_dot_dot 
     ];