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

1 vue (au cours des 30 derniers jours)
chiara garavelli
chiara garavelli le 24 Jan 2020
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

Réponses (1)

Nikhil Sonavane
Nikhil Sonavane le 22 Déc 2020
You may refer to the following documentation to understand solving non-linear ODEs-

Catégories

En savoir plus sur Programming dans Help Center et File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by