>

BDH1-7.mws

Bifurcations

Section 1-7 - Blanchard, Devaney, Hall

26 January 2000

Douglas Meade

Introduction

The purpose of this worksheet is to reinforce our in-class discussion of bifurcations (section 1.7 of the text). Definitions are not repeated in this worksheet. There are a number of questions for you to answer as you work through this worksheet. Take the time to do what is asked, and be sure to ask questions as they arise.

>

Getting started

Every Maple worksheet should begin by re-initializing the Maple "kernel" and loading the additional packages that we are most likely to use.

> restart;

> with( plots ):

> with( DEtools ):

>

One-Parameter Family with One Bifucation

The one-parameter family of differential equations introduced in the text uses

> f := y^2 - 2*y + mu;

in the differential equation

> MODEL := diff( y(t), t ) = subs( y=y(t), f );

>

To set this problem up for using DEplot , we have to define the variables and viewing window:

> VARS := { y(t) };

> DOMAIN := t=0..5;

> RANGE := y=-10..10;

>

Animated Slope Fields and Solution Curves

>

We will create an animated version of Figure 1.78 (p. 96). This animation will consist of the slope field and solution curves for several values of the parameter. The list of parameters will be created using the seq command. This command is essentially a `` do '' loop; the outer square brackets tell Maple to save this as a list .

> PARAM := [ seq( i/2, i=-8..8 ) ]; # change to -2..2 for testing

>

To maximize the impact of this animation we want to be sure that every equilibrium solution is displayed. One way to ensure this is to find all equilibria for each value of the parameter.

> IC := { seq( [0,i], i=-5..5 ) }: # start with some general IC

> for mu in PARAM do

> yROOT := solve( f=0, y ); # all equilibria for this mu

> yROOT2 := remove( has, {yROOT}, I ); # discard any complex roots

> IC := IC union { seq( [0,y0], y0 = yROOT2 ) }; # add new IC to set of ICs

> od;

>

Maple Question

Note that the IC are being kept as a set , not a list . Why do you thing I did this? (Consult the on-line help for information about sets and lists.)

>

>

The for .. do .. od command is used for the creation of the animation. (This could be done with the seq command, but it puts a lot into one command and can be rather imposing.)

>

> PLOTseq := NULL: # initialize list to NULL list

> for mu in PARAM do # create plot for each parameter value

> PLOTseq := PLOTseq, DEplot( MODEL, VARS, t=0..5, RANGE, IC, arrows=MEDIUM );

> od: # do NOT change to a semi-colon!

> display( PLOTseq, insequence=true ); # display list of plots as an animation

>

One-line version of animation (optional)

Here is the one-line version of the above animation created with seq .

> display(

> seq( DEplot( MODEL, VARS, t=0..5, RANGE, IC, arrows=MEDIUM ),

> mu=PARAM ),

> insequence=true );

>

>

Bifurcation Diagram

A bifurcation diagram can be quite difficult to create by hand. It's not easily created using Maple's plot command, but a reasonable job is done by the implicitplot command.

Cleaning Up

Prior to creating the bifurcation diagram, note that the definition of f no longer contains the parameter.

> f;

The problem is that the loops used above set the value of mu to a specific value.

> mu;

To restore mu to its status as a parameter, it has to be unassigned.

> mu := 'mu';

Now, the function is as it should be:

> f;

>

>

Recall the function that defines the RHS of the differential equation:

> f;

> implicitplot( f=0, mu=-4..4, y=-4..4, style=POINT );

Note that the plot ends with mu=1. This is because there are no equilibria for mu>1. If you really want to force Maple to create a plot with larger values of mu, use the view= argument to plot:

> implicitplot( f=0, mu=-4..4, y=-4..4, view=[-4..4, -2..4], style=POINT );

>

The pitchfork diagram in Figure 1.84 can be created as follows:

> f := y^3 - alpha*y;

>

> implicitplot( f=0, alpha=-4..4, y=-2..2, style=POINT );

>

Analytic Determination of Bifurcation Points

The discussion in the text leads to the conclusion that the only possible locations of bifurcation points of are values of and such that both and . We illustrate using the example from the text; this sequence of commands is easily modified for other examples. Note: In many cases, these computations are easy enough to do by hand.

> f := y*(1-y)^2 + mu;

> EQ1 := f = 0;

> EQ2 := diff( f, y ) = 0;

> EQ2 := factor( EQ2 );

> BIFpt := solve( { EQ1, EQ2 }, { mu, y } );

> BIFpt2 := [ [-4/27,1/3], [0,1] ];

>

> BIFdiag := implicitplot( f=0, mu=-2..2, y=-2..4, style=POINT ):

> BIFptPLOT := plot( BIFpt2, style=POINT, color=BLUE, symbol=CROSS ):

> display( BIFdiag, BIFptPLOT );

The axes get in the way of this plot. Use the context sensitive menu to change the axes (click on the graph, hold the right mouse button down, select Axes, then Boxed or Framed).

>

>

> restart;