Plotting vector fields, curves, and phase portraits

>	with(plots);with(DEtools);   # load packages

Here is the simplest way to plot a vector field. The field plotted below is V(,x,y) = ( y, 1-x^2).
After some experimenting I choice the box x = -2 to 2 and y = -1.3 to 1.3.

>	fieldplot([y,1-x^2],x=-2..2,y=-1.3..1.3);

Next we shall plot a curve. It is not a solution curve, just some curve.

>	implicitplot(x^2+y^2-2*x-y+1=0,x=-2.5..2.5,y=-2.5..2.5,color=red);

That looks a little rough. I will use an optional feature to get a better graph.

>	implicitplot(x^2+y^2-2*x-y+1=0,x=-2.5..2.5,y=-2.5..2.5,color=red,numpoints=10000);

Next we combine the two. For this we save the vector field and the curve and then display them together.
Notice I have ended the first two commends with a colon (:) instead of a semicolon (;).

>	vectorfield := fieldplot([y,1-x^2],x=-2..2,y=-1.3..1.3):

>	curve := implicitplot(x^2+y^2-2*x-1*y+1=0,x=-2.5..2.5,y=-2.5..2.5,color=red,numpoints=10000):

>	display(vectorfield,curve);

Notice that the vector field has a critical point on the circle. Let's use a different circle.

>	curve2 := implicitplot(x^2+y^2-2*x+1/2=0,x=-2.5..2.5,y=-2.5..2.5,color=green,thickness=3,numpoints=10000):

>	display(vectorfield,curve2);

Now you can compute the winding number of curve2 with respect to the vector field. In this case it is the index on the critical point (1,0).

Another command for getting a vector field is "phaseportrait". With it we can also graph solution curves.
The syntax is harder. Below we plot the same vector field. We have included two initial conditions, (-1.0, 0.2) and (-0.2, -0.2).
We have run time from -3 to 10 units.

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

Next we do an animation showing how the two solution curves change with time.

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

Below is an animation of the phase portrait of an example from Henle's textbook. See pages 61-64.
The equation is V(,x,y)=(2xy,Y^2-x^2-K^2). The parameter K goes from 1 to 0. Two solution curves
are shown. At K=0 there is a Single critical point at the origin. It is a dipole. For all other values of K
there are two critical points, (K,0) & (-K,0). Both are nodes, one repelling, the other attracting. This
qualitative change in the phase portrait is called a bifurcation. This example illustrates that a dipole
fix point is not structurally stable. Under even very small perturbations it is replaced by two critical
points that are linearizable.  

>	animate(phaseportrait,[[D(x)(t)=2*y(t)*x(t),D(y)(t)=y(t)^2-x(t)^2-K^2],[x(t),y(t)],t=-2..8, [[x(0)=0.01,y(0)=K+0.005],[x(0)=-0.1,y(0)=0.2]],x=-2..2,y=-2.0..2.0,stepsize=0.1,color=blue,linecolor=[yellow,pink]],K=1..0);

>