function ydot=RHS_SIPDB(t,y)


%  RHS_SI
%
%    define the RHS of the model equations. 
%        variables are 
%  y(1): S, susceptible, uninfected individuals [number]
%  y(2): I, infected individuals [number]
%  y(3): D, dead individuals [number]
%  y(4): P, infectious particles close 
%            to the susceptible population [number]
%   array indexes of variables

nVar=4;iS=1;iI=2;iD=3;iP=4;

global PAR 

% most rates are 1/day

ydot=zeros(nVar,1);


% EQUATIONS 

% Time changes of susceptable individuals,

ydot(iS)= PAR.birth*y(iS)- (PAR.Pinfect * y(iP)- PAR.Iinfect * y(iI)- PAR.Dinfect * y(iD)) * y(iS)-PAR.Bmort*y(iS);


% Time changes of infected individuals,

ydot(iI)=  (PAR.Pinfect * y(iP)- PAR.Iinfect * y(iI)- PAR.Dinfect * y(iD)) * y(iS)- PAR.Imort * y(iI);


% Time changes of dead individuals,

ydot(iD)=  PAR.Imort * y(iI)- PAR.Ddecay * y(iD);


% Time changes of susceptable individuals,

ydot(iP)=  PAR.Pburden *PAR.Prelease * y(iD)- PAR.Premove * y(iP);