wave equation, nonzero initial displacement

The basic wave equation is  (a^2) (Subscript[y, x x]) = (Subscript[y, t t])   for 0 < x < L, t > 0   
with boundary conditions 
y(0, t) = y(L, t) = 0 for t > 0 and
y(x, 0) = f(x),  (Subscript[y, t])(x,0) = 0 for  0 < x < L.

This represents the vertical displacement y along the length of a tight string of length L fixed at both ends for changing time.  It has an initial nonzero vertical displacement f.  Time begins as the string is temporarily held in the position so its initial velocity is 0.


notebook created by T. Shelton, Southwestern University
reference  Elementary Differential Equations and Boundary Value Problems, Boyce & DiPrima, 8th Edition, John Wiley & Sons, Inc.; 2005. 
reference  Differential Equations: Theory, Technique, and Practice, Simmons & Krantz, McGraw-Hill Co, Inc. 2007

Set-up

You may set the parameter a  for the wave equation and the length of the string, L.

Clear[a, L];
a=2;
L = 30;

Initial displacement of the string, f(x)

You may change the function definition f .

Clear[f];
f[x_] :=Piecewise[{{ (3/L)x , 0 ≤  x < L/3.},{ (3/(2 L))(L-x), L/3. ≤  x ≤ L}}]

Plot[f[x],{x,0,L},
PlotStyle->{Red,Thickness[0.01]},
 PlotLabel->"GRAPH 1: f, the initial shape of the string"]

Fit Fourier Series to f

You may change the number of terms to be calculated in any series and the number of periods to plot.

Clear[ numTerms,numPeriods];
numTerms = 10;
numPeriods = 3;

Use these automatic calculations or define your own.

Clear[minX, maxX, tMax,tInc];
maxX = numPeriods * 2*L; 
minX = - (numPeriods-1) *2* L;
tMax = numPeriods*L/a;
tInc = tMax/20.;

Fourier approximation to the initial displacement of the string

i.e. the partial sum for the Fourier sine series of f(x) = y(x, t=0)

Clear[b];
Table[b[n]=2/L NIntegrate[f[x]Sin[n π x/L],{x,0,L}],{n,1,numTerms}];

Clear[fPartialSum,t,n]
Print["A partial sum for the Fourier sine series"]
fPartialSum[x_]=  Sum[b[n] Sin[n π x/L],{n,1,numTerms}]
Plot[{fPartialSum[x],f[x]},{x,0,L},PlotStyle->{Blue,{Red}},AxesLabel->{" x "},
 PlotLabel->"GRAPH 2: f is red; partial Fourier is blue"]

approximate solution to the wave equation

The solution is a function whose space domain is -∞ < x < ∞, well beyond the interval 0 < x < L on which the string physically exists.  We see that the original displacement described by f(x) on [0, L] matches the solution y(x, t=0).

Clear[y];
Print["Some terms in the wave eqn solution y[x,t]"]
y[x_,t_] = Sum[b[n] Sin[n π x/L] Cos[n π a t /L],{n,1,numTerms}]

Show[Plot[y[x,0],{x,minX,maxX},PlotStyle->{Blue,Dashed},AxesLabel->{" x "}, PlotLabel->"GRAPH 3: f red; wave eqn approx at t=0 blue"],
Plot[f[x],{x,0,L},PlotStyle->{Red,Thickness[0.01]}]]

some graphs of the solution to the wave equation problem

3D graph

The solution is a two dimensional surface, plotted below, because the string displacement changes both along the string and in time.  Here we show it for the physical string (0 < x < L) and a time interval.

Plot3D[y[x,t],{x,0, L}, {t, 0,tMax},AxesLabel->{"x","t","y"}, 
PlotLabel->"GRAPH 4: y[x,t]"]

additional graph settings

You may change 
the max and min values for the vertical displacement of the string for subsequent graphs. 

These values should be a bit beyond the vertical values on the 3D graph

Clear[minY, maxY];
minY= -1.5;
 maxY = 1.5;

interesting slices of the 3D graph

Each image is a view of the string frozen in time.  See how the wave equation solution along the length of the string changes over  time by sliding the value of time. The boundary conditions guarantee that the ends always stay on the x-axis.
Each frame in GRAPH 5 is a slice of this surface for a particular t.
Each frame in GRAPH 6 is a slice of this surface through the x-value.

