Koda en pandemi

The following situation will be studied. The SIR model will be solved with a set of parameters 'par1' and a specified initial condition until the number of infected individuals potentially reaches a critical threshold 'upper_bd'. After this time, measures will be taken to decrease the number of infected (lockdown), i.e. another set of parameters 'par2' (with smaller infection rate) will be used. Once the number of infected individuals falls below the threshold 'lower_bd' the measures will be relaxed, i.e. the equation returns to the choice 'par1'Clarification:

This task asks to implement one lockdown only, i.e. the parameters will not be changed again after relaxing them after the lockdown (even if the upper threshold is crossed once more). *)

You will have access to working copies of the functions SIRmodel, check_upper and a corresponding function check_lower working in the same way as check_upper but for the lower threshold.Two outputs are required: a vector t corresponding to the time variable and a vector x containing the corresponding data points of S, I and R

Funktionerna check_upper, check_lower och SIRmodel ser ut som följande:

function [crossing, index] = check_upper(x, upper_bd) %output of the function

index=find(x>upper_bd, 1);

if x(index)>upper_bd

crossing=1;

else crossing=0;

end

end

function [crossing, index] = check_lower(x, lower_bd) %output of the function

index=find(x<lower_bd, 1);

if x(index)<lower_bd

crossing=1;

else crossing=0;

end

end

function y = SIRmodel(t,x,par)

a=par(1);

b=par(2);

dS=-a*x(1)*x(2);

dI=a*x(1)*x(2)-b*x(2);

dR=b*x(2);

y=[dS; dI; dR];

end

Vi bör också använda ode45 för att lösa S I R från SIRmodel

En början på koden är också given:

%initial set of parameters. Don't change these!

par1=[.0009,.07]; %initial parameters

par2=[.00002,.07]; %'lockdown' parameters

lower_bd=10; %lower threshold

upper_bd=100; %upper threshold

IC=[1000;5;0]; %initial values

T=500; %final time

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

t=[]; %solution time vector

x=[]; %solution vector containing S,I and R

%If you want to plot the solution

plot(t,x(:,1),'-',t,x(:,2),':',t,x(:,3),'--','LineWidth',2);

xlabel('time','FontSize',18)

ylabel('S,I,R', 'FontSize',18)

legend('S','I','R')