ME 406The Logistic MapsysidMathematica 6.0.3, DynPac 11.02, 4ê22ê2009plotreset; imsize = 250;intreset;‡1. IntroductionThe logistic map is discussed in many references. A very complete and readable discussion is given inChapter 10 of Nonlinear Dynamics and Chaos by Steven Strogatz, Addison-Wesley, 1994. Many of theinteresting properties of the map were discovered by the mathematical biologist Robert May ("Simple Mathemati-cal Models with Very Complicated Dynamics," Nature 261, 459, 1976.) The basic form of the map isxn+1 = rxn(1 - xn) .In many applications, the map is a model for the dynamics of a population, and xnis the population of the nthgeneration. As the work of May and others has shown, this map exhibits an astonishing range of behavior as thegrowth rate r is varied. We will use the range 0 < r b 4. For r = 4, the interval [0,1] is mapped onto itself; for 0< r < 4, it is mapped into itself. ,Many of the basic features of this map can be established analytically, and someof this analyis may be seen in section 10.3 of Strogatz. In this notebook, our approach will be primarily numeri-cal, using the functions built-in to DynPac. We begin by defining the system for DynPac. We start with thecommand setmap, which tells DynPac that this is a mapping and not a differential equation.setmap;Now we define the state variable, the mapping function, which is assigned to slopevec, and the optional systemname.setstate@8x<D; setparm@8r<D;slopevec = 8r * x * H1 - xL<; sysname = "Logistic Map";We check our definitions with a sysreport.sysreportSYSTEM DEFINITION H11.02LSystem name: sysname = Logistic MapState vector: statevec = 8x<State units: stateunits = 8<Slope vector: slopevec = 8r H1 - xL x<Parameter vector: parmvec = 8r<Parameter values: parmval = 8r<Parameter units vector: parmunits = 8<Time unit: timeunit =System Type: sysmode = mappingWe could have used this same function as the slope for a differential equation, in which case we woulduse the command setde to tell DynPac that we are working with a differential equation. The primary differencebe-tween the mapping and differential equation modes is just the actual stepping algorithm used in constructingsolutions -- a Runge-Kutta step for a differential equation, and a map iteration for the mapping. All of the othersupporting code is essentially identical in the two cases.‡2. Equilibrium PointsWe start by viewing the map for several different values of the parameter r.parmval = 80.5<;2 logistmap.nbviewmap@D;0.00.20.40.60.81.00.00.20.40.60.81.0xr H1 - xL x, 8r< = 8 0.50<parmval = 82.0<;viewmap@D;0.00.20.40.60.81.00.00.20.40.60.81.0xr H1 - xL x, 8r< = 8 2.00<For r = 0.5, there is a single equilibrium point at r = 0. For r = 2, there are two equilibrium points -- one at x = 0and the other at x = 0.5. For a general value of r, the fixed points arefindpolyfix@D:80<, :-1 + rr>>Thus there are always two equilibria, but the second one is relevant (in the range [0,1]) only for 1 b r b 4. Thestability of the equilibria is determined by the derivative of the mapping. We have D@slopevec, xD ê. x Ø 08r<logistmap.nb 3Simplify@D@slopevec, xD ê. x Ø H1 - 1 ê rLD82 - r<Thus the equilibrium at the origin is stable for 0 < r < 1. The other equilibrium is not relevant for 0 < r < 1. Forr > 1 it is relevant. It is stable as long as »2 - r» < 1, hence for 1 < r < 3. At r = 1, the two equilibria coincide at x= 0, and we recognize a transcritical bifurcation at r =1. We could also have determined the stability of the equilibria by using the function classifymap. We dothat now for three different values of r.eq0 = 80<; eq1 = 81 - 1 ê r<;parmval = 80.5<;classifymap@eq0Dstrictly stableparmval = 82.0<;classifymap@eq0Dunstableclassifymap@eq1Dstrictly stableparmval = 83.5<;classifymap@eq0Dunstableclassifymap@eq1DunstableNo surprises in the answers, although they do raise the question of the nature of the attractor when r > 3.We finish our discussion of equilibrium by looking at orbits approaching stable equilibria. We do thisfirst for r = 0.5, then for r = 1.5, and finally for r = 2.5. parmval = 80.5<;The only stable equilibrium for this value of r is the origin. We now use the functioniterate[init,initime,niter,ntoss,ncomp] to construct a sequence of iterations. The first four arguments of iterate arerequired. The fifth argument is optional. The starting value for the iteration is init, the starting time is initime,the total number of iterations performed is niter + ntoss, and the first ntoss of these are thrown away (useful insome cases to eliminate transient approaches to attractors). The optional fifth argument is the level of functionalcomposition for the map. If the fifth argument is omitted, the default is 1, which means the map itself is iterated.In looking for periodic solutions later, we will see uses for higher levels of composition. We try this now,carrying out 15 iterations, starting at x = 0.8 and time = 0, and throwing away none of the iterates. 4 [email protected], 0.0, 15, 0D:80., 0.8<, 81., 0.08<, 82., 0.0368<, 83., 0.0177229<,84., 0.00870439<, 85., 0.00431431<, 86., 0.00214785<,87., 0.00107162<, 88., 0.000535235<, 89., 0.000267474<,810., 0.000133701<, 811., 0.0000668417<, 812., 0.0000334186<,813., 0.0000167088<, :14., 8.35424 µ 10-6>, :15., 4.17708 µ 10-6>>We see that the output is a list of pairs. Each pair has the form {t , x }. We also see that the iterations areconverging well to the stable fixed point at the origin. We may visualize this process by constructing a cobwebplot. The function that does that is cobweb[init,niter,ntoss,ncomp]. The significance of the arguments is thesame as for iterate. We try it. [email protected], 15, 0D;0.00.20.40.60.81.00.00.20.40.60.81.0xr H1 - xL x, 8r< = 8 0.50<This gives us a good geometric view of the monotonic convergence to the equilibrium at the origin. Now we repeat these calculations for r = 1.5.parmval = 81.5<;The equilibria are now ateqstateval@eq0D80<logistmap.nb [email protected]<We can check directly that these are fixed points by evaluating the map at those values.mapval@eq0D80<[email protected]<From our calculations earlier, we know that eq0 is unstable and eq1 is stable. We check this.classifymap@eq0Dunstableclassifymap@eq1Dstrictly stableWe look at an iteration starting near x = [email protected], 0.0, 30,