MATLAB Codes for computing trapped mountain waves MATLAB
codes that compute linear modes that represent trapped waves in the
lee of mountain barriers can be
search for solutions to the nonlinear eigenvalue problem described in
of Durran et al
(2014) using the numerical method in the appendix of that paper.
Any questions or comments should be directed to Peter Blossey (email@example.com).
An example from Durran et al (2014)
By running the script TrappedWavesExample in MATLAB, the method can be
demonstrated on a number of the basic state profiles in Durran et al
(2014). For example, running with icase = 1 yields:
Searching for trapped or leaky trapped wave modes for:
Durran et al (2014) Prof. 1: No Stratosphere
Note that we are searching a wide variety of parameter space for
solutions so that error messages may appear below ...
Many error messages are produced because we are trying many potential
solutions for the nonlinear eigenvalue problem and not all of them
The few (if any) that satisfy the eigenvalue
problem are returned below.
Found 2 solutions for Durran et al (2014) Prof. 1: No Stratosphere
Mode 1: eigenvalue = 6.7671e-04 + 1i* 0.0000e+00, lambda x =
9.28 km, convergence residual = 1.2584e-12
Mode 2: eigenvalue = 4.8429e-04 + 1i* 1.1216e-03, lambda x =
12.97 km, convergence residual = 9.3655e-11
Here, two solutions have been found. The first one, a purely trapped
lee wave (Mode 1), has a wavelength of 9.28km and is familiar from
figure 2 of Durran et al (2014).
The second mode is a leaky trapped wave that decays strongly
downstream (losing all of its amplitude within one wavelength). The
vertical structure of base state and the two trapped wave solutions are shown in figures:
Because the second mode (at right) decays so strongly downstream, the linear trapped
wave response would be dominated by Mode 1 if that mode is excited by
the flow over a mountain.
The MATLAB script
TrappedWaveExample.m supplies ten
ready-made base states and find trapped wave solutions for each,
sometimes only one mode, and other times more than one. Try them all
by changing the limits on the
for icase = loop at the top
of the script.
Revisiting Corby & Sawyer (1958)
Even better, look for trapped wave solutions for a different basic
state. You just need to specify a profile that is piecewise linear in
U, piecewise constant in Brunt-Vaisala frequency N and has uniform
velocity over a layer at the top of the model. The last restriction
makes possible the application of our wave radiation condition at the
top of the profile.
As an example, let's
revisit Corby and
Sawyer (1958) and see if we can approximately reproduce the
wavelenths of the three modes in their figure 4. Let's build a
three-layer approximation to their four layer atmosphere -- our
radiation condition should let us cut off the profile in the third
layer. Let's use the following heights with vertically-uniform wind
and the Brunt-Vaisala frequencies (N^2 = l^2 * U^2, where l is the Scorer parameter):
H = [0 1e3 10e3 20e3];
Ubar = [10 10 10 10];
Nbrunt = sqrt(100*[4e-6 0.3e-6 6e-6]);
Note that we only need to specify three values of N because it is
uniform betwee the four heights specified in H. Also, note that Corby
and Sawyer only specified values for the Scorer parameter, so that we can
freely choose the (vertically-uniform) wind speed as long as l = N/U.
This base state is available under
TrappedWaveExample.m. Trying this yields:
Found 2 solutions for Corby & Sawyer (1958): four layer profile
Mode 1: eigenvalue = 7.6142e-04 + 1i* 1.4088e-08, lambda x =
8.25 km, convergence residual = 8.4203e-13
Mode 2: eigenvalue = 3.7669e-04 + 1i* 2.1257e-05, lambda x =
16.68 km, convergence residual = 2.3472e-12
so that we have found two modes. The first corresponds almost exactly
with the first mode in Corby and Sawyer (1958), matching all four
given digits of the wavenumber k (which they specify as k=0.7614 km^(-1)). The form of the second mode is very
similar to that of their third mode, but our wavenumber is slightly
larger, perhaps because we have truncated the base state at 20km
height in the stratosphere and applied a radiation condition there.
The two modes, shown below may be compared with figure 4a in Corby and
Note that our method only predicts the form, wavelength and downstream
decay rate of the trapped wave modes. A different method would be
required to predict their amplitude.