Discrete Dynamical Systems

by
K. D.  Stroyan
University of Iowa

copyright 1997 by Academic Press, Inc. - All rights reserved.

Use of  cobWeb[ f, initP, nrIter, pmax ] is described in a sample.

Needs["Graphics`Colors`"];

Off[General::spell1];
BeginPackage["Cobwebs`"];

cobWeb::usage = "cobWeb[f, p1, itnum, pmin,pmax, opns] calculates itnum
iterates of the discrete dynamical system p[i] = p[i-1] + f[ p[i-1] ]
beginning at p1 and draws a graph of the resulting 'cobweb'.  pmin, pmax
and opns are parameters that specify a range for the p-axis and
any options to use in drawing the plots.  cobWeb returns a list of
the iterates it calculates.";

Begin["`Private`"];

maxf[ g_,pmin_ , pmax_ ] :=
Module[ {prange, funcvalues},
		prange = Range[pmin, pmax, pmax/100] ;
		funcvalues = Map[ g[#]&, prange ] ;
		Max[ funcvalues, pmax ]
	]
minf[ g_, pmin_ , pmax_ ] :=
Module[ {prange, funcvalues},
		prange = Range[pmin, pmax, pmax/100] ;
		funcvalues = Map[ g[#]&, prange ] ;
		Min[ funcvalues, pmin ]
	]

makeWeb[ iterates_List, p1_ ] :=
Module[ {linepoints, funcpoints, pathpoints, depth},
		linepoints = Drop[ Transpose[ {iterates, iterates} ],
								1 ] ;
		funcpoints = Partition[ iterates, 2, 1 ] ;
		pathpoints = { } ;
		depth = Length[funcpoints] ;
		Do[
			AppendTo[ pathpoints,
				{funcpoints[[i]], linepoints[[i]]} ],
			{i, 1, depth}
		] ;
		pathpoints = Flatten[pathpoints, 1] ;
		pathpoints = PrependTo[pathpoints, {p1, 0} ] ;
		Line[Drop[pathpoints, -1]]
	]
	
cobWeb[ f_, p1_, itnum_, pmin_, pmax_,color___ ] :=
Module[ {iterates, path,ymin,ymax,gPlot},
	g[p_] := p + f[p] ;
	ymin = minf[ g, pmin,pmax ] ;
	ymax = maxf[ g, pmin,pmax ] ;
	gPlot=	Plot[ {p, g[p]}, {p, pmin, pmax}, 
			PlotRange -> {{pmin, pmax}, {ymin, ymax}},
			AxesOrigin -> {0, 0}, AspectRatio->Automatic,
			DisplayFunction->Identity ] ;
					
		iterates = NestList[ g, p1, itnum ] ;
		path = makeWeb[ iterates, p1 ] ;
Return[
{Show[gPlot,Graphics[{color,path}],DisplayFunction->$DisplayFunction],
		iterates}]
	]
	
End[ ];

EndPackage[ ] ;
On[General::spell1];

Notebook Overview

This notebook will allow you to experiment with discrete dynamical systems.  The "change form" of such systems is:

p[0] = P0
	p[t+1] = p[t] + f[p[t]]

where you define the change function f[p_] := ??? as in the samples.  The change in  p[.]  as  t  goes from  t  to  t+1  is p[t+1] - p[t], on one hand, and given by the function  f[.], on the other,

p[t+1] - p[t] = f[ p[t] ]   ⟺   p[t+1] = p[t] + f[p[t]]

The function iteration form that is used with "cobwebs" is:

p[0] = P0
	p[t+1] = g[p[t]]

where g[p_] := p + f[p]  or  f[p_] := g[p] - p.  Both of these forms are useful.

Comparison to Euler's Approximation

The change form is more like Euer's approximation to solution of differential equations:

(dx/dt)[t]= f[x[t]]

((x[t + δt] - x[t])/(δt))≈ f[x[t]]

x[t+δt] ≈ x[t] + f[x[t]]·δt   vs.   p[t+1] = p[t] + f[p[t]]

References to the Text

This NoteBook goes with Chapter 20 of the text.  Also see the FirstDscrDySy NoteBook.

Cobwebs and the Iteration Function  g[.]

The iteration approach to solving discrete dynamical systems can be viewed as a "cobweb" winding around the graph  q = g[p]  and  the graph  q = p.  This is explained in detail in Chapter 20.

The change form of a dynamical system is:

p[0] = P0
	p[t+1] = p[t] + f[p[t]]

where the function  f[p]  is known.  From the value of p[0] and the change law, we compute

p[1] = p[0] + f[p[0]]
	p[2] = p[1] + f[p[1]]
	p[3] = p[2] + f[p[2]]
	p[4] = p[3] + f[p[3]], etc.

If we define the iteration function  g[p] = p + f[p], this can be written

p[1] = g[p[0]]
	p[2] = g[p[1]] = g[ g[p[0] ]
	p[3] = g[p[2]] = g[ g[ p[1]] ] = g[ g[ g[p[0] ] ] ]
	p[4] = g[ g[ g[ g[p[0] ] ] ] ], etc.,
	p[n] = g[g[ ...n-times ...p[0]...]]

The cobWeb[.] function requires you to define  f[p]  in the form

f[p_] := 2.58 p - 3.58 (p^2);

The package itself defines the iteration function  g[p] = p + f[p].

If you have a formula for  g[p]  instead, say

g[p_]:=3.58 p (1-p);

you can simply define  f[p_] := 3.58 p (1  - p) - p, since  f[p] = g[p] - p, and let the package re-define the iteration function.

The "Cob Web" process needs the initial value of  p, p[0] = P0 as a specific number.  This is the
second variable in the cobWeb[f,initP,.] function.  The third variable is the number of times you wish to iterate the process.  The fourth variable is the maximum value of  p  in your graph.  To solve  15  steps of the problem

p[0] = (1/2)
	p[t + 1] = p[t] + 2.58 p[t] - 3.58 p[t]^2

we simply <Enter> the instructions in the next input cell (provided  cobWeb[·] was initialized above.):

\(Clear[f, p]; 
f[p_] := 2.58\ p - 3.58\ p\^2; 
\(initP = 0.5; \)\)nrIter = 15; pMin = 0; pMax = 1; cobWeb[f, initP, nrIter, pMin, pMax];

This sample does not save the sequence of terms as a list.  The list is useful when you want to see an explicit graph of the solution, rather than only the cobweb.  The list can be saved as follows:

\(Clear[f, p]; 
f[p_] := 2.58\ p - 3.58\ p\^2; 
\(initP = 0.05; \)\)nrIter = 15; pMax = 1; soln = cobWeb[f, initP, nrIter, pMin, pMax]; ListPlot[soln[[2]], 
      PlotJoined -> True];

A Nonlinear System p[t+1] = p[t] + f[ p[t] ]

for a nonlinear  f[.].  The function in the next cell has an attracting equilibrium point.  <Enter> the computation and watch the convergence.

\(Clear[f, p]; 
f[p_] := 1.5\ p\ \((1 - p\/50)\); 
\(initialPrice = 75; 
    \)\)numberIterations = 15; pMin = 0; pMax = 100; iterate1 = 
      cobWeb[f, initialPrice, numberIterations, pMin, pMax, Red]; ListPlot[iterate1[[2]], 
      PlotJoined -> True, PlotRange -> {0, 100}]; lastnrs = 
      Take[iterate1[[2]], -6]

Now let's try a different initial price:

initialPrice = 25; iterate2 = 
      cobWeb[f, initialPrice, numberIterations, pMin, pMax, Green]; ListPlot[iterate1[[2]], 
      PlotJoined -> True, PlotRange -> {0, 100}]; lastnrs = 
      Take[iterate1[[2]], -6]

Show[iterate1[[1]], 
      iterate2[[1]]];

Practice Exercise

1)  Change the initial price above from to  50 and the color to Blue and re-run the computation.

2)  Show the Red, Green and Blue cobwebs all together.

Chaos

p[t+1] = p[t] + f[p[t]]  when an equilibrium point is not attracting locally it need not be repelling either.  Points far away from equilibrium can jump back near equilibrium occasionally.  (This is called "chaotic behavior.)  Here is an example.

\(Clear[f, p]; 
f[p_] := 2.58\ p - 3.58\ p\^2; 
\(initialPrice = 0.95; 
    \)\)numberIterations = 100; pMin = 0; pMax = 1; iterates = 
      cobWeb[f, initialPrice, numberIterations, pMin, pMax, Red]; grf = 
      ListPlot[iterates[[2]], 
        PlotJoined -> True, PlotRange -> {0, 1}, 
        PlotStyle -> {Red}]; initialPrice = 0.1; iterate2 = 
      cobWeb[f, initialPrice, numberIterations, pMin, pMax, Blue]; grf2 = 
      ListPlot[iterate2[[2]], 
        PlotJoined -> True, PlotRange -> {0, 1}, 
        PlotStyle -> {Blue}];

Some cobWeb[.]  Experiments with The Linear System 
p[t+1] = p[t] + α p[t]	
or	
f[p_]  :=  α * p

The following cells experiment with different values of  alpha  to try to see when we have stablility in the simple system where the change is linear  f[p] = alpha p.  We will keep the maximum and minimum ranges the same in all our examples and always have 15 iterations.  <Enter> the definition of  f[.]  now.

\(Clear[f, p, α]; 
f[p_] := α\ p; 
f[p]\)numberIterations = 15; pMax = 10; pMin = -pMax;

Run the iterations below with various values of  alpha.  Notice that the change in  p is linear:

p[t+1] - p[t]  = f[p[t]] = α  p[t]

and change is only zero when  p[t] = 0.  Zero is the only equilibrium point.

You are primarily trying to observe whether points tend toward zero or tend away from zero.  We plot the solutions for 4 different starting values: p[0] = 1,-2,5,-6.

Sample computations with  α = -0.2

α = -0.2;

initialPrice = 1 ;
iterate1 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Red ];
ListPlot[ iterate1[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate1[[2]],-6]

initialPrice = -2 ;
iterate2 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Green ];
ListPlot[ iterate2[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate2[[2]],-6]

initialPrice = 5 ;
iterate3 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Blue ];
ListPlot[ iterate3[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate3[[2]],-6]

initialPrice = -6 ;
iterate4 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Maroon ];
ListPlot[ iterate4[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate4[[2]],-6]

Show[iterate1[[1]],iterate2[[1]],iterate3[[1]],iterate4[[1]]];

SUMMARY FOR THE EXPERIMENTS WHEN  α = -0.2:
In each case of initial value, the successive values of the solution sequence get smaller and smaller, tending toward zero.  The value  α =-0.2  makes the equilibrium of the linear discrete system attractive.

α = +0.2

α = +0.2;

initialPrice = 1 ;
iterate1 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Red ];
ListPlot[ iterate1[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate1[[2]],-6]

initialPrice = -2 ;
iterate2 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Green ];
ListPlot[ iterate2[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate2[[2]],-6]

initialPrice = 5 ;
iterate3 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Blue ];
ListPlot[ iterate3[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate3[[2]],-6]

initialPrice = -6 ;
iterate4 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Maroon ];
ListPlot[ iterate4[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate4[[2]],-6]

Show[iterate1[[1]],iterate2[[1]],iterate3[[1]],iterate4[[1]]];

SUMMARY FOR THE EXPERIMENTS WHEN  α = 0.2:
In each case of initial value, the successive values of the solution sequence get larger and larger, tending toward plus or minus Infinity.  The value  α = 0.2  makes the equilibrium of the linear discrete system repulsive.

α = -1.0

α = -1.0;

initialPrice = 1 ;
iterate1 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Red ];
ListPlot[ iterate1[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate1[[2]],-6]

initialPrice = -2 ;
iterate2 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Green ];
ListPlot[ iterate2[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate2[[2]],-6]

initialPrice = 5 ;
iterate3 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Blue ];
ListPlot[ iterate3[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate3[[2]],-6]

initialPrice = -6 ;
iterate4 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Maroon ];
ListPlot[ iterate4[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate4[[2]],-6]

Show[iterate1[[1]],iterate2[[1]],iterate3[[1]],iterate4[[1]]];

YOUR SUMMARY FOR THE EXPERIMENTS WHEN  α = -1.0:

α = -1.5

α = -1.5;

initialPrice = 1 ;
iterate1 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Red ];
ListPlot[ iterate1[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate1[[2]],-6]

initialPrice = -2 ;
iterate2 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Green ];
ListPlot[ iterate2[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate2[[2]],-6]

initialPrice = 5 ;
iterate3 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Blue ];
ListPlot[ iterate3[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate3[[2]],-6]

initialPrice = -6 ;
iterate4 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Maroon ];
ListPlot[ iterate4[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate4[[2]],-6]

Show[iterate1[[1]],iterate2[[1]],iterate3[[1]],iterate4[[1]]];

YOUR SUMMARY FOR THE EXPERIMENTS WHEN  α = -1.5:

α = -2.0

α = -2.0;

initialPrice = 1 ;
iterate1 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Red ];
ListPlot[ iterate1[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate1[[2]],-6]

initialPrice = -2 ;
iterate2 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Green ];
ListPlot[ iterate2[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate2[[2]],-6]

initialPrice = 5 ;
iterate3 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Blue ];
ListPlot[ iterate3[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate3[[2]],-6]

initialPrice = -6 ;
iterate4 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Maroon ];
ListPlot[ iterate4[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate4[[2]],-6]

Show[iterate1[[1]],iterate2[[1]],iterate3[[1]],iterate4[[1]]];

YOUR SUMMARY FOR THE EXPERIMENTS WHEN  α = -2.0:

α = -2.2

α = -2.2;

initialPrice = 1 ;
iterate1 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Red ];
ListPlot[ iterate1[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate1[[2]],-6]

initialPrice = -2 ;
iterate2 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Green ];
ListPlot[ iterate2[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate2[[2]],-6]

initialPrice = 5 ;
iterate3 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Blue ];
ListPlot[ iterate3[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate3[[2]],-6]

initialPrice = -6 ;
iterate4 = cobWeb[ f, initialPrice, numberIterations, pMin,pMax,Maroon ];
ListPlot[ iterate4[[2]], PlotJoined -> True] ;
lastnrs = Take[iterate4[[2]],-6]

Show[iterate1[[1]],iterate2[[1]],iterate3[[1]],iterate4[[1]]];

YOUR SUMMARY FOR THE EXPERIMENTS WHEN  α = -2.2:

Exercise: Analytical Support

A little theory might help you organize your observations on the above examples.  Let's try to solve this problem analytically: p[0] = P0,  p[t+1] = p[t] + alpha p[t] gives:

p[0] = P0
	p[1] = p[0] + α p[0] = p[0] (1 + α)
	p[2] = p[1] + α p[1] =  p[1] (1 + α) =  p[0] (1 + α) (1 + α) = p[0] (1 + α)^2
	p[3] = p[2] + α p[2] =  p[2] (1 + α) =  p[0] (1 + α)^??

Show that

p[n] = p[0] (1 + alpha)^n

How does

Limit[ p[0] (1+ alpha)^n , n->Infinity ] =  p[0] Limit[ (1+ alpha)^n , n->Infinity ] = ?

depend on alpha?  What can happen with this limit?

How does your analytical criterion compare with the experiments above?

α = -0.2

P0 Limit[(0.8^n), n -> ∞]

α = 0.2

P0 Limit[(1.2^n), n -> ∞]

α = -1.0

\(P0\ Limit[?? \^n\(, n -> ∞\)]\)

α = -1.5

\(P0\ Limit[?? \^n\(, n -> ∞\)]\)

α = -2.0

\(P0\ Limit[?? \^n\(, n -> ∞\)]\)

α = -2.2

\(P0\ Limit[?? \^n\(, n -> ∞\)]\)

Your Summary

FINALLY, GIVE YOUR SUMMARY:

The linear system

p[0] = P0
	p[t+1] = p[t] (1 + α)   ⟺   p[t] = P0 (((1 + α))^t)

1)  Tends toward zero for any initial  p[0]  provided  ??

2)  Tends toward +∞ for positive p[0] and -∞ for negative p[0]  provided ??

3)  Oscillates ???

4)  Oscillates more and more ???

