Engr 3: Introduction to Programming Fall 2020 Final Project Due Wednesday, December 16, 11:55PM (California time) It has been argued by at least one zombie expert that “engineers will make significant contributions during a zombie apocalypse” [1]. One such contribution is the computer simulation of mathematical models for the interactions between zombies and humans, with a goal of predicting if zombies or humans will prevail, and at what cost. To be ready, let’s practice performing such simulations using Matlab. A simple mathematical model for the interactions between zombies and humans is the SZR model [2]. In this model, the population is divided into three groups: (S)usceptibles (humans who have not been bitten and are trying to terminate zombies), (Z)ombies (the undead who are trying to bite humans, which turns them into zombies), and (R)emoved (zombies terminated by humans). There are two possible transitions: a human can become infected if they are bitten by a zombie, and a zombie can be destroyed by direct action by a human. The parameters that describe these transitions are β (the bite parameter that determines the rate at which a zombie will bite and infect a human), and κ (the kill parameter that determines the rate at which humans terminate the zombies). The transitions are illustrated in the following figure: RZS κβ The SZR model is a set of coupled ordinary differential equations: S˙ = −βSZ, Z˙ = (β − κ)SZ, R˙ = κSZ, where S, Z, and R denote the number of individuals in the Susceptible, Zombie, and Removed groups, respectively. We note that the interactions between humans and zombies are analogous to the “law of mass action”from chemistry, which states that a reaction A + B → C occurs at the rate k[A][B], where [ ] denotes concentration of the respective chemical, and k is the rate constant. In the SZR model, for example, the rate at which humans transition to zombies is βSZ, that is, it is directly proportional to the number of susceptibles times the number of zombies in the population, with proportionality constant given by the bite parameter β. Also, note that the variables in the SZR model are continuous, so they can take non- integer values. While this may seem strange - half a zombie? - it is the nature of such 1 models. We will consider a model in (c) which only allows integer values for the numbers of humans and zombies. For all programs, assume that the initial number of individuals in the Removed group is zero, that is, R0 = 0. a) (40 points) Write a Matlab program which simulates the SZR model for initial condi- tions S0 = 190, Z0 = 10, and parameters β = 0.01, κ = 0.008, given in units of (days)−1. The total simulation time should be 28 days. Curves for S(t), Z(t), and R(t) should all be shown on the same plot, and should be clearly labeled, for example by using the legend command. Please submit a file called szr.m for the main program, and a file called szrf.m for the function which encodes the righthand side of the SZR model, and which is called by the Matlab function ode23 in the main program. b) (20 points) It is shown in [2] that if β > κ, the zombies always win (S → 0 as t→∞), for all initial conditions S0 and Z0. On the other hand, if β < κ it is possible that either the humans or the zombies could win, depending on the initial conditions. Let’s consider β = 0.008 and κ = 0.01, and suppose that S0 + Z0 = 200, where S0 could be as high as 199 or as low as 1. Write a Matlab program which determines the smallest integer value that S0 can take for which S(t = 28) > 1. As a hint, start with S0 = 1 and Z0 = 199, and run the simulation. Is S(t = 28) > 1? If not, try S0 = 2 and Z0 = 198. Is S(t = 28) > 1? If not, try S0 = 3 and Z0 = 197, etc. Eventually you’ll find the smallest integer value of S0 for which S(t = 28) > 1; your program should determine and output this value. Please submit your file szr_ic.m which solves this problem. This should call your function szrf from (a). c) (40 points) While the results of the simulations of the SZR model might seem depressing, particularly when β > κ, there is another model called the Stochastic SZR (Stoch-SZR) model which, as we’ll see, is more hopeful. Quoting [2], this model suggests that “even a ferociously virulent zombie infestation might be killed early on by happy accident.” For the Stoch-SZR model, the variables S, Z, and R again denote the number of individuals in the Susceptible, Zombie, and Removed groups, respectively, but they only take non-negative integer values. There are two possible transitions: Susceptible + Zombie βSZ −−→ Zombie + Zombie, (1) Susceptible + Zombie κSZ −−→ Susceptible + Removed. (2) The interpretation of transition (1) is that a susceptible human comes into contact with a zombie, and is bitten and becomes a zombie; in this interaction, S decreases by 1, Z 2 increases by 1, and R stays the same. On the other hand, the interpretation of transition (2) is that a susceptible human comes into contact with a zombie, and kills the zombie; in this interaction, Z decreases by 1, R increases by 1, and S stays the same. The Stoch-SZR model is called a stochastic model because it uses random numbers to determine its behavior. In particular, let dt = 0.001 be a short time interval. The probability that (1) occurs in the interval dt is βSZdt, and the probability that (2) occurs in the interval dt is κSZdt. There is also a probability that neither of these transitions occurs during the interval dt, and in this case S, Z, and R do not change. We will assume that dt is short enough that at most one of transitions (1) and (2) can occur in a given interval. In our simulation, we’ll initialize S0 and Z0, then we’ll determine which of these possi- bilities occurs during the time interval from 0 to dt and update S, Z, and R appropriately, then we’ll determine which of these possibilities occurs during the time interval from dt to 2dt and update S, Z, and R appropriately, etc. Note that all of the time intervals are of size dt. To determine which possibility occurs during the time interval, we’ll draw a uniformly distributed random number r between 0 and 1 using the rand command, and apply the following rule: • If r < βSZdt, then (1) occurs. • If r > 1− κSZdt, then (2) occurs. • Otherwise, no transition occurs. Here S and Z are the current values of the variables. Notice that the minimum value that r can take is 0, and the length of the interval from 0 to βSZdt is βSZdt, which is the probability that (1) occurs; moreover, the maximum value that r can take is 1, and the length of the interval from 1−κSZdt to 1 is κSZdt, which is the probability that (2) occurs. This probability rule is illustrated in the following figure. To interpret this figure, recall that r takes some random value between 0 and 1: 0 11− κSZdt κSZdt βSZdt βSZdt S → S − 1 Z → Z + 1 R→ R S → S Z → Z − 1 R→ R + 1 S → S Z → Z R→ R We will run the simulation until either S = 0 (the zombies win) or Z = 0 (the humans win). 3 It turns out that for the initial conditions and parameters from (a), namely S0 = 190, Z0 = 10, β = 0.01, κ = 0.008, but now used for the Stoch-SZR model, it is possible for the humans to win! Write a Monte Carlo program which simulates 1000 trials to estimate the probability that the humans win. Here each trial is a simulation of the Stoch-SZR model which is executed until either the humans or the zombies win. Please submit your file szr_mc.m which solves this problem. Good luck - the fate of humanity is depending on this! And, please be sure to include comments in your programs! [1] J. Moehlis. Day of the Engineer: Engineering and the Zombie Apocalypse, in But if a Zombie Apocalypse Did Occur: Essays on Medical, Military, Governmental, Ethical, Economic, and Other Implications, ed. A. L. Thompson and A. S. Thompson, McFarland & Company, Inc., Jefferson, North Carolina, 122-137, 2015. [2] A. A. Alemi, M. Bierbaum, C. R. Myers, and J. P. Sethna. You can run, you can hide: The epidemiology and statistical mechanics of zombies. Physical Review E 92:052801, 2015. 4 Academic Honesty Agreement This project is open book, open notes. However, all work submitted must be your own. By submitting these programs, you are asserting that all work on this project is yours alone, and that you did not and will not provide any information to or discuss this project with anyone else doing the project, nor did you receive help from anyone in the class or anyone not in the class, including people on the internet. You are not allowed to post any information about this project on the internet at any time, either before or after the due date. Discussing any aspect of this project with anyone else before the due date constitutes a violation of the academic honesty agreement. A violation can lead to getting an F in the course. Please note that you may ask questions about this project to the Professor and the TAs, but we may answer in a limited way because we want to assess your coding ability, not ours. 5
欢迎咨询51作业君