## Linear algebra approaches to solving systems of constant coefficient ODEs

const_coeff_coupled

## Linear algebra approaches to solving systems of constant coefficient ODEs

John Kitchin

Today we consider how to solve a system of first order, constant coefficient ordinary differential equations using linear algebra. These equations could be solved numerically, but in this case there are analytical solutions that can be derived. The equations we will solve are:  We can express this set of equations in matrix form as: The general solution to this set of equations is where is a diagonal matrix of the eigenvalues of the constant coefficient matrix, is a matrix of eigenvectors where the column corresponds to the eigenvector of the eigenvalue, and is a matrix determined by the initial conditions.

In this example, we evaluate the solution using linear algebra. The initial conditions we will consider are and .

close all
A = [-0.02 0.02;
0.02 -0.02];


## get the eigenvalues and eigenvectors of the A matrix

[V,D] = eig(A)

V =

0.7071    0.7071
-0.7071    0.7071

D =

-0.0400         0
0         0



V is a matrix of the eigenvectors, with each column being an eigenvector. Note that the eigenvectors are normalized so their length equals one. D is a diagonal matrix of eigenvalues.

## Compute the matrix

V*c = Y0

We can solve this equation easily with the \ operator. As written, c will be a vector, but we need it in the form of a diagonal matrix, so we wrap the solution in the diag command.

Y0 = [0; 150]; % initial conditions
c = diag(V\Y0)

c =

-106.0660         0
0  106.0660



## Constructing the solution

We will create a vector of time values, and stack them for each solution, and .

t = linspace(0,100);
T = [t;t];

y = V*c*exp(D*T);
plot(t,y)
xlabel('t')
ylabel('y')
legend 'y_1' 'y_2' ## Verifying the solution

One way we can verify the solution is to analytically compute the derivatives from the differential equations, and then compare those derivatives to numerically evaluated derivatives from the solutions.

yp = A*y; % analytical derivatives

% numerically evaluated derivatves
yp1 = cmu.der.derc(t,y(1,:));
yp2 = cmu.der.derc(t,y(2,:));

plot(t,yp,t,yp1,'g--',t,yp2,'r--')
legend 'y''_1 analytical' 'y''_2 analytical' 'y''_1 numerical' 'y''_2 numerical'
xlabel('t')
ylabel('dy/dt') What you can see here is that the analytical derivatives are visually equivalent to the numerical derivatives of the solution we derived.

'done'
% categories: ODEs, linear algebra
% tags: math

ans =

done