Index exceeds the number of array elements. Index must not exceed 3. HELP

2 次查看（过去 30 天）
Elizabeth Pickford
Elizabeth Pickford 2024-4-26
回答： Sam Chak 2024-4-26
model code:
clear all;
close all;
S = 99;
I = 1;
R = 0;
N = 100; %Total population
beta= 0.1; % birth rate
alpha= 0.1; % infection person to person rate
lambda= 0.3; % infection by water rate
vac= 0.05; % recovery by vaccination rate
d= 0.03; % death rate
gamma= 0.8; % recovery rate
c= 0.9; % rate of contamination
m= 0.4; % rate of decay of V. cholerae
B= 0.0; % initial concentration of V. cholerae
t_f = 500; %Ending time of simulation
Q = [beta alpha lambda vac d gamma c m B];
[t,y]=ode45('cholera_de',[0:t_f/100:t_f],[S I R]',[],Q);
figure(1)
plot(t,y(:,1),'k-',t,y(:,2),'r--',t,y(:,3),'b:');
xlabel('\bf Time (days)');
ylabel('\bf Number of People by Category');
legend('S','I','R');
z=y(end,:)'
SN=y(:,1)/N;
IN=y(:,2)/N;
figure(2)
plot(IN,SN);
xlabel('\bf I/N');
ylabel('\bf S/N');
r0=beta/beta
[maxIN,y_maxtime]=max(y(:,2)/N);
maxIN
maxtime=y(y_maxtime)
eqIN=y(100,2)/N;
eqIN
function code:
function dy=cholera_de(t,Y,flag,Q)
beta= Q(1);
alpha= Q(2);
lambda= Q(3);
vac= Q(4);
d= Q(5);
gamma= Q(6);
c= Q(7);
m= Q(8);
S= Y(1);
I= Y(2);
R= Y(3);
B=Y(4);
N= S+I+R;
dy(1,1)= beta - alpha*I - lambda*B - vac*S - d;
dy(2,1)= alpha*I + lambda*B - d - gamma*I;
dy(3,1)= gamma*I + vac*S - d;
dy(4,1)= c*I - m*B;
Gives me this error:
Index exceeds the number of array elements. Index must not exceed 3.
Error in cholera_de (line 15)
B=Y(4);
Error in odearguments (line 92)
f0 = ode(t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode45 (line 104)
odearguments(odeIsFuncHandle,odeTreatAsMFile, solver_name, ode, tspan, y0, options, varargin);
Error in cholera_model (line 23)
[t,y]=ode45('cholera_de',[0:t_f/100:t_f],[S I R]',[],Q);
  1 个评论
Chuguang Pan
Chuguang Pan 2024-4-26
编辑：Chuguang Pan 2024-4-26
It seems that your initial condition is [S I R]', which has three elements. However, the cholera_de function need four states Y(1), Y(2), Y(3), Y(4). You should keep number of initial states the same as states transfered to cholera_de function.

请先登录，再进行评论。

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

回答（1 个）

Sam Chak
Sam Chak 2024-4-26
Ran in:
Hi @Elizabeth Pickford
Few places require fixes. Check out the annotated lines.
S0 = 99;
I0 = 1;
R0 = 0;
N0 = S0 + I0 + R0; % Total population % <-- fix it here
beta = 0.1; % birth rate
alpha = 0.1; % infection person to person rate
lambda = 0.3; % infection by water rate
vac = 0.05; % recovery by vaccination rate
d = 0.03; % death rate
gamma = 0.8; % recovery rate
c = 0.9; % rate of contamination
m = 0.4; % rate of decay of V. cholerae
B = 0.0; % initial concentration of V. cholerae
Q = [beta alpha lambda vac d gamma c m B];
t_f = 500; % Ending time of simulation
tspan = [0:t_f/100:t_f]; % <-- fix it here
y0 = [S0; I0; R0; N0]; % <-- fix it here
[t, y] = ode45(@(t, y) cholera_de(t, y, Q), tspan, y0); % <-- fix it here
figure(1)
plot(t,y(:,1),'k-',t,y(:,2),'r--',t,y(:,3),'b:'); grid on
xlabel('\bf Time (days)');
ylabel('\bf Number of People by Category');
legend('S', 'I', 'R');
z = y(end,:)'
z = 4x1
14.2061 -0.9011 91.6949 -2.0108
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
N = y(:,1) + y(:,2) + y(:,3); % <-- fix it here
SN = y(:,1)./N; % <-- fix it here
IN = y(:,2)./N; % <-- fix it here
figure(2)
plot(IN, SN); grid on
xlabel('\bf I/N');
ylabel('\bf S/N');
r0 = beta/beta;
[maxIN, y_maxtime] = max(y(:,2)/N);
maxIN;
maxtime = y(y_maxtime);
eqIN = y(100,2)/N;
eqIN;
function dy = cholera_de(t, y, Q) % <-- fix it here
%% parameters
beta = Q(1);
alpha = Q(2);
lambda = Q(3);
vac = Q(4);
d = Q(5);
gamma = Q(6);
c = Q(7);
m = Q(8);
%% definitions
S = y(1);
I = y(2);
R = y(3);
B = y(4);
N = S + I + R;
%% differential equations
dy(1,1) = beta - alpha*I - lambda*B - vac*S - d;
dy(2,1) = alpha*I + lambda*B - gamma*I - d;
dy(3,1) = gamma*I + vac*S - d;
dy(4,1) = c*I - m*B;
end

类别

在 Help Center 和 File Exchange 中查找有关 Numerical Integration and Differential Equations 的更多信息

标签

产品


版本

R2024a

Community Treasure Hunt

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

Start Hunting!

Translated by