Nomenclature
Latin Symbols
 Ar
nondimensional rate of angleofattack
 a
axis of rotation for a pitch and plunge aerofoil
 $a_{on}$

ONERA empirical parameter
 b

aerofoil semichord
 [C]

damping matrix
 C(k)

theodorsen function
 $C_{Am}$

Gangwani moment stall empirical parameter
 $C_{At}$

Gangwani trailing egde empirical parameter
 $C_{L_0}$

zeroangle lift coefficient
 $C_{L1}$

ONERA linear aerodynamics lift coefficient
 $C_{L2}$

ONERA nonlinear aerodynamics lift coefficient
 $C_{L_c}$

circulatory lift coefficient
 $C_{L\ell}$

static lift coefficient
 $C_{L_{nc}}$

noncirculatory lift coefficient
 $C_M$

pitching moment coefficient
 $C_{M_0}$

zeroangle pitching moment coefficient
 $C_{M1}$

ONERA linear aerodynamics pitching moment coefficient
 $C_{M2}$

ONERA nonlinear aerodynamics pitching moment coefficient
 $C_{M_c}$

circulatory pitching moment coefficient
 $C_{M\ell}$

static pitching moment coefficient
 $C_{M_{nc}}$

noncirculatory pitching moment coefficient
 $C_N$

normal force coefficient
 $C_{wR}$

Gangwani reattachment empirical parameter
 $C_{Wm}$

Gangwani moment stall empirical parameter
 c

chord
 $c_{\ell}$

lift coefficient
 $c_{l_{\alpha}}$

lift curve slope
 $c_{m,0.25c}$

pitching moment coefficient at quarter chord point
 E

Young’s modulus
 EF

stress excitation factor
 $E_{on}$

ONERA empirical parameter
 e

distance between mass and elastic axes
 GJ

torsional stiffness
 g

structural damping coefficient
 $H_n^{(2)}$

Hankel function of the second kind
 h

aerofoil plunging motion
 $I_1$

chordwise moment of inertia
 $I_y$

mass moment of inerta
 [K]

stiffness matrix
 $K_1$

Gormont pitch rate coefficient
 k

reduced frequency
 $k_0$

nondimensional aerodynamic centre offset from the quarter chord
 $k_{1/2}$

LeishmanBeddoes semiempirical coefficients
 $k_y$

blade area radii of gyration
 $k_{xx},k_{zz}$

blade mass radii of gyration
 L

aerodynamic lift
 LF

KirchoffHelmholtz unsteady lift factor
 M

aerodynamic pitching moment
 ${\rm M}_a$

Mach number
 m

section mass
 $m_i$

generalised mass for ith mode
 n

number of modes
 PF

KirchoffHelmholtz unsteady pitching moment factor
 Q

number of blades
 $Q_i$

generalised force for ith mode
 $q_i$

modal coordinate for ith mode
 R

blade radius
 $R_{hub}$

hub radius
 Re

Reynolds number
 r

nondimensional blade radius
 $r_{on}$

ONERA empirical parameter
 s

aerodynamic time variable, Wagner function
 $s_{1/2}$

LeishmanBeddoes semiempirical coefficients
 $s_{{1/2}_n}$

negative angle LeishmanBeddoes semiempirical coefficients
 $s_{mt}$

time at which the vortex detaches from the trailing edge
 $s_{on}$

ONERA empirical parameter
 Ta

blade axial tension
 t

time
 $t/c$

thickness/Chord ratio
 U

freestream velocity
 $V_{eq}$

equivalent airspeed
 $x_a$

distance from axis of rotation to aerodynamic centre
 $x_s$

nondimensional static upper surface separation point
 $x_{s_E}$

LeishmanBeddoes effective separation point
 $y_{CG}$

centre of gravity offset
 $y_{EA}$

elastic axis offset
 z

blade flapping motion
Greek symbols
 $\alpha$

angleofattack
 $\alpha_0$

initial angleofattack
 $\alpha_{Dm}$

pitching moment stall angle
 $\alpha_E$

effective angleofattack
 $\alpha_{QS}$

quasisteady angleofattack
 $\alpha_s$

staticstall angleofattack
 $\alpha_{s_n}$

negative staticstall angleofattack
 $\alpha_W$

unsteady decay parameter
 $\alpha_{Wm}$

unsteady decay parameter at moment stall
 $\beta$

PrandtlGlauert compressibility correction
 $\gamma$

blade taper ratio
 $\gamma_g$

Gormont nondimensional empirical parameter
 $\Delta \alpha_{DS}$

incremental dynamic stall angle
 $\Delta C_L$

difference in lift coefficient
 $\Delta C_M$

difference in pitching moment coefficient
 $\epsilon$

Gangwani empirical parameter
 $\zeta_i$

critical damping for the ith mode
 $\theta$

aerofoil pitching motion
 $\Lambda$

blade sweep
 $\lambda_{1/2}$

Dowell’s Theodorsen function coefficients
 $\lambda_{on}$

ONERA empirical parameter
 $\rho$

air density
 $\sigma$

rotor solidity
 $\sigma_{on}$

ONERA empirical parameter
 $\tau$

nondimensional time
 $\tau_1$

semiempirical relaxation time constant
 $\tau_2$

LeishmanBeddoes timedelay parameter
 $\phi_c(s,M_a)$

compressibility corrected Wagner function
 $\psi$

phasing between blades
 $\psi_i$

freestream inflow angle
 $\Omega$

rotational velocity
 $\omega$

natural frequency
 $\omega_i$

natural frequency for the ith mode
Acronyms
 BEMT

Blade Element Momentum Theory
 CFD

Computational Fluid Dynamics
 CG

Centre of gravity
 eVTOL

electrical Vertical TakeOff/Landing
 HMB3

Helicopter MultiBlock 3
 NACA

National Advisory Committee for Aeronautics
 NASA

National Aeronautics and Space Administration
 SAS

ScaleAdaptive Simulation
 SST

Shear Stress Transport
 URANS

Unsteady ReynoldsAveraged NavierStokes
 UTRC

