Dispersion of Gaussian Wave Packets

As you know, the dispersion relation for deep water waves is
		(ω^2)= gk.  
The phase speed is
			(c^2)=g/k.
Thus long waves (small k) move faster than short waves.
Here we consider a "wavetrain", or "periodic group" of deep
water waves.  The system can be scaled to be dimensionless so that period of the wave is  2π. and the dispersion relation is ω(k)=(Sqrt[k]) .  We make some pretty plots of:
		(Subscript[(Underoverscript[Σ, 1, k]), m])(Subscript[a, k])cos (kx-ω(k)t)

Some stuff

cs[k_,nm_] gives a list of the values of cos(kx) at nm_ points of x between -π and π

cs[k_,nm_]:=Module[{del,j},
		del=k*2*Pi/nm//N;
		j=(nm-1)/2;
		Table[Cos[del*i],{i,-j,j}]]

sn[k_,nm_] gives a list of the values of sin(kx) at nm points of x between -π and π

sn[k_,nm_]:=Module[{del,j},
		del=k*2*Pi/nm//N;
		j=(nm-1)/2;
		Table[Sin[del*i],{i,-j,j}]]

The next gives a list of lists of cos(kx) and sin(kx) for k=0 to k=(Subscript[k, m])

manycs[km_,nm_]:=Table[cs[k,nm],{k,0,km}]

manysn[km_,nm_]:=Table[sn[k,nm],{k,0,km}]

gauss[c_,k0_,km_]  gives a list  of the Fourier coefficients (Subscript[a, k])in
															(Subscript[(Underoverscript[Σ, 1, k]), m])(Subscript[a, k])cos (kx)				  
		that make a  pretty good approximation to the Gaussian wavepacket :
																		 (e^(-((x^2)/(4 c))))  cos((Subscript[k, 0]) x).   
Note that k=0 gets half the magnitude of the normal formula.

gauss[c_,k0_,km_]:={
		(Sqrt[c]/Sqrt[Pi]/2)*(Exp[-c*(0-k0)^2]+Exp[-c*(-0-k0)^2])//N,Table[(Sqrt[c]/Sqrt[Pi])*(Exp[-c*(k-k0)^2]+Exp[-c*(-k-k0)^2])//N,{k,1,km}]}//Flatten

setup is a module that sets up the lists of cosine and sine functions.
km should first be defined as an integer (of course)

setup:=Module[{},
		nm=2*km+1;
		ca=manycs[km,nm];
		sa=manysn[km,nm];
		]

Superimposes the waves (Subscript[(Underoverscript[Σ, 1, k]), m])(Subscript[a, k])cos (kx-ωt) to make a
group.  Uses identity cos(kx-ωt)=cos(kx)cos(ωt)+sin(kx)sin(ωt)

grp[t_]:=Module[{pc,ps},
		pc=Table[Cos[ω[k]t],{k,0,km}];
			ps=Table[Sin[ω[k]t],{k,0,km}];
		Sum[ak[[k]]*(pc[[k]]*ca[[k]]+ps[[k]]*sa[[k]]),
		{k,1,km+1}]]

plotgrp[t_] plots the group at time t_

plotgrp[t_]:=ListPlot[grp[t],Joined->True,Axes->False,PlotRange->{-1.1`,1.1`}]

Why it works:

Note(∫_(-∞)^∞)(E^(-c (((k-k0))^2))) Cos[k x]dk=(∫_(-∞)^∞)(E^(-c (((j))^2))) Cos[(j+k0) x]dj
For the latter integral, Mathematica (see below) tells us that:
(∫_(-∞)^∞)(E^(-c (((j))^2))) Cos[(j+k0) x]dj=(((E^(-((x^2)/(4 c)))) (Sqrt[π]) Cos[k0 x])/(Sqrt[c]))
Now we also have 
(∫_(-∞)^∞)(E^(-c (((k-k0))^2))) Cos[k x]dk
=(∫_0^∞)(E^(-c (((k-k0))^2))) Cos[k x]dk +(∫_(-∞)^0)(E^(-c (((k-k0))^2))) Cos[k x]dk 
=(∫_0^∞)(E^(-c (((k-k0))^2))) Cos[k x]dk+(∫_0^∞)(E^(-c (((-k-k0))^2))) Cos[k x]dk
Our Fourier series over the positive k approximates the above integrals.

(* Integrate[Exp[-c j^2]Cos[(j+k0) x],{j,-Infinity,Infinity}] *)

Now lets use it:

Now let's use it! Here is a dispersion relation appropriate for deep
water waves:

ω[k_]:=Sqrt[k]//N

Let's use km=40 wavenumbers.  c=.01 gives a nice, isolated hump, but not too isolated to abuse our truncated Fourier series.

km=40;

setup

Now let's make a wave packet.  Here it will be quite clear that the
group velocity is about half the phase velocity of the dominant wavenumber k0=10.

ak=gauss[.01,10,km];

Here is the group at t=0:

plotgrp[0]

Animate[plotgrp[t],{t,0,8,.2}]

Here are come other cool things to do.  We can plot the Fourier coefficents:

ListPlot[Table[{k-1,ak[[k]]},{k,1,km+1}],PlotRange->{Min[ak],Max[ak]},PlotStyle->{RGBColor[1,0,0]}]

We can plot a Fourier component (just to make sure it really is there):

ListPlot[ca[[3]]]