Beginning Mathematica

D.C.M. Burbulla and C.T.J. Dodson

University of Toronto

This NoteBook is intended to illustrate some basic operations in Mathematica and can be personalized for saving in your home directory for reference purposes. Note the use 
of square brackets [ ] to enclose arguments of functions, and note that Mathematica 
functions are usually standard names but start with capital letters. You can get help 
about a function by using ? before the function name. Click on the side bars to open 
interactive cells.

ELEMENTARY ALGEBRA CALCULUS AND PLOTS

?Expand

Expand[(x^2-3x)^3]

?Solve

Solve[x^2+6x+8==0,x]      (* Solves equation for x *)

NSolve[x^4 + x^3 -1==0,x]   (* NUMERICAL POLYNOMIAL SOLUTION *)

Plot[Sin[x], {x,-Pi,Pi}];

f[x_]:= x^4+3x^2-5x+1      (* Defines new function f *)

Factor[f[x]]

Plot[f[x],{x,-3,3}];

3-D PLOTS OF SURFACES

Plot3D[(-x^2-y^2),{x,-1,1},{y,-1,1}]; (* MAXIMUM AT (0,0) *)

Plot3D[(-x^2+y^2),{x,-1,1},{y,-1,1}];   (* SADDLE AT (0,0) *)

Plot3D[(x^2+y^2),{x,-1,1},{y,-1,1}];   (* MINIMUM AT (0,0) *)

ParametricPlot[{2 Cos[t],3 Sin[t]},{t,0,2Pi},
AspectRatio->Automatic];

There is a special package for polar plots:

<<Graphics`Graphics`

PolarPlot[Abs[Sin[t]],{t,0,2Pi}];

There is a special package for 3D parametric Plots

<<Graphics`ParametricPlot3D`

CylindricalPlot3D[Cos[t]Sin[t],{r,0.1,1},{t,0,2Pi}];

DIFFERENTIATION

D[x^n,x]         (* Differentiate x^n *)

D[Sin[x],x]

D[Sin[x y],x]

D[Sin[x y],{x,2}]       (* GIVES SECOND PARTIAL WRT x *)

FIND EXTREMA ANALYTICALLY AND NUMERICALLY USING f'[x]=0

?f

D[f[x],x]        (* Differentiate f at x *)

f'[x]            (* Evaluate derivative of f at x *)

