2simplephaseplots.mw

Chapter 9 examples: With mistakes Included!

> with(DEtools):with(plots):

Warning, the name changecoords has been redefined

Number 5, Section 9.3

x' = (2+x)(y-x) & y' = (4-x)(y+x)

Critical points are: (-2,2) & (4,4)

> phaseportrait([D(x)(t)=(2+x(t))*(y(t)-x(t)),
              D(y)(t)=(4-x(t))*(y(t)-x(t))],

[x(t),y(t)],t=0..1,

[[x(0)=1,y(0)=2]],

x=-4..6,y=-4..6,color=blue,linecolor=red);

[Plot]

But this does not look right! What's with the weird gap along x=y? For x=y we should have x'=0. So, at (1,1) we should have x'=0 and y'=6. But, that is not what we see. What! There is a typo in the y' equation: the second minus should be a plus! That explains it because the computer think x' and y' are both zero along x=y. It is vital that you learn how to spot errors. We all make typos. Let's try again.

> phaseportrait([D(x)(t)=(2+x(t))*(y(t)-x(t)),
              D(y)(t)=(4-x(t))*(y(t)+x(t))],

[x(t),y(t)],t=0..1,

[[x(0)=1,y(0)=2]],

x=-4..6,y=-4..6,color=blue,linecolor=red);


[Plot]

Cool! But it looks like I missed a critical point at (0,0). Plugging in we can see at (0,0) is indeed a critical point. Now I am going to add in a bunch more solution curves. But, I am not picking points at random. I will choose points (curves) the illustrate interesting features. It looks like (0,0) is a saddle, (-2,2) has only one real eigenvalue and (4,4) is a spiral attractor.

> phaseportrait([D(x)(t)=(2+x(t))*(y(t)-x(t)),
              D(y)(t)=(4-x(t))*(y(t)+x(t))],

[x(t),y(t)],t=-1..1,

[[x(0)=1,y(0)=2],[x(0)=-2,y(0)=3],[x(0)=-2,y(0)=-1],[x(0)=-2.1,y(0)=2],[x(0)=3,y(0)=-0.4]],

x=-4..6,y=-4..6,color=blue,linecolor=red);

[Plot]

I think I can do better than that. Notice one curve is a bit jagged. This is caused by numerical estimation errors. We'll need to change the step size.

> phaseportrait([D(x)(t)=(2+x(t))*(y(t)-x(t)),
              D(y)(t)=(4-x(t))*(y(t)+x(t))],

[x(t),y(t)],t=-1..2,

[[x(0)=1,y(0)=2],[x(0)=-2,y(0)=3],[x(0)=-2,y(0)=-1],[x(0)=-1.9,y(0)=1.8],[x(0)=-1.9,y(0)=1.784],[x(0)=-2.1,y(0)=2.3],[x(0)=3,y(0)=-1]],

x=-4..6,y=-4..6,color=blue,linecolor=red,stepsize=0.01);

[Plot]

I took me about ten tries before I got a picture I liked.

Example 1, Section 9.7. This example has a limit cycle.

x' = x+y - x(x^2 + y^2) & y' = -x + y - y(x^2 + y^2)

Think about what happens if we are on the unit circle: x^2 + y^2 = 1. Thus, x'=y & y' =-x. So, at a point (x,y) the velocity vector is (y,-x). Hence the velocity  and is perpendicular to the position vector. That implies circular motion as you undoubtedly remember from your physics or dynamics course. I turns out that for x(0) = 1 , y=(0) = 0, the parametric equation <x(t),y(t)> = <cos t, -sin t> is the solution. Check this by direct substitution.

> phaseportrait([D(x)(t)=y(t)+x(t)-x(t)*(x(t)^2+y(t)^2),D(y)(t)=-x(t)+y(t)-y(t)*(x(t)^2+y(t)^2)],
[x(t),y(t)],t=0..10,[[x(0)=1,y(0)=2],[x(0)=0,y(0)=0.2]],x=-2..2,y=-2..2,color=blue,linecolor=red,stepsize=0.1);

[Plot]

>