Computing an ODE in Matlab

Question

Given a system of the form y' = A*y(t) with solution y(t) = e^(tA)*y(0), where e^A is the matrix exponential (i.e. sum from n=0 to infinity of A^n/n!), how would I use matlab to compute the solution given the values of matrix A and the initial values for y?

That is, given A = [-2.1, 1.6; -3.1, 2.6], y(0) = [1;2], how would I solve for y(t) = [y1; y2] on t = [0:5] in matlab?

I try to use something like

t = 0:5
[y1; y2] = expm(A.*t).*[1;2]

and I'm finding errors in computing the multiplication due to dimensions not agreeing.

Answer 1

Please note that matrix exponential is defined for square matrices. Your attempt to multiply the attenuation coefs with the time vector doesn't give you what you'd want (which should be a 3D matrix that should be exponentiated slice by slice).

One of the simple ways would be this:

A = [-2.1, 1.6; -3.1, 2.6];

t = 0:5;
n = numel(t);  %'number of samples'

y = NaN(2, n);
y(:,1) = [1;2];


for k =2:n
    y(:,k) = expm(t(k)*A) * y(:,1);
end;

figure();
plot(t, y(1,:), t, y(2,:));

Please note that in MATLAB array are indexed from 1.