1 Background The background information is from the book Mathematical Biology by J.D. Murray [2]. The spruce budworm is a species of pests which can defoliate the balsam fir. That is, they feed on the needles of the balsam fir and if there are enough of them, they can eat all of the needles. The population dynamics was proposed to be modeled [1] by: dN = TBN(1 dt - N -) - p(N). KB Here, rB is the linear growth rate of the species (that is, the growth rate when there are no predators and abundant resources). The parameter KB is the carrying capacity. This is related to the density of needles on the trees. Finally, the term p(N) represents predation. Birds are the major predator. When the number of budworms is below a certain threshold, the birds don't really eat the budworms and look for food elsewhere. When the number of budworms is above that threshold, p(N) ~ p*. Thus, p(N) is essentially a "switch". When N is smaller than some threshold, p(N) ~0 and when N is bigger than this threshold, p(N) ~ p*. A good form for p(n) is: p(N): = BN2 A + N2 1 The parameter A is the "threshold" and B is the value p*. Thus, for values of N less than A, p(N) 0 and for values of N bigger than A, p(N) ~ B. (It might be good to graph this for some particular values of A and B. This function is an example of something called a sigmoid function. If you want, read more about them here: https://en.wikipedia.org/ wiki/Sigmoid_function.) So, the ODE becomes: dN dt N = B(1 - -)N KB - BN2 A + N2 Thus, there are four parameters. The parameters A and KB have the same dimensions as N (that is, number of budworms). The parameter rg has dimensions of (time) -1 and B has the dimensions of (dimensions of N) (time) -. Thus, overall, the dimensions of BN (1- A2+N2 are number of budworms - BN2 time - which is what we want. 2 Kill All The Birds We are first going to assume that the trees on which the budworms feed are located in an area where there are no birds. Thus, A = or B = 0 whichever way you want to view it, the point is that p(N) = 0 and so we have the new equation: dN dt == TB (1 - N -)N. KB This is an autonomous equation. Recall that this particular one is called the logistic equation. This can be solved using methods similar to the ones we have used so far. The first part of the project is to analyze this equation. You must do the following: 2 (1) What are the equilibrium solutions? Are they stable or unstable? (2) Let's say you are in charge of the trees on which the budworms feed. What are some steps you can take to limit the number of budworms? How would these changes affect the differential equation (e.g., what happens if you change some of the parameters). (3) Let B .25 and KB = 4. Use dsolve to find and plot several solutions to this ODE with various initial population values. In this problem, there are two equilibrium solutions; I'll call them no and n. Find two solutions with starting values in the (no, n) interval and two solutions with starting values in the (n,) range. (4) Now, using the euler method, find numerical solutions to the ODE. You are going to do this with step sizes .5, .25, .1. For each step size, plot them solutions on the same graph as the solutions from (3). So, you will end up with three graphs: each graph will contain solutions plotted from part (3) and numerical solutions with a given step size. If you are using Matlab, use the script template "budworm_part1.m. (Of course, you will also have to fill out the euler_308.m function that I provided). 3 Bring The Birds Back Now, we are going to go back to the original equation: dN = TB (1 dt N - -)N- KB BN2 A + N2 3 For various reasons, it is good to make this equation "dimensionless". So we introduce the parameters: = N A' r = B KB Bt q= T = A A Thus, we can re-write the equation as: du dT | = ru( 1 1-4) u - = f(u;r,q). 1+u In addition to making the equation dimensionless, we have reduced the number of parameters. We will study this equation. We are going to set r = .25 and q = 16. This gives the ODE: du dT = 42 (1 - = f(u). 16 1+u You will analyze this ODE and answer the following questions: (1) Find the equilibrium solutions (round them to 4 decimal places). You can use either the formulas for roots of quartic or cubic equation, or you can use the solve function in Matlab (or something else). (2) You can't use dsolve to solve this ODE it can't find an explicit solution. So you will use your Euler Method function to approximate solutions. Let's say that the equilibrium solutions are up, u, U2, u3. Use Euler's method to plot numerical solutions. As in the previous part, you will want to plot two solutions with initial values in the interval (uo, u), two solutions in (u, u2) and so on. As above, you will also do this with t step values of of .25, .1, .01. (3) Use your graph above to determine which equilibrium solutions are stable or unstable. 4 If you are using Matlab, use the script template "budworm_part2.m". (Of course, you will also have to fill out the euler 308.m function that I provided). 4 Deliverables You will answer the questions asked above. You must turn in the following files: (1) Your euler 308.m file. (2) Your budworm_part1.m file. (3) Your budworm_part2.m file. (4) A pdf with your answers from Part 1 and 2 and with the graphs on the PDFs. This should be easy for me to read. If you are using something other than matlab, you need to turn in the equivalent of what I asked for in parts 1-3. If you are not using Matlab, you need to use Java, C, C++, python, or Octave. If you want to use something else, you need to check with me first. References [1] D. Ludwig, D.D. Jones, and C.S. Holling, J. Anim. Ecol 47, 315-332. [2] J. D. Murray, Mathematical biology, Biomathematics, vol. 19, Springer-Verlag, Berlin, 1989. 5