Problems replacing for-loops

Hello,

it has been a while since I needed to touch matlab seriously, or any programming language with heavy array-operation, for that matter.

We've been tasked to simplify and vectorise this function ('s for-loop). Aaaaand, I am stuck. Stuck trying to wrap my head around it, it just doesn't happen. This is supposed to be an introductory task, but this topic has always baffled my mind, and so far I didn't have the actual need to stop using for-loops. I know that you should not use them, but I was always having problems removing them. THey are the quick-and-dirty way to get stuff done, and my usual workloads are not big enough that I would care about saving the marginal amount of time I am saving each iteration.

Could someone give an explanation as to how this is supposed to be done?

Thank you.

Sincerely,

Claudius

clear all
close
clc
% set start param
S0=20;
I0=15;
R0=10;
a=3;
d=0.1;
lambda=0.3;
beta=0.1;
m=0.01;
alpha1=1;
alpha2=0.1;
r=0.2;
tmax=200;


SIR_01(S0,I0,R0,a,d,lambda,beta,m,alpha1, alpha2,r,tmax)



function SIR_01(S0,I0,R0,a,d,lambda,beta,m,alpha1, alpha2,r,tmax)
    S=zeros(1,tmax);
    I=zeros(1,tmax);
    R=zeros(1,tmax);
    S(1)=S0;
    I(1)=I0;
    R(1)=R0;
    T(1)=0;
    sMat=S
    tMat=0:1:tmax
    for t=1:1:tmax
        
        T(t+1)=t;
        S(t+1)=S(t)+(a-d*S(t)-(lambda*S(t)*I(t)/(1+alpha1*I(t)+ alpha2*I(t)^2))+beta*R(t));
        I(t+1)=I(t)+((lambda*S(t)*I(t)/(1+alpha1*I(t)+alpha2*I(t)^2))-(d+m+r)*I(t));
        R(t+1)=R(t)+(m*I(t)-(d+beta)*R(t)+r*I(t));
    end
%     plot(T,S,T,I,T,R); % deactivated the plot stuff to skip that while
%     working. I am currently not interested in the plot itself.
    IsequalNewT=isequal(T,tMat) 

%       title('SIR Model'); 
%       legend('S(t)','I(t)','R(t)');
%       xlabel('time, t'); 
%       ylabel('Population');
    
end

Answer 1

[David Hill]

There is nothing wrong with for-loops. The only thing you can vectorize is T

T=0:tmax;
for t=1:1:tmax
    S(t+1)=S(t)+(a-d*S(t)-(lambda*S(t)*I(t)/(1+alpha1*I(t)+ alpha2*I(t)^2))+beta*R(t));
    I(t+1)=I(t)+((lambda*S(t)*I(t)/(1+alpha1*I(t)+alpha2*I(t)^2))-(d+m+r)*I(t));
    R(t+1)=R(t)+(m*I(t)-(d+beta)*R(t)+r*I(t));
end