To obtain the solution and its graph for other values of the parameter or other initial values and , just re-execute the last three commands, changing only the first one as needed. For further details on plotting, especially on simultaneous plotting of several graphs, see the handout " More on Basic Maple Commands".

2. Plotting Direction Fields and Numerical Solutions. The basic Maple command for the graphical representation of direction fields and approximate solution curves is DEplot ; it is part of a package of routines called DEtools , which must be explicitly loaded with the command

> with(DEtools):

To plot a direction field for a first-order d.e. of the form ' = , in the rectangle , , use the DEplot command in the form

DEplot(D(x)(t)=f(t,x(t)),x(t),t=a..b,x=c..d) .

The following command, for example, generates a direction field for ' = in the rectangle , :

> DEplot(D(x)(t)=x(t)^2,x(t),t=-3..3,x=-6..6);

The default number of grid points used is 20 in either coordinate direction. You can change that, if needed, by adding the option

dirgrid=[m,n]

inside the DEplot command, where m and n denote the desired number of grid points in the horizontal and vertical direction, respectively. Example:

> DEplot(D(x)(t)=x(t)^2,x(t),t=-3..3,x=-6..6,dirgrid=[15,30]);

In addition to a direction field, Maple will compute and display one or more approximate solution curves of a given d.e. if suitable initial conditions are specified in the form

{[x(t0)=x0],[x(t1)=x1],[x(t2)=x2],...} ,

immediately following the range specifications in the DEplot command. Note that the set of initial conditions is enclosed in curly braces, while each individual condition appears in brackets. Both, brackets and braces, must be present even if only a single initial condition is given.

For example, the following command produces a direction field for ' = , as before, along with the (approximate) graph of a particular solution, satisfying the initial condition :

> DEplot(D(x)(t)=x(t)^2,x(t),t=-3..3,x=-6..6,{[x(0)=1]});

Now we will add a second solution curve, satisfying the initial condition :

> DEplot(D(x)(t)=x(t)^2,x(t),t=-3..3,x=-6..6,{[x(0)=1],[x(0)=-1]});

In general, the given initial points do not have to lie inside the plot window.

If initial conditions are given, specification of the vertical range ( -range) is optional; if omitted, DEplot will choose the -range automatically, depending on where the solution curves "go". This is sometimes useful, but not always. Let's try it with the above example (same command as before, but without -range specification):

> DEplot(D(x)(t)=x(t)^2,x(t),t=-3..3,x=-6..6,{[x(0)=1],[x(0)=-1]});

Can you explain the poor outcome? (Hint: The exact solution with initial condition has a vertical asymptote at .)

If only the solution curves are of interest, display of the direction field arrows can be suppressed by adding the option

arrows=none

inside the DEplot command, following the specification of initial conditions. Example:

> DEplot(D(x)(t)=x(t)^2,x(t),t=-3..3,x=-6..6,{[x(0)=1],[x(0)=-1]},arrows=none);

Many further options are available. In fact, most of the options available with the basic plot command will also work with DEplot . Ver y useful for printing purposes is the option

linecolor=black ,

which causes the solution curves to be drawn in black rather than yellow (the default, which looks pitiful on a black-and-white printout). The color and appearance of the direction field arrows can be customized, too, but this will hardly ever be necessary. For details, consult the relevant help pages (by typing ?DEplot or ?plot,options).

Note: Options may be added to the DEplot command in any order, but must be preceded by the range specification(s) and initial condition(s) (if any). Example:

> DEplot(D(x)(t)=x(t)^2,x(t),t=-3..3,x=-6..6,{[x(0)=1],[x(0)=-1]},
dirgrid=[10,10],linecolor=black,color=green,title=`dx/dt=x^2`);

The default method that DEplot uses to compute approximate solutions of the given d.e. is "rk4", a so-called 4th order Runge-Kutta method, much more sophisticated than the classical Euler method. You can force Maple to use the Euler method by specifying the option

method=classical .

Both "rk4" and "Euler" employ a fixed stepsize; the default stepsize is ( )/20, where is the horizontal range ( -range) of the plot. Frequently, it will be necessary to reduce the stepsize, for improved numerical accuracy. This can be done by specifying the option

stepsize=h ,

where h denotes the desired stepsize, smaller than ( )/20 . Let's try the following example, where we approximate the same solutions of ' = as before, but with a larger vertical range ( ) and using the classical Euler method (instead of "rk4"):

> DEplot(D(x)(t)=x(t)^2,x(t),t=-3..3,x=-9..9,{[x(0)=1],[x(0)=-1]},method=classical,linecolor=black);

The horizontal part of the lower solution curve is clearly erroneous -- it certainly doesn't match the direction field! The (default) stepsize used in the plot is .3 (namely, ( )/20 with and ), and apparently this is not accurate enough. Reducing the stepsize should remedy the problem:

> DEplot(D(x)(t)=x(t)^2,x(t),t=-3..3,x=-9..9,{[x(0)=1],[x(0)=-1]},method=classical,stepsize=.1,linecolor=black);

To accurately track a rapidly growing or highly oscillatory solution, it may be necessary to employ a very small stepsize and compute millions of data points. Not only will this consume a lot of time -- storing millions of data points may also tax Maple's memory. Fortunately, it is not nessary to store all the computed data points (the number of points used in creating the plot is anyway limited by the resolution of computer screen and laser printer). With the option

iterations=n ,

you can instruct Maple to store and plot only every n th of the computed data points. This option, if applied judiciously, can save a lot of time and memory, with no loss of numerical accuracy. As a drastic example, we will compute and display an approximate solution of the i.v.p.

, ,

on the interval , using "rk4" (the default method) with a stepsize of .001 and 100 iterations between plotted points. (Thus, Maple will compute 15,000 data points, but store and plot only 150.)

> DEplot(D(x)(t)=1-exp(t)*x(t)^2,x(t),t=0..15,{[x(0)=0]},
stepsize=.001,iterations=100,arrows=none,linecolor=black);

Generating this plot still takes quite some time, but without the option iterations=100 , it would take about ten times as long, without any visible improvement! On the other hand, if you tried to speed things up by using a stepsize significantly greater than .001, the plot would not be reliable. Check it out!

Beyond "rk4" (the default method) and "Euler" ( method=classical ), you have a choice of several, much more sophisticated d.e. solvers. One of them is "rkf45", a so-called 4th/5th order Runge-Kutta-Fehlberg method, which can be invoked by specifying

method=rkf45

in the DEplot command. While "rk4" and "Euler" use a fixed stepsize, "rkf45" employs a variable stepsize that's automatically adjusted in accordance with the "stiffness" of the problem. As an adaptive (self-correcting) method, "rkf45" is more reliable than "rk4" or "Euler", but usually slower. But in "stiff" problems (like the previous example), "rkf45" offers a double advantage: You don't have to worry about the stepsize, and it is faster than "rk4" or "Euler"! (In general, controlling and adjusting the stepsize as needed is more efficient than using a very small stepsize throughout.)

Since "rkf45" controls the (numerical) stepsize automatically, specifying a stepsize when using this method has no effect on the numerical accuracy of the procedure; it does, however, affect the number of points plotted and, thereby, the apparent smoothness of the graph. More precisely, in conjunction with "rkf45", the specified stepsize becomes simply the horizontal distance between consecutive points plotted. For example, if we solve the previous problem with "rkf45", we should specify a stepsize of .1 to obtain the same number of points plotted as before, when we used "rk4" with a (numerical) stepsize of .001 and 100 iterations between points plotted. Let's try it:

> DEplot(D(x)(t)=1-exp(t)*x(t)^2,x(t),t=0..15,{[x(0)=0]},
method=rkf45,stepsize=.1,arrows=none,linecolor=black);

DEplot will not only solve first-order , but also higher-order i.v.p.s (although no direction fields can be generated for the latter). The following command, for example, produces an approximate solution of the i.v.p. '' + , , x'(0) = -3, whose symbolic solution is :

> DEplot(D(D(x))(t)+x(t)=0,x(t),t=-10..20,{[x(0)=2,D(x)(0)=-3]},stepsize=.2,linecolor=black);

Note that the two initial conditions appear inside a single pair of brackets (because they determine only one solution curve ). The stepsize specification was necessary, due to the oscillatory nature of the graph and the length of the plot interval; the default stepsize would have been 30/20 = 1.5 ! (Check out what you would get without stepsize specification.)

Basically, you have the same options as in the first-order case (except for those relating to direction fields, of course), and just as in the first-order case, several approximate solutions can be generated simultaneously, by specifying several sets of initial data in the form

{[x(t0)=x0,D(x)(t0)=v0],[x(t1)=x1,D(x)(t1)=v1],[x(t3)=x3,D(x)(t3)=v3],...} .

Note that, in analogy to the first-order case, the entire set of initial conditions is enclosed in curly braces, while the initial conditions for each individual solution curve appear in brackets. Example:

> DEplot(D(D(x))(t)+x(t)=0,x(t),t=-10..20,{[x(0)=2,D(x)(0)=-3],[x(0)=-3,D(x)(0)=2]},linecolor=black);

Pretty lousy! Obviously, we should reduce the stepsize or use a more flexible solver. Let's do both and also add some perks:

> DEplot(D(D(x))(t)+x(t)=0,x(t),t=-10..20,{[x(0)=2,D(x)(0)=-3],[x(0)=-3,D(x)(0)=2]},
method=rkf45,stepsize=.2,linecolor=[blue,green],
title=`x(0)=2, x'(0)=-3 (blue); x(0)=-3, x'(0)=2 (green)`,
titlefont=[TIMES,ROMAN,12]);