BACK
Lab 2: Maps and Orbits

Section 1: Introduction

1.1 Flows.  Broadly speaking, a dynamical system is any system whose evolution in time is described by deterministic laws.  An example could be Newton's laws of motion.  Given the velocities and positions of particles at some instant of time, these laws enable the prediction of the velocities and positions at all future times.  In mathematical terms, such laws can be best represented by a set of ordinary differential equations:

dx/dt = F(x) ,
where x is a vector whose components are the positions and velocities of the particles. Since time t is a continuous quantity (a real number), such dynamical systems are called continuous, or flows. The solutions of the differential equations can be though of as lines in the space of  the vector x  (called the phase space) and are called orbits. Orbits starting from different points in phase space indeed resemble the flow of a fluid (the orbits cannot intersect each other), hence the name.

1.2. Maps.  Another type of dynamical systems can be imagined, where the evolution law is given by a function f(x) that maps the system from one state, say x0, to another state, say x1, as follows:  x1= f(x0).  Next, x1 is mapped into x2 according to the same rule: x2= f(x1), and so on, so that, in general

xn+1 = f(xn) .
In this case, the states of the system are numbered by an integer n which can be though of as a "time" variable.  Such dynamical system is called discrete, or simply a map.  An orbit of the discrete dynamical system is an infinite sequence of points {x0,  x1,  x2,  x3, ...}.

Section 2: Plotting orbits

During this lab session, we shall be examining the orbits generated by a simple quadratic map 

yn+1 = yn2- c

with a particular focus on how the properties of the orbits depend on the parameter c.  To get started, we will write a simple Matlab program to study the behaviour of orbits generated by this map.

Below is a program which computes an orbit of length N for a quadratic map with a given parameter c starting from a specified initial condition.  Copy this program into a file and then run it in the Matlab command window (by typing the name of the file on the command line).  You should be able to copy and paste the lines from the Explorer window into the Matlab editor without re-typing them.

Program 1

c = 0.5;  N = 50;
y = zeros(N,1);
y(1) = 0;
for n = 1:N-1,
  y(n+1) = y(n).^2 - c;
end,
plot(1:N, y, 'o-');
title('Orbit of the quadratic map');
xlabel('n'); ylabel('y_n');

Comments to the program:

 
Q2.1.  Use the program above to explore the behaviour of the orbits for different values of c and different initial conditions.  (Remember to save the file each time you change it, before running the program).  Find values of parameter c for which the orbits a) tend to the same constant value independent of the initial condition; b) oscillate between two values (period-2 orbit); c) oscillate among four values (period-4 orbit); d) escape to infinity; e) exhibit complicated (chaotic?) behaviour. Produce a figure for each case, indicating the value of parameter c and briefly characterising the orbit behaviour.


As a consequence of this experiment, we see that even a very simple dynamical system - as simple as a quadratic function - can exhibit fairly complicated behaviour.  Moreover, the behaviour changes as the parameter c changes.

Next, we shall try to identify values of the parameter c for which the behaviour of orbits changes from one type to another.  Before doing that, note that it is not very convenient to have to save the program each time you make changes to the parameters.  It would be much more convenient if the parameters could be entered directly into the program from the command line.  To this end,
 

Q2.2.  Transform the program into a Matlab function, so that the same results can be achieved by typing, for example,

>> quadratic(0.5,50,0.0);

where 'quadratic' is the name of the file containing the function.  The result should be the same as that of the program above.
 

Q2.3.   Use the function from Q2.2. or the program above to find (approximately) the values of the parameter c for which all orbits
 a) stop escaping to infinity and begin converging to a fixed point;
 b) stop converging to a fixed point and start converging to a period-2 orbit;
 c) stop converging to a period-2 orbit and start converging to a period-4 orbit.


Section 3: Graphical analysis

3.1. Graph of a function.  The goal of this section is to understand the behaviour of the orbits of a given map without having to compute the orbits, but simply by looking at the graph of the function f(x).   To plot the graph of a function in Matlab, we can use the command 'fplot'.  For instance, to plot the function f(x) = 2.5 x (1-x) in the range of x values from 0 to 1, we type

>> fplot('2.5*x*(1-x)',[0 1],'r-');

(The string 'r-' specifies that the graph should be a red line, the default is a blue line).

Alternatively, we can plot the function using 'plot' command: first specify the range of x values:

>> x = 0:0.01:1;

Next, compute the values of the function and plot the result:

>> f = 2.5.*x.*(1-x);
>> plot(x, f, 'r-');

(Note that we use '.*' for multiplication instead of just '*'.  The reason is that '*' is interpreted by Matlab as a vector multiplication.   The '.*' symbol tells Matlab to multiply vectors element by element.)

3.2. The cobweb diagram.  There is a simple geometric procedure for describing the behaviour of orbits using only the graph of f(x).  It is sometimes referred to as the cobweb diagram.  Let's say we want to show on the function graph an orbit starting from x0.  To help us do that, we first draw the diagonal line y = x, which makes a 45o angle with x- and y-axes.

>> hold on; plot([0 1],[0 1],'k-');

(The command 'hold on' holds the current plot of the function, so that the black line from (0, 0) to (1, 1) is added to it).

The next point on the orbit is the number f(x0).  The graph y = f(x) allows us to read off this point, since (x0, f(x0)) is the point on the function graph directly over x0.  To visualise this,  we can draw a vertical line from the point (x0, x0) on the diagonal to the point (x0, f(x0)) on the graph of the function f(x).

>> x0 = 0.1;
>> x = x0; f = 2.5.*x.*(1-x);
>> plot([x x], [x f]);

Next, we need to assign the number f(x0) to the next value of x:  x1 = f(x0) .  In order to do it on the graph, we draw a horizontal line from the point  (x0, f(x0))  till it meets the diagonal line x = y precisely at the point (f(x0),  f(x0)) .

>> plot([x f], [f f]);

We can continue graphing the orbit by assigning the value of f to x:

>> x = f; f = 2.5.*x.*(1-x);

and repeating the above steps again.  All this can be arranged in a loop and written as a Matlab program.  As before, copy the program into a file and run it in Matlab.  Explore different values of r, N and x0.

Program 2

% Plotting the graph of the function and the diagonal
r = 2.5;  N = 10;  x0 = 0.1;
x = 0:0.01:1;
f = r.*x.*(1-x);
clf;
plot(x,f,'r-',[0 1],[0 1],'k-');
hold on;
% Plotting the 'cobweb' diagram
f = x0;
for i = 1:N,
  x = f;
  f = r.*x.*(1-x);
  plot([x x], [x f]);
  pause;
  plot([x f], [f f]);
  pause;
end,

Comments to the program:

In the diagram, each successive point of the orbit appears on the diagonal y = x after one 'reflection' off of the function graph.
 
Q3. Modify the above program to produce the cobweb diagrams for the following maps:
a) xn+1 = xn2 - c;  c = 0.6;  x0 = 0,  in the range [-1 1];
b) xn+1 = a sin(xn);  a = 2.0;   x0 = 0.1;  in the range [0 pi];
You may need to use Matlab command 'axis([xmin xmax ymin ymax])' to set the axis range.  For example,
axis([-2 1 -4 3]) sets the axis range from -2 to 1 for the x-axis and from -4 to 3 for the y-axis.
 

For the Lab Report No. 2, give answers to all the questions above, providing texts of the programs and figures wherever necessary.

BACK