United Technologies Research Centre
1.0 Introduction
There are at least three known types of propeller flutter: bendingtorsion flutter, stall flutter and whirl flutter. Bendingtorsion flutter, or classical flutter, is driven by the coupling and excitement of selected modes of vibration. For this type of flutter, the aerodynamics can often operate in the linear regime with nonlinearity introduced via the coupling of modes due to the combined effects of structural and aerodynamic forces. To understand the coupling procedure, a stepbystep approach can be taken. Any change in effective angleofattack, whether this be driven by internal or external effects such as inflow fluctuations, results in changes to the blade lift. This results in a flapping of the blade which triggers the torsional mode. Depending on the current freestream velocity, this excitation can grow, dampen or selfsustain. As a result of the linear aerodynamics, lower fidelity, fast models can be utilised which predict the critical flutter boundaries accurately. Although a good level of accuracy is found when assessing for the critical flutter points, the models ability to predict the blade response outside of this critical flutter point loses accuracy. Moreover, such tools can be used to understand and, potentially, adjust such points to ensure that the propeller operation avoids such excitement.
Although several sources of nonlinearity can be found within a propeller blade when looking at classical flutter, such as the the four sources of deflection arising from nonlinear inertial loading, a type of flutter which is inherently more nonlinear than classical flutter is that of whirl flutter. This is due to the combination of bladewingnacelle structure with the propeller gyroscopic and aerodynamic effects adding to the sources of nonlinearity already found within a standard propeller. Such a response can be found in helicopters, turboprop aircraft and tiltrotors, with turboprop aircraft being the principle test case. Whirl flutter for tiltrotors is more complex than for standard turboprops due to the ability of the tiltrotor to operate in different flight regimes, including helicopter, transition and propeller modes. As a result of the aerodynamics associated with propeller systems, nonlinear aerodynamics becomes more critical and, hence, the modelling becomes more complex. Therefore, the combination of more complex structural models with nonlinear aerodynamics significantly increases the sophistication of the modelling requirement. Simpler models that assume rigid propellers, for example, are available but these often produce conservative flutter boundaries with respect to the experiment. Although conservative boundaries, from an operational pointofview and ensuring efficiency of the flight path, is not terminal as whirl flutter is typically a catastrophic oscillation and therefore, avoidance of the boundary is required at all costs.
Stall flutter takes a combination of the blade structural nonlinearity with the unsteady aerodynamics associated with detached flow. Stall flutter is typically a one degree of freedom oscillation, which is triggered by detached flow. It is often found in the torsional mode for propellers due to the introduction of a pitching moment as a result of a change in the local centre of pressure. It is the detached aerodynamics that introduces significant nonlinearity to this phenomenon, and the trigger of it can be a result of the blade design or flight conditions. For example, typical conditions for a propeller to undergo stall flutter is during takeoff. This is due to the requirement of high thrust/pitch during the acceleration phase of the aircraft. Due to the highly dynamic and, sometimes, sudden nature of stall, significant reductions to the overall propeller flutter boundary can be found due to this phenomenon. Such results are discussed and presented in this paper.
Although stall flutter is typically a transient and shortlived response, the development of modern propeller blades with high sweep/taper and thin aerofoil sections is likely to result in changes to the flutter boundary. Furthermore, a new generation of propellers is needed to power electrically driven aircraft, both fixed and rotarywinged. In the latter category, there is a great variety of multirotor configurations over a range of scales and Reynolds numbers (including eVTOL), operating with variable angular speeds. As a result, a reinvestigation of the stall flutter is required to ensure safety of the new designs and to also push the boundaries of the design.
A high number of stall flutter studies focus on aerofoil dynamic stall. In this contribution, a description is presented on dynamic stall; links to aerodynamic damping and propeller stall flutter boundaries are demonstrated. Details are then shown on the propeller stall flutter. This includes historic experiments and modern numerical studies. Stall flutter associated with turbomachinary and helicopter rotors is not included in this review, except within the context of the general stall flutter description.
2.0 Stall flutter description
One of the characteristics associated with helicopter rotors in skewed flight is the stall effects on the retreating side due to the combination of the angular and freestream velocities. Similar stall effects can also be found in a propeller during a high thrust/pitch setting such as takeoff. If the blade exhibits elastic behaviour, there is a potential for further blade oscillations. This type of flutter is classed as stall flutter and is due to the combination of structural dynamics and nonlinear aerodynamics. It remains a challenging problem to this day.
Stall on a rotor blade due to the control pitch oscillations is termed as dynamic stall and can be studied as a twodimensional aerofoil problem. The inclusion of dynamic stall on an aerofoil introduces a unique set of characteristics not found in static aerodynamic stall. The use of quasistatic lookup tables for lift, drag and moment allows the calculation of rotor loads for hover, vertical flight or lowspeed level flight. This method has been verified with rotor load studies and is justified due to the relatively low angles of attack. However, at high speeds such quasistatic tables cannot be used; investigations such as those of Kufield and Bousman [Reference Kufield and Bousman1] pointed out the difficulty when the airloads of the UH60A helicopter rotor were measured in flight. The airloads indicate a significant deviation from the twodimensional aerofoil results. To understand this discrepancy, an examination of both the staticstall and dynamic stall effects on an aerofoil is presented.
2.1 Aerofoil staticstall
The range of angleofattack is significant in rotorcraft performance. The lift, drag and moment coefficients are examined from attached flow to deep stall; the typical trend for each case is presented in Fig. 1. This data is adapted from the experiments of Critzos et al. who studied the NACA 0012 static coefficients between $0^{\circ}$ and $180^{\circ}$ at Reynolds numbers $Re = 1.8\times10^6$ and $0.5\times10^6$ [Reference Critzos, Heyson and Boswinke2].
A linear increase in lift is found during the initial stages and this remains until the static stall angle ( $\alpha_s$ ). This angle is defined based upon the location of the maximum lift coefficient ( $c_{\ell_{max}}$ ) and can often be accompanied by detached flow. Following this static stall angle, the lift coefficient decreases in value for a short range before gradually increasing towards a secondary peak. This secondary peak is typically under the presence of deep stall where threedimensional detached flow effects introduce pressure reductions causing the increase in lift.
Up to the point of static stall, the drag and moment coefficients remain fairly constant. Following this point, both the drag and moment increase sharply; the pitching moment experiences a sharp nosedown effect resulting from a rearward shift in the centre of pressure. Both coefficients continue to increase with the angleofattack, with the drag reaching a maximum at $90^{\circ}$ , when the aerofoil is perpendicular to the flow.
2.2 Aerofoil dynamicstall
A common feature among aerofoils undergoing light and deep dynamic stall is the presence of hysteresis, which is caused by a time lag between forcing and aerodynamic response. This hysteresis depends on several parameters, including frequency, angleofattack amplitude, Mach number, turbulence, etc. Typical trends are shown in Fig. 2 for the lift and moment coefficients. One of the effects of a dynamically oscillating aerofoil is the increase in maximum lift coefficient and stall angle, Fig. 2(a). This gives an indication as to the difference in airloads between the lowangleofattack lookup tables and what is observed in the true blade. An increase in maximum lift is potentially beneficial in terms of the blade operation, as momentary increases can be obtained at lower angles of attack. However, an increased load can potentially lead to an increase in blade torsion, which pushes the blade towards a potential stall flutter boundary due to the introduction of a larger pitching moment.
Likewise, a pitching moment hysteresis loop can be obtained, Fig. 2(b). A similar increase in moment above the static stall is found; this is related to the increase in lift. The number of loops included within the profile depends on the angle at which the aerofoil is oscillating around and the amount of detached flow associated with the condition. These loops can categorise the stability of an aerofoil section via the determination of aerodynamic damping in pitch. Aerodynamic damping is calculated as the summation of areas of each individual loop, with a clockwise loop indicating negative damping and anticlockwise loop producing positive damping.
2.3 Aerodynamic damping
Aerodynamic damping is the result of aerodynamic forces and moments exerted on to a structure. Aerodynamic damping often opposes structural damping, and can potentially result in aeroelastic instabilities. Aerodynamic damping is critical to the flutter characteristics of a structure. Above the flutter velocity, the work of the given fluid on a structure is said to be negatively damped. Thus, the structure’s oscillatory motions tend to increase with time unless an external factor such as a cyclic input adjusts the given pitch angle.
Stall flutter originates from separated flow and is found to be present in helicopter rotors, propellers and other rotating wings. A stall flutter instability can only be corrected via positive structural damping, change of stiffness or a change in the aerodynamic conditions. As a result, an investigation into stall flutter can begin from the aerodynamic damping of a system. Damping estimation is often performed for aerofoils and full threedimensional calculations are rare.
Early investigations of stall flutter were twodimensional and experimental [Reference McAlister, Pucci, McCroskey and Carr4–Reference McCroskey6]. These investigations focused on the determination of the aerodynamic coefficients during dynamic stall to improve helicopter rotor performance in forward flight. Such 2D investigations of oscillating aerofoils highlighted the trends seen during differing stall regimes. These regimes are highlighted in Fig. 3, where the pitching moment coefficient is shown for a pitching NACA0012 aerofoil.
No stall : The aerofoil motion remains below the static stall angle and the use of quasisteady aerodynamics is sufficient enough to predict the aerofoil loading. Both the lift and pitching moment coefficients are found to circle anticlockwise with no crossing of the downstroke and upstroke profiles. Therefore, aerodynamic pitch damping is found to be positive.
Stall onset : The aerofoil motion reaches the static stall angle. There is often found a slight reduction within the area of the anticlockwise loop, however, no crossing of the profiles are seen and therefore quasisteady aerodynamics can be used to estimate the loads [Reference Leishman and Beddoes7, 8].
Light stall : It is within this regime that dynamic stall vortices are present. These vortices are a combination of shedding and passage events that occur across the upper surface of the aerofoil. A typical, and clear, pattern often exists where there is a buildup of detached flow at the leading edge following an attached period. This then passes over the upper surface before shedding as it reaches the trailing edge. Often, secondary events can occur under the right conditions. The cycle then returns to the attached state. During light stall, the aerofoil motion reaches values higher than the static stall angle, with the aerofoil loads characterised by a hysteresis effect and it is the dynamic stall vortices which introduce the large hysteresis. The development of separated flow regions are found to be sensitive to the aerofoil geometry, freestream Reynolds number and Mach number, reduced frequency of the aerofoil oscillation, and the mean and harmonic angles of attack. It is also within the regime that there is the highest tendency towards negative aerodynamic damping.
Deep stall : The aerofoil motion is often found to pitch entirely beyond the static stall angle and it is where the strongest effects of the dynamic stall vortex are seen. The aerofoil loading is characterised by a strong hysteresis effect, with significantly larger peak load coefficients and multiple damping loops.
Following the initial propeller stall flutter studies, a transition was made from propeller driven aircraft to turbofan engines. This resulted in a reduction in propeller research with attention paid to the more immediate issue of helicopter dynamic stall. It is from this research that a significant amount of understanding of aerodynamic damping was gleaned. The above description of dynamic stall is focused on helicopter applications, and as such, involves large amplitude variations in pitch due to the cyclic inputs of a helicopter in forward flight. In propeller dynamic stall, significantly smaller oscillations in pitch are found around the static stall angle due to the lack of cyclic inputs, and therefore the mechanism which triggers the stall flutter phenomena must be taken from the context of aerodynamic damping.
Following some of the initial experimental propeller stall flutter investigations [Reference Baker9], attempts were made to better understand the stall flutter boundaries with respect to aerodynamic damping. This included the experiments of Fanti et al. [Reference Fanti, Carta and Pitt10] & Lemnios [Reference Lemnios11]. However, as described by Carta & Lorber [Reference Carta and Lorber12], such estimations were seen to be inadequate for design and deeper understanding of the mechanism due to insensitive instruments. Carta & Lorber conducted an experimental investigation into the variation in aerodynamic damping due to small amplitude oscillations in pitch of a propeller aerofoil section [Reference Carta and Lorber12]. This study followed the successful demonstration of the damping extraction for a cascade environment [Reference Carta13].
Experiments were performed in the United Technologies Research Centre (UTRC) wind tunnel. Sinusoidal oscillations of the Sikorsky SC1095 aerofoil section were performed over a range of amplitudes, mean angles and reduced frequencies. Pitch amplitudes ranged from $0.5^{\circ}$ to $4^{\circ}$ , with the mean angle between $8^{\circ}$ and $10.5^{\circ}$ and reduced frequencies between 0.05 and 0.16. From initial static results, it was shown that the static stall angle was $\sim9.5^{\circ}$ for this section. Pressure transducers were installed in the section in order to determine surface pressure and load profiles. The results indicate a significant influence of these parameters. The static derivative results were shown to be sensitive to the mean value when oscillating at low angles of attack for a range of reduced frequencies. The derivatives were found to vary significantly once the mean angle approached the static stall. This variation was seen to be random, with this random nature also present within the load time histories at low oscillating amplitudes [Reference Carta and Lorber12]. However, at higher oscillation values such as those seen in helicopter rotors, reduced frequency becomes more influential with almost no trend captured due to changes in mean angle. Therefore, in the interest of propeller stall flutter which takes into account smaller oscillation values, the mean value sensitivity is critical to the test case in question.
Pitching moment hysteresis loops similar to Fig. 3 were presented, which determine the aerodynamic damping estimations. For this experiment, lowamplitude oscillations were found to induce a dominant clockwise loop when the mean angle was selected above the reference static stall angle. The trend of aerodynamic damping from these given profiles is presented in Fig. 4. For the highest oscillation value, the mean angleofattack is seen to have very little influence on the aerodynamic damping value with an almost constant negative value produced. However, the greatest variation in aerodynamic damping comes from the lower values with the destabilising range being dependent on the mean angle. An overall negative damping result can be mitigated via the blade design. However, the variation in damping seen from the $0.5^{\circ}$ and $2^{\circ}$ plots, is more difficult to design for; therefore, such lower amplitude profiles become more critical.
Using the data gained from experimental investigations, such as Carta & Lorber [Reference Carta and Lorber12] and McCroskey [Reference McCroskey6], validation of aerodynamic models could be undertaken with respect to dynamic stall and aerodynamic damping. This involved theoretical and semiempirical models, with a review conducted by Ericsson & Reding [Reference Ericsson and Reding14]. However, with the need for experimental parameters to ensure the accuracy of semiempirical formulations, a review was conducted by Ekaterinaris & Platzer [Reference Ekaterinaris and Platzer15], which focused on the use of the NavierStokes equations for aerofoil aerodynamic damping extraction. This was brought on by the development of the numerical methods and increased computational power to allow for NavierStokes methods to be used in dynamic stall modelling.
The use of NavierStokesbased modelling can not only assist in the extraction of damping estimations but also determine steady airloads for the semiempirical formulations, should no experimental data exist. Such a procedure was conducted by Beedy et al. [Reference Beedy, Barakos, Badcock and Richards16] who utilised the ONERA aerodynamic model developed by Tran and Petot [Reference Tran and Petot17] and Dat and Tran [Reference Dat and Tran18] in the 1980s, with CFD coefficients for the study of aerofoil stall flutter. Such a procedure was found to qualitively correlate well to experimental data [Reference Dunn and Dugundji19] in terms of flutter velocity and frequency, with a comparison shown in Fig. 5. As observed, a rapid increase in frequency is found around 30m/s, correlating with the experiment. At this point, further work was needed to assess the ability of the turbulence model to predict the required ONERA coefficients; however the results were deemed to be satisfactory and thus highlighted the ability to use NavierStokes CFD with semiempirical models should typical experimental values be unavailable.
Ekaterinaris & Platzer [Reference Ekaterinaris and Platzer15] discussed the effects of numerical schemes, turbulence modelling and the effect of transition, with the turbulence and transition modelling seen to have a significant effect on the experimental correlation. They concluded that compressibility, transition, and flow reattachment were key flow effects limiting the correlation between the experiments and CFD results. Compressibility effects are driven by the partially supersonic profile observed across the aerofoil section. This is typically seen on the advancing side of the rotor due to the combination of freestream and rotational velocities combining. The use of transition models was seen to have a significant effect on the leading edge separation and, with the transition models of the time utilising empirical coefficients, an improvement was required to enhance the correlation.
The final factor was related to flow reattachment, and this was seen during the aerofoil downstroke where poor agreement was shown. This is often found in the light and deep dynamic stall regimes where the influence of the threedimensional detached flow effects become important. Twodimensional CFD cannot capture this effect, and therefore, an increase to the computational expense is required. This comes in the form of threedimensional calculations and/or resolved flow modelling. Resolved flow modelling allows for the shedding process associated with dynamic stall process to be captured, potentially improving the downstroke performance.
In order to improve the correlation and mitigate the effects of flow reattachment, an investigation was conducted by Higgins et al. who utilised resolved flow modelling across aerofoil sections undergoing deep dynamic stall [Reference Higgins, Barakos and Jinks20]. Quasi3D simulations of the NACA 0012, $70\%R$ and $90\%R$ Commander blade sections were conducted. The simulations were termed quasi3D due to the use of periodic boundary conditions where typical twodimensional planes are applied. This allowed for the use of scaleresolving techniques to be used, and for this investigation the ScaleAdaptive Simulation (SAS) method was employed. The SAS model is a hybrid URANSLarge Eddy Simulation method used to resolve the stalled flow [Reference Menter and Egorov21]. For the SAS formulation, an additional source term is added to the $k\omega$ Shear Stress Transport (SST) equations, which allows for the local adjustment to the vonKarman length scale and balances the contributions from resolved and statistical components. Using this method, an improvement on the NACA 0012 experimental load correlation and aerodynamic damping estimation was found against the twodimensional URANS results. This was driven by the ability of the method of capture the shedding of vorticity during the heavily detached downstroke/reattachment phase. This emphasises the potential need for highfidelity aerodynamic modelling during detached flow conditions.
Furthermore, many studies have been conducted in recent times which focused on twodimensional aerofoils. Such investigations made use of pitch and plunge rigid models [Reference Verstraelen, Habib, Kerschen and Dimitriadis22, Reference Razak, Andrianne and Dimitriadis23, Reference Yamasaki, Isogai, Uchida and Yukimura24], cyberphysical models (where the structural response of the aerofoil is given based upon a known set of derived equations of motion) [Reference Culler and Farnsowrth25, Reference Culler, Fagley, Seidel, McLaughlin and Farnsworth26], or if numerical, used dynamic stall based aerodynamic models [Reference Bryant, Gomez and Garcia27, Reference Arena, Lacarbonara and Marzocca28]. The use of such studies can provide fundamental analysis into aerofoil stall response, with conservative boundaries when applied to threedimensional test cases.
3.0 Propeller stall flutter
The topic of propeller flutter is split into three components, and this includes bendingtorsion, stall and whirl flutter. These types of flutter, particularly whirl and stall, require detailed modelling of the aerodynamics and structural response of a propeller. For stall flutter, the nonlinear aerodynamics is a result of the detached flowfield, which triggers the aeroelastic excitation. Successful capture of such aerodynamic interactions allows for increased accuracy in surface loads, and increased accuracy in predicting the resultant flutter boundary.
3.1 Experimental studies
Propeller stall flutter research began during World War II and was primarily motivated by the use of propellerdriven aircraft such as the Supermarine Spitfire and North American P51 Mustang. This research began with two studies: one from the USA and another by the UK. The UK study was conducted by Sterne and focused on the flutter analysis of the 4bladed Spitfire propeller [Reference Sterne29]. The USA investigation was conducted by Theodorsen & Regier and used a range of windtunnel designed propellers [Reference Theodorsen and Regier30].
Focusing on the UK study [Reference Sterne29], static experiments were conducted to determine the flutter boundary of the propeller, with the boundary presented in the form of the blade pitch versus blade rotational velocity, Fig. 7. The blade operated at 19 different pitch angles ranging from $8^{\circ}$ to $32^{\circ}$ (the pitch was measured at the $70\%R$ station). Note that a reference pitch of $8^{\circ}$ at $70\%R$ results in a $0^{\circ}$ pitch angle at the blade tip. As observed from Figure 7(a), the experiment found that the blade would flutter at a much lower velocity within a specific pitch region. This region ranged from $20^{\circ}$ to $28^{\circ}$ with a 24% reduction in flutter velocity found.
Sterne concluded that the distinct decrease in velocity between $20^{\circ}$ and $28^{\circ}$ was the result of the presence of stall, hence stall flutter oscillations. The remaining pitch angles were judged to be a result of classical bendingtorsion flutter. This observation was also found by the work of Stüder who found a similar drop in velocity across a model wing due to stall and theorised it to be a result of negative aerodynamic damping [Reference Stüder31].
In addition to this spin test on the Spitfire propeller, Sterne conducted an experiment to investigate the influence of the number of blades on the flutter characteristics. Those tested included a 2, 3 and 4bladed Firefly propeller. The Fairey Firefly blade was constructed in the same manner as the Spitfire, using compressed wood, with a slightly larger diameter of 13 feet compared to 10 feet 9 in. A comparison of the thickness/chord ratio and the solidity for the 4bladed cases is presented in Fig. 8.
From the results in Fig. 7 it was implied that the number of blades have an influence on the flutter boundary. In the classical flutter boundary range, the 2 and 4bladed propellers produced similar results with the 3bladed propeller operating at an angular velocity $\sim90 \mbox{rpm} $ higher. One difference in the experimental setup between the 2 or 4bladed propeller and the 3bladed was a slightly different hub. It was expected that such a difference would have an influence on the flutter boundary and, hence, the increase within the classical range could be explained via this difference. In terms of the stall flutter range, a significant reduction in stall flutter boundary was observed in the 3bladed propeller compared to the 2 and 4bladed. Sterne concluded that the reduction in stall velocity could be a result of the greater influence of the preceding blade stall upon the oncoming blade. With the experiment conducted in static conditions, and no thrust data known, the blade wake could very well interact with a subsequent blade. Additionally, the influence of the propeller induced velocity is unknown due to the lack of thrust data and therefore, this propeller could be operating on the edge of a wake interaction that is felt for the 3bladed propeller but not the 4bladed due to the expected increase in induced velocity (assuming all other parameters remain the same).
The experimental investigation of Theodorsen & Regier (the USA wartime report on propeller flutter) [Reference Theodorsen and Regier30] included a range of experiments on model propellers for use, and examination, in the windtunnel. This investigation was conducted in connection with the construction of the Langley, Ames and Cleveland blades for use within the respective tunnels. In this study, two propeller blades were presented. As in the case of the Spitfire/Firefly blades, both were made of laminated wood (spruce) with flatbottom Clark Y aerofoil sections. The first propeller (Propeller A) was designed as a 6bladed propeller with a 45inch diameter. The second propeller (Propeller B) was designed as single bladed with a reduction in the chord length and thickness/chord ratio of one third. This reduction in size was to reduce the flutter velocity in order to examine the results.
The experiments were conducted in an open wind tunnel in static conditions, with only the induced velocity present at the tunnel exit. For Propeller B, this induced velocity was enhanced by a booster fan at the rear of the motor due to being of single blade design. A stroboscope was used at the blade tip in order to measure the bending and torsion, with strain gauges also used. With the blade pitch fixed, the rotational velocity was increased until flutter was determined for a range of propeller lift coefficients. These lift coefficients were adjusted through the change in the wind tunnel exit area, i.e. essentially compressing or expanding the wake. The change in lift coefficient was confirmed via pitot tubes.
From the experiment it was observed that at the lowest flutter velocities, the oscillations were found to be pure torsional fluctuations. This was determined from the strain gauges of Propeller B. It was concluded from this experiment, as in the Sterne experiments, that the stall flutter velocities are much lower than the classical boundary. In order to mitigate this effect, Theodorsen & Regier concluded that the design blade angle should be such that the lift coefficient produces zero blade twist. Equations are presented to determine such angles, with the effect of addition twist seen as inducing further detached flow. Approaching the classical flutter boundary, less and less detached flow was required to excite the blade.
Although both studies provide a good insight into propeller stall flutter, particularly the influence of stall on the propeller flutter region, the construction of the blades make it difficult to study numerically. In comparison to modern manufacturing processes, the compressed wood construction has the potential to introduce manufacturing inconsistencies and, therefore, make the understanding of the blade structural properties potentially difficult. From this basis, the Spitfire, Firefly and NASA blades could not be used for a numerical study.
Following this observation of classical and stalled flutter within the working range of the wooden Spitfire and Firefly propellers, Ewing et al. conducted a experiment on the duralumin blades of the Hawker Tempest 5bladed aircraft [Reference Ewing, Kettlewell and Gaukroger32]. This experiment was conducted in the spinning tower of the Royal Aircraft Establishment, with strain gauges applied to the blades to measure the vibratory stress. Each propeller was pitched below and above the stalling region with only stall flutter found during the runs.
The experiments found a heavy dependence of the stress results on the phase differences between subsequent blades. This phase difference ( $\psi$ ) was given by $\psi = 2 \pi p / Q$ , where Q is the number of blades and $p = 0, 1, 2 \ldots$ . As p changed value, there was a change in the phase relationship and waveform. The experiments found the value of p to increase with blade rotational velocity at a given pitch angle. Due to this influence of phase, an exact comparison of the stress values for the number of blades could not be made, hence general stresses are compared. For the 3, 4 and 5bladed propellers, the greatest stress values occurred at $24^{\circ}$ . This reduced to $20^{\circ}$ for the 2bladed case. Either side of $24^{\circ}$ , the 4 and 5bladed propeller stress values drop rapidly with a progressively lower gradient for the 3 and 2bladed propellers. From this experiment, it was observed that the examination of the number of blades of the same design could only provide a rough indication of stress level due to the influence of phase. For example, the velocity at which the peak stress occurred may coincide with a phase change and thus be at a lower stress level than the other.
In 1955, [Reference Baker9] Baker conducted a static experimental investigation into the effects of various parameters on the flutter boundary of a single model propeller blade to determine the minimum stall flutter condition. These parameters included the effects of the blade structure, the blade geometry, and the freestream flow conditions, with the key list and range presented in Table 1. Due to the change in certain parameters, a range of blades were manufactured. For the construction, a 16series NACA aerofoil was used with the material selected as aluminum alloy, steel or maple, depending on the model and test case. All propeller blades were tested across the full pitch angle range.
The experiments were conducted in the Langley vacuum sphere to allow for the use of Freon12 to achieve more realistic Reynolds numbers during the density and Mach number studies. Strain gauges were attached to the blades to measure the bending and torsion oscillations with the rotational velocity gradually increased until flutter was observed.
The experimental results showed very little difference in the boundaries for the torsional stiffness, blade taper, blade twist, length/chord ratio and density studies. Substantial changes were found for the centreofgravity, sweep, thickness and Mach number runs with the nondimensional flutter velocity results presented in Fig. 9 for each. The models used for the results shown are described in Table 2.
^{*} At tip, polynomial sweep
The most significant change in the stall flutter boundary was found in the centreofgravity results (Fig. 9(a)), where a 14.5% shift in CG results in a 44% reduction in the minimum flutter boundary. Little difference in the thickness ratio results (Fig. 9(b)) was found between the 3% and 6% plots with the 9% model resulting in a $\sim25\%$ increase in the boundary. Although a benefit was found using a thicker aerofoil section, propeller performance dictates thin aerofoils.
The introduction of sweep (Fig. 9(c)) to the blade design has a positive influence on the stall flutter boundary, with all sweep configurations seeing an increase in the flutter velocity of greater than 25%. Examining the results from a pitch angle of $10^{\circ}$ to $35^{\circ}$ , the nonlinear sweep has the greatest continuous increase in the boundary across this region. The final observation from the investigation was the increase in flutter boundary for Mach ${\rm M}_a > 0.6$ (Fig. 9(d)). The Mach number influence was tested across two models at different pitch angles. The second test (Model 6), was seen to have a more gradual increase in the boundary, with the first (Model 5) having a steeper gradient above 0.6.
Baker did note that the resultant flutter boundaries were not within the normal operating cruise conditions for typical propellers, and it was only during takeoff like conditions that the blades were found to suffer from stall flutter. Following the experiment of Baker, greater insight into several aerodynamic and structural influences on the stall flutter boundary were found due to this investigation. Although some details are available for the blade used (Table 2), greater detail would be required to map the blade inertias along the radius for a numerical study.
In 1956, Hubbard et al. [Reference Hubbard, Burges and Sylvester33] followed on from the work conducted by Baker by focusing on three specific parameters. The primary aim was to determine the effect of Mach number, structural damping and builtin twist on the flutter boundary. For these parameters, four singlebladed model propellers were constructed (one each for the Mach number and twist study with the final two for the structural damping). Each blade was defined in the same way as the Baker propellers [Reference Baker9] with Table 3 presenting the blade descriptions. Model 3B was constructed using laminated steel in order to change the internal structural damping. In the same manner as the Baker experiments, the propeller models were tested in the Langley vacuum sphere in static conditions. During the experiment, strain gauges were used to measure the vibratory stress levels with runs conducted at a fixed pitch angle, only varying the rotational velocity.
Figure 10 shows the experimental results for the twisted blade, Mach number and structural damping studies. Focusing on the twisted blade result (Fig. 10(a)) and using the strain gauges, the type of blade oscillations across the examined pitch range was determined. Through the low pitch range ( $8^{\circ}$ to $8^{\circ}$ ), the blade oscillations were dominated by bendingtorsion or pure bending oscillations. The highest flutter velocities were found during this region with an average value 228% above the minimum nondimensional velocity of 1.26. The reduction of the pitch angle to highermagnitude negative angles resulted in the activation of the torsional mode, reducing the flutter velocity. At high positive angles, the traditional shape expected from the Baker [Reference Baker9] and Sterne [Reference Sterne29] experiments was found. Due to the detached flow at these angles, the flutter velocity was reduced significantly compared to the classical flutter bounds due to the active torsional mode.
The effect of Mach number (Fig. 10(b)), where the effect of Mach number is depicted by values of soundspeed coefficient which is defined as $b \omega_{\alpha} / a$ with b being the semichord, $\omega_{\alpha}$ the natural torsional circular frequency and a the speed of sound, on the propeller flutter velocity mirrored that found by Baker [Reference Baker9]. As the speed of sound was decreased, thus increasing the Mach number, the extent of the flutter region was reduced. This is likely down to a compressibility effect which limits maximum lift and delays the stall onset.
In terms of structural damping (Fig. 10(c)), an increase in the internal damping was seen to increase the flutter velocity. This was to be expected and therefore any oscillations due to a negative aerodynamic damping were mitigated by the increase in structural damping.
The combined works of Baker [Reference Baker9] and Hubbard et al. [Reference Hubbard, Burges and Sylvester33] allowed for a greater understanding of propeller stall flutter due to the volume of parameters investigated. However, the blades used in these experiments do not represent a realistic modern blade and this has a significant effect on the observed flutter boundary.
Similarly, Allis and Swihart conducted a static experimental investigation to determine the effect of bladesection camber [Reference Allis and Swihart34]. This work followed on from the experiments of Wood & Swihart [Reference Wood and Swihart35] who found an increase in stall flutter boundary due to an increase in camber at one specific pitch angle. This pitch angle was extended to test over a range from $16^{\circ}$ to $38^{\circ}$ . Three 2bladed NACA propellers, of varying camber, were tested in the Langley dynamometer with the blades fitted with strain gauges. Figure 11 presents the installed blades on the dynamometer. At a fixed pitch angle for the given blade, the propeller rotational velocity was slowly increased until flutter was encountered.
Presented in Fig. 12 is the nondimensional flutter velocity and thrust coefficient results for the range of cambered propellers. Focusing on the flutter boundary result (Fig. 12(a)), it was shown that at $16^{\circ}$ an increase in flutter velocity of 45% was achieved with the highest cambered propeller compared to the symmetrical blades. This difference slowly reduces with increased pitch angle until $30^{\circ}$ where the flutter velocity values are approximately the same. Above $30^{\circ}$ , an increase in the flutter velocity was seen by all, with the symmetrical propeller having the highest value at the largest pitch angle. Allis and Swihart observed that at $16^{\circ}$ , only the inboard radial stages were stalled and this increased to the entire blade at $30^{\circ}$ .
One of the key observations from the experiment was the significant difference in the flutter boundary at $16^{\circ}$ . The flutter boundary was not only driven by the combination of detached flow and structural damping, but also the total airloads applied to the blade. At higher thrust, higher fluctuations in blade deformation was expected. This is what was observed at the $16^{\circ}$ pitch angle, with this confirmed in Fig. 12(b) where the highest camber blade was shown to have a 36% increase in thrust in comparison to the symmetrical blade. With greater thrust comes greater induced velocity and therefore a reduced effective angleofattack at a given rotational velocity. In order to assess the true effect of bladesection camber, a trimming process would be required to ensure the same thrust was generated from the blade.
Rogallo & Yaggy conducted an experimental campaign on a 3bladed 10foot propeller to determine the effect of positive/negative thrust and thrust axis inclination on the stall flutter characteristics [Reference Rogallo and Yaggy36]. The experiments were conducted in the 40by80 foot NASA Ames wind tunnel facility, with the installed blades presented in Fig. 13(a). The blades were instrumented with strain gauges to determine the torsional and bending stresses. The positioning of the strain gauges, and definition of the blade profiles, shown in Fig. 13(b). These gauges were accompanied by flowfield rakes upstream and downstream of the propeller, with the rakes shown in Fig. 13(a). The use of the NASA Ames wind tunnel facility allowed for the assessment of freestream velocity, the first published dataset to take into account its effect.
In order to capture positive and negative thrust conditions, the blade pitch angle was varied throughout the experiment from $17.5 ^{\circ}$ to $+42^{\circ}$ . The exact range was dependent on the desired condition and also the freestream velocity. The freestream velocity was varied from zero to 183.5 ft/s. This maximum velocity represents 50% of the maximum achieved tip velocity. For the majority of the tests, the rotational velocity was fixed with the blade pitch angle varied until stall flutter was encountered. When examining the effect of the thrust axis inclination (yaw), the positive thrust pitch range was examined along with the inclusion of a freestream velocity. The axis angle was varied by $10^{\circ}$ and $15^{\circ}$ , with only the $15^{\circ}$ yaw at the highest velocity presented.
Presented in Fig. 14 is a selection of the stall flutter results for the Rogallo & Yaggy experimental investigation [Reference Rogallo and Yaggy36]. A selection using different formats is presented to highlight the key observations of Rogallo & Yaggy. For the positive thrust tests, similar trends were found in the torsional stress results with varying freestream velocity. To highlight the sharp increase in stress observed at the flutter boundary, the results in static conditions are shown in Fig. 14(a). With increasing rotational velocity, lower pitch angles are required before entering stall flutter. For rotational velocities between 1,800 and 1,200rpm, a sharp rise in stress was captured at a given pitch with very few additional angles tested before reaching the stain gauge limits. This sharp rise trend transitioned at 1,000rpm, in static conditions, to a more gradual gradient with the lower velocity runs able to reach the full pitch angle range.
This transition velocity, from very high to lower stress gradients, reduced with increasing freestream. The overall effect of increasing freestream is presented in Fig. 14(b) where the stall flutter boundary in terms of pitch and rotational velocity is shown for the positive thrust cases. An increasing freestream velocity pushed the stall flutter boundary further, i.e. a greater pitch and rotational velocity was required before encountering stall flutter.
In terms of the negative thrust study, similar trends and values were found for all tested freestream velocities. Unlike the positive thrust cases, no coupling between the pitch angle and rotational velocity was found at a given freestream, with increases in stress observed at similar pitch angles. With the sharp rise in shear stress occurring around $0^{\circ}$ , the stress values remained fairly stable before slowly decreasing at the lowest pitch angles.
One of the final investigations was the effect of axis yaw on the stall flutter boundary. This effect is shown in Fig. 14(c) where the highest freestream velocity results are presented for a yaw angle of $15^{\circ}$ over a range of rotational velocities. The introduction of a yaw angle had little effect on the shear stress results. To understand, the thrust and torque coefficients were obtained along the blade radius from the rake data. The thrust values indicated only the outer section of the blade (20%) to be stalled with a small magnitude found in the change of angleofattack. With the change in angleofattack, driven from the axis yaw, too small, the detached flow dominated resulting in similar stress results.
In order to understand the trend between generated thrust and the stall flutter boundary, Rogallo & Yaggy extracted the thrust coefficient results along the blade radius and compared them to the values indicated from the shear stress plots when above $\pm 1,000\ psi$ . In Fig. 14(d) are the thrust coefficient results for the $90\%R$ radial station at the highest freestream velocity in positive thrust. It was shown that the stall flutter boundary appears to coalesce with the dropoff in thrust coefficient, with this observation clearest at $90\%R$ compared to inner stations. This gives further indication to the dominance of the separated flow.
A similar investigation was conducted by Dowty Propellers in 1979 at the spinning tower of the Royal Aircraft Establishment, Farnborough [Reference Burton37]. The aim was to determine the torsional stress levels, and subsequently the flutter characteristics, of the Aero Twin Commander threebladed propeller when spun at fixed pitch over a range of rotational velocities. The torsional stress levels were measured via strain gauges placed along the propeller radius.
During the first stage of the test, high levels of torsional stress were seen for a given pitch angle, with the test terminated due to excessive oscillations. Upon examination of the propeller apparatus, the propeller crosshead had failed and hence a redesign was required. Following, overspeed and 30 minute power run tests were completed at lower blade pitch angles, without further attempts to probe the stall flutter boundary and therefore further trends in terms of the stall flutter boundary could not be obtained. From the torsional stress results, a clear stall flutter boundary was observed with a sharp increase in stress seen around 1,650rpm, in a similar manner to the plots of Rogallo & Yaggy (Fig. 14(a)).
Smith conducted a static experimental investigation into three propfan model propellers [Reference Smith38]. These models were designated SR2, SR3 and SR5, with the blades featuring increasing levels of sweep from the unswept, SR2, design to the highly swept, SR5. Additionally, different blade numbers are found across these test cases with the SR2 and SR3 having eight blades and ten for the SR5. No attempt was made to understand the effect of these changes on the flutter response. As in the Dowty tests, for a fixed pitch angle, the propeller rotational velocity was increased to maximum before returning to its baseline value. The stress levels were measured via strain gauges. The positioning of the strain gauges and the blade profiles are presented in Fig. 15. This study was motivated by the transition to thin and highly swept blade designs, conclusions that were drawn from past studies, with the aim to ensure the structural stability of the blades.
Each propeller was installed on an isolated nacelle in the UTRC large subsonic wind tunnel where previous inhouse tests had been conducted [Reference Goepner39]. The SR2 propeller was an 8bladed rotor constructed of steel. The SR3 and SR5 were both made using aluminum with the SR5 having an additional two blades. All configurations were derived from the full scale designs at an intended operating Mach ${\rm M}_a = 0.8$ . A reduction to $1/8^{th}$ scale was required due to the wind tunnel constraints.
The total vibratory stress results were compared against previous experiments of Smith in the NASALewis Research Center [Reference Smith40]. The additional tests of Smith [Reference Smith40] not only looked into static conditions, but also introduced a forward velocity and yaw angle to the freestream to determine its effect. Focusing on the static results, plots of the total vibratory shear stress are presented in Fig. 16 for each blade. Outputs were also presented for the mid and tip blade bending, with shear only presented here due to stall flutter typically being active in the torsional mode.
For each propeller blade, results are presented for a pitch angle of $32^{\circ}$ and a subsequent higher value. This higher value was selected as $40^{\circ}$ , $60^{\circ}$ and $50^{\circ}$ for the SR2, SR3, and SR5 propellers, respectively. Focusing on the comparison between the experiments (Fig. 16(a–c)), it was shown by Smith [Reference Smith38] that good correlation was found between the UTRC and NASALewis investigations. The UTRC results at the critical speeds were found to be higher in all blades with this theorised to be the effect of turbulence.
Figure 16(d) presents the results for the $32^{\circ}$ pitch runs only for each propeller. It was observed that during the lower velocity range, similar trend and gradients were seen between the blades. The results start to deviate around 6,000rpm. At this stage, the SR2 propeller was stopped due to excessive vibration with the maximum speed of 9,000rpm not reached. At this point, the SR2 propeller was seen to have the highest stress values. As the stress levels are lower during this stage for the SR3 and SR5 propellers, the runs continue towards the maximum velocity. Upon reaching the maximum velocity, the SR3 propeller was seen to have a greater stress value, with the SR5 results dropping following 7,000rpm. This overall trend mirrors previous experiments of Smith [Reference Smith40] in that an increase in blade sweep results in a decrease in the stall flutter boundary. However, a trend which was not previously captured was the excitation of a bending mode during a stall flutter response and this was found, particularly, for the SR3 propeller. It was also concluded that the SR5 oscillated in the second bending mode, however, the excitement is found to be limited and therefore any such conclusion may be difficult to correlate. This response is different from expected responses in that it is usually a torsional mode that is triggered and this was found for the SR2. Upon examination of the results, Smith found that the high stress regions were found around a pitch angle of $50^{\circ}$ with a rotational velocity of $6000 $ rpm ( $195.4{\rm m/s}$ tip velocity) with these values way beyond stall. Hence, the response at these angles may be that of a buffeting response.
Along with the static experiments, Smith conducted a similar range of tests on the SR2, SR3 and SR5 propellers in low forwardflight conditions in axial and yawed conditions [Reference Smith40]. A maximum forwardflight ${\rm M}_a = 0.35$ was selected with a yaw range of $2^{\circ}$ to $+15^{\circ}$ . In Fig. 17 is the total stress results for the Mach number and yawed inflow effects on the SR3 propeller. In Fig. 17(a) the total vibratory stress results are presented with Fig. 17(b) showing the total stress per excitation factor (termed stress sensitivity). The excitation factor was used by Smith due to its linear dependence on yaw and was defined as:
where $\psi_i$ is the inflow angle and $V_{eq}$ the equivalent airspeed. Focusing on the inflow effect (Fig. 17(a)) and the results presented are for a pitch angle of $35^{\circ}$ over a range of inflow velocities. A significant decrease in the stress values was seen with the introduction of an inflow velocity. The reduction in stress was mainly seen for the lowest forwardflight speed ( ${\rm M}_a = 0.1$ ), with similar trends observed during the higher velocities. A comparable profile was also found for the SR2 and SR5 blades with these profiles presented in the experimental article. Along with these trends, an increase in stress was observed with an increase in blade pitch. This result pushed the propeller further towards detached flow, and was observed across all tested Mach numbers.
The secondary factor analysed by Smith was the introduction of an angle to the inflow velocity. Furthermore, Smith presented the results for the $15^{\circ}$ yawed inflow for all blades, with the SR3 at ${\rm M}_a = 0.1$ , across a range of blade pitch angles, presented in Fig. 17(b). The introduction of a yawed inflow resulted in a 1/revolution fluctuation in the stress. This was seen as the main component and therefore resulted in a linear increase in stress with yawed inflow angle. A resonancelike effect was found when the 1/revolution frequency equaled the rotational velocity, resulting in a spike in stress. Although not directly related to stall flutter, this highlights the influence of the yaw effect on the blade vibrations. This spike was mitigated with the introduction of greater sweep. As observed from Fig. 17(b), the stress sensitivity increased with increasing blade pitch and rotational velocity.
Due to the amount of data available and the publicly available blade definitions, they remain a viable option for future stall flutter analysis. The structural properties would need to be defined with a clear flutter boundary not present within the test range. Smith also acknowledged the effect of the highlystiff titanium model propellers due to being of single material construction used for this study. The full discrepancy between model and fullscale would need to be examined.
Presented in Table 4 is a summary of the experimental datasets currently available. The table includes some of the pros and cons highlighted within this review along with some of the conditions examined. It is clear from this section, and table, that no ideal test case is currently available and a factor which may have to be revisited, especially with the introduction of the eVTOL market.
3.2 Numerical studies
In addition to the experiments, Smith [Reference Smith38, Reference Smith40] conducted numerical simulations of the SR propellers with a comparison made to the wind tunnel results. Both investigations used the F203 stability analysis code developed by Turnberg [Reference Turnberg41]. Although initially derived for classical flutter, through the aerodynamic update it was used to analyse for stall flutter. This code required the supply of the blade mode shapes and frequencies, along with the aerodynamic loads, to determine an eigenvalue solution.
The mode shapes and frequencies were derived via a finite element analysis, with a beam based method used for the SR2 propeller and NASTRAN for the SR3 and SR5 blades. The beam method was limited to long, slender, isotropic blades. Due to the sweep, offset and large chord lengths of the SR3 and SR5 blades, the beam models would fail and, hence, finite element grids were derived in NASTRAN to utilise its nonlinear static analysis method. The grids for these blades used a combination of quadrilateral and triangular elements. These results were then transformed into the required coordinate system. The transformation splits the displacements into radial sections with two translational and one rotational degree of freedom. This was done in order to linearise the aerodynamics during the stability analysis. For the unsteady aerodynamics, a Goldsteintype [Reference Goldstein42] performance strip analysis was used to determine the steadystatic sectional lift and moment curves.
Figure 18(a) shows the estimated stall flutter predictions of the static experiments of Smith [Reference Smith38], with the prediction assumed to be at the point where the damping goes negative for the third mode of type torsion, where the only instability was found. As shown, the stall flutter predictions mirror the experiments in that the introduction of greater sweep was seen to increase the stall flutter boundary or reduce torsional stress. At this current condition (7,000rpm), the highly swept SR5 propeller does not suffer from stall flutter and requires the velocity to be increased towards the maximum value to induce the instability.
The forwardflight results (Fig. 18(b)) predicted the same trend as the experiments [Reference Smith40] where the difference of Mach number indicates a change to the freestream velocity with constant speed of sound. With increasing freestream velocity, higher blade pitch angles and velocities were required to induce stall flutter. Although the trend was captured through the numerical method, a poorer correlation was found in terms of the exact boundary. This will be primarily a result of the aerodynamic modelling and, hence, greater effort is required to model the unsteady flow.
Based upon the experimental results of Smith [Reference Smith38], Reddy & Kaza conducted a numerical investigation into the SR2 propeller [Reference Reddy and Kaza43]. They sought to improve on the numerical methods of Smith [Reference Smith38, Reference Smith40], in particular the aerodynamic modelling. As was described in Section 2.3, the instability of stall flutter can be derived from the understanding of dynamic stall and such a study was conducted by Reddy & Kaza in 1987 [Reference Reddy and Kaza44]. This study compared three semiempirical dynamic stall models in terms of their lift and moment hysteresis loop prediction for a sinusoidally oscillating aerofoil.
The three semiempirical models include the Gormont model [Reference Gormont45, Reference Gross and Harris46], the Gangwani model [Reference Gangwani47, Reference Gangwani48] and the ONERA model [Reference Tran and Falchero49]. The Gormont model required only one empirical parameter, with a number of stall parameters required for the Gangwani and ONERA models. A short description is given for each modelling method:
Gormont model : This model uses a similar technique as per the LeishmanBeddoes model [Reference Leishman and Beddoes50] in the form of an effective angleofattack. The lift and moment coefficient are then obtained from static aerofoil data. This correction for the effective angleofattack is a function of the rate and is based on oscillating aerofoil test data. The correlation between the true and effective angles is given by:
where $\Delta \alpha_{DS}$ is an incremental dynamic stall angle given by:
Here U is the resultant velocity, c the aerofoil chord, $\dot{\alpha}$ the rate of angleofattack and $\gamma_g$ the nondimensional empirical function derived from experimental tests. This empirical parameter is determined by:
where, Ar is the nondimensional rate of angleofattack. This parameter is a function of Mach number and thickness/chord. The final parameter $K_1$ is given by:
Gangwani model : This synthesised approach was development by Gangwani [Reference Gangwani47, Reference Gangwani48]. The model utilises semiempirical expressions representing the qualitative physical features of dynamic stall. The model is seen to capture the formation and streamwise travel of the dynamic stall vortex.
In this model, three stages of dynamic stall are defined: the stall onset; the trailing edge vortex; and the reattachment. Each of these stages are related empirically by the static stall angle, $\alpha_s$ , the nondimensional rate of angleofattack, Ar, and the unsteady decay parameter, $\alpha_W$ . This decay parameter was previously defined as part of the findings of Beddoes [Reference Beddoes51] using a compressibility correction for the approximated Wagner function [Reference Fung52] (Equation (6)).
Stall onset: The angle at which the dynamic stall vortex is initiated is given based upon the instance of moment stall (since moment stall occurs before lift) and this is represented as:
where the empirical parameters $\epsilon$ , $C_{Am}$ and $C_{Wm}$ are obtained from curve fitted experimental data.
Trailing edge vortex: Following the stall onset, the moment coefficient increases significantly in magnitude and this increase is a result of the vortex travel across the aerofoil. The point at which the moment coefficient reaches its maximum value is seen as the point at which the vortex reaches the trailing edge of the aerofoil. The time at which the vortex detaches from the trailing edge is given by:
where $C_{At}$ and $C_{\alpha t}$ are the empirical parameters defined using experiments.
Reattachment: The point at which the flow reattaches to the aerofoil depends on the freestream Mach number. For ${\rm M}_a < 0.4$ , the reattachment occurs at an angle which is less than the static stall angle. When ${\rm M}_a > 0.4$ , the reattachment angle occurs at an angle greater than the static stall. Using the empirical parameters $\epsilon$ , $C_{AR}$ and $C_{wR}$ , this reattachment angle ( $\alpha_{re}$ ) is defined as:
Using these equations for the dynamic stall events, the lift and moment coefficients can be derived alongside additional empirical parameters.
ONERA model : The ONERA model defines the lift and moment coefficients as a set of ordinary differential equations [Reference Tran and Falchero49]. When attached, the equations are first order with this increasing to third order for the detached flow regime. A total of 11 empirical parameters are required for the differential equations (6 for lift and 5 for moment) with both lift and moment containing a linear ( $C_1$ ) and nonlinear ( $C_2$ ) expression resulting in the following simple equations:
where each linear and nonlinear expression is given as follows:
The total aerodynamic angleofattack ( $\alpha$ ) is a summation of the angles due to the pitching ( $\theta$ ) and plunging ( $\dot{h}/b$ ) motions. $C_{L\ell}$ and $C_{M\ell}$ represent the static lift and moment coefficients in the linear region, with $\Delta C_L$ and $\Delta C_M$ the difference between the extended linear curve and the true static curve. The operators $\dot{()}$ and $\ddot{()}$ represent derivatives with respect to the nondimensional time, $\tau = U t / b$ . This leaves the empirical parameters $\lambda_{on}$ , $s_{on}$ , $a_{on}$ , $\sigma_{on}$ , $r_{on}$ and $E_{on}$ which are derived from wind tunnel data.
In addition to the pitching aerofoil study of Reddy & Kaza [Reference Reddy and Kaza44], this model has been used to obtain the dynamic response of a typical helicopter blade section [Reference Rogers53] and the entire helicopter blade [Reference Peters54]. For the dynamic stall study [Reference Reddy and Kaza44], all models were used across a range of test cases for the NACA 0012 aerofoil. For the lift coefficient, all models correlate well across the range of test conditions to the experimental data. The moment coefficient was not presented for the NACA 0012 test cases using the ONERA model, however, it was seen to predict the moment coefficient well for a reference OA212 aerofoil section. The moment coefficient was not predicted well for the Gormont model, with very good agreement for the Gangwani. The results from the Gangwani model were expected due to the amount of empirical data taken from the experimental results.
Reddy & Kaza [Reference Reddy and Kaza44] concluded from this initial dynamic stall investigation that the Gormont and ONERA models can be used for propeller flows in typical operating conditions due to the blade being primarily within the light stall regime. The ONERA model involves fewer experimental parameters than the Gangwani, with the Gormont only requiring one. This could potentially be a factor in the choosing of an aerodynamic model as it depends on the availability of static aerofoil data of the blade sections. This situation was found in the propeller stall flutter study [Reference Reddy and Kaza43] and, hence, the Gormont model was used due to the lack of blade sectional data. This model was used to provide the stall effects with quasisteady, incompressible twodimensional strip theory used to calculate the generalised aerodynamic forces across a range of sections along the blade.
The structural modelling used by Reddy & Kaza [Reference Reddy and Kaza43] involved the derivation of the staticdisplacement and normal modes and frequencies from a nonlinear finite element analysis which was then supplied to the standard modal equations of motion for a propeller blade [Reference Kaza, Mehmed, Narayanan and Murthy55]. The final equations were derived as:
where $m_i$ was the generalised mass, $\omega_i$ the natural frequency, $\zeta_i$ the critical damping, $Q_i$ the generalised force terms and $q_i$ the normal coordinate for the $i^{th}$ mode where n modes were supplied. These equations were solved within the timedomain using the Wilson $\theta$ method [Reference Wilson and Bathe56]. This method is an implicit integration method that assumes a linear variation in acceleration between timesteps. The static and modal displacements were derived using NASTRAN, with four normal modes supplied. The aeroelastic flowchart used by Reddy & Kaza is presented in Fig. 19.
The experimental results of Smith [Reference Smith38] presented the stall flutter boundaries in terms of stress across the blade. Using the Gormont dynamic stall model, qualitative correlation was found between the experiment and numerical study, however, the quantitative comparison of the damping plots showed a conservative estimation of the boundary using this model. This conservative estimation may be acceptable from a design point of view if the boundary is outside the normal operating range of the propeller.
Due to the fact that only a qualitative comparison could be made to the experiments, a numerical comparison of the modal amplitude results was conducted against the study of Bielawa [Reference Bielawa, Johnson, Chi and Gangwani57]. The investigation of Bielawa used the Gangwani dynamic stall model with a nonlinear beam model for the structural modelling. The use of a nonlinear beam model for the unswept SR2 is acceptable due to the design of the blade. However, for a complex structure, such as the swept SR5, the general finite element analysis of Reddy & Kaza [Reference Reddy and Kaza43] is more suitable. The modal amplitude predictions of both the simulations agreed well. A slight improvement on the dampingboundary plots was found in the Bielawa study whilst using the Gangwani model with the predictions shown in Fig. 20. This slight improvement is shown via the closer approximation of the zero damping pitch angle to the single experimental point.
In 2007, DelamoreSutcliffe [Reference DelamoreSutcliffe58] conducted a numerical investigation into propeller stall flutter with a comparison made to the experiments of Baker [Reference Baker9]. The investigation involved the development of an unsteady aerodynamic model, coupled with the Brookes and Houbolt equations of motion for a rotating beam, to determine the stall flutter boundaries. The main objective of this work was to transition from the empirical models still used by industry to a more advanced formulation without penalising the fast prediction characteristics of such methods.
DelamoreSutcliffe builtup the aerodynamic model from the linear to nonlinear regime in twodimensions, with experimental tests conducted to ensure a validated model was derived. Three fundamental models were tested for the attached flow (linear) regime, with this involving two established models and a third derived in the timedomain. The experiments used to validate the attached flow model were those of Piziali in 1994 [Reference Piziali59]. Details such as the pressure distribution over a semispan wing and oscillating aerofoil were presented, thus ensuring an extensive range of data from which to compare.
The first model was defined as the indicial response model and considers the instantaneous changes between two steady states. A review of this method was conducted by Beddoes [Reference Beddoes60]. With an indicial response function defined for the circulatory and noncirculatory terms as a function of semiempirical coefficients, the total force coefficient is seen as a sum of the subsequent parts. These semiempirical coefficients have been discussed by Leishman in which a method for obtaining [Reference Leishman61] and validating [Reference Leishman62] was presented. Using a finite difference approximation to Duhamel’s integral [Reference Bisplinghoff, Ashley and Halfman63], it was possible to determine the response to an arbitrary timehistory. A full derivation was presented by Leishman [Reference Leishman64].
The second model used was the Theodorsen model [Reference Theodorsen65, Reference Theodorsen and Garrick66] where a pitching and plunging aerofoil is modelled in static conditions. Taking the fixedwing expressions for lift (L) and moment (M) on an aerofoil with unit span:
where C(k) is the Theodorsen function. This function defines the unsteady aerodynamic component. It is complex and typically defined as:
where $H_n^{(2)}$ are Hankel functions of the second kind, which are Bessel functions of the first and second kinds. With assumptions made for the rotorcraft equivalent, taking the differential across the blade radius allows the equivalent lift and moment equations are to be derived (Equation (15)). Although such a theory would appear restrictive to the Theodorsen function, by substituting different liftdeficiency functions the unsteady airloads for other assumed conditions can be expressed.
Due to the fact that the Theodorsen model was constructed within the frequency domain, an ideal coupling to a timemarching aeroelastic model remains tricky. This leads on to the third model, which was a timedomain approximation of the Theodorsen model derived by Dowell [Reference Dowell67]. Here, it was found that a link between the frequency and timedomain can be obtained via a Fourier transform pair. This results in an updated Theodorsen function with coefficients $\lambda_1$ and $\lambda_2$ , which can be substituted into the full lift and moment equations (Equation (13)) before converting to the final timedomain differentials by using the inverse of the differential rule for Laplace transforms.
All attached flow models were tested and provided satisfactory correlation with the experimental results [Reference Piziali59]. The third model was chosen as the basis for the full aerodynamic model due to the lack of empirical data and formulation within the timedomain.
Following the derivation of the linear aerodynamic model, an extension was required to include the dynamic stall effects. Using a statespace approach, a nondimensional parameter was introduced in order to provide the detached flow mechanism. This parameter was the static upper surface separation point ( $x_s$ ) and this was combined with the Kirchoff and Helmholtz theory [Reference Woods68] to determine equations for the normal force and moment coefficients. The separation point was defined as $x_s \in [0,1]$ , where $x_s = 0$ corresponds to detached flow at the aerofoil leading edge and $x_s = 1$ to a fully attached flow. The normal force coefficient ( $C_N$ ) was given by:
where $c_{l_{\alpha}}$ is the lift curve slope and $\alpha$ the aerofoil angleofattack. For the moment coefficient, a reference to the work of Khrabrov and Ol was made as they produced a theory which combined thin aerofoil theory and Kirchoff concepts to determine a closed solution for the pitching moment coefficient at high angles of attack [Reference Khrabrov and Ol69]. However, to overcome the small angle assumption, the work of Leishman and Beddoes is cited [Reference Leishman and Beddoes50] for the ratio equation between the normal force and moment coefficients ( $C_M$ ), due to the fact it used the separation point. The ratio was given as:
where $k_0 = (0.25x_{ac})$ is the nondimensional aerodynamic centre offset from the quarter chord. The remaining two parameters, $k_1$ and $k_2$ , are semiempirical coefficients derived from experiments. In order to estimate the separation point, further semiempirical coefficients were required depending on the current state of the flow and the static stall angle. Leishman and Beddoes stated that the static stall angle ( $\alpha_s$ ) corresponds closely to the point where $x_s=0.7$ . This fact was used alongside experimental data to determine the static stall angle. The separation point was then given by:
where $s_1$ and $s_2$ were semiempirical coefficients chosen to define the shape of the state variable curve pre and poststall. In order to model the unsteady effects of detached flow, an effective angleofattack was introduced:
where $\tau_2$ defines the time delay associated with quasisteady effects. The model was completed by modelling the transient aerodynamics associated with the relaxation process following detached flow. These effects are determined via a first order differential equation for the effective separation point $x_{s_E}$ :
where $\tau_1$ was the semiempirical relaxation time constant.
This dynamic stall model was validated against the experiments of Jumper et al. [Reference Jumper, Schreck and Dimmick70]. Jumper et al. conducted wind tunnel tests of the pitching NACA 0015 aerofoil over a range of constant pitch rates. The aerodynamic model was tested over a range of the experimental rates with the semiempirical parameters $s_1$ , $s_2$ and $\alpha_s$ determined from the data. Initial estimates of the two dynamic coefficients ( $\tau_1$ and $\tau_2$ ) were sourced from the studies of Goman and Khrabrov [Reference Goman and Khrabrov71] and Khrabrov and Ol [Reference Khrabrov and Ol69]. These estimates were then adjusted to provide improved correlation. Excellent correlation was found across all pitch rates.
In order to combine the attached and detached flow effects, a similar summation as the linear model was used. This involved the addition of the zeroangle, circulatory and noncirculatory lift coefficients, with the high angleofattack effects included via an unsteady factor (LF).
This unsteady factor was multiplied by the circulatory lift and was derived in the Kirchoff/Helmholtz theory [Reference Woods68] as:
A similar process was found for the pitching moment term. An additional unsteady pitching moment factor (PF) was introduced and this was multiplied by the unsteady lift factor and circulatory lift. This term was then included in the summation of the zeroangle, circulatory and noncirculatory moment coefficients.
The full aerodynamic model was verified against the pitching aerofoil experiments of Piziali [Reference Piziali59], with prerequisite testing used to determine the semiempirical coefficients. Good correlation was found in static and dynamic tests, with some small deficiencies at high static angles of attack and the dynamic overshoot of the lift and pitching moments.
In order to verify the aerodynamic model for the stall flutter phenomena, DelamoreSutcliffe conducted a range of wind tunnel experiments of a twodimensional aerofoil with a single degree of freedom in pitch. This verification was conducted via a timemarching and eigenvalue stability analysis through the coupling of the aerodynamic and structural models. The structural model was derived using Lagrange’s equations of motion [Reference BismarckNasr72]. Due to the motion of the aerofoil, the potential remains for the aerofoil to stall at a negative angleofattack. As a result, the equations for the separation point needed updating to include the negative range with respect to the negative staticstall angle ( $\alpha_{s_n}$ ):
where semiempirical parameters $s_{1_n}$ and $s_{2_n}$ were updated for the negative range. The equations were coupled in a second order matrix equation including the mass, damping and stiffness matrices. A RungeKutta scheme was used for the timemarching formulation and a linearised solution in the form of an exponential was assumed in order to obtain eigenvalues for the stability analysis. Further details were presented by DelamoreSutcliffe et al. [Reference DelamoreSutcliffe, Whiting and Greenwell73]. A comparison of the experimental and numerical results showed good correlation, with the numerical results from both the stability and timemarching methods slightly below the experiment.
In preparation for the transition to the propeller test cases, the aeroelastic model developed for the twodimensional aerofoil was extended to look at a threedimensional cantilever wing. The wing was structurally modelled using strip theory and assumed modes with two degrees of freedom in pitch and heave. A key benefit to using this method was the limitation of twodimensional aerofoil characteristics and the computational solution time. Disadvantages included the lack of threedimensional crossflow modelling and the requirement for empirical data, a common theme in rotor aerodynamic modelling.
As per the structural model, the aerodynamic model was updated to compute the spanwise lift and moment coefficients in both attached and detached conditions. With the use of strip theory, the aerodynamic model was applied at each spanwise location with a stationbystation update required for the semiempirical coefficients and lift curve slope. Experiments [Reference Piziali59] were again used for the coefficients, along with curve fitting, to determine the spanwise distribution of the values. Theoretical expressions used to convert the twodimensional slope to threedimensions with a quintic power distribution over the span was assumed for the lift curve slope with a value of zero applied at the tip. Furthermore, the effect of cantilever wing twist was taken into account via liftingline theory and Prandtl’s integral equation for circulation [Reference Houghton, Carpenter, Collicott and Valentine74].
These expressions for aerodynamics and structure were combined in a similar manner to the twodimensional aerofoil equations. The structural equations were generalised with respect to the assumed modes and then integrated over the span. Two bending and a single torsional mode were included, with higher modes having a minimal effect on the accuracy. A similar process was conducted for the aerodynamic terms including the circulatory and noncirculatory coefficients. The solution to these equations were obtained via a timemarching and eigenvalue analysis. Due to the increase in matrix terms, the timemarching calculations were limited due to the increase in computational cost. These were subsequently compared to the eigenvalue analysis which was conducted across the entire test case range.
This threedimensional wing method was compared against the experiment of Dunn and Dugundji [Reference Dunn and Dugundji19]. Good correlation between the experimental and numerical calculations were found. For the first test case examined, the timemarching method appeared to produce an improved correlation of the flutter frequency in comparison to the stability analysis. However, the flutter velocity was underpredicted at high initial pitch angles. Only the stability analysis was used for the second test case with an underprediction of the flutter velocity at low and high initial angles, and underprediction of the frequency across the full range.
Following the verification of the threedimensional wing aeroelastic model, the aerodynamic model was extended and used with the Brookes and Houbolt equations of motion for a rotating beam [Reference Houbolt and Brookes75] to conduct a propeller stall flutter investigation. For the structural model derived by Brookes and Houbolt, several assumptions were made to reduce the complexity of the equations without losing accuracy. These were:

