(* Mathematica notebook, Physics 237, Frank L. H. Wolfs.  
This notebook plots the probability densitry distribution of a free particle at different times.  It is assumed that the wavefunction describing the free particle is the sum of a series of wavefunctions with different wavenumbers k.  This calculation allows us to study the motion of the particle and the spread of the wavefunction as function of time. *)


(* Define the wavefunction Psi.  Psi depends on x, t, and k.  Energy and omega depend on k, hbar, and mass. *)
Psi[x_,t_,k_, hbar_:1,mass_:1] := 
Block[{energy, omega},energy=hbar^2*k^2/(2 mass);omega = energy/hbar;
Exp[I (k x - omega t)]];


(* Define the sum of a set of wavefunctions.  The sum of wavefunctions is weighted with a Gaussian weight factor.  The wavenumbers of the waves being added cover the range between k-0.01*nk and k+0.01*nk. *)
PsiSum[x_,t_,k_, nk_, hbar_:1,mass_:1] :=
 (∑_(n=-nk)^nk)Exp[-((0.01*n)^2)/(2*0.5^2)]Psi[x,t, k+0.01*n,hbar,mass];


(* Create an array of plots to show the real part of Psi at different times. The wavefunction is evaluated between x = -150 and x = +150 at six different times between t = 0.0 and t = 5.0. *)
figures=Table[Plot[Evaluate[Conjugate[PsiSum[x,t,20,100]]* PsiSum[x,t,20,100]],{x,-25,110},PlotRange->Full, AxesLabel->{"x (a.u.)","(SuperStar[Ψ])Ψ (a.u.)"},PlotLabel->"t = "<>ToString[t],Prolog->Thickness[0.001]],{t,0.0,5.0,1}];


(* Plot the wave packets. *)
figure=Show[GraphicsArray[Partition[figures,2]]]