% Aph161 - PSet5, Spring 2009
% ODE Solver Tutorial - Harmonic Oscillator
function dXdt = HO_ODE(T,X,w_0)
% This function will define the set of ODEs for the harmonic oscillator.
% This is the syntax that is necessary.
% dXdt holds the values of the ODEs for a given timepoint as the ode45
% solver evaluates the ODEs numerically (dXdt is ONLY defined inside of the
% function and CANNOT be called in the main script)
% HO_ODE is the name of the function. Note: YOU MUST NAME THIS function
% after what you put after the equals sign OR ELSE IT WILL NOT WORK (i.e.
% in this case, save the function as 'HO_ODE.m' without the quotation
% marks)
% The parameters inside HO_ODE refer to time points, the position/velocity 
% vector, and the parameter w_0 for the harmonic oscillator DEs, 
% respectively.

% First, initialize the dXdt vector as a 2x1 vector of zeros (Note: the
% vector MUST be a column vector, so saying zeros(1,2) is not allowed).
% This is just good coding practice.
 dXdt = zeros(2,1);
% Now, define the position and velocity ODEs
 dXdt(1) = X(2);     % Derivative of position is velocity!
 dXdt(2) = -w_0^2*X(1); % Derivative of velocity is acceleration!
