function ode_5
% ode_5.m  solve  d^2x/dt^2 + 5 dx/dt - 4 x(t) = sin(10 t)
% initial conditions x(0)=0, x'(0)=0
  t = 0:0.001:3.142; % time scale
  initial_x = 0;
  initial_dxdt = 0;
  [t,x] = ode45(@rhs, t, [initial_x, initial_dxdt]);
  figure(1)
  plot(t,x(:,1));
  title('x(t)');
  xlabel('t');
  ylabel('x');

  figure(2)
  plot(t,x(:,2));
  title('dx/dt');
  xlabel('t');
  ylabel('dx');

  figure(3)
  n = size(t)
  for i=1:n
    ddx(i) = -5*x(i,2) + 4*x(i,1) + sin(10*t(i));
  end % i
  plot(t,ddx);
  title('dxx/dtt');
  xlabel('t');
  ylabel('dxx');

  figure(4)
  for i=1:n
    s10t(i) = sin(10*t(i));
  end % i
  plot(t,s10t);
  title('sin(10*t)');
  xlabel('t');
  ylabel('s10t');

  function dxdt=rhs(t,x) % move everything past =
    dxdt_1 = x(2);
    dxdt_2 = -5*x(2) + 4*x(1) + sin(10*t);
    dxdt   = [dxdt_1; dxdt_2];
  end % rhs
end % ode_5
