遅延微分方程式（MATLABのdde23）

Question

私はdde23を使って遅延微分方程式を解くことを試みましたが、正しく理解できなかったようですので、私が書いた関数にエラーがあり、出力は正しいです。

私はこのシステムを解決したい：

私はプログラムに最後の5つの方程式を追加する方法を理解できませんでした。私は、画像内の1を図ようとして出力を得るために表1パラメータを使用してこのシステムを解決しなければなりません。

これは私が書いたコードは次のとおりです。

clear all; 
clc; 
lags=1; 
sol=dde23(@eq24,lags,@eqh,[0 80]); 
plot(sol.x,sol.y) 


function v=eqh(t) 
v=zeros(6,1); 


function v=eq24(t,s,Ia,Is,R,N) 
Alfa=0.1; 
beta1=0.09; 
beta2=0.1; 
sigma1=0.3; 
sigma2=0.4; 
mu=0.01; 
alfa=0.2; 
rho=0.4; 
r1=0.4; 
r2=0.2; 
d1=0.2; 
d2=0.15; 
k=0.1; 
p=0.8; 
tau=1; 
T=4; 
dsdt=Alfa-beta1*((s*Ia)/(1+sigma1*s))-beta2*((s*Is)/(1+sigma2*s))-mu*s+alfa*R; 
dIadt=rho.*exp(-mu*tau).*s(((beta1.*Ia)/(1+sigma1.*s))+((beta2.*Is)/(1+sigma2.*s)))-(r1+d1+mu).*Ia; 
dIsdt=(1-rho).*exp(-mu.*tau).*s(((beta1.*Ia)/(1+sigma1.*s))+((beta2.*Is)/(1+sigma2.*s)))+(1-k).*r1.*Ia-(r2+d2+mu).*Is; 
dRdt=k.*r1.*Ia+r2.*Is-mu.*R-alfa.*R; 
dNdt=Alfa-mu.*N-d1.*Ia-d2.*Is; 

Answer 1

私は、実行中のコードを思い付くために管理しています。しかし、私はまだあなたの質問のt⁺の部分について混乱しています。私のコードの結果は、あなたが提示した結果とほぼ同じです。私はこれがあなたが何をしているのかに役立つことを願っています。

function sol = eq242 

clf 

global lambda beta1 beta2 sigma1 sigma2 mu alpha rho r1 r2 d1 d2 k tau 

lambda=0.1;beta1=0.09;beta2=0.1;sigma1=0.3;sigma2=0.4;mu=0.01;alpha=0.2; 
rho=0.4;r1=0.4;r2=0.2;d1=0.2;d2=0.15;k=0.1; 

opts = ddeset('RelTol',1e-5,'AbsTol',1e-8); 
sol = dde23(@eq24,tau,[1, 1, 1, 1, 1],[0, 50],opts); 

figure(1) 
plot(sol(1).x,sol(1).y(1,:),sol(1).x,sol(1).y(2,:),sol(1).x,sol(1).y(3,:),sol(1).x,sol(1).y(4,:),sol(1).x,sol(1).y(5,:)) 


function dydt=eq24(t,y,Z) 

global lambda beta1 beta2 sigma1 sigma2 mu alpha rho r1 r2 d1 d2 k tau 
s = y(1); 
Ia = y(2); 
Is = y(3); 
R = y(4); 
N = y(5); 
slag1 = Z(1,1);  
ialag2 = Z(2,1); 
islag3 = Z(3,1); 

dsdt=lambda-beta1*((s*Ia)/(1+sigma1*s))-beta2*((s*Is)/(1+sigma2*s))-mu*s+alpha*R; 
dIadt=rho*exp(-mu*tau)*slag1*(((beta1*ialag2)/(1+sigma1*slag1))+((beta2*islag3)/(1+sigma2*slag1)))-(r1+d1+mu)*Ia; 
dIsdt=(1-rho)*exp(-mu*tau)*slag1*(((beta1*ialag2)/(1+sigma1*slag1))+((beta2*islag3)/(1+sigma2*slag1)))+(1-k)*r1*Ia-(r2+d2+mu)*Is; 
dRdt=k*r1*Ia+r2*Is-mu*R-alpha*R; 
dNdt=lambda-mu*N-d1*Ia-d2*Is; 
dydt= [dsdt; dIadt; dIsdt; dRdt; dNdt];