(1) The chordwise bending was assumed to be very small, and hence negligible, in comparison to the flapwise bending. This was due to the fact that the chordwise bending stiffness is much greater than the flapwise stiffness when the thickness/chord ratio is less than 15%. This is acceptable with modern propellers due to their need to utilise thin aerofoil sections for performance. This assumption may fail for tiltrotor blades where thicker sections are used.

(2) The blade was assumed to have no set twist distribution in the spanwise direction. This was applicable within this study due to the fact that the experimental blades have zero twist. For more complex, highly twisted blades, this assumption would fail and the terms associated would need to be included.

(3) The final assumption was that the tensile and rotational axis of the blade align with the elastic axis. This was, again, driven by the experimental data were no information was supplied regarding these offsets from the elastic/mass axis.
These assumptions resulted in the following structural equations for the propeller blade:
Where the axial tension was defined as $Ta = \int_y^L \Omega^2 m (y + R_{hub}) dy$ and e the distance between the elastic and mass axes was given by $e = (x_{ea}x_{cg})c$ . An important point to note is the reduction of the above equations from propeller to wing configurations is possible via the setting of the rotational velocity ( $\Omega$ ) and hub radius ( $R_{hub}$ ) to zero. This provides continuity in the aeroelastic modelling.
In terms of the aerodynamic modelling, the strip theory used in the cantilever wing was updated to include the spanwise variation in induced velocity. One of the key physics of rotorcraft aerodynamic modelling is the influence of the wake on the preceding blade. This effect can be captured by methods such as Loewy [Reference Loewy76] who took into account this effect via the modification to the liftdeficiency function from Theodorsen’s theory. Although such a method aims to capture a greater amount of the physics, it was assumed that such an effect is minimal in comparison to the overall load results obtained via Theodorsen’s method. This is true in high thrust/inflow condition were the propeller wake is induced far away from the disc.
In order to update the aerodynamic model for the rotorcraft environment, the effect of the induced incidence must be taken into account. The induced incidence is the effect of the blade rotation and inflow velocity, and was modelled via Blade Element Momentum Theory (BEMT) [Reference McCormick77]. This model was used for each section with an iterative process conducted to obtain the induced velocity effect.
Both the structural and aerodynamic terms were combined as the previous versions, in statespace form of stiffness, damping and mass matrices. In the cantilever wing model, the structural damping was set to zero. This was not the case for the propeller model due to the fact it was taken into account in the experimental studies [Reference Baker9] and, hence, must be included within the equations. To include structural damping, the theory of Theodorsen [Reference Theodorsen and Garrick66] was applied where the structural damping was seen as a function of the generalised stiffness matrix, i.e.:
where, g is the structural damping coefficient. In order to solve the equations of motion, the stability analysis conducted for the aerofoil and cantilever wing was used. Due to the increased computational cost associated with the timemarching method, it was subsequently not used for this method. Further details regarding the aeroelastic modelling was presented by DelamoreSutcliffe and Greenwell [Reference DelamoreSutcliffe and Greenwell78].
A comparison to the experimental data of Baker [Reference Baker9] was made with three test cases selected and presented in Fig. 21. A good correlation between the experiment and numerical analysis was found. The trend of lowering velocity with increased pitch angle was seen in both results with the minimum flutter velocity well captured. Differences at higher angles of attack were seen within the comparison and these were determined to be an effect of the high angleofattack detached flow modelling. However, the pitch angles where this was seen are outside the blades operational regime.
Ognev in 2011 investigated different unsteady aerodynamic models and their influence on the determined flutter boundary [Reference Ognev and Rosen79]. These unsteady aerodynamic models included Theodorsen’s model [Reference Theodorsen80], three cascade models [Reference Whitehead81, Reference Smith82, Reference Jones and Moore83] and a detailed threedimensional model [Reference Iosilevskii, Iosilevskii and Rosen84]. The cascade models took into account the effect of the returning wake (something which is neglected in the Theodorsen model). However, they were limited by their twodimensionality with the effect only applied from equivalent cross sections. This was mitigated by the use of the threedimensional model by Isoilevskii et al. [Reference Iosilevskii, Iosilevskii and Rosen84] where bound, trailing and shed vortices were included. For all models, the effect of viscosity is neglected and this resulted in the loss of detached flow modelling, not ideal for a stall flutter analysis, hence the comparison between the simulations and experiments [Reference Hubbard, Burges and Sylvester33] were only fair within the classical bendingtorsion flutter regime.
In 2019, Higgins et al. conducted a numerical investigation into the Aero Twin Commander propeller blade [Reference Higgins, JimenezGarcia, Barakos and Bown85]. This investigation used the HMB3 CFD solver and compared the results of URANS and SAS [Reference Menter and Egorov21] calculations against experimental data. The use of the SAS method for aerodynamic modelling had proven beneficial in a typical twodimensional form via the study of pitching aerofoils [Reference Higgins, Barakos and Jinks20]. Such a method allows for a balance between the resolved and statistical content to be found, and improves the accuracy without significantly increasing the computational cost.
For the structural modelling, a loosely coupled timemarching modal method was used. This method, again, allows for a balance to be made between modelling accuracy and the heavy computational cost associated with coupled methods. For the modal method, the mode shapes and frequencies of the Commander propeller blade were required, and for this, a NASTRAN nonlinear beam model was built. A nonlinear beam model was selected due to its ability in previous helicopter rotor aeroelastic studies [Reference Dehaeze, Barakos, Garipova, Kusyumov and Mikhailov86], with good validation found between the experimental and numerical mode shapes and frequencies for the propeller blade.
As in the experimental investigation [Reference Burton37], the blade was accelerated from 1,400 to 1,750rpm at a fixed pitch angle of $27^{\circ} 40^{\prime}$ . Good correlation was found in terms of the trend of torsional stress across the blade radius between test and simulations (Fig. 22(b)). The SAS results were found to be slightly higher than the URANS, indicating the need for scaleresolving. This was confirmed with the modal amplitude results where a qualitative agreement to the experiment was found in terms of the third mode amplitudes. As the blade accelerates, a significant increase in amplitude was found which mirrors the trend in stress from the test results (Fig. 22(a)).
This investigation highlighted the ability and requirement for scaleresolving CFD methods when studying stall flutter. Advances in computing power has enabled the use of such methods, thus allowing for closer validation which can lead to greater advances to be made in design. Expanding the boundary which the blade can operate.
Presented in Table 5 is a summary of the numerical investigations conducted for propeller stall flutter. A list of the techniques used, the range of conditions and the date is presented. From the table, and this section, it can be found that a limited number of aeroelastic numerical investigations have been conducted. This is likely down to the small range of potential trigger points during a typical flight path. Thus, allowing for fast, eigenvalue type studies to be dominant within the field as the ultimate aim of avoidance could be achieved from conservative boundaries. However, again with the introduction of the eVTOL market, propeller stall flutter is becoming a more active area of research due to shift in propeller flight path and the redesign of the blade for greater efficiency.
4.0 Concluding remarks
4.1 Discussion
A range of experimental investigations have been carried out over the years. Due to the heavy use of propeller driven aircraft between 1945 and 1960, a high number of experiments were conducted to understand blade flutter. As these studies were not specifically focused on stall flutter, merely propeller flutter in general, correlation with more recent studies is difficult as only a few fundamental parameters were examined. This included rotational velocity, pitch and number of blades. However, these studies by Sterne [Reference Sterne29] and Theodorsen & Regier [Reference Theodorsen and Regier30] provided some of the initial trends and influences stall flutter can have on propellers, with this being the significantly reduced boundary as a result of a propeller entering stall flutter.
One of the additional factors studied by Sterne [Reference Sterne29] was the effect of number of blades. A different number of blades was studied by Theodorsen & Regier, but due to the fact that one of the two propellers was singlebladed, a trend with blades numbers would be invalid due to lack of influence on the single blade from a proceeding blade wake. However, the study by Ewing et al. [Reference Ewing, Kettlewell and Gaukroger32] did examine different blade numbers. Ewing et al. found that the effect of number of blades had an influence on the vibratory phase between blades and, as a result, any comparison in terms of stress levels between propellers of different blade numbers would also have to take into account the phase. Such a comparison remains difficult due to the fact that a peak stress level for one propeller may aline with another which induces a phase change and thus the stress levels are lower. This is potentially one of the influences that was found in the study by Sterne in that a different phasing was introduced by the 3bladed propeller that resulted in higher stress levels, and thus lower stall flutter boundaries, than the 2 and 4bladed alternatives.
In the 1950s, a range of model blades were designed for wind tunnel investigations [Reference Baker9, Reference Ewing, Kettlewell and Gaukroger32, Reference Hubbard, Burges and Sylvester33, Reference Allis and Swihart34, Reference Rogallo and Yaggy36]. Such blades were made of a modern alloy to allow for a reduction in the variability of the manufacturing process. A range of parameters could be investigated and the greatest insight into propeller stall flutter was found during this period. One of the parameter investigated during this period was the effect of blade twist, with this parameter being the only factor examined across the different studies. This was primarily conducted by Baker [Reference Baker9] and Hubbard et al. [Reference Hubbard, Burges and Sylvester33]. However, during the tests, only two data points were examined. In the case of Baker [Reference Baker9], an untwisted and $17^{\circ}$ twisted blade were compared with very little difference found in terms of minimum flutter speed. For the investigation by Hubbard et al. [Reference Hubbard, Burges and Sylvester33], again an untwisted and a $31^{\circ}$ twisted blade were tested but the untwisted blades were conducted at different conditions. Hence, a true sidebyside comparison could not be conducted.
With the transition from propellerdriven aircraft to turbofan engines, interest was lost in propeller stall flutter; therefore, very few studies were conducted after the 1950s. A privatesector test was conducted by Dowty Propellers in 1979 using an inservice blade [Reference Burton37]. Although the primary aim of the test was not specifically for stall flutter understanding, stall flutter was encountered due to the pitch angle range. Furthermore, interest within propfans grew in the 1980s and this resulted in new studies to understand the blade flutter boundaries [Reference Smith38, Reference Smith40]. These investigations into the propfans by Smith looked into a range of factors previously examined. This included blade sweep previously studied by Baker [Reference Baker9] and freestream velocity and yaw which was previously examined by Rogallo & Yaggy [Reference Rogallo and Yaggy36].
Focusing on the sweep trends produced by Smith [Reference Smith38, Reference Smith40] and Baker [Reference Baker9], Baker initially found the positive benefit of sweep. An unswept blade was compared to three swept blades, two of which had a linear profile and the other a nonlinear. Although the nonlinear profile had the greatest impact across the stall pitch region, both the linear swept blades produced a benefit in terms of the flutter velocity. The comparisons of sweep by Smith looked at the unswept SR2 against the nonlinear swept SR3 and SR5. Both the SR3 and SR5 reduced the levels of stress acting on the blade, confirming the trend found by Baker. One of the benefits to the introduction of sweep is the delay in separation onset with minimal performance losses. This correlates with the stall flutter results as increasing sweep is seen to delay the onset of high stress.
The majority of the experimental studies have been conducted in static/hover conditions. This is mainly due to the fact that propeller stall flutter is likely to occur during low inflow conditions, such as takeoff. However, the experiments of Rogallo & Yaggy [Reference Rogallo and Yaggy36] and Smith [Reference Smith40] looked into the effect of added freestream velocity and found that the overall stress acting on the blade was significantly reduced and thus an improvement in the stall flutter boundary was achieved. The improvement being the requirement for higher rotational velocity and pitch angles before stall flutter is induced compared to the hover trigger values. Such a response is expected when the effective angleofattack is taken into account. Although the geometric pitch remains the same, the additional inflow vector reduces the overall effective angle and thus, stall is delayed. Such a trend is not only observed in propellers but also wings and aerofoils.
The final crossover factor examined experimentally was the effect of yaw and this was again studied by Rogallo & Yaggy [Reference Rogallo and Yaggy36] and Smith [Reference Smith40]. In order to examine for the effect of yaw, a freestream velocity is required. For the Rogallo & Yaggy test, only the highest freestream condition was reported and therefore, any yaw influence on the stall flutter boundary is minimised due to the dominating reduction in boundary of the freestream. This is essentially confirmed due to the fact that almost no change in stress was found as a result of yaw. Smith did present results for the yawed cases across a range of freestream velocities, however, redefined the stress values presented by taking into a yaw factor. As a result, no direct comparison to a nonyawed case could be made with only observations presented by Smith available. At the lowest freestream velocities, Smith presented a linear trend of total maximum stress with yaw angle and this is induced due to the oneperrevolution fluctuation in loads due to the variation of velocity across the blade disc. This additional fluctuation could have a significant impact on the stall flutter boundary, however, the dominating effect of the freestream must be minimised in order to gain a true picture of its effects.
With the minimal development in aeroelastic models during the 1950s, interest in numerical stall flutter investigations only began with the development/study of the SR propfans. These studies, conducted in the 1980s, involved both eigenvalue and timemarching modal methods with semiempirical aerodynamic models. The use of the eigenvalue analysis was selected to ensure a fast calculation and this trend continued into the modern studies. With the availability of faster computing power, more complex finite element models could be used and a greater amount of timemarching analysis was conducted. However, the focus still remained on the use of the semiempirical dynamic stall models and therefore conservative boundaries were achieved.
4.2 Conclusion
From the initial studies of Sterne [Reference Sterne29] and Theodorsen & Regier [Reference Theodorsen and Regier30], wooden propellers are used. One of the benefits of using modern alloys is the ability to refine the manufacturing process to reduce any inconsistencies in structural properties. Modal analysis can indicate frequency differences in blades of the same design. Therefore, to ensure no additional inconsistencies are applied on top of numerical errors, metal alloy propellers are preferable, making wooden blades nonideal.
It is clear from the studies in the 1950s [Reference Baker9, Reference Ewing, Kettlewell and Gaukroger32, Reference Hubbard, Burges and Sylvester33, Reference Allis and Swihart34, Reference Rogallo and Yaggy36] that a great deal of data and insight was gleaned. However, due to the purposebuilt nature for each of the studies, the blades were not fully defined within the reports. As a result, such investigations remain nonideal for numerical investigation. Additionally, no significant crossover of parameters was conducted and therefore any trends, regardless of propeller type, remain sought after.
For each experimental investigation described, the structural properties are unknown and unreported. These properties, therefore, have to be estimated in order to conduct numerical investigations.
From the propeller stall flutter review it is concluded that there is a lack of definitive propeller stall flutter data. The experimental investigations are historic and often lacking in clear boundary assessment. This minimises the range of validation that can be carried out with modern numerical methods.
4.3 Recommendations
For any future study, the influence of the freestream velocity must remain minimal with ideal test conditions involve zero to low speed conditions. The exact boundary of the low speed condition is dependent upon the specific design. However, from the studies of Rogallo & Yaggy [Reference Rogallo and Yaggy36] and Smith [Reference Smith40], even the lowest freestream increment had an effect during low rotational velocity conditions. In the Rogallo & Yaggy study, the lowest ratio of freestream to tip speed was $\sim5\%$ and this $\sim5\%$ introduction of freestream found an increase in pitch of $5^{\circ}$ , or 22%, with respect to the hover boundary. This highlights the significant influence freestream velocity has and hence the recommendation.
In order to improve the numerical techniques, and additionally obtain a greater insight into the phenomena via nonintrusive tools, a new set of experiments is required to extract greater datasets. Not only with the current datasets being minimal and historic, it is further emphasised that modern blades for novel aircraft configurations are likely to be considerably different from wooden propellers of 80 years ago, and operate in completely different regimes (from eVTOL to tiltrotor aircraft), and are made of completely new materials, such as composites, which have different aeroelastic response. For these new systems, there are virtually no data.
A range of numerical methods are available from which to study propeller stall flutter, with each model’s application dependent upon the overall goals. Should a significant range of calculations be required, an eigenvalue analysis with semiempirical aerodynamic model would be of use as this would allow for a full conservative boundary to be formed. Should greater detail about a specific condition be required, a timemarching method with highfidelity CFD aerodynamics could be used as this would allow for an assessment of the structural and aerodynamic response in time. With the latter more computationally expensive, improvements in computing power are continually being made and therefore the likelihood of using a timemarching method is increasing.
Acknowledgments
The support provided by The Engineering and Physical Sciences Research Council (EPSRC), the U.K. Vertical Lift Network (UKVLN), the U.K. Chapter of the Vertical Flight Society and Dowty Propellers is gratefully acknowledged.
The support of EPSRC via the MENtOR project (EP/S013814/1 for Glasgow and EP/S010092/1 for Manchester) is gratefully acknowledged.