Mathematica for physics tutorial
Mathematica has multiple types of cells. This line is a "text line" (Alt 7), It will be ignored when Mathematica goes through an executes your commands. It is great for comments.
This line is an "input line" (Alt 9). Mathematica will try to execute these commands and will give a red error box to the right, because this is not a mathematical input. This is the default type of line created when you just start typing.
To have mathematica calculate a line, put your cursor on that line and hit Shift+Enter to run.
x=5
x+9
The previous input should be interpreted as "x=5"->Define x as 5 for the rest of the mathematica session.
You should check that basic math works correctly:
3x
a=x!
b=x^x; (* Notice that the semicolon supresses the output, You don't always want to see what Mathematica has calculated.*)
c=Log[x]
N[c] (* N[ ] is a great function for getting a numerical estimate for something*)
x+y
Notice: x has changed from blue font to black, this means that x is defined, while y is still an undefined variable.
Functions and commands in Mathematica are case sensitive and their arguments need to be square brackets
Integrate[y^3,y]
One of the great powers Mathematica bestows upon its user is Simplify and FullSimplify. Use it well
Simplify[(1-Cos[θ]^2)(1+Cot[θ]^2)]
Before we move on with the document we want mathematica to "forget" the definitions we gave it above.
Clear[x] (* This clears x only, notice that after excuting the symbol x should appear blue again. *)
ClearAll/@Names["Global`*"] ; (*this should clear everything*)
Calculus
This is how you define a function (Notice both the Underscore after the variable x_ and the "colon equals" := are both important.. The := operator ensures that definition of the function is evaluated every time )
f[x_]:=a1 x^3
g[y_]:=f[y]+2
g[z]
Derivatives and Integrals
Command for derivatives:
D[f[x],x]
Higher derivatives:
D[f[x],{x,2}]
Command for indefinite integrals:
Integrate[f[x],x]
Command for definite integrals:
Integrate[f[x],{x,0,1}]
Warm up exercises
1) Create a function h[x] that is quadratic in x. Use numbers rather than variables for the coeffiecents.
2) Have mathematica calculate dh/dx for you.
3) Use the Solve[ ] function to find the roots of your quadratic equation. To get to the right documentation quickly, type out Solve and then hit F1 while the cursor is in that key word.
Taylor Series
ClearAll/@Names["Global`*"] ;
Let us now expand the function about the point x=0, keeping up to 5th order terms
Series[Cos[x],{x,0,5}]
You can get rid of the O[x]^2 symbol via Normal[ ]
% calls the result of the previous line, %n calls the result of the n^th line
Normal[%]
If you want to call just the forth order coefficient use:
SeriesCoefficient[Cos[x],{x,0,4}]
So what is the difference between a function and it' s Taylor series? Let's look at an example.
f[x_]:=Cos[x];
g[x]=Normal[Series[f[x],{x,0,7}]];
Plot[{f[x],g[x]},{x,-5,10}, PlotLabels->"Expressions"]
Warm up exercises
1) Create a function h[x] that is quadratic in x. Use numbers rather than variables for the coeffiecents.
2) Have mathematica calculate dh/dx for you.
3) Use the Solve[ ] function to find the roots of your quadratic equation. To get to the right documentation quickly type out Solve and then hit F1 while the cursor is in that key word.
Differential equations
ClearAll/@Names["Global`*"] ;
This is an example of solving a Differential equation with specified boundary conditions
soln = Assuming[ k>0 , DSolve[{ ψ''[x] + k ψ[x] == 0, ψ[0] == 2.0, ψ[4.0] == 0.0}, ψ , x ]] ;
kind of ugly. but if we want to plot the solution, we must give a value to k and extract the solution into
ψ[x]=FullSimplify[ψ[x]/.soln[[1]]]
to plot the function we need to define k so all variables can be evaluated numerically
k=7;
Plot[ψ[x],{x,0,10}]
k=.
A more difficult example: Lets compare the heat equation to the wave equation
ClearAll/@Names["Global`*"] ;
initialCondition[x]=Abs[x-1]-1;Plot[initialCondition[x],{x,0,2}] (*here is the initial condition, a rod that is cold in the middle and has the endpoints at a fixed temperature OR a string bent into a V but fixed on both ends*)
(* The Heat equation states that the Concavity is related to the velocity of each point*)
HEAT=D[D[T[x,t],x],x]== D[T[x,t],t]
(* The wave equation states that the concavity is related to the acceleration of each point*)
WAVE=D[D[y[x,t],x],x]==D[D[y[x,t],t],t]
solHEAT=NDSolve[{HEAT,T[x,0]==initialCondition[x],
(*fix the end points*) T[0,t]==0,T[2,t]==0 },
T[x,t],{x,0,2},{t,0,4}];
solWAVE=NDSolve[{WAVE,y[x,0]==initialCondition[x],
(*fix the end points*) y[0,t]==0,y[2,t]==0 },
y[x,t],{x,0,2},{t,0,4}];
T1HEAT[x_,t_]=T[x,t]/.solHEAT[[1]]; (* This is an annoying formatting step that turns the array/list of solutions into a function we can plot later*)
y1WAVE[x_,t_]=y[x,t]/.solWAVE[[1]];
(*This manipulate shows how the temperature of the rod changes as time goes on*)
(*Remeber, concavity controls velocity*)
Manipulate[Plot[T1HEAT[x,t], {x,0,2}, PlotRange->{-1,1}],{t,0,1}]
(*This manipulate shows how the wave propogates in the string as time goes on*)
(*Remember, concavity controls acceleration*)
Manipulate[Plot[y1WAVE[x,t], {x,0,2}, PlotRange->{-1,1}],{t,0,1}]
(* Now both at the same time*)
Manipulate[Plot[{T1HEAT[x,t],y1WAVE[x,t]}, {x,0,2}, PlotRange->{-1,1}],{t,0,4}]
A totally different example: Wave function of Hydrogen
ClearAll/@Names["Global`*"] ;
(* The wave function of hydrogen needs 6 inputs. Three of them are normal spherical coordinates, (r, θ, ϕ) and three of them are special integers called quanum numbers (n,l,m) *)
(* These functions are famous enough that they specialized names that Mathematica has predefined as Laguerre and SphericalHarmonicY*)
ψ[n_,l_,m_][r_,θ_,ϕ_]:=Sqrt[(n-l-1)!/(n+l)!] E^(-(r/n)) ((2 r)/n)^l 2/n^2 LaguerreL[n-l-1,2 l+1,(2 r)/n] SphericalHarmonicY[l,m,θ,ϕ];
(*To see what these functions look like we can just pick a few n, l, m values that are integers
This should work any time n is an positive integer, l{Sqrt[x^2+y^2+z^2],ArcCos[z/Sqrt[x^2+y^2+z^2]],Arg[x+I y]}];
n=2;l=1;m=1;
Manipulate[RegionPlot3D[Abs[ψ[n,l,m][r,θ,ϕ]/.sphericalToCartesian]>Sqrt[density]/n^3,{x,-n^2/zoom,n^2/zoom},{y,-n^2/zoom,n^2/zoom},{z,-n^2/zoom,n^2/zoom},AxesLabel->{x,y,z},PlotStyle->Directive[Opacity[0.3`]],Mesh->None,PlotPoints->20,ColorFunction->Function[{x,y,z},Hue[Rescale[Arg[ψ[n,l,m][r,θ,ϕ]/.sphericalToCartesian],{-Pi,Pi}]]],ColorFunctionScaling->False],{{density,.02},.01,.1},{{zoom,.5},.1,2}]