# Math 21a Fall 2019

## Multivariable Calculus

# Mathematica Laboratory

- The Mathematica project is available
here. - Strict deadline: Wednesday, December 4. Midnight. (The email has to have a December 4 time stamp).
- Submit by email to knill@math.harvard.edu Please remove all graphics first from the notebook: (Menu: Cell -> Delete all Output) and attach the .nb file to your email. Do not use canvas nor any dropbox/google drive/cloud type service to submit. Just attach the .nb file to a regular email.

# Mathematica

You can get Mathematica from here. During installation you will be prompted for an Activation Key. To get it, try this direct link.## Double integrals

f[x_]:=If[0 < x < 1,1,0]; Integrate[ f[x] f[y] Abs[x-y],{x,-Infinity,Infinity},{y,-Infinity,Infinity}]

## Lagrange problems

Solve a Lagrange problem with 2 variablesf=2x^2+4 x y; g=x^2 y; Solve[{D[f,x]==L*D[g,x],D[f,y]==L*D[g,y],g==1},{x,y,L}]With 3 variables:

f=2x^2+4 x*y+z; g=x^2 y + z; c=1; Solve[{D[f,x]==L*D[g,x], D[f,y]==L*D[g,y], D[f,z]==L*D[g,z], g==c},{x,y,z,L}]With 3 variables and two constraints

f=z; g=z^2-x^2-y^2; h=4x-3y+8z; c=0; d=5; Solve[{D[f,x]==L*D[g,x] + M D[h,x], D[f,y]==L*D[g,y] + M D[h,y], D[f,z]==L*D[g,z] + M D[h,z], g==c, h==d}, {x,y,z,L,M}]

## Classifying critical points

Here is example code on how to compute the gradient and the discriminantf=x^3 y^3- x y^2; Grad[f,{x,y}] D[f,{x,2}]*D[f,{y,2}]-D[f,{x,1},{y,1}]^2

## Classifying critical points

This produces a nice table.f=4 x y - x^3 y - x y^3; ClassifyCriticalPoints[f_,{x_,y_}]:=Module[{X,P,H,g,d,S}, X={x,y}; P=Sort[Solve[Thread[D[f,#] & /@ X==0],X]]; H=Outer[D[f,#1,#2]&,X,X];g=H[[1,1]];d=Det[H]; S[d_,g_]:=If[d<0,"saddle",If[g>0,"minimum","maximum"]]; TableForm[{x,y,d,g,S[d,g],f} /.P,TableHeadings->{None,{x,y,"D","f_xx","Type","f"}}]] ClassifyCriticalPoints[f,{x,y}]If the numerical solution does not work, replace Solve with NSolve as follows:

f=4 x y - x^3 y - x y^3; ClassifyCriticalPoints[f_,{x_,y_}]:=Module[{X,P,H,g,d,S}, X={x,y}; P=Sort[NSolve[Thread[D[f,#] & /@ X==0],X]]; H=Outer[D[f,#1,#2]&,X,X];g=H[[1,1]];d=Det[H]; S[d_,g_]:=If[d<0,"saddle",If[g>0,"minimum","maximum"]]; TableForm[{x,y,d,g,S[d,g],f} /.P,TableHeadings->{None,{x,y,"D","f_xx","Type","f"}}]] ClassifyCriticalPoints[f,{x,y}]

#### Surfaces:

Plot3D[Sin[x y],{x,-3,3},{y,-3,3}] ContourPlot3D[ x^2+y^2-z^2==0,{x,-2,2},{y,-2,2},{z,-2,2}] RevolutionPlot3D[t, {r, 0, 1}, {t, 0, 2 Pi}] SphericalPlot3D[Sin[5 s], {s, 0, Pi}, {t, 0, 2 Pi}] ParametricPlot[{Cos[3t],Sin[5t]},{t,0,2Pi}] ParametricPlot3D[{Sin[s] Cos[t],Sin[s] Sin[t],Cos[s]},{s,0,Pi},{t,0,2Pi}]In spherical plot the first variable is the phi variable, the second the theta variable. You can name the variables as you like, for example

SphericalPlot3D[Sin[5 phi], {phi, 0, Pi}, {theta, 0, 2 Pi}]Plot two surfaces at once, like plotting the tangent plane to a surface

f = x^2 + y^2; L = 2 + 2 (x - 1) + 2 (y - 1); Plot3D[{f, L}, {x, -2, 2}, {y, -2, 2}, AspectRatio -> 1]

#### Curves:

ContourPlot[ Sin[x y],{x,-3,3},{y,-3,3}] ParametricPlot3D[{t Cos[10t],t Sin[10t],t},{t,0,2Pi}] PolarPlot[ Cos[7 t],{t,0,2Pi}]

#### Integration

Integrate[Sin[x],x] Integrate[Sin[x]^10,{x,0,Pi}] ArcLength[{Cos[t], Sin[t], t}, {t, 0, 9}] Integrate[Norm[D[{Cos[t], Sin[t], t}, t]], {t, 0, 9}] NIntegrate[Sin[x^10 Sin[x],{x,0,Pi}]NIntegrate does a numerical integration. One can use it if there is no closed form answer.

#### Dfferentiation

D[Sin[x],x] r= {Cos[t],Sin[t],t^5}; D[r,t] (* velocity *) D[r,{t,2}] (* acceleration *)

#### Checking differential equations

Here is a check that a function f satisfies the heat equation. In HW 12, you have to verify such a thing using technology but for a bit more complicated function and also for a bit more complicated PDE:f=(1/Sqrt[t])*Exp[-x^2/(4t)]; Simplify[ D[f,t] == D[f,{x,2}]]

#### Partial differential equations

This solves the heat equation and plots the function f(t,x).s=NDSolve[{D[f[t,x],t]==D[f[t,x],x,x],f[0,x]==Exp[-x^2]},f,{t,0,3},{x,-5,5}]; Plot3D[Evaluate[f[t,x] /. s],{x,-5,5},{t,0,3}]

## Roller Coaster

In one Data project, we have used this notebook about roller coaster to learn how to run a notebook. At the end of the course, we will have a Mathematica project. There is no need to buy the software. Harvard has a site license.## Elevation heights

A = Reverse[ Normal[GeoElevationData[ Entity["City", {"Cambridge", "Massachusetts", "UnitedStates"}], GeoProjection -> Automatic]]]; ReliefPlot[A] ListDensityPlot[A] ListPlot3D[A] ListContourPlot[A]