Hi,
To recreate the phase space for the system described in the above equations, you can use the "ode45" function. This will help you to integrate the above equations over time. To change the "y" units to be "s/K" instead of just "s", I scaled the "s" values by "1/K" in the plotting function. You can refer to the example code below for better understanding:
function phase_space_plot
% Parameters
alpha = 1;
beta = 1;
k = 0.1;
K = 1;
mu = 0.1;
u0 = 0.75;
s0 = 0.01 * K; % Initial condition for s
% Time span
tspan = [0 100]; % Adjust as necessary
% Initial conditions
y0 = [u0, s0];
% Solve the ODE system
[~, Y] = ode45(@(t, y) systemODEs(t, y, k, K, alpha, beta, mu), tspan, y0);
% Extract solutions
u = Y(:, 1);
s = Y(:, 2) / K; % Scale s by K for plotting
% Plotting
figure;
hold on;
% Plotting the trajectory
plot(u, s, 'b', 'LineWidth', 2);
% Generating a grid of points for the direction field
[U, S] = meshgrid(0:0.05:1, 0:0.05:1.2); % Adjust grid density and range as needed
Udot = -k.*U.*(1-U).*(1-S/K);
Sdot = -k.*S.*U.*(1-S/K) + alpha.*(1+beta.*U).*S.*(1-S/K) - mu.*S;
% Normalizing vectors for uniform arrow size
Vnorm = sqrt(Udot.^2 + Sdot.^2);
Udot = Udot ./ Vnorm;
Sdot = Sdot ./ Vnorm;
% Adding the directional arrows to the plot
quiver(U, S, Udot, Sdot, 0.3, 'k'); % Adjust scaling factor as necessary
xlabel('u');
ylabel('s/K');
title('Combined Phase Space and Directional Plot');
grid on;
hold off;
end
function dydt = systemODEs(~, y, k, K, alpha, beta, mu)
u = y(1);
s = y(2);
% System of equations
du_dt = -k*u*(1-u)*(1-s/K);
ds_dt = -k*s*u*(1-s/K) + alpha*(1+beta*u)*s*(1-s/K) - mu*s;
dydt = [du_dt; ds_dt];
end
phase_space_plot
I defined a function "phase_space_plot" that solves the system of the given differential equations using "ode45" and plots the phase space trajectory with directional arrows, where "u" is on the x-axis and "s/K" on the y-axis. You can adjust the code to your requirements.
For more information on the "ode45" function, you can refer to below documentation:
