Geostrophic Adjustment

Consider an initial displacement η(x,0) in a layer of shallow water.  The displacement is shaped like a "tophat" , with value of 1 for -L <x< L and a value 0 elsewhere.  The initial displacent is the superposition of a steady geostrophic solution (Subscript[η, g])(x,0) and gravity waves that will propagate away.  We find the  (Subscript[η, g])(x,0).

First find the Fourier transform of η(x,0):

tran=Integrate[Cos[k x],{x,-L,L}]

Check to see if the inverse transform gives back the tophat function.  (Is (1/(2 ))(Sign[L-x]+Sign[L+x]) the tophat function?).

invtran=Integrate[tran*Cos[k x],{k,-Infinity,Infinity}]/(2 Pi)

Recall the amplitude of geostrophic part of a sine wave disturbance in η is:

geopart=1/(1+(LR*k)^2)

We now find (Subscript[η, g])(x,0) using an inverse Fourier transform method:

geotran=Integrate[tran*geopart*Cos[k x],{k,-Infinity,Infinity}]/(2 Pi)

All that stuff in the middle of the If statement is what we are looking for.  We grab it with the mouse and clean it up with TrigToExp:

trythis =TrigToExp[ 
  1/2*(Sign[L - x]*(-Cosh[Sqrt[(L - x)^2/LR^2]] + 
        Sinh[Sqrt[(L - x)^2/LR^2]] + 1) + 
     Sign[L + x]*(-Cosh[Sqrt[(L + x)^2/LR^2]] + 
        Sinh[Sqrt[(L + x)^2/LR^2]] + 1))]

Make a function so that we can plot it out.  Notice I used = and not := to define the function.

h[x_, LR_, L_]=trythis

Let's see if the function works...well it does give back a number...that's nice...

h[1.9,.1,1.]

For comparison, we need to be able to plot η(x,0):

h0[x_,L_]:=1/2*(Sign[L- x]+ 
     Sign[L+ x])

We can nondimensionalize distance in terms of L, so L=1 in the argument list.  We plot the case for (Subscript[L, R])=0, which is the same as the initial condition, and the solution for three other values of (Subscript[L, R]).

fourplots=Plot[{h0[x,1],h[x,0.25,1],h[x,1,1],h[x,4,1]},{x,-10,10},PlotRange->{-0.2,1},AxesOrigin->{-10,0},PlotStyle->{{RGBColor[1,0,0]},{RGBColor[0,0,0]},{RGBColor[0,0,0]},{RGBColor[0,0,0]},{RGBColor[0,0,0]}},AxesLabel->{x,h},Exclusions->None]

Export["fourplots.eps",fourplots,"EPS"];