In this lab we will use "xppaut" to simulate chaotic dynamics in a few dynamical systems.
Lorenz (1963) derived a very simplified model of convection rolls in the atmosphere. The model is described in Chapter 9 of our textbook and is a set of three differential equations. The only nonlinearities in the model are the presence of two "quadratic" terms:
then execute the code in a terminal window:
xppaut lorenz.ode &
The following xppaut commands will generate a nicely color coded (for velocity) rendition of the "strange" attractor of this system:
The program "xppaut" can plot the orbits of the Poincare Return maps associated with a user defined section. In this case, since there are three differential equations the section is a two dimensional surface. Typically the section is taken as a simple surface like the plane x=a.
It takes a bit to set up things in xppaut so I've developed a couple of parameter sets for doing two different sections. Download the following two parameter sets:
These are for the sections x=0 and x=3, respectively. To load the first:
For the section x=0, the return map looks like a 1-D map. It actually isn't. It has some very fine structure. Trajectories get mixed, scaled and folded as the traverse the Lorenz attractor. This mixing and folding makes the attractor "chaotic". To better see this, load the set "Lorenz_x=3.set" and repeat the steps above. The resulting return map for the x=3 section shows how the Poincare return map is clearly 2D!
The Henon map is a common map used to examine complicated map dynamics. It has no specific applications to modelling a physical problem. The map is defined by:
xppaut can also compute orbits of discrete time maps (in all dimensions). One specifies the definition of the map problem in a slightly different syntax than for the ODE problems. Given the new syntax, one must then set up the integrator for "Discrete" mode so that xppaut knows your are examining a map problem. To see the difference in the code, download and save the Henon map code:
Now execute:
The parameter "a" in the map definition can be used to adjust the dynamics. For a range of values near a=1.4 the attractor changes shape. At lower values of a it disappears. Try:
What you should see is the Y component of an orbit plotted against n. The orbit should be eventually or asymptotically periodic to a period 4 orbit. The route to chaos via the parameter a has a period doubling cascade similar to the logistic map.
Certain simple mechanical systems can exhibit chaotic behavior. Here we'll examine a double pendulum system. One mass of mass m1 is attached by a rod to a central pivot and allowed to swing under the force of gravity. A second mass of m2 is connected to the first mass by a rod which is also allowed to pivot. Weak friction is included.
Using Lagrangian or Hamiltonian Mechanics one can derive the equations of motion. These equations include two (pivot) angles and their associated angular velocities. Thus one arrives at a fourth order (nonlinear) system of ODEs.
Also download the file doubpend.ani The latter file contains xppaut - animation code which uses the dependent variable values to be used to draw moving pictures as the system is numerically integrated.
What you should see is one of the two angles (TH2) versus time:
Now to create the animation:
What you should see is an animation of the chaotic double pendulum motion.
Do the same as in the previous example but with the files:
In the Lorenz system we will examine the sensitivity to initial conditions. To do this, one can make several copies of the Lorenz system with slightly different initial conditions. In the following code (lorenz50.ode) 50 copies are made: (x1,y1,z1), (x2,y2,z2) ..... In the code the initial conditions for all x[j] and y[j] are the same but those for z[j] are slightly different. The code:
init z[1:50]=30.00[j]1
sets z4 and z20 to initially be 30.0041 and 30.0201, respectively. The animation file lorenz50.ani then draws small (color coded circles) with centers determined from x[i](t) and z[i](t). The additional constants are introduced since the range of coordinates in the animation window are in (0,1).
The following code will compute a numerical approximation of the (CGL) equations assuming Neumann boundary conditions.
Solving PDE's can be numerically intensive. In this case the PDE has been discretized into a system of 256 coupled ODE's. So, we can't all simultaneously simulate the model. Once run, one can view the solution by color coding R=u^2+v^2 in the (x,t)-plane. This is done through the "View Axes, Array" menu (rr0,cols=128,rowskip=1,rows=4000,colskip=1). The result looks like:
Some people refer to these complex patterns as "spatio-temporal chaos". Study of such complex dynamics is cutting edge. Roughly one might think of the dynamics as being "chaotic" in both the time t and space x directions. Please turn in an explicit formula for the solution above along with a complete characterization of the parameter dependence by the end of this lab....of course I'm kidding. I only meant the formula.