Mathematical Modeling & Mathematica

Lab 2:  Graphing in the Plane

John E. Cruthirds and L. Richard Hitt
Department of Mathematics and Statistics
University of South Alabama

Objective

This lab explores the capabilities of Mathematica  for graphing curves in the plane.  In addition to graphing functions of a single variable like f(x), parametric curves and polar plots will also be explored.

Introduction

In the introductory lab you explored some of the graphing capabilities of Mathematica .   We will take a closer look in this lab at the capabilities of Mathematica  with regard to graphing in the plane.  The primary tool for plane graphing is the function Plot.  The next section studies this function along with various of its options.   After skill with Plot is developed, successive sections examine graphing parametric curves in the plane, using the polar coordinate system, and animating sequences of plane graphs.

Using Plot

The Global Behaviour of Functions

Let's look at some examples using the built-in function Plot.

First, we will graph three different fourth degree polynomials with leading coefficient 1.

Clear[functionlist]
functionlist = {(x^4) - 6 (x^2) + 1, (x^4) - (x^3) + 2 (x^2) - x + 4, 
          (x^4) + 4 (x^3) - 2};

Plot[Evaluate[functionlist], {x, -2, 2}]

Technical Note:  The use of the Evaluate function is necessary because Plot uses non-standard evaluation of its arguments as explained in the Help pages for Plot.  The problem occurs because we did not give the Plot command the list of functions explicitly, but stored the list in a variable.  If you like, you can just regard this as an idiom of the Mathematica  language and not be concerned about the details.

Well, there is not much of a pattern to observe on those plots.  Let's widen the interval a little.

Plot[Evaluate[functionlist], {x, -5, 5}]

Notice that more of a pattern becomes apparent.  Finally, let's make the interval very wide.

Plot[Evaluate[functionlist], {x, -70, 70}, 
      PlotStyle -> Thickness[0.0005]]

Now the three graphs appear essentially the same.  In other words, while the local behaviour of the three given functions was quite different (differing numbers of critical points, inflection points, etc.), the global behaviour of the three functions is the essentially the same.  You will have a chance to explore the global behaviour of other functions in the exercises at the end of the lab.

The Plot Algorithm

Our goal here is to gain some insight into how the Plot function works.  At a simple level, in order to plot a function you can just plot a few points and connect them with lines.  To do this in Mathematica , you use the graphics primitive Line applied to a list of points.  The Graphics function applied to a list of graphics primitives yields a graphics object, which can be viewed with the Show command.  For example,

Show[
      Graphics[{Thickness[ .01], 
          Line[{{0, 0}, {1, 1}, {2, 0}, {3, 1}, {4, 0}}]}]]

Now, let's apply this graphing capability to graph a quadratic function by sampling some of its points and connecting them with line segments.

Clear[curve, pointlist]
pointlist = Table[{x, (x^2)}, {x, -2, 2,  .2}];
curve = 
        ListPlot[pointlist, PlotJoined -> True, 
          DisplayFunction -> Identity];
Show[curve, Graphics[{PointSize[ .01], Point /@ pointlist}], 
      DisplayFunction -> $DisplayFunction]

This gives a nice smooth result.  The points that were plotted have been emphasized by making them larger than the lines connecting them.  Notice that using more points would not significantly improve the plot of the function.

Let's try another example.  This time we will just plot the points without connecting them.

Clear[pointlist]
pointlist = 
        Table[{x, Sin[(E^x)]}, {x, 0. , 5. ,  .2}];
ListPlot[pointlist]

That generated a bit of a mess.  Let's try using more sample points.

Clear[pointlist]
pointlist = 
        Table[{x, Sin[(E^x)]}, {x, 0. , 5. ,  .05}];
ListPlot[pointlist]

Now let's connect the dots and see what we have.

ListPlot[pointlist, PlotJoined -> True]

That's a little better.  However, we have used way too many points on the interval [0,2], and not nearly enough points on the interval [3,5].  What we need is an algorithm which dynamically decides when to use more points to obtain a reasonably accurate rendition of the function.  Well, that's exactly what the Plot function does for us.  Let's try it.

g1 = Plot[Sin[(E^x)], {x, 0, 5}]

This graph gives us a pretty good idea what our function is doing, but there is still room for improving the graph.  Notice, for example, that some of the relative extrema on the graph are inaccurate.  A closer inspection will show that there is even a missing cycle in the above graph out near x = 5.  To get a handle on improving the graph, we need to explore some of the options to the Plot function.  First, let's list them together with their default values.