Manipulate[
Plot[y[x,t],{x,0, L}, PlotRange->{minY,maxY},AxesLabel->{" x "}, 
PlotLabel->"GRAPH 5: time-slice of wave eqn approx along the string"],
{{t,tInc,"t"},0,tMax,tInc}]

Each image focuses on a single point on the string, showing how that point on the string changes in time.  Slide the space variable x to see the time changes for other points on the string.

Manipulate[
Plot[y[x,t],{t,0,tMax}, PlotRange->{minY,maxY},AxesLabel->{" t "}, 
PlotLabel->"GRAPH 6: space-slice of wave eqn approx in time"],
{{x,0.5,"x"},0,L,.05}]

extended solution to the wave equation

The solution actually exists on a "virtual string" that is infinite in both directions, as we saw in GRAPH 3 above.  For each position along the virtual string, see how that point changes in time.  Slide along for different points in the extended space variable.

Manipulate[
Plot[y[x,t],{t,0,tMax}, PlotRange->{minY,maxY},AxesLabel->{" t "}, 
PlotLabel->"GRAPH 7: space-slice of wave eqn approx in time; virtual string"],
{{x,minX,"x"},minX,maxX,.05}]

Plot3D[y[x,t],{x,minX, maxX},  {t, 0,tMax},AxesLabel->{"x","t","y"},PlotLabel->"GRAPH 8: y[x,t] for a virtual string"]

a traveling wave approach

d'Alembert and Euler described the solutions as waves generated by the initial displacement that radiate in both directions.  
The waves travel with velocity "a" given by the parameter in the wave equation.

Manipulate[
Plot[{f[x], f[x - a*t],f[x + a*t]},{x,minX, maxX}, 
PlotStyle->{{Red,Thickness[0.01]},{Blue },{Purple}},
PlotRange->{minY,maxY},AxesLabel->{" x "}, PlotLabel->"GRAPH 9: f red; right-moving wave f[x + a*t] blue; left purple"],
{{t,0,"t"},0,tMax,tInc}]

See how the waves that travel from the original f on [0, L] meet the solution at t=0.

Manipulate[
Plot[{f[x], f[x - a*t],f[x + a*t],y[x,0]},{x,minX, maxX}, 
PlotStyle->{{Red,Thickness[0.01]},{Purple,Thickness[0.005] },{Purple,Thickness[0.005]},
{Brown}},
PlotRange->{minY,maxY},AxesLabel->{" x "}, PlotLabel->"GRAPH 10: f red; traveling waves purple; wave eqn approx at t=0 brown"],
{{t,0,"t"},0,tMax,tInc}]

The solution y[x, t] is the average of the traveling waves for the odd periodic extension of the initial displacement function.  Below we plot the average traveling wave dashed so the Fourier solution shows through.  We see these are the same function.

Clear[fOddPeriodic,x];
fOddPeriodic[x_]:= f[x] /;(0≤x<L)
fOddPeriodic[x_]:=-f[-x] /;(-L≤x < 0)
fOddPeriodic[x_]:= fOddPeriodic[x-2*L * Sign[x]]/; ((x <-L)||(x >L ))

Manipulate[
Plot[{( fOddPeriodic[x - a*t]+fOddPeriodic[x + a*t])/2,y[x,0]},{x,minX, maxX}, 
PlotStyle->{{Purple,Thickness[0.005], Dashing[Large]},
{Black}},
PlotRange->{minY,maxY},AxesLabel->{" x "}, PlotLabel->"GRAPH 11: traveling waves purple and dashed; wave eqn approx black"],
{{t,0,"t"},0,tMax,tInc}]

Manipulate[
Plot[{fOddPeriodic[x - a*t],fOddPeriodic[x + a*t],y[x,t]},{x,0, L}, 
PlotStyle->{{Purple },{Green,Dashed},{Brown,Thick}},
PlotRange->{minY,maxY},AxesLabel->{" x "}, PlotLabel->"GRAPH 12: traveling waves purple, green; average brown; physical string"],
{{t,tInc,"t"},0,tMax,tInc}]

Manipulate[
Plot[{f[x], fOddPeriodic[x - a*t],fOddPeriodic[x + a*t],y[x,t]},{x,0, L}, 
PlotStyle->{{Red,Thickness[0.01]},{Purple,Thickness[0.005] },{Purple,Thickness[0.005]},
{Brown,Thick}},
PlotRange->{minY,maxY},AxesLabel->{" x "}, PlotLabel->"GRAPH 13: f red; traveling waves purple; wave eqn approx brown"],
{{t,0,"t"},0,tMax,tInc}]