Math 454: Ordinary Differential Equations and Stability Theory (Intro to Dynamical Systems)
Welcome to the Math 454 homepage. Here you will find information on homework, projects, and computer resources.
Suggested Homework
These homework problems are recommended. Every 2-3 classes, a quiz will be given, based upon problems corresponding to the most recent lectures. Note that the order of the problems listed may not be the same as the lecture order.- Chapter 2: 2.1.1-4; 2.2.4-6,10; 2.3.2-3; 2.4.3-6,9; 2.5.3-4; 2.7.1,6,7;
- Chapter 3: 3.1.1,2,4,5; 3.2.2-4; 3.4.1,2,5,7,9,11; 3.5.2,3,7,8; 3.6.1-3;
- Chapter 4: 4.1.5-7; 4.3.1,7,10 ; 4.5.1;
- Chapter 5: 5.1.4, 5, 9; 5.2.3,5,7,9,12
- Chapter 6: 6.1.2,4,6; 6.2.1; 6.3.1-4,10(a-c); 6.4.4(c-e),7; 6.5.1,2,8,9,19(c-d); 6.8.3,4,7,12
- Chapter 7: 7.1.1,3; 7.2.9,12; 7.3.1,7,10; 7.5.1,3;
(you might find pplane to be useful for some of these) - Chapter 8: 8.1.1,7,8,11; 8.2.1,5, 13; 8.6.2(a)
(you might find pplane to be useful for some of these) - Chapter 10 and 8.7: 10.1.6,10,11; 8.7.2,9; 10.3.3,4,7; 10.5.1-3;
- Chapter 9: 9.2.1,3,6; 9.3.1,8;
Projects
Due dates for these will be posted as necessary. Project #1: MATLAB introduction and numerical solutions
(DUE Sept 8)
This project is meant to be a gentle introduction to some of the features
of MATLAB, as well as getting some experience with numerical solutions
of dynamical systems.
First, look at this Matlab Runge-Kutta
implementation and save it as a file "rk.m". Read the code carefully
to see what's going on.Notice that many things are done for you. You will need to edit two things:
(1) The loop over initial conditions. You need to select a sensible range which demonstrates features like fixed point solutions.
(2) The function f(x) in the ODE dx/dt = f(x). Once you understand what the code is doing, launch MATLAB, and type 'rk' to run the script.
Do problems 2.8.2 (a-d) in Strogatz computationally instead of by hand. For each problem, identify the fixed points and stability. Write up your results with a little discussion of each problem following the guidelines explained in the syllabus.
Finally, look at the MATLAB tutorial on this web site to get a more in-depth look at what MATLAB does. We will not, of course, use most of the features, but it's nice to see what can be done.
Project #2: Finding bifurcations in a fishery model (DUE Oct 1)
This project is based on Exercise 3.7.4, which involves a model that has logistic growth (for example in a fish population) together with a harvesting term. Obtain the MATLAB function bifurcation.m , which has already been set up for the model with a = 1/2 and "h" as the bifurcation parameter. This program plots the bifurcation digram for you (!) and lets you interactively see how this corresponds to the phase portrait. You may find this is useful to develop intuition, but your answers to the following also need to be backed up by hand-written analysis.
The parts of this project are as follows:
(1) Include exercise parts (a)-(e) in your write-up. (part a hint: consider what happens to the harvesting term if N is much bigger or smaller than a). These involve some analytic work, and you need to explain and justify your answers.
(2) If the parameter "a" is changed to 1, the two bifurcations happen to merge, creating a different bifurcation type. What is this type? Show a picture to justify.
(3) Finally, a little science-based policy-making: What rate of harvesting should be allowed to ensure a positive, stable equilibrium? What rate should be allowed to ensure that this equilibrium is attracting for ALL positive initial conditions?
Project #3: Love Affairs (DUE Oct. 13 )
Read section 5.3, which is an amusing application of linear systems.
(1) For fun, begin with problem 5.3.1: invent nomenclature for various models.
(2) Do problems 5.3.3-5.3.6. Justify your answers with written analysis and explanations and include a little discussion of the ``real-world" implications.
Project #4: Chemical oscillators (DUE Nov 5)
There are a number of chemical reactions which are known to produce oscillations.
The oldest and most famous is the Belousov-Zhabotinsky reaction. You can view a
laboratory demonstration of it here.
This project deals with the oscillating chemical system known as the "Brusselator",
which has the form
dx/dt = 1 - (b+1) x + a x^2 y
dy/dt = bx - a x^2 y
Here x,y,a,b are all >0.
Part A. Find all fixed points and classify them using the Jacobian.
Part B. (You'll want to use PPLANE for this and part D) Choose some reasonable values for the parameters (this might take some trial and error). Graph the nullclines and the vector field, and draw in a "trapping region" by hand.
Part C. Find and explain analytically that a Hopf bifurcation occurs at some critical value of b (which will depend on the parameter a). Does the limit cycle exist for b larger or smaller than this critical value?
Part D. Plot some sample trajectories for b larger and smaller than the critical value. Identify the limit cycle using PPLANE when it exists.
Part E. Approximate (analytically) the period of the limit cycle when b is near its critical value.
Project #5: Maps, period doubling, chaos, and Lyapunov exponents (Due Dec. 1)
Here is Matlab code for plotting that groovy bifurcation diagram associated with the logistic map.A1. Change the code to compute the sine map (see text problem 10.5.6 ), and BE SURE TO ALSO CHANGE THE RANGE of the parameter r,
so that the results look similar to the logistic map. Show labeled plot(s) of your results.
A2. Explain why things look different when r > 1 (does it map to [0,1]?)
A3. Find the values of the parameter $r$ where the period doubling bifurcations occur (use the zoom feature to get these accurately). Compute the ratio of their successive differences (see text) to see how close you can get to the Feigenbaum constant.
A4. Find the point where chaos sets in, and other features, such as a window where 3-cycles live. Label these on your plot(s).
Part B: Here is Matlab code for finding the Lyapunov exponent of the attractor. It plots out a graph whose LIMIT tends toward the exponent for the logistic map.
B1. Set r=3.6 and determine the Lyapunov exponent. Show that you get the same result independent of the initial data that you use. Why is the result positive?
B2. Determine the exponent for a stable 2-cycle and 4-cycle, indicating the value of r you are using. Why are the exponents negative?
Project # 6: Complex mappings, basins of attraction, and the Gallery of Chaos (due on final exam day)
We will meet in the usual room at the final exam time (Thursday, December 17, 2009, 2:00 - 4:00 p.m.) where everyone will share their results for part E of the project (see below).This project uses the MATLAB code complexmap.m , which displays the basin of attraction for a particular complex mapping z_(n+1) = f(z_n). The mapping and parameters can be specified along with the viewing window specified by (xmin,xmax,ymin,ymax). The code iterates the mapping N times, and plots points only where the initial condition leads to a sequence for which z_N is fairly small.
A. Consider the complex mapping f(z) = z^2 +c. Argue that if z is large enough, the modulus of z_n will always tend toward infinity as n gets large. This means not every point will be in the basin of attraction.
B1. Run the code for the complex mapping f(z) = z^2 +c, where c = -.12256117 + .74486177 i. You should (after a few seconds) see the famous "Douady's Rabbit". Plot this.
B2. Change (xmin,xmax,ymin,ymax) so that you are looking at a region of radius ~0.1 near the point (-.5, .3). You should see quasi-self similarity. Plot this as well.
C1. Set c=i and make the viewing window large enough to see the whole basin of attraction, which looks like lightning. With N=40, you'll probably see nothing. This is because the basin of attraction is very "thin", and almost every initial condition will explode. Plot using N = 20 to capture initial conditions that are close enough to the basin of attraction that their iterates don't get too large.
C2. Determine (using, say, the MATLAB command line) the first several iterates starting at z=0. Is z=0 in the basin of attraction?
C3. Zoom in on zero by changing (xmin,xmax,ymin,ymax) and adjust N so that you capture a basin structure that looks self-similar.
D. There appears to be symmetry in the basins of attraction for parts A and B. Prove the following: if z_0 is in the basin of attraction (so {z_n} stays bounded) so is -z_0. (hint: do these two initial points create the same sequence?)
E. In class presentation: create a one-page peice of artwork, showing one of the following (include captions to explain them):
1. Your own complex mapping, showing the basin of attraction (either from experimentation of from a reference book),
2. Another fractal or a chaos related picture.
We will hang up our artwork during the final exam period and allow the class to view them.
MATLAB resources
Matlab Primer (PDF) This will help you get started using MATLAB.Simple phase portrait m-file
This is a script which plots 2-d phase portraits and sample trajectories. It is small and easy to customize, but requires some MATLAB knowledge.
PPLANE7 (courtesy John Polking) This is a sophisticated tool for generating 2-D phase diagrams which can identify fixed points and plot trajectories. It has a user-friendly interface, but cannot be easily customized. You'll need two files: pplane7.m and dfield7.m . You only need to run pplane7. There is also a stand alone version (google pplane to find it).
Bifurcation tool (m-file)
This script will plot a bifurcation diagram for one dimensional systems. It requires the user to edit the function definition and provide bounds for the parameter and dependent variable.