I have mathematical model with 7 equation . and its difficult to find endemic equilibrium point with analitical , i want to solve with numerical method

6 次查看（过去 30 天）

Latifah 2023-4-12

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1945654-i-have-mathematical-model-with-7-equation-and-its-difficult-to-find-endemic-equilibrium-point-with

评论： Walter Roberson 2023-4-12

can you help me to find the endemic equilibrium point , with all the parameter known

close all; clear all; clc;
Pg=0.0068;
Ps=0.012;
Pr=0.006;
pe=0.002;
beta1=0.00702;
beta2=0.00702;
beta3=0.00702;
d4=1.7143;
u=0.05;
a1=1.10;
a2=4.6205;
a4=4.6666;
a5=1.10;
c1=0.0002;
c2=0.032;
c4=0.032;
c5=0.0002;
d10=1.2*10^-7;
d11=4.2*10^-8;
d12=1.0*10^-7;
d20=0.2051;
d21=0.00431;
d22=19.4872;
d50=1.2*10^-7;
d51=4.2*10^-8;
d52=1.0*10^-7;
tau=0.1615;
mu=0.00371;
psi=0.01813;
delta=2.4*10^-4;
gamma=0.136;
alpha=2.0;
%t = linspace(0,0.1,100)';
%tspan = [0 1000];
X0 = [1.0 0.01 0.0 1.0 0.0 0.01 0.01];
t0 = 0;
tf = 1000;
options = [];
%options=odeset('Abstol', 1e-6, 'Reltol', 1e-6);
[t,y1]=ode45(@glioma1,[t0:0.1:tf],X0,options,Pg,Ps,pe, beta1, beta2,a1, a2, a4, a5, c1,c2,c4,c5,u,d10,d11,d12,d20,d21,d22,d50,d51,d52,d4,tau,mu, beta3,psi, delta, gamma,alpha, Pr);
function hdot = glioma1(t,x,Pg,Ps,pe, beta1, beta2,a1, a2, a4, a5, c1,c2,c4,c5,u,d10,d11,d12,d20,d21,d22,d50,d51,d52,d4,tau,mu, beta3,psi, delta, gamma,alpha, Pr);
% The time-dependent term is A * sin(w0 * t - theta)
%H = 1*sin(-1*t);
%G = x(1)
hdot=zeros(7,1);
%H= heaviside(x(5));
v=0.6;
Phi=3.3*10^-3;
    H=0;
if x(6)>0 || x(7)>0
    H = 1;
%else%if (x(5)<=0);
%   H = 0;
end
hdot(1)=Pg*x(1)*(1-x(1))-beta1*x(1)*(x(2)+x(3))-(d10+d11*x(4)+d12*x(7))*(x(1)*x(6))/(a1+x(1));
hdot(2)=Ps*x(2)*(1-(x(2)+x(3))/1+tau*x(4))-beta2*x(1)*x(2)-u*x(2)*H-(d20+d21*x(4)+d22*x(7))*(x(2)*x(6))/(a2+x(2));
hdot(3)=Pr*x(3)*(1-(x(2)+x(3))/(1+tau*x(4)))-beta3*x(1)*x(3)+u*x(2)*H-v*H*x(3);    
hdot(4)=mu*(x(2)+x(3))+pe*x(4)*(1-x(4))-d4*((x(4)*x(7))/a4+x(4));
if (hdot(1)<0)
hdot(5)=alpha*hdot(1)*x(5)-(d50+d51*x(4)+d52*x(7))*(x(5)*x(6))/a5+x(5);    
    
elseif (hdot(1)>0)
hdot(5)=-(d50+d51*x(4)+d52*x(7))*(x(5)*x(6))/a5+x(5);
end
hdot(6)=Phi-(psi+c1*x(1)/(a1+x(1))+c2*x(2)/(a2+x(2))+c5*x(5)/(a5+x(5)))*x(6);
hdot(7)=delta-(gamma+c4*x(4)/(a4+x(4)))*x(7);
%hdot = hdot';
%hdot=[hdot(:);H,h1]; 
% To make xdot a column
% End of FUN1.M
end

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Walter Roberson 2023-4-12

The mathematics of ode45 is not valid when you have discontinuities in the first or second deriviatives of the equations. You need to construct event functions to detect each condition under which you are changing your equations -- detect x(6) or x(7) changing between negative and positive, detect hdot(1) changing between negative and positive.

Question: what should your hdot(5) be if hdot(1) is exactly 0?

请先登录，再进行评论。

请先登录，再回答此问题。