The purpose of this lab is to introduce you to some of the things you can do with Mathematica. Once you sit down at the machine and get the program running you get a window where you can type commands and see the answers.  You can also type text commentary like this.
	To get started try typing the command below and then executing it.  After you are done typing it you execute it by doing a Shift-Return or Shift-Enter, that is, hold the Shift key down while pressing the Return or Enter key.   You should get a nice picture.

Plot[Sin[1/x], {x, 0.01, 0.5}, PlotPoints -> 100]

You can expand the size of graphics using the mouse.  Click on the graphics object and drag the corner to expand the picture.  Another useful feature is that you can go back to a command you have previously executed, change it, and execute it again.  To get a better look at the interesting part of the graph above, try changing the range of the plot to {x,0.01,0.05}.
	Next we will plot some algebraic varieties.  Here is the four-leaved rose from the homework assignment.  First we will do a parametric plot.  Note that to type multiline command like the one below, you just use the Return or Enter key without holding the Shift key down.

r  =  Sin[2*t];
x =  r*Cos[t];
y  =  r*Sin[t];
ParametricPlot[{x, y}, {t, 0, 2*Pi}, AspectRatio -> 1, 
      PlotPoints -> 100]

Note:  a semicolon after a commands prevents the output from being printed. Inside the commands you set various options by using the arrow -> .  To get some help on a command, you can do the following:

?ParametricPlot

There is also a help menu,  To get a list of options for a command do this:

Options[ParametricPlot]

After executing the commands above, the variables r, x, y which were previously treated as just symbols have now been assigned values. For example, try typing this.

f  =  x^2 + 1

You can set a new value for x by just typing x = (whatever).  If you want to unassign x so that it again becomes a pure symbol, you type x = . as in the next cell.  Then you can enter an expression using the symbols again.

x =. ; y =. ; r =. ;
f  =  (x^2 + y^2)^3 - 4*x^2*y^2

The four-leaved rose is the algebraic variety V(f). To plot this variety directly from its equation you use a ContourPlot.

ContourPlot[f, {x, -1, 1}, {y, -1, 1}, 
      Contours -> {0, 0.01, -0.01}, ContourShading -> False, 
      PlotPoints -> 100, AspectRatio -> 1]

The program has difficulty near the origin.  You can try increasing the number of PlotPoints but the result is not really satisfying.  One other option is to put in more contours, that is to plot the curves where f = c as well as where f=0.  Try  Contours->{0, 0.01, -0.01}.  The curve changes a lot when you make a small change in c.

For another example we will make a movie of some Cassini curves.  These are obtained by setting the product of the distances to two fixed points to a constant.  Using the two points (1,0) and (-1,0) and setting the product of the distances equal to 2a give (after squaring) a polynomial f(x,y):

f  =  ((x + 1)^2 + y^2) ((x - 1)^2 + y^2) - 
        4*a^2

You can make the computer multiply this out for you by using the command Expand[ ].  We will reset the value of f to be the expanded version of itself.

f =  Expand[f]

To plot it you need to set a value for the variable a.  There is a way to set the value of a variable in an expression like f and still keep the variable around as a symbol.  This is called making a substitution and it is done as in  the command below.

f1  =  f /. a -> 1

Now we can plot f1 and see the curve.  Instead of typing this command you could go back up to the ContourPlot cell above and "cut and paste" using the mouse and the Edit menu.  You still have to edit the command to change the old variable f to f1 and fix up the plot ranges.

ContourPlot[f1, {x, -2, 2}, {y, -2, 2}, Contours -> {0}, 
      ContourShading -> False, PlotPoints -> 100, 
      AspectRatio -> 1]

Here is a sequence of ContourPlots for twenty different values of a.  You do it by putting the plot command inside a Table[   ] command. The Table command makes the variable a run over the range from 0 to 1 in steps of 0.05.

Table[
      ContourPlot[f, {x, -2, 2}, {y, -2, 2}, Contours -> {0}, 
        ContourShading -> False, PlotPoints -> 100, 
        AspectRatio -> 1], {a, 0, 1, 0.05}]

Now you can make a movie out of this as follows.  Use the mouse to click on the outer cell bracket surrounding the pictures. This will select them all.  Now go to the Cell menu and choose the command Animate Selected Graphics.   You will get "Cassini -- the Movie".

Here is a parametric surface from the text.  You need to specify the three Cartesian coordinates as functions of two parameters, say t and u.

x  =  t*(u^2 - t^2)
y = u
z  =  u^2 - t^2
ParametricPlot3D[{x, y, z}, {t, -1, 1}, {u, -1, 1}]

Now we turn to some of the more algebraic features of Mathematica.  One of the most basic is just the ability to do exact computations with integers and rational numbers.  Here are some large integers.  It is possible to get an overflow error if you push too far.  The largest number you can handle is $MaxNumber below.

11^1111
2^(2^16) - 1
$MaxNumber

Rational numbers are automatically reduced.

12345678910/111213141516

The Euclidean algorithm and its extended version allow you to find d = GCD[m,n] and to express it as an integer combination of m and n.

m  =  20978;
n  =  1479;
GCD[m, n]
ExtendedGCD[m, n]

The answer for the extended GCD is a list  {d, {a, b}} where d  = a*m+b*n.  To get at the entries of this list you can do a multiple assignment of variables.

{d, {a, b}}  =  ExtendedGCD[m, n]
d
a
b

Now you can test whether it worked

a*m + b*n

You can also find the GCD  or  ExtendedGCD of more than two numbers.  For example, to understand the ideal of integers  < 6, 8, 9 >  you could do the following.

ExtendedGCD[6, 8, 9]