Exercises

Exercise CD20.1.1    Stability of p[t+1] = p[t] + a p[t] + b

Experiment with your choice of calculations in this NoteBook to formulate a conjecture about the role of a and b in determining both the stability and limiting value of the general affine linear dynamical system.  Consider the cases:

(1)  a = -0.2, 	b = 100,	P0 = 100, 500, 1000
	(2)  a = +0.2,	b = 100,	P0 = 100, 500, 1000
	(3)  a = -0.2,	b = 50,		P0 = 100, 500, 1000
	(4)  a = +0.2,	b = 50,		P0 = 100, 500, 1000
	(5)  a = -0.1,	b = 100,	P0 = 100, 500, 1000
	(6)  a = +0.1,	b = 100,	P0 = 100, 500, 1000
	(7)  a = -1.1,	b = 110,	P0 = 100, 500, 1000
	(8)  a = -2.0,	b = 200,	P0 = 100, 500, 1000
	(9)  a = -2.2,	b = 220,	P0 = 100, 500, 1000

HINT:  When is the change  f[p]  = 0, so that you are at an equilibrium?  When was the linear system in the sample comptation above stable?

Apply the results of your solution to make predictions about the Bytie and Corn economy models of the text, that is, re-formulate your conjectures in terms of the parameters k, aD, aS, bD, bS, in both cases.  (These can be prhased in economic terms once you apply your analytical criterion.)

