Solution of Practice Problem for Exam 3

clearvars

Consider the ODE system

syms x y
f = 4*x+y+x*y;
g = x+4*y+y^2;

Find all stationary points:

Find

with

[xs,ys] = solve([f;g]==[0;0],[x;y])

xs =

ys =

m = length(xs) % number of solutions we found

m = 3

We found three stationary points (0,0), (3,-3), (-5,-5).

Find Jacobian matrix J:

J = jacobian([f;g],[x;y])

J =

For each stationary point : find matrix

Find eigenvalues and eigenvectors

For LINEARIZED PROBLEM: find type and stability of stationary point, sketch phase portrait

Result:

(0,0) is nodal source (unstable, repelling)
(3,-3) is saddle point (unstable, not repelling)
(-5,-5) is nodal sink (stable, attracting)

for k=1:m

stationary_point = [xs(k),ys(k)]

A = subs(J,{x,y},{xs(k),ys(k)}) % evaluate Jacobian at stationary point

[eigenvectors,eigenvalues] = eig(A)

figure

phaseportrait(A);

title(sprintf('linearized problem for stationary point (%g,%g)',stationary_point))

end

stationary_point =

A =

eigenvectors =

eigenvalues =

stationary_point =

A =

eigenvectors =

eigenvalues =

stationary_point =

A =

eigenvectors =

eigenvalues =

For nonlinear problem: What can you say about type and stability of the stationary points?

What does the phase portrait for the nonlinear problem look like?

For the nonlinear problem we will have the same types and stabilities as for the linearized problem (only in the case of a center things can be different for the nonlinear problem).

We will have the same tangent directions of the trajectories at the stationary point.

We use

in place of

and define the function

as follows:

f = @(t,y) [ 4*y(1)+y(2)+y(1)*y(2) ; y(1)+4*y(2)+y(2)^2 ];

For many initial points

: use ode45 and

solve ODE with initial condition for t going from 0 to 5, plot trajectory in phase plane
solve ODE with initial condition for t going from 0 to , plot trajectory in phase plane

for a1 = -6:5

for a2 = -6:1

[ts,ys] = ode45(f,[0,5],[a1;a2]); plot(ys(:,1),ys(:,2),'c'); hold on

[ts,ys] = ode45(f,[0,-5],[a1;a2]); plot(ys(:,1),ys(:,2),'c')

end

vectorfield(f,-6:5,-6:1);

axis([-6 5 -6 1]); hold off

We can see

nodal source at (0,0)
saddle at (3,-3)
nodal sink at (-5,-5)

and the phase portrait close to each stationary points corresponds to the phase portrait of the linearized problem.

Function for drawing phase portrait

function phaseportrait(A)
syms y(t) [2 1]
syms C1 C2
sol = dsolve(diff(y)==A*y);             % solve ODE 
ysol = subs(y(t),sol);
f = @(t,y) A*y;                         % define function f(t,y) for vectorfield
vectorfield(f,-5:5,-5:5); hold on
for c1=[-1 0 1]
  for c2=[-1 0 1]
    ys = subs(ysol,{C1,C2},{c1,c2});    % use C1=c1, C2=c2
    fplot(ys(1),ys(2),[-5 5]);          % plot trajectory for t=-5...5
  end
end
hold off; axis equal; axis([-5 5 -5 5]); xlabel('y_1'); ylabel('y_2')
end