%% Solving Linear Systems with the Laplace Transform
%% Kurt Bryan and SIMIODE
% This worksheet illustrates using the Laplace transform to solve linear system of ODEs (and also just using Matlab's dsolve command).
% Consider the linear constant-coefficient nonhomogeneous system

syms x1(t) x2(t); %Declare functions
de1 = diff(x1(t),t) == -5*x1(t) + 6*x2(t)-2*sin(t)+10*cos(t)-6*exp(t)
de2 = diff(x2(t),t) == -3*x1(t) + x2(t)+6*cos(t)

% with initial data x1(0) = 2 and x2(0) = 1.
% To solve, Laplace transform both sides of both equations

syms s;
de1lap = laplace(de1, t, s);
de2lap = laplace(de2, t, s);

% Substitute in the initial conditions, and (for convenience) let X1 = laplace(x1(t),t,s) and X2 = laplace(x2(t),t,s):

syms X1 X2;
de1lap2 = subs(de1lap,[x1(0),x2(0),laplace(x1(t),t,s),laplace(x2(t),t,s)],[2,1,X1,X2])
de2lap2 = subs(de2lap,[x1(0),x2(0),laplace(x1(t),t,s),laplace(x2(t),t,s)],[2,1,X1,X2])

% Solve for the transforms X1 and X2

Xsols = solve([de1lap2,de2lap2],[X1,X2])

% Let X1sol and X2sol denote the transforms

X1sol = Xsols.X1
X2sol = Xsols.X2

% Inverse transform to find the solutions

x1sol(t) = ilaplace(X1sol,s,t)
x2sol(t) = ilaplace(X2sol,s,t)

% A quick check using the dsolve command:

sol = dsolve([de1, de2],[x1(0)==2,x2(0)==1]);
sol.x1
sol.x2

%
%