Exercise CD20.2.4

Sketch the first three iterates of the 'cobweb' for the dynamical systems by hand and machine.  Then sketch 30 by machine.

(1)	p[0] = 10
	p[t+1]  =  p[t] + 2 p[t] ( 1 - p[t] / 50 )

(2)	p[0] = 70
	p[t+1]  =  p[t] + 2 p[t] ( 1- p[t] / 50 )

Exercise CD20.3.1

Show that the only equilibrium of the system   p[t+1] = p[t] + a p[t]   is  Pe = 0.  Draw the cobwebs and write the closed formula solutions in the cases

(1)  a = -1/2	
	(2)  a = +1/2	
	(3)  a = -1	
	(4)  a = +1

Exercise CD20.4.2

(1)  Find the two equilibrium values of the Logistic Growth Model

p[0] = P0
p[t+1] = p[t] + α p[t](1 - p[t]/c)

(2)  Formulate conjectures about the stability of these two equilibria by doing experiments with this nobebook in various cases, for example:

α = 1.5		C = 50		P0 = 10,  60,  70
	α = 2.5		C = 50		P0 = 10,  60,  70

Exercise CD20.5.2

Find the equilibria of the discrete dynamical system

Make cobweb diagrams for a variety of initial conditions above and below these equilibria and decide by observation which are stable attractors and which are not.  Apply the Local Stability Theorem to the equilibria.

a

p[0] = P0
	p[t + 1] = p[t] + (1/3) p[t] (1 - p[t]^2)( has equilibria at  p = -1, p = 0 and p = +1)

b

p[0] = P0
	p[t + 1] = (1/3) (p[t]^2 + 5  p[t] - 3)

c

p[0] = P0
	p[t + 1] = p[t] + (1/3) (p[t]^3 - p[t])

d

p[0] = P0
	p[t + 1] = p[t] + (1/3) (p[t] - p[t]^3)

e

p[0] = P0
	p[t + 1] = 2 Sin[p[t]]

Exercise CD20.5.4

Show that the following system fails the local stability criterion for both of its equilibria.  Experiment with its behavior by using the cobWeb[.] function.

p[0]  =  P0
	p[t + 1] = p[t] + 2.58 p[t] - 3.58 p[t]^2
		=3.58 p[t] (1-p[t])