Here is an alternate derivation of reducing Laplace's tidal equations along the equator. Recall that the behavior of QBO (and ENSO) is isolated to the equator.
For a fluid sheet of average thickness , the vertical tidal elevation , as well as the horizontal velocity components and (in the latitude and longitude directions).
This is the set of Laplace's tidal equations (wikipedia). Along the equator, for at zero we can reduce this.
where is the angular frequency of the planet's rotation, is the planet's gravitational acceleration at the mean ocean surface, is the planetary radius, and is the external gravitational tidal-forcing potential.
The main candidates for removal due to the small-angle approximation along the equator are the second terms in the second and third equation. The plan is to then substitute the isolated and terms into the first equation, after taking another derivative of that equation with respect to .
Taking another derivative of the first equation:
Next, on the bracketed pair we invert the order of derivatives (and pull out , which may not be right for ENSO where the thermocline depth varies!)
Notice now that the bracketed terms can be replaced by the 2nd and 3rd of Laplace's equations
The terms are in longitude so that we can use a separation of variables approach and create a spatial standing wave for QBO, where is a wavenumber.
The next bit is the connection between a change in latitudinal forcing with a temporal change:
to describe the external gravitational forcing terms, then the solution is the following:
where is an aggregate of the constants of the differential equation.
So we can eventually get to a fit that looks like the actual QBO, for a fit from 1954 to 2015 for the 30 hPa data.
Besides the excursion-limiting behavior imposed by a modulation, this formulation can also show amplitude folding at the positive and negative extremes. In other words, if the amplitude is too large, the outer modulation starts to shrink the excursion, instead of just limiting it. Note around the 2016 time frame, the negative peak starts folding inward.
One caveat in this derivation is that the QBO may actually be the term - the horizontal longitudinal velocity of the fluid, the wind in other words - which can be derived from the above by applying the solution to Laplace's third tidal equation in simplified form above.
EDIT (8/22/2016) — The caveat in the paragraph above is resolved in this post.
What I have derived is a first step, good enough for modeling the QBO for the time being. But the dynamics of ENSO substantially differ in character from the QBO, showing more extreme excursions. How could that come about?
Consider if instead of being a constant, it varies in a similar sinusoidal fashion, then that term will go to the RHS as a forcing, and the above impulse response function will need to be convolved with that forcing . The same goes for the factor which I foreshadowed earlier. Having these vary with time changes the character of the solution. So instead of having a relatively constant amplitude envelope characteristic of QBO, the waveform amplitude will become more erratic, and more similar to the dynamics of ENSO, where the peak excursions vary a considerable amount.
Modifying the factor in fact makes the solution closer to a Mathieu equation formulation. The stratosphere is constant in depth (i.e. stratified) so that's why we set that to a constant. Yet remember that with ENSO, it is the thermocline depth that is sloshing wildly. The density differences of water above and below the thermocline is what is sensitive to the lunar gravitational forcing. That's why with a strict biennial modulation applied to the term allows an ENSO model that shows a behavior that matches the data so well:
This is a fit for ENSO from 1880 to 1950, with the validation to the right
Note how well the amplitudes vary in prediction.
In contrast for QBO, a training run from 1953 to 1983 shows this (this may be sign inverted, I don't care since the positive direction of QBO is defined arbitrarily):
The validated model to the right of 1983 shows much less variation with amplitude, in keeping with the excursion-limiting behavior of the sin-of-a-sin formulation (i.e. ).
These kinds of approximations are very common in the physics that I was taught, especially in the hairy world of solid-state physics. Simplifying the starting set of equations will not only reduce complexity but it often brings a level of understanding not normally accessible from looking at the original formulation.