This shows that < 6, 8, 9 > = <1>  which means that every integer can be written as an integer combination  a*6 + b*8 + c*9.

Next we will try some polynomial algebra.  The are basic commands are Expand and Factor.  Factor finds a factorization of a polynomial with rational number coefficients.

x =. ; y =. ;
f  =  Expand[(x + y)^20]
g =  Factor[4  x^5 + 4  x^4 + 5  x^3 + 8  x^2 + 5  x + 1]

You can differentiate polynomials and other functions using the D[ ] command.  Here are the partial derivatives of f = (x+y)^20.

D[f, x]
D[f, y] // Factor

In the last command we used the so-called "Postfix"  way of applying a command.  It gives the same result as Factor[D[f,y]] would have done (but with fewer square brackets).

You can read off the coefficients of a polynomial as follows.

Coefficient[f, x^18  y^2]
Coefficient[f, x^18]
CoefficientList[f, x]

Next we turn to division of one-variable polynomial and the polynomial GCD.  First we enter a couple of polynomials to play with.

f  =  x^5 + 4  x^2 - 3
g =  x^3 + 1

The division algorithm finds polynomials q and r such that  f = g*q + r.

q =  PolynomialQuotient[f, g, x]
r = PolynomialRemainder[f, g, x]

To test is out try:

g*q  +  r  // Expand

The theory discussed in class showed that the ideal < f, g > generated by f and g in the polynomial ring Q[x] has a single generator:  < f, g > = < h >  where h is a greatest common divisor of f and g.  There is a built-in command for finding this in Mathematica.

h  =  PolynomialGCD[f, g]

This polynomial h should be a common divisor of f and g of highest possible degree.  To see that it is a common divisor you could try taking the PolynomialRemainder or else try Factoring.

f // Factor
g // Factor

Since < f, g> = < h > you should be able to express h as a polynomial combination of f and g. This is accomplished the PolynomialExtendedGCD[ ] command.  Unfortunately this command is not part of the standard Mathematica function library. You have to load in a special package by using the following command.  Note the slanted single quote marks at both ends.

<< Algebra`PolynomialExtendedGCD`
?PolynomialExtendedGCD

Now you can use the command to find both the GCD and the coefficients such that h = a*f + b*g.

{h, {a, b}}  =   PolynomialExtendedGCD[f, 
        g]
h
a
b

And to test whether it worked:

a*f + b*g // Expand

A variety in a one-dimensional space consists of finitely many points -- the common roots of the polynomials defining the variety.  These are just the roots of the GCD of the polynomials.  Here is an example.  Define a variety V(f, g) where f and g are defined below.  Of course you can probably see before you start what the answer is going to be but the computation will illustrate the theory.

f  =  (x - 1) (x - 2)^2 (x^2 + 1) (x + 7) // 
        Expand
g  =  x (x - 2)^2 (x^2 + 1)^2 // Expand

h = PolynomialGCD[f, g]
Factor[h]

The variety depends on what field you are using.  Over the rational numbers Q the factorization above reveals all the roots.  Namely, there is a single rational root  x = 2 and the variety is just that one point.  Over the complex numbers there are more roots.  These can be found approximately as follows.  Note the double equals sign.

NSolve[h == 0, x]

So over the complex number field, the variety contains three points V(f, g) = { 2, i, -i }.  Now to find the ideal of the variety you can use some theory which will be covered in the second homework assignment.   Since the polynomial GCD  h generates the ideal < f, g > we know that the ideal  < h >  is a subset of I(V).  But in this case I(V) is larger than < h >.  The problem is the repeated factor (x-2)^2 in the factorization of h.   The theory tells us that  I(V) is generated by the so-called square-free part of h , that is, h with all repeated factors removed.  Here we could get this by just dividing by one factor of (x-2) but there is a general method which does not require knowing the factorization first.  It involves the GCD of h and its derivative.

h
hprime  =  D[h, x]
PolynomialGCD[h, hprime]

Note that the polynomial GCD of h and hprime is exactly the factor we want to eliminate. The theory says that this will always be the case, so by dividing h by the GCD we get a generator for I(V).

h1 =  
      PolynomialQuotient[h, PolynomialGCD[h, hprime], 
        x]
h // Factor
h1 // Factor

If you get this far, here are a few problems (some of which are from the text) to try.

1.  Find a way to express 1 (one) as an integer combination of  30*a+ 42*b+ 70*c +105*d.  Extract the values of a, b, c, d from the computer's answer and test your solution.

2. Without using factoring, find the square-free part of the polynomial 
f  = x^11- x^10 +2x^8 -4x^7 +3x^5 -3x^4 +x^3 +3x^2 -x -1
Test your answer by factoring both f and its square-free part.

3. Find a basis for the ideal I(V) of the variety  V( f, g) where
f = x^5-2x^4+2x^2-x
g = x^5-x^4-2x^3+2x^2+x-1

4.  To find a single generator of an ideal < f1, ... , fm > generated by more than two polynomials you need to find the PolynomialGCD of all of them.  This can be done two at a time.  For example, if h12 = Polynomial[GCD[ f1, f2] then 
PolynomialGCD[f1, f2, f3] = PolynomialGCD[ h12, f3 ]
Test this equation for the polynomial 
f1= x^6-1		f2 = x^12-1 		f3 = (x^2+x+1)^3
by computing both sides of the equation using Mathematica.

5. Unfortunately, the PolynomialExtendedGCD command only works for two polynomials at a time.  See if you can figure out how to express  h = PolynomialGCD[ f1, f2, f3] as a polynomial combination a1*f1 + a2*f2 + a3*f3  where fi are the polynomials in problem 4.  Test your answer.