Math Modeling Lecture 9: Discrete Dynamical Systems Page 14452 Mathematical Modeling Lecture 9: Discrete Dynamical Sys-temsMany processes in the world are not continuous, but discrete. Compound interest, population growth, andfeedback control all may be more suitably modeled as a discrete time event. Discrete dynamical systems ordifference equations leads to recurrence relations which must be analyzed. Recurrence relations are easy toset up and solve, using a calculator, a spreadsheet, or a computer. Although you may be able to get a closedform for the solution, you can always just iterate the recurrence relation to determine the behaviour of thesystem.Recurrence Relation BasicsA first order recurrence relation has the form:yn+1= f(n, yn), n = 0, 1, 2, 3, . . . ,which is sometimes written asy(xn+1) = f(n, y(xn)), n = 0, 1, 2, 3, . . . .This is first order since the value of yn+1depends explicitly on the previous yn, but not earlier valueslike yn−1, yn−2. It is linear if f is a linear function of yn, otherwise it is nonlinear. It is autonomous iff(n, yn) = f(yn). An initial condition y0= α is required to initialize the difference equation, and thesolution is the sequence of numbers y0, y1, y2, . . ..Aside: A Difference Equation would look like:yn+1− yn= f(n, yn).Example Let’s analyze the recurrence relationyn+1= −12yn+ 6. (1)After choosing an initial condition y0= α, you could just iterate this equation and plot it and try anddetermine the b ehaviour that way. However, Mathematica has a package which is very similar to DSolvethat can solve recurrence relations for you. It is RSolve, and it has to be loaded into the kernel before youcan use it. Using RSolve, we find that the solution isyn=−4 (−1)n+ 22+n+ (−1)nc2nThis is a family of sequences, where the c can b e identified as the analogue of the constant of integration ina differential equation. We can use the initial condition to determine c:y0= α = cMath Modeling Lecture 9: Discrete Dynamical Systems Page 2yn=−4 (−1)n+ 22+n+ (−1)nα2nThe solution for α = 2 is shown in Fig. 1. The sequence of solutions looks like this:{{0, 2}, {1, 5}, {2, 3.5}, {3, 4.25}, {4, 3.88}, {5, 4.06}, {6, 3.97}, {7, 4.02}, {8, 3.99}, {9, 4.00}, . . .}We can see an oscillation, and that the solution appears to be approaching 4 as n → ∞ (∞ = 10, in thiscase...).0 24 6810n22.533.544.55y_n024 6810n22.533.544.55y_nFigure 1: The solution to Eq. (1) when y0= 2. Typically, any figures you find in the literature will have thepoints joined for ease of viewing, but really the discrete function is not defined except at the points.The equilibrium solution is found by solving yn= yn+1= f (yn) for yn. We guessed that the equilibriumsolution would be 4 based on our graph, and we can prove that numerically now.Another way of describing the series is sometimes used for autonomous equations. This method is related tothe idea of nullclines we saw in differential equations. Here you let yn= x, plot the function f(x) = −12x + 6and the straight line yn= x. The intersection of these two lines will be at the equilibrium point. I will drawthis representation for this case in class, but I will include the graph here in Fig. 2.12 3456123456Figure 2: The solution to Eq. (1) when y0= 2.Math Modeling Lecture 9: Discrete Dynamical Systems Page 3Programming NotesBefore we get to chaos, I should say a few words about programming with smarts. Here is the obvious wayto program a recursion relation in Mathematica:y[0] := 1/3.y[n_] := 4 y[n - 1] ( 1 - y[n - 1])data2 = Table[{n,y[n]}, {n, 0, 40}]This works fine logically, and is alright for small n. However, think about what this is doing. What ishappening with this construct is Mathematica climbs down the “tree” all the way to the base to get eachnew yn! It is not storing values, so if it wants y5it has to recompute y4, y3, y2, and y1(y0is simply definedin this case). Then for y6it has to recompute all those values! It would be much better (read: faster) if westored our values as we went, and then Mathematica could just look ups each value each time rather thanrecomputing. Storing each value can be accomplished by the following Mathematica commands:y[0] := 1/3.y[n_] := y[n] = 4 y[n - 1] ( 1 - y[n - 1])data2 = Table[{n,y[n]}, {n, 0, 40}]ChaosVirtually any book on chaos begins with simple recursion relations, and chaos becomes readily apparent.Chaos can be thought of as a system which is highly sensitive to its initial conditions. Over short timeperiods, small changes in the initial conditions do not affect the solution much at all, but the cumulativeeffect of small differences at each time step results in two solutions which are greatly different over longertime intervals. This obviously has some relation to the sensitivity analysis you have been performing.Why this chaotic nature is readily apparent in discrete systems is due to the time delay that is built intothe system. Figure 3 shows two solutions to the discrete dynamical systemyn= 3.65yn−1(1 − yn−1), y0= α, n = 1, 2, 3, . . . (2)for different initial values α. This is a logistic equation which can be used to model population growth. Youcan see that for short times the solutions remain similar, but for longer times they are different. It turnsout that the solution is chaotic because of the choice of 3.65 in our model. It is difficult to determine thebehaviour in an analytic sense for this problem, since it is nonlinear. Chaos would tell us that if the growthrate (which is what 3.65 represents) is too large it will be impossible to make long range predictions on thebehaviour of the system [1].Comparing Continuous and Discrete ModelsWhat about the continuous analogue of Eq. (2)? You might think that it would be the differential equationdydt= 3.65y(1 − y), y0= α,Math Modeling Lecture 9: Discrete Dynamical Systems Page 4but this isn’t quite right. Let’s take a moment to examine what is meant by the relationship between acontinuous and discrete model. It is easiest to begin from the continuous, and work towards the discrete.dydt= f(y)∆y∆t= f(yn−1)yn− yn−1∆t= f(yn−1)yn= f(yn−1)∆t + yn−1, n = 1, 2, 3, . . .and we see that our original guess was wrong! We really want to be working in the frame of differenceequations, not recurrence relations, when we compare discrete and continuous models.So if we have the discrete model which is given by a recurrence relationyn= 3.65yn−1(1 − yn−1), y0= α, n = 1, 2, 3, . . . (3)the