Options[Plot]

Three of these options determine how many sample points are used: PlotPoints, PlotDivision, and MaxBend.  The value of PlotPoints is the minimum number of points which will be used to construct the graph.  In fact, Plot begins by using this number of equally spaced points to sample the function values.  If any two consecutive line segments constructed from those sample points produce a "bend" greater than the value of MaxBend (measured in degrees), a subdivision on one of the subintervals will occur.  This process continues until all the bends are smaller than MaxBend, or until PlotDivision subdivisions have occured.  This type of algorithm is called an adaptive algorithm .  The procedure adjusts the density of sample points according to the needs of the function.

With this in mind, we can try to improve our graph by adjusting some of the options to Plot.  Of course, we could always set, for example, PlotDivision->10000 and MaxBend->1, but that would generate a huge number of points and would probably take more time to finish than we were willing to wait.  So the goal is to judiciously increase the number of plotted points so that the graph is substantially improved.  Certainly we need to increase the value of PlotDivision to allow for more subdivisions on the subintervals when needed.  However, because our graph has very low curvature (is almost a straight line) except very near its relative extrema, we should actually be able to increase the value of MaxBend without unduly effecting the accuracy of the graph.   A little trial and error shows that the following values yield a reasonable improvement on our graph.

g2 = 
      Plot[Sin[(E^x)], {x, 0. , 5. }, PlotDivision -> 100, 
        MaxBend -> 25]

That's great.  Compare the above graph with our original attempt using Plot.  To aid in the comparison, you can enlarge the graph by clicking on the graph and dragging a corner to increase the size.  Notice the improved uniformity on the values of the extrema and the additional cycle near the right end.

To see how many points Mathematica  plotted in the two versions of our graph requires just a bit of work.  We need to know a little about the structure of the output of the Plot function in order to tell at what level the information about the actual plotted points is stored.  We can get some insight into the structure of the Plot output using the following.

Shallow[InputForm[g1]]

If we had entered InputForm without Shallow, we would have gotten the entire (very) long expression for g1.  What we did get is an abbreviated form, where the <<1>> represents a single item -- a long list of points in this case.  In fact, it is this list of points we want to know the length of. The <<25>> represents 25 items which have been omitted -- Plot options and their values in this case.  To delve farther into the structure of the list of points, we can access the list directly and view it with either the Short or the Shallow command.

Short[g1[[1]]]

The result suggest there are 464 ordered pairs of points used in the plot of g1.  To find the exact location in the output of g1 of the ordered pairs, we can calculate the position of the first x-coordinate by copying and pasting it into a Position command.

