• Click anywhere in the chart to reset the operating point
• Throttle the compressor towards instability by increasing \(k\)
• Try throttling the compressor with different values of \(B\)
  • Low values of \(B\) give stall: small, even oscillations
  • High values of \(B\) yield surge: larger, cyclic changes

Introduction

Surge and stall are two problems that may occur in a jet engine. This video gives a good layman’s explanation of the fluid mechanics involved (first minute is sufficient):

Stall is usually recoverable, with smaller oscillations in mass flow. Surge, on the other hand, results in complete flow reversal, with extreme forces on the machine entailing damage. It is useful to know which of stall and surge will occur in a given compression system; in the past, the distinction was not well understood.

Part of the difficulty was that simulating surge in a laboratory experiment was difficult and potentially dangerous, as suggested by the video. Greitzer (1976) developed an elegant theoretical model for the dynamic behaviour of a compressor, and showed that surge was a system instability, which occurs due to coupling between fluid inertia and storage effects in all parts of the machine.

The Greitzer model is neat because it is a good use of lumped parameter modelling and computers. The governing equations, four first-order ordinary differential equations, are simple enough in concept to be understood by a calculus beginner. We can only make things so simple by judicious application of fluid mechanics to lumped control volumes. As the equations are non-linear, they cannot be solved by hand, but the computational cost is so low we can solve them numerically right in a web browser!

Governing equations

Jet engines are complicated machines. However, we can make progress by analysing a simpler system that exhibits the same physical phenomena.

Here, we have a compressor, denoted by subscript \(\mathrm{c}\), that provides pressure rise \(\Delta p_\mathrm{c}(t)\) in an inlet duct of length \(\mathcal{L}_\mathrm{c}\) and area \(\mathcal{A}_\mathrm{c}\). The compressor discharges to a plenum with volume \(\mathcal{V}\) at a pressure \(p_\mathrm{p}(t)\). To drop the pressure back down to atmosphere at \(p_\mathrm{a}\), the flow passes through a throttle of strength \(K\) (defined later) in a duct denoted by a subscript \(\mathrm{t}\). The lengths and areas are selected so that the fluid inertia in the ducts is comparable to the total in the real machine; similarly, the plenum volume is selected to yield comparable mass storage capacity to the real machine.

Now, we need to lay out a few assumptions:

1. Low Mach number flow throughout the compression system;
2. Flow is incompressible with constant density everywhere apart from in the plenum;
3. Geometry length scales are much smaller than acoustic wavelengths;

