Mathematica: Beyond the Basics

Workshop notes

Basic animation: Turning a secant line into a tangent line

Let's say we want to animate a secant line to a curve at a point turning into a tangent line as Δx -> 0.  Now, you could make a really fancy version of this that does everything automatically for you (try it once you get the hang of things), but let's set this up step by step.

Solution

First, let's define our curve:

f[x_] := Sin[(x^2)]

Our point of tangency:

xPt = 0

So, for a given Δx, the slope of the secant line is:

slope := ((f[xPt + Δx] - f[xPt])/Δ)x

This works great until Δx = 0, at which time bad things happen...  So, a way to avoid this is to use the following trick:

slope := ((f[xPt + Δx] - f[xPt])/Δ)x /; Δx ≠ 0

slope := f'[xPt] /; Δx == 0

Thus, the equation of the secant line is:

secLine[x_] := slope (x - xPt) + f[xPt]

Table[
      Plot[{f[x], secLine[x]}, {x, -1.5, 1.5}], {Δx, 1, 
        0, -0.1}]

This gets the basic idea across, but it jumps around a lot.  The way to fix this is to use the PlotRange option:

Table[
      Plot[{f[x], secLine[x]}, {x, -1.5, 1.5}, 
        PlotRange -> {-0.5, 2}], {Δx, 1, 
        0, -0.1}]

Let's spice it up a little more with different colors:

<< Graphics`

graphs = 
      Table[Plot[{f[x], secLine[x]}, {x, -1.5, 1.5}, 
          PlotStyle -> {Red, Green}, 
          PlotRange -> {-0.5, 2}], {Δx, 1, 
          0, -0.1}]

Another cool trick is to plot all of these together on the same axes:

Show[graphs]

Fancier animation: Stewart Calculus problem with limits

Consider a fixed circle (Subscript[C, 1])with equation (((x - 1))^2) + (y^2) = 1 and a shrinking circle (Subscript[C, 2]) with radius r and center at the origin.  P is the point (0, r), Q is the upper point of intersection of the two circles, and R is the point of intersection of the line PQ and the x-axis.  What happens to R as (Subscript[C, 2]) shrinks, that is as \(TraditionalForm\`r -> \(0\^+\)\)?  (Stewart, James - Calculus Concepts and Contexts, Section 2.3, problem 44, pg. 119)

(See How I made the diagram below to see how this figure was generated.  It was a bit of a pain, but does make for a nice illustration.)

Solution

<< Graphics`

First, find the point of intersection (Q):

intersections = 
      Solve[{(((x - 1))^2) + (y^2) == 1, 
          (x^2) + (y^2) == (r^2)}, {x, y}]

pointQ = {x, y} /. intersections[[2]]

We know point P by inspection:

pointP = {0, r}

So, the equation for the line is:

slope = ((pointQ[[2]] - pointP[[2]])/(pointQ[[1]] - pointP[[1]]))

lineEq[x_] = 
      slope (x - pointP[[1]]) + pointP[[2]]

pointR = {x, 0} /. Solve[lineEq[x] == 0, x] // 
        FullSimplify

This has an extra set of curly brackets (and that will matter later), but we can fix that by using Flatten:

pointR = pointR // Flatten

So, now for the grand finale...  As r -> 0,

Limit[pointR, r -> 0]

I have to admit that this didn't fit in with my initial guess about what would happen (why 4?).  So, to help visualize this, I set up an animation:

exampleAnimation = 
      Table[ImplicitPlot[
          Evaluate[{(((x - 1))^2) + (y^2) == 1, 
              (x^2) + (y^2) == (r^2), 
              y == lineEq[x]}], {x, -1.5, 4.5}, {y, -1.2, 1.2}, 
          PlotRange -> {-1.4, 1.4}], {r, 1. , 
          0, -0.2}]

