function [dX] = BeamDynamics( X, T )
% BeamDynamics(X, T);
% ECE 463 Lecture #8
% m = 1 kg
% J = 0.2 kg m^2

r = X(1);
q = X(2);
dr = X(3);
dq = X(4);
g = 9.8;
m = 1.0;


M = [1.4*m, 0; 0, 0.2 + m*r*r];
B1 = m*r*dq*dq - m*g*sin(q);
B2 = T - 2*m*r*dr*dq - m*g*r*cos(q);

ddX = inv(M)*[B1; B2];
dX = [dr; dq; ddX];

end