A set of three ordinary differential equations, derived from conservation laws, governs the state of the model system. Applying Newton’s Second Law to the fluid inside the compressor duct, $$ \renewcommand{\cs}{_\mathrm{c}} \renewcommand{\ts}{_\mathrm{t}} \renewcommand{\ps}{_\mathrm{p}} \renewcommand{\Lc}{\mathcal{L}\cs} \renewcommand{\Lt}{\mathcal{L}\ts} \renewcommand{\Ac}{\mathcal{A}\cs} \renewcommand{\At}{\mathcal{A}\ts} \renewcommand{\AR}{\mathit{AR}} \renewcommand{\Vc}{V\cs} \renewcommand{\Vt}{V\ts} \renewcommand{\Vref}{U} \renewcommand{\rc}{\rho} \renewcommand{\Vol}{\mathcal{V}} \renewcommand{\pa}{p_\mathrm{a}} \renewcommand{\pp}{p\ps} \renewcommand{\Dpc}{\Delta p\cs } \renewcommand{\Dpt}{\Delta p\ts } \renewcommand{\df}{\mathrm{d}} \renewcommand{\rp}{\rho\ps} \renewcommand{\mdot}{\dot{m}} \renewcommand{\mc}{\mdot\cs} \renewcommand{\mt}{\mdot\ts} \renewcommand{\dt}[1]{\frac{\df #1}{\df t}} \renewcommand{\dT}[1]{\frac{\df #1}{\df \tau}} \Ac \left(\pa + \Dpc - \pp \right) = \rc \Lc \Ac \dt{\Vc} , $$ then putting \(\mc = \rc \Vc \Ac \) and simplifying gives, $$ \dt{\mc} = \frac{\Ac}{\Lc} \left (\pa + \Dpc - \pp\right) . \tag{I} $$ Similarly for the exit duct, $$ \At \left(\pp + \Dpt - \pa \right) = \rho \Lt \At \dt{\Vt} . $$ A loss coefficient \(K\), characterises the pressure drop across the throttle, \(\Dpt = -\rho \Vt^2 K /2\). Rewriting in terms of mass flow, $$ \dt{\mt} = \frac{\At}{\Lt} \left (\pp - \frac{1}{2}\frac{ \mt^2 K }{\rho \At^2} - \pp\right) . \tag{II} $$ Applying conservation of mass to the plenum, $$ \mc - \mt = \Vol\dt{\rp} , $$ where the density inside the plenum \(\rp\) is not constant, and allowed to vary in time. For an isentropic compression process, \(\df \rho = \df p/a^2\), and so we substitute for \(\rp\), $$ \dt{\pp} = \frac{a^2}{\Vol}\left(\mc - \mt\right) . \tag{III} $$

To clean things up a bit, we make Eqns. (I–III) non-dimensional by using, $$ \phi = \frac{\mdot}{\rho \Vref \Ac} , \quad \psi = \frac{\pp - \pa}{\frac{1}{2} \rho \Vref^2} , \quad \tau = t \omega ,\quad \omega = a \sqrt{{\Ac}/{\Vol \Lc}} , $$ where \(\Vref\) is the compressor blade speed, and \(\omega\) is the angular frequency of a Helmholtz resonator formed by the plenum and duct. After substitution and simplification, the final governing equations are, $$ \dT{\psi} = \frac{1}{B} \left(\phi\cs - \phi\ts\right) \tag{IV.i} , $$ $$ \dT{\phi\cs} = B \left(\Psi - \psi\right) \tag{IV.ii} , $$ $$ \dT{\phi\ts} = gB \left(\psi - k \phi\ts^2\right) \tag{IV.iii} , $$ involving the parameters, $$ B = \frac{\Vref}{2\omega\Lc} ,\quad \Psi = \frac{\Dpc}{\frac{1}{2}\rho\Vref^2} ,\quad g = \frac{\At\Lc}{\Ac\Lt} , \quad k = K \left(\frac{\Ac}{\At}\right)^2 . $$ The physical meanings of these are,

• \(B\) is a reduced frequency, proportional to the ratio of the natural period of the Helmholtz resonator, and the convection time through the compressor duct. \(B\) takes higher values for large plenums, where there is more mass storage capacity;
• \(\Psi\) is the usual non-dimensional form of a compressor pressure rise characteristic;
• \(g\) is a geometry parameter, relating the inertias of fluid contained within the compressor and throttle ducts;
• \(k\) quantifies pressure loss through the throttle, at a fixed value of inlet mass flow.

Compressor characteristic

As variations in compressor operating point are not quasi-steady, the pressure rise at a given instant in time is not equal to that expected in a steady flow. To model this, we introduce a (purely empirical) lag equation to be solve with Eqns. (IV), $$ \dT{\Psi} = \frac{1}{T}\left(\hat{\Psi}(\phi\cs) - \Psi \right) , \tag{V} $$ where \(T\) is a time constant, and \(\hat{\Psi}(\phi\cs)\) is the pressure rise characteristic of the compressor in steady state as a function of flow coefficient. This equation implies that the response of the compressor to a step change in mass flow is an exponential approach to the new steady-state value of pressure rise.

The steady-state pressure rise characteristic is a function of compressor design, and must in general be found from experimental measurements. To proceed, we can prescribe a polynomial form for the characteristic after Moore and Greitzer (1986), $$ \hat{\Psi}(\phi\cs) = \hat{\Psi}_0 + h \left[1 + \frac{3}{2}\left(\frac{\phi\cs}{w}-1\right) - \frac{1}{2}\left(\frac{\phi\cs}{w}-1\right)^3 \right] , \tag{VI} $$ where the parameters \(\hat{\Psi}_0\), \(h\), and \(w\) are chosen to fit experimental data.

Numerical solution method

We can solve Eqns. (IV) and (V) as follows:

1. Choose an initial state of the system, that is, values for \(\Psi , \psi , \phi\cs , \phi\ts \) at time \(\tau=0\);
2. Evaluate the right-hand sides to get new values for the time derivatives;
3. Update the state of the system at time \(\tau + \Delta \tau\) using a numerical integration scheme, namely a fourth-order accurate Runge–Kutta method.

In this way, we march the system forward in time and can observe its behaviour.

References

• Greitzer, E. M. (1976). Surge and Rotating Stall in Axial Flow Compressors: Part I—Theoretical Compression System Model. J. Eng. Power, 98:190–198.
  https://doi.org/10.1115/1.3446138
• Moore, F. K., and Greitzer, E. M. (1986). A Theory of Post-Stall Transients in Axial Compression Systems: Part 1—Development of Equations. J. Eng. Gas Turbines Power, 108:68–76. https://doi.org/10.1115/1.3239887