Plotting time series of ODE for several initial values in subplot
3 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
I am checking the stability of Equilibrium points for a system of 3 ODE numerically via ODE54 and need to plot the result for several initial conditions for example 4 initial condition for each equilibrium point and plot the result in subplots. Is it possible ? is it possible to change the initial value in the title accordingly?
Function
function [dxdt] = Subsystem_NonDim_LFR(t,x,s,epsilon,omega,beta,g12,gamma12)
%UNTITLED2 Summary of this function goes here
% Detailed explanation goes here
dxdt(1,:)=x(1)*(1-x(1)-x(2)-x(3))+s;
dxdt(2,:)=x(2)*(-epsilon(1)*(1+omega(1)*x(2)+omega(2)*x(3))-g12*x(3)+beta(1)*x(1));
dxdt(3,:)=x(3)*(-epsilon(2)*(1+omega(3)*x(2)+omega(4)*x(3))+gamma12*x(2)+beta(2)*x(1));
end
Script:
s=0
omega=[1,2,3,4]
g12=1
gamma12=2
beta=[1,2]
epsilon=[4,3]
s=0
x0=[0, 200 ,100];
t=linspace(0,30,30000);
[t,x]=ode45(@(t,x) Subsystem_NonDim_LFR(t,x,s,epsilon,omega,beta,g12,gamma12),t,x0);
figure
p=plot(t, x(:,1),'r--',t, x(:,2),'b:',t,x(:,3),'g-.')
p(1).LineWidth = 2.5;
p(2).LineWidth = 3;
p(3).LineWidth = 3.5;
xlabel('Time (day)')
ylabel('Populations (num/ml)')
legend('Bacteria','protozoa','Daphnia')
grid
ylim([-5,10 ])
title(' $$x_{0}=[0.00000001,200,100]\ ,\beta_{1}<\epsilon_{1} and\ \beta_{2}<\epsilon_{2}$$','interpreter','latex')
2 Kommentare
Akzeptierte Antwort
darova
am 29 Jan. 2020
SOmething like that i think
x0=[0, 200 ,100
3 2 1
3 1 2
31 2 1]; % some initial conditions (random)
t=linspace(0,30,30000);
for i = 1:4
[t,x]=ode45(@(t,x) Subsystem_NonDim_LFR(t,x,s,epsilon,omega,beta,g12,gamma12),t,x0(i,:));
subplot(2,2,i)
p=plot(t, x(:,1),'r--',t, x(:,2),'b:',t,x(:,3),'g-.')
p(1).LineWidth = 2.5;
p(2).LineWidth = 3;
p(3).LineWidth = 3.5;
xlabel('Time (day)')
ylabel('Populations (num/ml)')
legend('Bacteria','protozoa','Daphnia')
grid
ylim([-5,10 ])
str = sprintf('$$x_{0}=[%d,%d,%d]',x0(i,:));
title([str '\ ,\beta_{1}<\epsilon_{1} and\ \beta_{2}<\epsilon_{2}$$'],'interpreter','latex')
end
Weitere Antworten (0)
Siehe auch
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!