Name:

Section 3.2 Maple commands for component curves, solution curves, orbits and time-state curves.

by Tom Linton.

While DEplot and its 3 dimensional counterpart DEplot3d can handle all of these plotting tasks, they require rewriting the second order IVP as a system for some plots in section 3.2. These commands are part of the DEtools package, and we'll utilize the spacecurve command in the plots package as well:

> map(with, {DEtools, plots}):

The normal form for a second order IVP is:

y'' = F(t, y, y'), y(t0) = y0, y'(t0) = v0.

This can be translated to the first order system (where both y and v are the unknown functions):

y' = v, y(t0) = y0,

v' = F(t, y, v), v(t0) = v0.

which is legal for input into more of Maple's plotting commands than is the original version. Let's explore the ODE of problem 1d on page 164 of the Borrelli and Coleman text.

> ode1:=diff(y(t),t,t)+0.2*diff(y(t),t)+10*y(t)+0.2*y(t)^3=-9.8;

Solution curves:

These can be handled in a variety of ways. To see a variety of solutions to various IVPs associated with the ODE, DEplot is handy, and can take the original form of the ODE as input. The initial value of y' uses the differential operator "D". It differentiates its input and returns a function for the derivative. This means that D(y) is a function, and can be evaluated at 2.0 by the command D(y)(2.0). It looks funny, but that is how you need to specify initial values of y'.

> DEplot(ode1,y(t),t=0..40,
[[y(0)=0,D(y)(0)=0],[y(0)=-4,D(y)(0)=0],[y(0)=1,D(y)(0)=.5]],
y=-5..3,linecolor=[red,blue,black],
stepsize=0.1);

Here is a handy way to zoom in on such a DEplot . The number 1 corresponds to the horizontal variable (t) and 2 to the vertical variable (y). The percent sign, %, means "the last executed output". Thus, we're asking Maple to zoom in the last plot from t = 5 to t = 10.

> zoom(%,1=5..10);

Notice (as is allowed since we are not looking in [t,y,y'] space) the solutions cross, but with different slopes (so the values of y' are different, and they wouldn't cross in [t,y,y'] space).

If you need to evaluate a solution at many points (and perhaps plot it too), using dsolve(..., numeric,...) and turning the solution into a function will work well. In this case, a single initial condition is used, inside the first parameter to dsolve .

> ans1:=dsolve({ode1,y(0)=-4,D(y)(0)=0},y(t),numeric):

Make ans1 into a function (or procedure).

> soln1:=u->subs(ans1(u),y(t)):

You can get a value just by asking for it:

> soln1(3.4);

Several y values as follows:

> map(soln1,[1.0,1.1,1.2,1.3,1.4]);

A single "triple":

> ans1(3.4);

Several triples via:

> map(ans1,[1.0,1.1,1.2,1.3,1.4]);

Notice that the values of y'(t) are given as well (from ans1). We can therefore get a procedure for y'(t) (numerically) in the same manner we make ans1 into a function for y(t):

> dsoln1:=u->subs(ans1(u),diff(y(t),t)):

Values of y'(t):

> dsoln1(1.2), dsoln1(5.3);

You can plot soln1 to get a solution curve (which is also a component curve ). Since we're plotting a procedure, we don't use "t =" in the range.

> plot(soln1,0..5,color=blue,thickness=2,
labels=[`t`,`y`]);

And, you can plot dsoln1 to get the other component curve . Labels are useful here!

> plot(dsoln1,0..5,color=green,thickness=2,
labels=[`t`,`y'`]);

You can get an orbit by parametrically plotting soln1 and dsoln1 over a given range (put the 2 functions and the range inside a list and use the standard plot command). Once again, there is no variable specified in the range.

> plot([soln1,dsoln1,0..5],
labels=[`y`, `y'`],title=`Orbit`);

To get a time-state curve , we need a procedure for z(t) = t, like this:

> tproc:=t->t:

And then (so long as the plots package is loaded) we feed spacecurve the three state variables in the order [y, y', t] (since the text plots them that way) and add a few options to get labels and axes.

> spacecurve([soln1, dsoln1, tproc],0..5,
thickness=2,axes=frame,
labels=[`y`, `y'`, `t`]);

Click on the plot above (to enclose it in an anchored frame) and then drag it around to get a better feel for what it looks like.

>

By translating the ODE into a system, you can use DEplot to get (several simultaneous) component curves and orbits. You can use DEplot3d to get (several) time-state curves.

> sys1:={diff(y(t),t)=v(t),
diff(v(t),t)=-0.2*v(t)-10*y(t)-0.2*y(t)^3-9.8};

Then, use the scene option to get components or orbits. Here are solution curves (or t-y components )

> DEplot(sys1,[y(t),v(t)],t=0..5,
[[y(0)=0,v(0)=0],[y(0)=-4,v(0)=0],[y(0)=1,v(0)=.5]],
y=-5..3,linecolor=[red,blue,black],
scene=[t,y],
stepsize=.05);

Now, t-y' components :

> DEplot(sys1,[y(t),v(t)],t=0..5,
[[y(0)=0,v(0)=0],[y(0)=-4,v(0)=0],[y(0)=1,v(0)=.5]],
y=-5..3,linecolor=[red,blue,black],
scene=[t,v],
stepsize=.05);

Orbits:

> DEplot(sys1,[y(t),v(t)],t=0..5,
[[y(0)=0,v(0)=0],[y(0)=-4,v(0)=0],[y(0)=1,v(0)=.5]],
y=-5..3,linecolor=[red,blue,black],
scene=[y,v],
stepsize=.05);

If we use DEplot3d , with its scene option, we can get time-space curves :

> DEplot3d(sys1,[y(t),v(t)],t=0..5,
[[y(0)=0,v(0)=0],[y(0)=-4,v(0)=0],[y(0)=1,v(0)=.5]],
y=-5..3,linecolor=[red,blue,black],
scene=[y,v,t],
stepsize=.05);

The homework assignment is numbers 1b,f, 2b and 3 from section 3.2 on pages 164 and 165.

>