function dydt = osc_td4(t,y,varargin)
% harmonic oscillator system
k=varargin(1);
m=varargin(2);
%%  we could use
% dydt(1) = y(2);
% dydt(2) = -(k/m)*y(1);
% dydt=dydt(:)  % this makes dydt a column vector as required by odexx
%% or we could just use 
dydt=[y(2); -(k/m)*y(1)];
end

