Phase Plane

MA 507

November 6, 2000

Richard Hitt

Introduction

This worksheet provides Maple implementations for some of the material in Chapter 7 of Greenberg's book. These qualitative methods are useful for studying non-linear differential equations when we cannot solve to equation exactly.

The Phase Plane for the Free Oscillator

We begin with a free oscillator given by . If we set then the xy-plane is called the phase plane . A solution plotted in the phase plane is called a phase trajectory , and a collection of such trajectories is called a phase portrait for the DE.

In order to get Maple to plot phase trajectories, we need to write the DE as a system of equations solved for each of the coordinate axis variables. For example, we can use the following Maple commands to accomplish this.

> sysde := [y(t)=diff(x(t),t), diff(y(t),t)\+x(t)=0];

> with(DEtools):

> phaseportrait(sysde,[x(t),y(t)],t=0..12,[[x(0)=1,y(0)=1]],scene=[x(t),y(t)],stepsize=.1);

>

Using DEplot instead of phaseportrait produces the same graph.

Exercise: Try the last Maple command without the stepsize option and explain what you see.

One can add in as many initial points as needed to obtain a phase portrait.

> phaseportrait(sysde,[x(t),y(t)],t=0..2*Pi,
[[x(0)=1,y(0)=1],[x(0)=2,y(0)=1],[x(0)=3,y(0)=1]],scene=[x(t),y(t)],stepsize=.1);

>

Exercise: The orbits (trajectories) in the above phase portrait appear to be ellipses. Verify this (using the fact that we can solve the DE exactly).

To illustrate the relationship between the phase plane and the xt-plane that we are more accustomed to, we can plot a solution curve in the xyt-space and project the curve into various of the coordinate planes.

> IC :=[[0,3,0],[0,2.5,0],[0,2,0],[0,1.5,0],[0,1,0],[0,.5,0]];

> DEplot3d(
{diff(x(t),t)=v(t),diff(v(t),t)+9.8*x(t)=0},
[x(t),v(t)], t=0..2, IC, stepsize=.05, orientation=[0,90], linecolor=blue);

>

Exercise: Left-click and drag the graphic above to rotate the curves in 3-space. You will be able to see the xt-plane graph as well as the vt-plane graph after suitable rotation.

Stiff Oscillator

If the restoring force for the oscillator increases as a cubic function of the displacement from equilibrium rather that as a linear function, we get models for some of the non-linear oscillators that arise in certain applications. Elastic springs (rubber bands, bungee cords, etc.) fall closer to this model than the Hooke's Law model in the previous section. As a first example, suppose our DE is where we have an added cubic term.

> sysde2 := [y(t)=diff(x(t),t), diff(y(t),t)+x(t)+(x(t))^3=0];

> phaseportrait(sysde2,[x(t),y(t)],t=0..2*Pi,[[0,.25,0],[0,.5,0],[0,1,0],[0,1.5,0],[0,2,0]],scene=[x(t),y(t)],stepsize=.1);

>

When x remains close to 0, the DE will be close to the simple oscillator in the previous section, so the curves in the phase plane will be approximately elliptical. As x moves away from 0, however, the curves deviate from ellispes and become more rectangluar.

Exercise: In the above Maple plot of the phase portrait, change the t-interval to go from 0 to Pi. Interpret the result. What does this tell you about the motion for different initial conditions?

In situations where elastic devices are self-limiting, we could model the restoring force with a function that is approximately linear around 0 but which has a vertical asymptote at a particular value of x (producing a restoring force that becomes unbounded as the displacement approaches that x-value). The tangent function is such a function.

> sysde3 := [y(t)=diff(x(t),t), diff(y(t),t)+tan(x)=0];

> phaseportrait(sysde3,[x(t),y(t)],t=0..4,[[0,.25,0],[0,.5,0],[0,.75,0],[0,1,0],[0,1.25,0],[0,1.5,0]],scene=[x(t),y(t)],stepsize=.1);

Soft Spring

This time, we use . Note that if then the solution have no physical meaning for an oscillator since that would correspond to a "restoring" force in the direction away from the equilibrium.

> sysde3 := [y(t)=diff(x(t),t), diff(y(t),t)+x(t) - (x(t))^3=0];

> phaseportrait(sysde3,[x(t),y(t)],t=0..2*Pi,[[0,.25,0],[0,.5,0],[0,.75,0],[0,1,0]],scene=[x(t),y(t)],stepsize=.1,arrows=none);

>

Exercise: Again, what happens to the period of the motion as the initial conditions take x farther away from the equilibrium position? What happened to the x=1 solution?