Solve[f'[x]==0,x]

N[Solve[f'[x]==0,x]]   (* Finds horizontal tangents *)

Plot[f[x],{x,0.6,0.7}];  (* Plots f[x] near real minimum *)

INTEGRATION

?Integrate

Integrate[x^2, x]

Integrate[Sin[x],x]

Integrate[Sin[x],{x,0,Pi/15}]

NUMERICAL INTEGRATION

NIntegrate[Sin[x],{x,0,Pi/15}]

NIntegrate[E^(-x^2),{x,0,1}]

er[y_]:= NIntegrate[E^(-x^2),{x,0,y}]  (* FUNCTION AS INTEGRAL *)

Plot[er[y],{y,0,2}];

erdata=Table[{y,er[y]},{y,0,2,0.1}]; TableForm[erdata]  (*TABULATE FUNCTION *)

MULTIPLE INTEGRALS

Integrate[E^y + Sin[x], {y,0,1},{x,0,Pi/2}]

f1[a_]:= Integrate[E^y + Sin[x], {y,0,a},{x,0,Pi/2}]  (*FUNCTION AS INTEGRAL*)

Plot[f1[a],{a,0,1}];    (* PLOT FUNCTION *)

t=Table[{a,N[f1[a]]},{a,0,1,0.2}]        (*TABULATE IN 0.2 STEPS FROM 0 TO 1*)

TableForm[t]      (* TABLE FORM FOR DATA *)

PLOT OF NUMERICAL SOLUTION OF ODE

NDSolve[{y'[x] == 1 - y[x]^2, y[0] == 4}, y, {x,0,1}];

Plot[Evaluate[y[x] /. %], {x,0,1}]  ;

ORDINARY DIFFERENTIAL EQUATIONS

You can think of a first order differential equation as a way to specify a vector field in the x-y plane. Thus, dy/dx is the size of the y-component when the x-component has value 1. Then 
a solution to the equation is a curve that follows the vector field flow.  Mathematica can show
this as follows. The equation dy/dx =x has as solution  (x^2)/2 + c  for any constant c. The 
vector field has components (1,x) and the solution curves are y(x) = (x^2)/2 + c. We show 
these as follows using Mathematica's special package

<<Graphics`PlotField` (*INPUTS PACKAGE WITH SPECIAL FUNCTIONS*)

pf = PlotVectorField[{1,x}, {x, -1, 1}, {y, -1, 1}];

pd = Plot[x^2 /2, {x,-1,1},
 PlotStyle -> {{Thickness[0.02],RGBColor[1,0,0]}}];

Show[pd,pf];

Mathematica has a special package to help solving differential equations.

<<Calculus`DSolve`  (* INPUTS PACKAGE WITH SPECIAL FUNCTIONS*)

DSolve[{y'[x] == x, y[0] == 0}, y[x], x]

DSolve[{y'[x] == 1 - y[x], y[0] == 4}, y[x], x]

The following example concerns the logistic growth equation which has differential equation

		DE:    x'(t) =  k x(t) (1 - x(t)/m)

In the following example we shall take   k = 1   and   m = 4.

?NDSolve

The following function is defined to produce a numerical solution to the above differential equation with initial condition, x = a    if   t = 0 :

curve[a_] := NDSolve[{x'[t] == x[t] (1 - x[t]/4), x[0] == a},
                             x, {t, 0, 6}]

curve[1]

Plot[ Evaluate[ x[t] /. curve[1] ], {t, 0, 6}];

graphs = Table[ Plot[Evaluate[ x[t] /. curve[a] ], {t, 0, 6}],
                        {a, 0, 6}]

Show[graphs];

Note that for all of the above initial conditions, except for  a = 0,  x always tends to  x=4  as   t  tends to infinity. That is, in the long run the value of  x  approaches  4. (If  x = 0, when t = 0 then   x(t) = 0 , for all values of t.) What happens if   x = 4  when  t = 0 ?

MORE ADVANCED TOPICS

VECTOR ANALYSIS

Note that Mathematica has a special package Calculus`VectorAnalysis`
which contains functions to compute Jacobians, change coordinates among
most of the standard systems, provide arclength functions,as well as 
the operator functions Div, Grad, Curl and Laplacian.
To input this package into your Mathematica session use

<<Calculus`VectorAnalysis` (*INPUTS PACKAGE WITH SPECIAL FUNCTIONS*)

Div[{x^4 y , y^3 , z x}, Cartesian[x,y,z]]

Curl[{x^4 y z , y^3 , z x},Cartesian[x,y,z]]

Grad[1/(x^2+y^2+z^2)^(1/2), Cartesian[x,y,z]]

Grad[1/r, Spherical]

Grad[(x^3+y^2-z^4)].{1/3,2/3,-2/3}(* DIRECTIONAL DERIVATIVE ALONG (1,2,-2) *)

CoordinatesToCartesian[{1,Pi/2,Pi/4},Spherical]

CoordinatesFromCartesian[{x,y,z},Spherical]

JacobianMatrix[Spherical]

JacobianDeterminant[Spherical]

Integrate[r^2 JacobianDeterminant[Spherical], {r,0,2},
         {theta,0,Pi},{phi,-Pi,Pi}]

FINDING EXTREMA USING Grad[H[x,y]==0]

H[x_,y_]:=1+4/3 x^3 + 4y^3 -x^4 -y^4

Plot3D[H[x,y],{x,-10,10},{y,-10,10}];

sol=Solve[Grad[H[x,y]]=={0,0,0},{x,y}]; (* SOLVE FOR Grad[H] = 0 *)

H[x,y] /. sol                         (* SUBSTITUTE SOLUTIONS *)

Max[H[x,y] /. sol]                    (* FIND MAXIMUM *)

PLOTTING VECTOR FIELDS

<<Graphics`PlotField` (* INPUTS PACKAGE THAT PLOTS 2-D FIELDS*)

PlotVectorField[{-y,x}, {x, -1, 1}, {y, -1, 1}];

PlotVectorField[{-y , x Sin[x]}, {x, -1, 1}, {y, -1, 1}];

PlotGradientField[1/(x^2 + y^2), {x, -1, 1}, {y, -1, 1},
      MaxArrowLength->600]

<<Graphics`PlotField3D` (* INPUTS PACKAGE THAT PLOTS 3-D FIELDS*)

PlotVectorField3D[{-y , x , z x}, {x, -1, 1}, {y, -1, 1}, {z, -1, 1},
      VectorHeads->True];

PlotGradientField3D[1/(x^2 + y^2 + z^2), 
{x, -1, 1}, {y, -1, 1}, {z, -1, 1}, VectorHeads->True, 
      MaxArrowLength->100];

SAVING DEFINITIONS AND GRAPHICS

(* SAVING DEFINITIONS OF FUNCTIONS INTO A FILE *)
Save["/yourpathname/myfile",f,g]          (* use your home directory path *)

(* READING IN A FILE *)
<< /yourpathname/myfile         (* use the full  directory path *)

(* TYPING CONTENTS OF A FILE WHILE IN MATHEMATICA *)
!! /yourpathname/myfile         (* use the full  directory path *)

(* EXECUTING A SHELL COMMAND WHILE IN MATHEMATICA *)
! ls -l 

(* SAVING A GRAPHIC OUTPUT TO A POSTSCRIPT FILE *)
Display["/yourpathname/graphic1.ps",Out[3]]  (* use your home directory path *)

(* YOU NEED TO CONVERT GRAPHIC graphic1.ps TO EPS  FORM USING* ) 

psfix /yourpathname/graphic1.ps > /yourpathname/graphic1.eps

(* NOW YOU CAN DELETE /yourpathname/graphic1.ps AND KEEP THE eps VERSION *)
(* USE ghostview TO VIEW eps FILES *)