Tools › Aerodynamics › Stability and Control › Flutter and Divergence › Method
Flutter and Divergence: Method
Open the Flutter and Divergence tool
Flutter and Divergence: the Typical Section, Theodorsen's Function, the U-g and p-k Solutions, Divergence, Reversal and the Fin Criterion
What the page is for
Flutter is a self-excited oscillation in which the air forces on a moving wing feed its bending and torsion; divergence is the static loss of torsional stiffness when the lift twists the wing faster than the structure resists; aileron reversal is the speed at which the twist produced by an aileron cancels the roll it was meant to give. The page gives first estimates of the three from the public reports of the National Advisory Committee for Aeronautics. It is a study and preliminary-design aid for light aircraft, models and fins: a real aircraft is cleared by a full analysis with measured modes, a ground vibration test and a flight flutter test.
The typical section
The wing is replaced by one section of semichord b (Theodorsen, NACA Report 496). Positions are measured in semichords from midchord, positive aft: the leading edge is at −1, the quarter chord at −½, the elastic axis at a, and the center of gravity at a + xα. The section has mass m per unit span and radius of gyration rα b about the elastic axis, and is held by a bending spring and a torsion spring that give the uncoupled frequencies ωh and ωα. The mass ratio μ = m / (π ρ b²) compares the wing with the cylinder of air standing on its chord; Theodorsen's κ is 1 / μ. With the plunge h positive down and the pitch α positive nose up, the equations of motion are
m ḧ + m xα b α̈ + Kh h = −L, m xα b ḧ + Iα α̈ + Kα α = M,
and for harmonic motion in incompressible flow the lift and the moment about the elastic axis are
L = π ρ b² [ḧ + V α̇ − b a α̈] + 2 π ρ V b C(k) [ḣ + V α + b (½ − a) α̇],
M = π ρ b² [b a ḧ − V b (½ − a) α̇ − b² (⅛ + a²) α̈] + 2 π ρ V b² (a + ½) C(k) [ḣ + V α + b (½ − a) α̇].
The first bracket of each is the reaction of the air that must move with the section (apparent mass); the second is the circulatory lift acting at the quarter chord, proportional to the downwash at the three-quarter chord and reduced and delayed by the wake through the circulation function.
Theodorsen's function
C(k) = F + i G = H1(2)(k) / [H1(2)(k) + i H0(2)(k)], with the reduced frequency k = ω b / V and the Hankel functions Hn(2) = Jn − i Yn. In real form (Report 496, above its Table II), F = [J1 (J1 + Y0) + Y1 (Y1 − J0)] / [(J1 + Y0)² + (Y1 − J0)²] and −G = [Y1 Y0 + J1 J0] / [the same denominator]. C is 1 at k = 0 (steady flow) and falls to ½ at high frequency. The page computes the Bessel functions from their ascending series below an argument of 14 and from the Hankel asymptotic expansions above. Three cells of the report's tables do not follow its own formulas and are not used: Y1(0.1) is printed −7.0317 for −6.4590 (which also spoils F and G of that row), F(0.5) is printed 0.6030 for 0.5979, and T4 at c = −½ in Table I is printed −1.6614 for −2.5274.
Flutter by the U-g method
Assume harmonic motion at a chosen k. Then V = ω b / k, every air force is proportional to ω², and the equations become K (1 + i g) q = ω² [M + A(k)] q, in which g is an artificial structural damping added to both springs. For each k the 2 × 2 eigenvalue λ = (1 + i g) / ω² gives a frequency, a speed and the damping that would be needed to keep the motion harmonic. Plotted against speed, a branch with g below zero is stable (the wing has more damping than it needs) and flutter is where g rises through zero, or through the structural damping actually present. The page scans k from 6 down to 0.02, follows the two branches by continuity, and refines each crossing by bisection. Speeds are reported as V / (b ωα), the flutter speed coefficient of the reports, and in physical units with the semichord and torsion frequency entered.
Flutter by the p-k method
The U-g curves are only true at g = 0. The p-k method (Hassig) looks for the actual roots s = (γ + i) ω of det[s² M + s B(k) + K(k)] = 0 at each speed, taking the aerodynamic matrices at the reduced frequency of the root itself, k = Im(s) b / V, and repeating until the two agree. The page starts at low speed from the roots of the section in still air (with its apparent mass), steps up in speed following each root, and plots 2 Re(s) / Im(s), which is comparable with g. The flutter speed is where a real part turns positive. Because both methods use the same harmonic aerodynamics at the neutral point, their flutter speeds must agree; the page prints the difference.
Divergence, reversal and the approximations
Divergence of the section follows from the static balance of the lift 2 π ρ V² b α at the quarter chord against the torsion spring: VD / (b ωα) = √[μ rα² / (1 + 2 a)] (NACA Report 685, page 17). Report 496 prints this relation with ½ + a in place of √(1 + 2 a); NASA/TP-2015-218765 lists that as an error. For aileron reversal the aileron adds lift CLδ δ and a nose-down moment Cmδ δ; the section twists until the net lift is zero when VR / (b ωα) = √[μ rα² CLδ / (4 (−Cmδ))], independent of the elastic-axis position. The page takes the thin-airfoil values CLδ = 2 T10 and Cmδ = −(T4 + T10) / 2 from Theodorsen's constants T4 = −cos⁻¹c + c √(1 − c²) and T10 = √(1 − c²) + cos⁻¹c, with the hinge at c = 1 − 2 E for an aileron of chord ratio E; real ailerons give less lift than this, so the estimate is for orientation.
Three approximations are printed beside the full solution. The quasi-steady result keeps every term but sets C = 1, that is, ignores the wake's lag; it is often far below the true flutter speed and for many sections predicts instability from zero speed, which is why it is not a design method. The steady result keeps only the static lift and finds the speed at which the bending and torsion frequencies merge. The approximate formula of Report 685, Vf / (b ωα) ≈ √[μ rα² ⋅ ½ / (½ + a + xα)], was offered for heavy wings (κ below 1/10) with a low frequency ratio.
Altitude
The mass of the wing is fixed, so μ rises as the air thins. The page repeats the solution through the standard atmosphere and lists the flutter speed as true and equivalent airspeed with the Mach number. The aerodynamics are incompressible; above a Mach number of about 0.6 the result is only an indication. A proper clearance uses aerodynamics computed for each Mach number and looks for the match point, where the flutter speed, the density and the Mach number all belong to the same altitude.
Wing divergence with sweep (NACA TN 1680)
Diederich and Budiansky treat one wing clamped at its effective root (where the elastic axis meets the fuselage side), with air forces by strip theory in planes perpendicular to the elastic axis. On a swept wing the bending slope changes the angle of attack of those strips: sweepback unloads the wing as it bends, sweepforward loads it. Their approximate formula is
qD = [(GJ)r / (me cr L³ cos² Λ)] ⋅ [L / (e1 cr)] ⋅ K1 / [1 − K2 ((GJ)r / (EI)r) (L / (e1 cr)) tan Λ],
with L the length of the wing along the elastic axis, cr the root chord normal to it, e1 the distance from the aerodynamic center aft to the elastic axis as a fraction of the chord, and (K1, K2) = (2.81, 0.614), (2.74, 0.497), (2.47, 0.390), (2.22, 0.326) for taper ratios 0.2, 0.5, 1.0 and 1.5 with stiffnesses varying as the fourth power of the chord; the page interpolates between them. For the uniform unswept wing the exact strip-theory value is π² G J / (4 me e1 c² L²), and 2.47 is π²/4 to three figures. The effective lift-curve slope is me = m0 A / (A + 4 cos Λ) at subsonic, subcritical speeds. A negative qD means the wing cannot diverge. Against an independent eigen-solution of the uniform swept wing the formula is within 2 percent. The divergence speed is the equivalent airspeed √(2 qD / ρ0).
Aileron reversal of a uniform wing
For a uniform unswept wing with a full-span aileron, torsion only, the twist under the aileron loads satisfies G J θ″ + q c² e1 me θ = −q c² (e1 CLδ + Cmδ) δ with θ(0) = 0 and θ′(L) = 0. Requiring the rolling moment ∫ (me θ + CLδ δ) y dy to vanish gives (e1 CLδ + Cmδ) φ(z) = Cmδ / 2 with z = q c² e1 me L² / (G J) and φ(z) = (sec √z − 1) / z, which the page solves for z. With e1 = 0 it reduces to qR = 2.4 CLδ G J / (me c² (−Cmδ) L²). When e1 CLδ + Cmδ is positive the aileron's own lift twists the wing nose up and there is no reversal below divergence. A swept or tapered wing is treated as uniform with its root properties and flagged.
Fin flutter criterion (NACA TN 4197)
Martin starts from the approximate formula above, evaluates the section at three-quarter span, inserts the torsional frequency of a solid tapered panel and an empirical aspect-ratio factor A / (A + 2), and arrives at
(Vf / a)² = GE / {[39.3 A³ / ((t/c)³ (A + 2))] ⋅ [(λ + 1) / 2] ⋅ (p / p0)},
with GE in lb/in², A the aspect ratio of one exposed panel (span over the chord at its midspan), t/c the streamwise thickness ratio, λ the taper ratio, and 39.3 lb/in² = 24 ε γ p0 / π for a center of gravity a quarter chord behind the quarter chord (ε = 0.25) and γ = 1.4. The note uses the bracket as the ordinate of a chart on which missiles that fluttered and missiles that did not are separated by a band, and calls the criterion crude. The page prints the note's X, the ordinate, Vf / a and the corresponding speed at the altitude entered. It covers bending-torsion flutter only, for solid sections of constant thickness ratio. The thickness ratio is taken on the midspan chord, which is this page's choice for a flat plate of constant thickness (the note assumes a constant ratio).
Comparison with the dive speed
14 CFR 25.629(b)(1) requires freedom from aeroelastic instability within the VD/MD envelope enlarged by 15 percent in equivalent airspeed. The page compares each critical equivalent airspeed with the entered factor times VD as a plain ratio. It is not a showing of compliance.
What is not covered
Control-surface degrees of freedom (aileron, tab), compressibility, three-dimensional unsteady aerodynamics, mode shapes of a real wing, whirl flutter, stall flutter, panel flutter, buzz, and the Regier-number estimate.