Summary
Our motion model is based on the spring-mass system shown in the figure on the right. The foil has two degrees of freedom: pitch (rotational) and plunge (vertical). The foil is connected to a linear plunge spring and a torsional pitch spring, and the zero for pitch and plunge is defined as the equilibrium point of the springs. The equations of motion for this system produce a system of two second-order ODEs (one for each degree of freedom). The characteristic equation of this system of ODEs can be used to analyze the system behavior and determine what configurations are stable and unstable.
We also programmed a Python script that calculates flutter parameters based on this motion model. An interactive calculator using this script can be found on the "calculator" subpage of this website.
The detailed explanation below is adapted from [4].
Resource [1]
Interlude: Notation
The following equations and derivations use a fair bit of aerodynamic notation, so here is a key for the meaning of some of the common symbols. A point worth noting is that for a given fluid (e.g. air) at a constant density, the dynamic pressure is a function of the airspeed. In the following text, it is assumed the reader understands this relationship, and the effects of dynamic pressure and airspeed are used interchangeably as appropriate.
Definition of Dynamic Pressure
Equations of Motion
The first step in analyzing the system is to write equations of motion by applying Newton's second law that the net vertical force is equal to the vertical acceleration, and the net torque is equal to the angular acceleration. This yields two equations (equations 1 and 2): one for each degree of freedom (i.e. pitch and plunge). It is of note that the effect of torque about the elastic axis due to gravity, and of the vertical displacement of the center of mass due to rotation, creates a coupling between the equations. This coupling is captured in the coupling inertia, which is defined in equation 3.
If we assume the solution will take the form of simple harmonic motion, the solutions for pitch and plunge will be as shown in equations 4 and 5, where both ω and the amplitudes may be complex numbers. Note that expanding these complex exponentials creates a waveform that is a damped oscillation with the real part of ω corresponding to the frequency, and the imaginary part of ω corresponding to the net damping term.
The assumption of simple harmonic motion also allows us to relate the spring constants, mass/moment of inertia, and natural frequencies of the structural modes alone (i.e. assuming no aerodynamic forces). These are shown in equations 6 and 7. This relationship is not used directly in the derivation of the characteristic equation of the system of ODEs, but allows us to see the limit behavior at very low airspeeds when the aerodynamic forces are negligible.
By plugging in the assumed form of solution into the equations of motion we obtain equations 8 and 9.
Aerodynamic Forcing
The equations of motion as presented in equations 8 and 9 treat lift as an unknown function of the angle of attack. Any further simplification requires treating lift as a known function of the angle of attack, so the next step is to get a function for lift. For our model, we assumed that torsion dominates, we neglected 2nd order aeroacoustic effects, we neglected compressibility effects, we assume that the airfoil is symmetric (i.e. has a coefficient of lift of zero at zero angle of attack), and we assumed the coefficient of lift varies linearly with the angle of attack. Starting with the lift equation (equation 10) [9], we can simplify to equation 11.
By plugging in this lift function, canceling the exponential term, and expressing the equations of motion in matrix form, we obtain equation 12, which relates the amplitude of the motion, the frequency, the geometry of the foil section, the mass properties, and the dynamic pressure.
Characteristic Equation
The number of unknowns in equation 12 is large enough that we can't plug in numbers and solve directly; however, we can use the structure of equation 12 to analyze the behavior of the system. In essence, equation 12 states that a matrix multiplied with the amplitude vector is equal to the zero vector. It can be clearly seen that there are two major cases for which this is true: one is when the amplitude vector is zero, the other is when the determinant of the matrix is zero.
The trivial solution is when the amplitude vector is zero. In this case, the system is "trivially stable" and nothing happens. This is the most common case (assuming a properly designed wing).
The other family of solutions is when the determinant of the matrix is zero. In this case, the amplitudes can be any number and the vector product will still be zero. This means that the amplitude is unbounded and the wing can oscillate with an amplitude that grows until the structure fails (or higher order effects take over).
By setting the determinant of the matrix to zero, we obtain equation 13, which is the characteristic equation for this system.
System Behavior
We can express equation 13 as a quadratic in terms of ω^2, with factors A, B, and C found by using algebra and shown below. Solving this quadratic for ω^2 using the quadratic formula yields equation 14. The determinant of this quadratic (B^2 - 4AC) can give a lot of information about the behavior of the system.
If the determinant is positive, there are two distinct, real roots. The fact that these roots are distinct means that the pitch and plunge modes will oscillate at different frequencies. Recall from equations 4 and 5 that the real part of ω describes the frequency, and the imaginary part of ω describes the net damping. The purely real ωs from a positive determinant represent pure oscillations with zero net damping. Because second-order effects are not taken into account by our model, these oscillations are damped in real life and will decay back to the trivially stable case.
As the dynamic pressure, q, increases, the determinant will decrease and eventually equal zero. When the determinant equals zero, the frequencies of the pitch and plunge modes are the same. This means the modes are synchronized and can feed energy into each other and cause flutter. The airspeed at which this happens is the critical flutter airspeed.
If the airspeed is increased any higher than the critical flutter airspeed, the determinant becomes negative. This means that ω becomes a complex number. It can also be seen that the imaginary component of this complex number is negative. From equations 4 and 5, we can see this means the oscillations have negative net damping and will grow in amplitude to infinity (in reality, the amplitude is limited by higher order effects or failure of the airframe).
What this all looks like in terms of real-world observation is that for low airspeeds, when given an initial disturbance, oscillations will decay to zero, and the frequencies of pitch and plunge will be different and close to their natural frequencies found in equations 6 and 7. For these low airspeeds, the system is underdamped. As the airspeed is increased, the frequencies of the pitch and plunge modes converge toward each other, eventually reaching a single value. Once they reach this single value, the modes synchronize and the oscillations grow in amplitude until the structure fails or the displacement is large enough that the linear spring and coefficient of lift assumptions break down. For these post-critical airspeeds, the system has negative damping
Divergence
By looking at equations 4 and 5, one may notice that if the real part of ω is zero, but the imaginary part is non-zero, the motion is a pure exponential. If the imaginary part in this case is negative, then the motion is an exponential growth to infinite displacement (in reality displacement is limited by second-order effects or failure of the airframe). This phenomenon when the displacement grows to structural failure without oscillation is called divergence. By setting ω to zero in the characteristic equation, we can see that the divergence condition occurs when equation 15 is satisfied.
Model Limitations
"All models are wrong, but some are useful" - George E. P. Box
Although very useful for analyzing system behavior (including non-intuitive relationships as explained on the "super cool discoveries" page), this simple model has several limitations.
One limitation particularly applicable to our test setup is that this model assumes a single mass for both pitch and plunge, which is the case for a real wing; however, the linear sliders and springs of our test rig have mass, so there is some mass in our system that is affected by plunge but not pitch. Splitting up the masses into a plunge mass and a pitch mass would make the model more accurate for our test setup.
Another limitation is that our model does not account for structural damping, nonlinear spring properties, nonlinear airfoil lift properties, or any other higher-order effects (e.g. friction, aeroacoustics, slop). For small displacements, these assumptions are accurate, but for large displacements, these effects can have a non-negligible effect. For instance, our model indicates that the amplitude of the oscillations grows to infinity for flutter above the critical airspeed. In the real world, these amplitudes are limited by effects not taken into account by our model.