Solving the ODE system

We are given a matrix

and want to find

such that

Example: Solve the initial value problem

with initial condition

Key idea: use the guess

We need to find a number r and a vector

such that

We obtain the equation

. With the identity matrix

we can write this as

Since

the matrix

must be singular:

This gives the characteristic equation

where

and

Solving this quadratic equation gives the eigenvalues

A = sym([1 2; 4 3])

A =

eig(A) % find eigenvalues of matrix A

ans =

For each eigenvalue r we need to find eigenvectors

satisfying

For the eigenvalue

we obtain the linear system

Note that the matrix is singular, so we can find nonzero solution vectors. We obtain

, so

is an eigenvector for

(multiplying this vector by any nonzero number would also be an eigenvector).

For the eigenvalue

we obtain the linear system

Note that the matrix is singular, so we can find nonzero solution vectors. We obtain

, so

is an eigenvector for

(multiplying this vector by any nonzero number would also be an eigenvector).

Find eigenvalues and eigenvectors with Matlab: The command [V,D]=eig(A) gives

a matrix containing the eigenvectors
a diagonal matrix containing the eigenvalues

[V,D] = eig(A) % find matrix V of eigenvectors, diagonal matrix D with eigenvalues

V =

D =

Matlab gives

the eigenvalue with eigenvector
the eigenvalue with eigenvector

The general solution of the ODE

is now given by

Here we obtain

For the initial condition we plug in

and obtain the equation

We can solve this linear system in Matlab as follows

c = V\[1;0] % solve linear system V c = [1;0]

c =

Hence the solution of the initial value problem is

Solving an ODE system with dsolve

We have to tell Matlab that y(t) is an array of size [2,1]:

syms y(t) [2 1]             % define vector valued function y(t) = [y1(t);y2(t)]
syms C1 C2

Find general solution

sol = dsolve( diff(y)==[1 2;4 3]*y);

ysol = subs(y(t),sol) % in y(t) substitute the values from sol

ysol =

This is the general solution, with

. We obtain

by substituting

Y1 = subs(ysol,{C1,C2},{1,0})

Y1 =

We obtain

by substituting

Y2 = subs(ysol,{C1,C2},{0,1})

Y2 =

Plot vector field and trajectories in the phase plane

You need to download the file vectorfield.m and put it in the directory where your mlx files are.

Plot the curves

and

for

f = @(t,y) A*y; % define function f(t,y) for vectorfield

vectorfield(f,-8:8,-8:8); hold on

fplot(Y1(1),Y1(2),[-3,3],'b'); % plot trajectory for Y1 in blue

fplot(-Y1(1),-Y1(2),[-3,3],'b:'); % plot trajectory for -Y1 in blue, dotted

fplot(Y2(1),Y2(2),[-3,3],'k'); % plot trajectory for Y2 in black

fplot(-Y2(1),-Y2(2),[-3,3],'k:'); % plot trajectory for -Y2 in black, dotted

axis equal; axis([-8 8 -8 8]); xlabel('y_1'); ylabel('y_2')

title('four trajectories in the phase plane')

Solve initial value problem

We have the initial condition

sol = dsolve( diff(y)==[1 2;4 3]*y, y(0)==[1;0] ) % solve initial value problem

sol = struct with fields:
    y2: (2*exp(5*t))/3 - (2*exp(-t))/3
    y1: (2*exp(-t))/3 + exp(5*t)/3

ysol = subs(y(t),sol) % in y(t) substitute the values from sol

ysol =

Plot vector field and trajectory in the phase plane

You need to download the file vectorfield.m and put it in the directory where your mlx files are.

We plot the curve

in the

plane for

. The inital point

for

is marked with a green dot.

f = @(t,y) A*y; % define function f(t,y) for vectorfield

vectorfield(f,-8:8,-8:8); hold on

plot(1,0,'g.','MarkerSize',20) % mark initial point as green dot

fplot(ysol(1),ysol(2),[-3,3],'g'); % plot trajectory for ysol as green curve

hold off; axis equal; axis([-8 8 -8 8]);

xlabel('y_1'); ylabel('y_2')

title('solution of initial value problem (green curve)')