% rk4pen.m
% solves driven-damped pendulum 
% using Runge-Kutta 4nd order

function [y1,y2] = rk4pen(t0,y0,dt,nSteps)

y1 = zeros(1,nSteps);
y2 = zeros(1,nSteps);

t = t0;
y = y0;

fdt = 1.0/6.0*dt;

for k=1:nSteps
	k1 = rkFunc(t,y);
	k2 = rkFunc(t+0.5*dt,y+0.5*dt*k1);
        k3 = rkFunc(t+0.5*dt,y+0.5*dt*k2);
        k4 = rkFunc(t+dt,y+dt*k3);
	y = y + fdt*(k1+2*k2+2*k3+k4);
        t = t+dt;

	tTrace(k) = t;
	y1(k) = y(1);
        y2(k) = y(2);
end

% System of ODEs 

function yPrime = rkFunc(t,y)

q = 2;
g = 1.0;
om = 2/3;

yPrime = zeros(1,length(y));

yPrime(1) =  y(2);
yPrime(2) = -y(2)/q-sin(y(1))+g*cos(om*t);
