%% Gnarly ODE
% BJ Furman
% 24SEP2015
clear all
close all
syms s t y(t) Y         % create symbolic variables
dy(t)=diff(y(t),t)      % create dy/dt
d2y(t)=diff(dy(t),t)    % create d(dy/dy)/dt
f(t)=1+5*dirac(t-5)     % create the forcing function (the LHS of the ODE)
eqn = d2y(t) + 2*dy(t) + 10*y(t) - f(t)  % the complete equation: LHS - RHS = 0
Leqn = laplace(eqn,t,s)     % take the Laplace transform of the complete equation
Leqn_N = subs(Leqn,{'y(0)','D(y)(0)'},{1,2})  % substitute numeric values for the y(0) and dy/dt(0)
Leqn_Y = subs(Leqn_N,'laplace(y(t),t,s)',Y)   % substitute Y for the operational description of the Laplace transform
Sol = solve(Leqn_Y,Y)       % solve for Y(s)
y = ilaplace(Sol,s,t)       % take the inverse Laplace transform of Y(s)
y_num = matlabFunction(y)   % turn the symbolic y(t) into a numeric function
t=linspace(0,10,200);       % create a time vector to plot the function
plot(t,y_num(t))            % plot the numeric function over the time vector
grid