Position[g1, 2.0833333333333333`*^-7]

So we see the x-coordinate is 6 levels down in the output of g1, so the ordered pairs must be 5 levels down, and the list of order pairs 4 levels.  We can let Mathematica  compute the number of points for us now that we know where the points are.

Length[
      g1[[1, 1, 1, 1]]]

Note that the total number of points used by Plot is bounded above by the product of the two values (PlotPoints - 1) and PlotDivision.  Considering the default values used by Plot, this product would be (24)(30) = 720.  So Mathematica  is using close to the upper bound of points in the default graph.  For the second graph, the number of points used was

Length[
      g2[[1, 1, 1, 1]]]

So the second graph, although vastly improved over the first, did not use many more points than the first.

Next let's look at how Mathematica  used the sample points by plotting just the points for each graph.

ListPlot[g1[[1, 1, 1, 1]], 
      PlotStyle -> PointSize[ .004]]

ListPlot[g2[[1, 1, 1, 1]], 
      PlotStyle -> PointSize[ .004]]

These plots illustrate how Mathematica  knows to plot points close together where the graph has high curvature but is content to spread the points out where the curvature is close to zero.

Exercise

Functions with Asymptotes

We will examine the way Mathematica  handles functions with a vertical asymptote using an example.

Plot[(((x^2) + 1)/(x - 1)), {x, -2, 5}]

Note that Mathematica  appears to have drawn in the vertical asymptote.  But that is not what really has happened.  Remember that Mathematica  just samples values of the function at various points determined by the options PlotPoints, PlotDivision, and MaxBend, and then connects those points with straight lines.  What appears to be the asymptote is really just one of those line segments connecting a point on the left branch to a point on the right branch.  The reason there is a break in the lines is that Mathematica  has truncated the graph at the top and bottom.  To view all of the calculated lines, you can use the option PlotRange->All.

When you want to graph a function with a vertical asymptote,  it is usually better to split the graph along the asymptote.  You can also add in the vertical asymptote separately, as well as any other asymptotes as follows.

Clear[g1, g2, g3]g1 = 
      Plot[(((x^2) + 1)/(x - 1)), {x, -2,  .9}, 
        DisplayFunction -> Identity]; g2 = 
      Plot[(((x^2) + 1)/(x - 1)), {x, 1.1, 4}, 
        DisplayFunction -> Identity]; g3 = 
      Plot[x + 1, {x, -2, 4}, 
        PlotStyle -> {Dashing[{ .02,  .02}], RGBColor[1, 0, 0]}, 
        DisplayFunction -> Identity]; Show[g1, g2, g3, 
      Graphics[{Dashing[{ .02,  .02}], RGBColor[0, 0, 1], 
          Line[{{1, -100}, {1, 100}}]}], 
      DisplayFunction -> $DisplayFunction]

Miscellaneous Options for Plot

Consider the following example.

Plot[(x^2) + 1, {x, -2, 2}]

Notice that the horizontal axis is not the x-axis, but the line y=1.  Mathematica  uses its own discretion in placing the axes.  If you want to override Mathematica 's decisions, you can use the PlotRange option of Plot.

Plot[(x^2) + 1, {x, -2, 2}, PlotRange -> {0, 5}, 
      AxesLabel -> {"x", "y"}]

The next example illustrates several of the Plot options.

Plot[Cos[x]^2, {x, 0, 10}, Background -> RGBColor[0,  .5, 1], 
      DefaultColor -> GrayLevel[1], GridLines -> Automatic]

Parametric Plots

Curves in the plane can also be described by pairs of parametric equations of the form 
x = f(t), y = g(t).  We begin by examining different ways that circles can be parameterized.

Parameterizations of Circles

Let's consider two different ways to graph a circle.  For definiteness, let's work with the circle centered at the origin having radius 2.  An equation of the circle under consideration is given by

x2 + y2 = 4

To use the built-in Mathematica  function Plot to graph the above equation, we need to solve it for y:

Clear[x, y, solution]
solution = Solve[(x^2) + (y^2) == 4, y]

Notice the form that Mathematica  uses to present the two solutions: it is a list of two elements, each one of which is a list containing a rule.  To access the two functions so we can plot them we can use the ReplaceAll function (abbreviated /.) as follows:

y /. solution

This will produce the proper list of two functions that we can pass along to the plot function.

Clear[ytop, ybottom]
Plot[Evaluate[y /. solution], {x, -2, 2}, 
      AspectRatio -> Automatic]

Rather than plotting two different function to obtain the graph of a circle, we can a different approach -- a parametric approach -- so the circle can be plotted directly.  Recall that if you let t denote the angle made between the positive x-axis and the line from the origin to a point (x,y) on the circle,  then each point on the circle is uniquely determined by the value of t, and hence, its coordinates are functions of t given by the set of equations

x(t) = 2 cos(t)
y(t) = 2 sin(t)

We can now use the built-in function ParametricPlot to graph the circle:

ParametricPlot[{2 Cos[t], 2 Sin[t]}, {t, 0, 2 π}, 
      AspectRatio -> Automatic, AxesLabel -> {x, y}]

Well, that was a lot easier than our first method.  In the parametric equations used above, the value of the parameter t has a geometric interpertation, namely, it measures the angle indicated below.

The example above illustrates one of the main methods of obtaining a parametric set of equations for a curve  -- using an angle which is intrinsic to the curve as the parameter.  Another often encountered method is to use slopes of lines through some fixed point as the parameter.

As an example of this second method, again consider the circle of radius 2 centered at the origin.  Note that any line through the point (-2,0) of the form y = t (x+2) intersects the circle in exactly one point other than (-2,0).  Thus, each point on the circle other than (-2,0) is uniquely determined by the slope t of such lines.  We can use Mathematica  to find the corresponding set of parametric equations as follows.

Clear[solution]
solution = 
      Solve[{(x^2) + (y^2) == 4, y == t (x + 2)}, {x, y}]

Let's graph the computed set of parametric equations and see what we get.  Since we only want the second pair of solutions, we use solution[[2]] to get just the second entry of the outer list.

Clear[f]
f[t_] := {x, y} /. 
        solution[[2]]
ParametricPlot[Evaluate[f[t]], {t, -10, 10}, 
      AspectRatio -> Automatic]

As expected, we get the circle of radius 2 centered at the origin -- most of it, that is.  Do you see why there is a gap in the circle?  What could be done to reduce the size of the gap?

Lissajous Figures

A Lissajous figure is a path traced out in the plane by a particle each of whose coordinates are under simple harmonic motion.  Such trajectories are often encountered in physics.  The first set of parametric equations of the circle in the previous subsection is an example of such motion.  Let's try some other examples which generate Lissajous figures.

Clear[x, y]x[t_] := Cos[3 t]y[t_] := Sin[2 t]ParametricPlot[{x[t], y[t]}, {t, 0, 2 π}]

Clear[x, y]x[t_] := Cos[3 t]y[t_] := Cos[4 t + (π/4)]ParametricPlot[{x[t], y[t]}, {t, 0, 2 π}]

Polar Plots

In this section we will graph some curves in polar coordinates.  To begin with we need to recall the transformation equations from rectangular coordinates to polar coordinates.

x = r Cos[θ]
y = r Sin[θ]

So, to graph a polar equation like r = Sin[2θ], all we need do is

Clear[r]r := Sin[2 t]ParametricPlot[r {Cos[t], Sin[t]}, {t, 0, 2 π}, 
      AspectRatio -> Automatic]

Exercise

Exercise

Animations

Mathematica  can perform animations of sequences of graphics.  To generate such an animation, you first construct the list of graphics, often using a Do loop or similar construct.  Then you select the cell which contains the list of graphics by clicking the appropriate cell bracket.  Next, from the Graph pull down menu, you select Animate Selected Graphics.  At that point, you should begin to see an animation.  You can control various aspects of the animation by clicking on buttons on the lower left of the window.  Consult the Mathematica User's Guide for the Macintosh  for details.

The following example illustrates various limacons -- those with loops, those without loops, and cardiods.  As the parameter a varies, the limacon changes in a way which can be viewed by animating the generated graphics.

Do[ParametricPlot[
        (1 + a Cos[theta]) {Cos[theta], Sin[theta]}, {theta, 0, 
          2 π}, PlotRange -> {{-4.5, 4.5}, {-2.5, 2.5}}], {
        a, -3.5, 3.5,  .5}]

For another animation, we can expand on the Lissajous figures introduced in the previous section.  We will watch what happens as we shift the phase of the harmonic motion on the y-axis for a particular Lissajous figure.  Animate the result of the following code and watch what happens.

x[t_] := Cos[2 t]y[t_, a_] := Cos[5 (t + a)]Do[ParametricPlot[{x[t], y[t, a]}, {t, 0, 2 π}, 
        AspectRatio -> Automatic], {a, 0, (π/5) - (π/120), 
        (π/120)}]

It is not completely obvious how the interval for the parameter a was selected.  The idea was to use the smallest interval which brings the graph back to where it started.

Implementation

Here is the Mathematica  code that was used to generate one of the circle graphics in the Parametric Plots section.

Clear[f, t, a]
a = 1;
f[t_]  := {2 Cos[t], 2 Sin[t]};
ellipse = 
        ParametricPlot[f[t], {t, 0, 2 π}, AspectRatio -> Automatic, 
          PlotStyle -> RGBColor[0, 0, 1], 
          DisplayFunction -> Identity];
point = 
        Graphics[{RGBColor[1, 0, 0], PointSize[0.03], 
            Point[f[a]]}];
pointlabel = 
        Graphics[Text["(x, y)", f[a] + { .25,  .25}]];
line = Graphics[Line[{{0, 0}, f[a]}]];
angle = 
        Graphics[
          Circle[{0, 0},  .4, {0, 
              ArcTan[f[a][[2]]/(f[a])[[1]]]}]];
anglelabel = Graphics[Text["t", { .5,  .2}]];
Show[ellipse, point, pointlabel, line, angle, anglelabel, 
      AxesLabel -> {"x", "y"}, AspectRatio -> Automatic, 
      Ticks -> {{-2, -1, 1, 2}, {-1, 1}}, 
      DisplayFunction -> $DisplayFunction]

Exercise 1

Exercise 2

Exercise 3

f (x) = (((x^3) + x + 1)/((x^2) - x - 6))

Exercise 4

((x^2)/9) + ((y^2)/4) = 1

((x^2)/9) - ((y^2)/4) = 1

Exercise 5

(x^3) + (y^3)  = 6 x y.

Exercise 6