plot the phase plane for the SIR model

I need to plot the phase plane for the SIR model (looks like the attached image). Anyone can tell me how to do it please? (assume the values of the variables are given including s(0) and i(0)) I searched out that I should use the meshgrid and quiver function. I have been trying so hard but still can't figure out how these two functions work, can someone help me please?

it means a lot to me, thanks.

Answer 1

[Star Strider]

You can select the phase plane output, or you can simply substitute the solved values for ‘t’ and ‘v’ (in my code here) back into the original ODE function to get the derivatives.

This code calls the 'odephas2' output function:

% MAPPING: s = v(1),    i = v(2)
ode = @(t, v, beta, mu, gamma) [-beta.*v(1).*v(2) + mu - mu.*v(1);  beta.*v(1).*v(2) - (gamma + mu).*v(2)]; 

beta  = rand;                                               % Choose The Correct Values
mu    = rand;
gamma = rand;

v_init = rand(2,1);
tspan  = linspace(1, 20, 50);

opts = odeset( 'OutputFcn','odephas2');                     % Set To Plot Phase Plane

[t, v] = ode45(@(t,v) ode(t, v, beta, mu, gamma), tspan, v_init, opts);

figure(2)
plot(v(:,1), v(:,2))                                        % Plot The Phase Plane Manually
grid

for k1 = 1:length(t)
    derivs(:,k1) = ode(t(k1), v(k1,:), beta, mu, gamma);      % Calculate The Derivatives
end

figure(3)
plot(derivs(1,:), derivs(2,:))                              % Plot Derivatives
grid