Warning: Failure at t=-1.769748e+00. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (3.552714e-15) at time t.
Solution Exam 3 (MATH246, Spring 2024)
Problem 1(a),(b)
syms x y
f = 1-2*y;
g = x^2-y^2;
Find stationary points:
[xs,ys] = solve([f;g]==[0;0],[x;y])
m = length(xs) % number of solutions we found
m = 2
We found two stationary points: 
and 
.
and .J = jacobian([f;g],[x;y])
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
Result:
- is saddle (unstable, not repelling)
- is spiral sink (stable, attracting)
Problem 1(c)
The stationary points for the nonlinear system have the same type and stability as in 1(b) (only for the center we may have differences for the nonlineary ODE system).
f = @(t,y) [ 1-2*y(2) ; y(1)^2-y(2)^2 ]; % we now use y1,y2 instead of x,y
for a1 = -1:.25:1
for a2 = 0:.125: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
end
vectorfield(f,-1:.25:1,0:.125:1);
hold off; axis equal; axis([-1 1 0 1]);
Problem 2
A = sym([6 2; 5 9])
[eigenvectors,eigenvalues] = eig(A)
We have two different positive eigenvalues, so this is a nodal source (unstable, repelling).
phaseportrait(A)
Problem 3
A = sym([0 -2; 8 0])
eigenvalues = eig(A)
We have two imaginary eigenvalues, so this is a center (stable, not attracting), counterclockwise.
phaseportrait(A)
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