How to solve non linear system of ODE about trophic levels?

2 views (last 30 days)
Hi,
I have to solve this system of ODE in Matlab, describing nitrogen exchanges between different trophic levels.
The problem is that I am not satisfied by the plot, I expected something different according to Lotka_Volterra equations.
I Include the script, the function, the plot and the text of the exercise.
Thaks for any advice,
C.
Cattura.PNG
global r K p1 p2 rho1 rho2 rho3 rho4 mu1 mu2 mu3 mu4 l1 l4 k50 F05y U F05t
r=0.1;
K=20000;
p1=0.1;
p2=0.05;
rho1=0.6;
rho2=0.6;
rho3=0.6;
rho4=0.6;
mu1=0.2;
mu2=0.4;
mu3=0.5;
mu4=0.6;
l1=0.05;
l4=0.1;
k50=0.3;
U=13*10.^12;
F05t=(linspace(0,20,200)); % external time dependent input for the fifth equation
F05y=zeros(1,200);
F05y(1:end)=U;
Y0=[200, 100, 50, 50, 500]; %condizioni iniziali
t0=0;
tf=20;
tRange=[t0,tf]; %intervallo di tempo
%ode45
[tSol,YSol]=ode45(@prede_predatoriCopia, tRange, Y0); %runge kutta 4-5 ordine
passiode45=length(tSol)
y1=YSol(:,1);
y2=YSol(:,2);
y3=YSol(:,3);
y4=YSol(:,4);
y5=YSol(:,5);
figure
plot(tSol,y1,'*'), hold on
plot(tSol,y2,'*')
plot(tSol,y3,'*')
plot(tSol,y4,'*')
plot(tSol,y5,'*'), hold off
legend('y1','y2','y3','y4','y5')
%%
function [dYdt] = prede_predatori(t,Y)
global r K p1 p2 rho1 rho2 rho3 rho4 mu1 mu2 mu3 mu4 l1 l4 k50 F05y F05t
y1=Y(1); %primar producers
y2=Y(2); %herbivorous
y3=Y(3); %carnivores
y4=Y(4); % decomposers
y5=Y(5); % nutritious
% F05=zeros(size(y1));
% F05(1000:end)=U;
F05=interp1(F05t,F05y,t);
dy1dt=r*y1*(1-y1/K)-p1*y1*y2-(rho1+mu1)*y1;
dy2dt=p1*y1*y2-(rho2+mu2)*y2-p2*y2*y3;
dy3dt=p2*y2*y3-(rho3+mu3)*y3;
dy4dt=mu1*y1+mu2*y2+mu3*y3-(rho4+mu4)*y4;
dy5dt=F05-l1*r*y1*(1-y1/K)+l4*mu4*y4-k50*y5;
dYdt=[dy1dt;
dy2dt;
dy3dt;
dy4dt;
dy5dt];
end

Answers (1)

Nikhil Sonavane
Nikhil Sonavane on 22 Dec 2020
You may refer to the following documentation to understand solving non-linear ODEs-

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!