The errors are due to round-off induced division by 0 (i.e., not really a problem).  However, this animation looks lousy.  Even after explicitly giving a PlotRange (usually this fixes the "jumpiness" here), it still jumps around all over the place.  (I consider this a bug in Mathematica; it isn't actually doing what I asked it to.)  However, you can fix this by using the command ShowAnimation (this is why I named the above animation):

ShowAnimation[exampleAnimation, 
      PlotRange -> {{-1.4, 4.5}, {-1.2, 1.2}}]

How I made the diagram

You must have already executed all the definitions above in the solution (pointP, pointQ, etc.) first.

firstframe = 
      ImplicitPlot[
        Evaluate[{(((x - 1))^2) + (y^2) == 1, 
              (x^2) + (y^2) == (r^2), y == lineEq[x]} /. 
            r -> 1], {x, -1.5, 4.5}, {y, -1.2, 1.2}, 
        PlotRange -> {-1.4, 1.4}]

Show[firstframe, 
      Graphics[Evaluate[{PointSize[ .02], Point[pointP], 
              Text["P", pointP + {0.1, 0.15}], Point[pointQ], 
              Text["Q", pointQ + {0, 0.15}], Point[pointR], 
              Text["R", pointR + {0, 0.15}]}] /. r -> 1]]

Building a program: Faking the data

When I was writing a lab on least-square curve fitting for Calculus III, I ran into a problem.  I needed a bunch of data that fit the curves I was working with (mostly exponential, limited, and logistic growth models) but looked at least moderately realistic.  Of course, one option would have been to look it up or to use some sort of computer based data collection technique to get "real" data.  On the other hand, we can think of this as an experiment in virtual reality, so let's generate our own data. 

What I want to design here is a program that will accept a generating function (provided by the user) and some sort of measure of "perturbation" and generate a list of data that could then be manipulated by the students.  Also, it would be nice to plot this data so I can see how closely it fits the generating curve (so I don't make it too far out or too obvious).  Ideally, it should be able to do this for any kind of function I choose.  Let's further assume that I want to ensure that all my data is positive (not necessary, but why not?).

Solution

There are certainly a lot of ways this could be done and the "solution" I give here could be expanded significantly if you so desire (or made more statistically rigourous and/or useful).  

When you want to make one function dependent on another function, there are different ways to do it.  Each of these has some advantages and disadvantages.  For this particular case, we are going to develop a one command function that does everything we want (i.e., you don't have to define your generating function in advance).  

First, let's define a function, perturbFunc, that takes your generating function, the name of your independent variable, plugs in a value for x, and then displaces that value by some (pseudo) random amount:

perturbFunc[funct_, varName_, valVar_, pertCoef_] := 
      Abs[(Random[Real, {-pertCoef, pertCoef}]) + funct /. 
          varName -> valVar]

Here, Random[Real, {-pertCoef, pertCoef}] gives you a pseudo-random real between -pertCoef and +pertCoef.  Also, notice that I take the absolute value of everything to make sure nothing becomes negative (we could have skipped this if negative values are acceptable).  For certain cases, this might really screw up your data (if perturbCoef is too big compared to your function values, for example).

Let's try this out for a few different functions and values. The third argument to this function (perturbCoef) is a way to control the maximum amount the value will vary from the value of your function.

perturbFunc[(x^2), x, 1, 1]

perturbFunc[(x^2), x, 1, 1]

perturbFunc[(x^2), x, 2, 5]

perturbFunc[(x^2), x, 2, 5]

Now, we need to generate a table (or "list" in Mathematica-speak) of data and graph it along with your original function:

makeData[function_, {variable_, startValue_, endValue_, howManyPts_}, 
        perturbCoef_] := 
      Table[{x, perturbFunc[function, variable, x, perturbCoef]}, {x, 
            startValue, 
            endValue, ((endValue - startValue)/(howManyPts - 1))}] // 
        N

dataList = makeData[(x^2), {x, 0, 10, 21}, 5]

We can plot this using ListPlot to see how it looks:

<< Graphics`

dataPlot = ListPlot[dataList, PlotStyle -> Red]

What about the original function?

curvePlot = 
      Plot[(x^2), {x, 0, 10}, PlotStyle -> LightBlue]

Now let's superimpose them:

Show[curvePlot, dataPlot]

Okay, now it's time to put all this together into one big function that will take your input and return the list of data (along with a convenient plot) that you can then work with.  There are a few changes that need to be made, but we will try to keep them as similar to the steps above as possible.  Comments to explain what is going on are enclosed in the code (* like this *).  

Notice that the entire "program" definition is enclosed in a Module[...] command and statements are separated by semicolons.  The three variables listed at the start of the module {dataList, dataPlot, curvePlot} are all "local" variables that we use inside our program; we won't use them again outside of the program, so by "declaring" them in the module like this, we make sure they don't cause us problems later (i.e., we don't have to Clear them).

\(fakeData[function_, {variable_, startValue_, endValue_, howManyPts_}, 
        perturbCoef_] := 
      Module[
 (*\ First, \ 
          we\ list\ the\ variables\ defined\ internal\ to\ this\ program, \ 
          called\ "local variables"\ in\ computer\ \(\(programming\)\(:\)\)\ *) 

{dataList, dataPlot, 
          curvePlot}, 

 (*\ 
          This\ seeds\ the\ random\ number\ generator\ so\ your\ results\ don' t\ repeat\ themselves\ *) 
SeedRandom[]; 

 (*\ Now, \ 
          we\ define\ our\ \(\(function\)\(:\)\)\ *) 

perturbFunc[funct_, varName_, valVar_, pertCoef_] := 
          Abs[\((Random[Real, {\(-pertCoef\), pertCoef}])\) + funct /. 
              varName -> 
                valVar]; 

 (*\ Now, \ 
          we\ start\ making\ use\ of\ these . \ \ Notice\ that\ I\ don' 
            t\ define\ the\ makeData\ function\ this\ time; \ 
          I\ just\ use\ the\ right - hand - 
            side\ of\ that\ function\ to\ generate\ the\ data\ \(\(below\)\(:\)\)\ *) 

dataList = 
          Table[{x, perturbFunc[function, variable, x, perturbCoef]}, {x, 
                startValue, 
                endValue, \(endValue - startValue\)\/\(howManyPts - 1\)}] // 
            N; 

 (*\ Plot\ the\ data, \ 
          but\ don' 
              t\ display\ it\ yet . \ 
                You\ do\ this\ with\ the\ DisplayFunction -> 
            Identity\ option\ *) 

dataPlot = 
          ListPlot[dataList, PlotStyle -> Red, 
            DisplayFunction -> 
              Identity]; 

curvePlot = 
          Plot[function, {variable, startValue, endValue}, 
            PlotStyle -> LightBlue, 
            DisplayFunction -> 
              Identity]; 

 (*\ Now, \ 
          put\ them\ both\ together\ on\ a\ graph\ and\ display\ it; \ 
          to\ do\ this, \ 
          you\ must\ use\ the\ option\ DisplayFunction -> \(\($DisplayFunction\)\(:\)\)\ *) 

Show[
          dataPlot, curvePlot, 
          DisplayFunction -> $DisplayFunction]; 

 (*\ 
          The\ actual\ value\ returned\ by\ this\ function\ is\ the\ last\ result\ output\ by\ it . \ \ By\ omitting\ the\ semicolon\ after\ the\ next\ command, \ 
          we\ ensure\ that\ the\ list\ of\ data\ we\ generated\ is\ the\ "value"\ returned\ by\ this\ \(\(function\)\(:\)\)\ *) 

dataList]\)

The following generates fake data for the function f(x) = (x^2) and generates 21 data points evenly spaced between x = 0 and x = 10.  It will randomly perturb the y values from their correct values by ±6:

fakeData[(x^2), {x, 0, 10, 21}, 6]

You can also give a name to the results for later manipulation:

joeData = 
      fakeData[9 - 4  (E^(-0.35 t)), {t, 0, 20, 51}, 
        0.5]

Let's list this last set of data in a nice table form:

joeData // TableForm

Or, let's graph the data with the origin where it should be and bigger dots:

ListPlot[joeData, PlotRange -> {0, 10}, 
      PlotStyle -> {Red, PointSize[0.01]}]

Exercises for further study (i.e., homework)

Some of these I have already done myself, others just seemed like a good idea at the time...

Easier problems:

Animate a Sin curve whose amplitude varies over time.

Animate the "beating" of two Sin curves (i.e., when their period differs by only a small amount, with one period held steady and the other varied).

Animate a tangent line moving along a Sin curve.

Trickier problems:

Write a function that will animate a tangent line moving along an arbitrary function (provided by the user).

Write a function that will animate the rotation of a two-dimensional graph around the X axis to produce a surface.  (Bonus: Try the same thing around the Y axis or, better yet, around an arbitrary vertical or horizontal line.)

Develop an animation to graphically demonstrate Newton's Method for an arbitrary function (supplied by the user).

Develop a comparison of the Binomial distribution, Poisson's distribution, and the Normal distribution, demonstrating when they are similar and when they differ significantly, using tables of data and overlaid graphs.  (This one isn't so difficult to program, but takes some thought about what exactly it is you want to do...)