Flow past a self-oscillating airfoil with two degrees of freedom: measurements and simulations

The paper focuses on investigation of the unsteady subsonic airflow past an elastically supported airfoil for subcritical flow velocities and during the onset of the flutter instability. A physical model of the NACA0015 airfoil has been designed and manufactured, allowing motion with two degrees of freedom: pitching (rotation about the elastic axis) and plunging (vertical motion). The structural mass and stiffness matrix can be tuned to certain extent, so that the natural frequencies of the two modes approach as needed. The model was placed in the measuring section of the wind tunnel in the aerodynamic laboratory of the Institute of Thermomechanics in Nový Knin, and subjected to low Mach number airflow up to the flow velocities when self-oscillation reach amplitudes dangerous for the structural integrity of the model. The motion of the airfoil was registered by a high-speed camera, with synchronous measurement of the mechanic vibration and discrete pressure sensors on the surface of the airfoil. The results of the measurements are presented together with numerical simulation results, based on a finite volume CFD model of airflow past a vibrating airfoil.


Introduction
In aerospace engineering, fluid-structure interaction can play a very important and potentially dangerous role: under certain circumstances, the coupling between flow and structure may lead to unstable exponentially increasing oscillations. The classical example is the flutter instability of airfoils, which occurs for systems with two degrees of freedom when the critical flow velocity is surpassed [1]. The consequent large-amplitude oscillations can result in destruction of the elastic part due to low-cycle fatigue [2], with possibly catastrophic consequences. Thus, the importance of theoretical and experimental analysis of the aeroelastic phenomena in these applications cannot be placed in doubt.
The purpose of this study is the research of the flowinduced vibration of airfoils using interferometry and CFD. The current paper presents the results of first windtunnel measurements of a newly reconstructed model of a NACA0015 airfoil with two degrees of freedom (DOF)pitching (rotation about the elastic axis) and plunging (vertical translational motion). The results are compared to numerical simulations of airflow past an airfoil, whose motion is prescribed based on the kinematic data identified from the measurements.

State-of-the-art
The majority of experimental works dealing with 2-DOF vibration of airfoils is still limited to the simplified case of forced (externally excited) vibration [3][4][5][6][7][8][9][10][11][12][13]. The standard method for the measurement of the flow field in the proximity of the vibrating profile is the Particle Image Velocimetry (PIV), complemented with local pressure measurements using pressure probes. The usage of interferometric or schlieren methods has not yet been published in this specific context, although these methods are well suited especially for the case of higher subsonic, transonic and supersonic velocities where the PIV method is hardly applicable. Moreover, the interferometric measurements can provide pressure field data for the calculation of the forces and moments acting on the airfoil. The drawback of the interferometric approach -the assumption of the isentropic flow in the calculation of the pressure field -can be minimized by concurrent pneumatic measurements used to calculate correction factors.
The interferometric images provide results, whose precision is proportional to the test section width and flow velocity (the number of interferometric fringes increases linearly with the light path length, and also with the Mach number). The PIV and interferometric methods can be seen as complementary, since with the increasing Mach number the difficulties of a typical PIV setup generally grow, but the applicability of the interferometric method is not influenced and the precision of the results evaluated from the fringes increases.

Physical design
The experimental setup, which can be seen in figure 1, builds on previous experience obtained in the Institute of Thermomechanics in the years 2006-2012, when selfoscillating models of the NACA0015 and DCA18% airfoils with two degrees of freedom were designed and the methodology of the measurements developed [14][15][16]. A model of a symmetric NACA0015 airfoil with a chord length of 65 mm and width of 80 mm was placed in a 80x210 mm test section of a vacuum wind tunnel in the laboratory of the Institute of Thermomechanics in Nový Knín. The supporting frame of the airfoil allows two degrees of freedom -vertical shift (plunge) and rotation about elastic axis (pitch), which is located at 1/3 of the chord length (see figure 2). The translation stiffness is realized by adjustable flat springs, the torsional stiffness comes from a thin torsion bar inside the profile. For higher flow velocities and pitching amplitudes, this bar is highly dynamically loaded and prone to destruction. During first preliminary measurements, the torsion bar underwent rupture due to excessively high pitching amplitudes, destroying also the cables and some of the sensors inside the airfoil. Certain parts of the model had to be remanufactured, and the guiding mechanism was also slightly improved. The total mass of the dynamic system can now be tuned by mounting additional masses onto the moving frame. In this way, the eigenfrequency of the plunging mode can be shifted apart from the torsional mode eigenfrequency. From the theoretical point of view, this should postpone the flutter onset to higher flow velocities.
Prior to wind tunnel measurements, the natural frequencies and damping ratios of the system with various values of additional masses attached was measured and identified in laboratory [17]. For the system with all three additional masses of m = 350 g, the eigenfrequency of the translational mode is 17.7 Hz, about 2-3 Hz lower than without the additional masses. The natural frequency of the rotational mode is 24.8 Hz and should not depend on the additional masses.

Sensors and wiring
The airfoil is equipped with six miniature pressure transducers Freescale MPXH6115 mounted flush with the airfoil surface -three on the top surface, three on the bottom. Due to damage of the model during previous measurements, only two sensors (V1h and V2h, both on the top surface) remained functional. The signal from the transducers, which have built-in preamplifiers, leads off the measuring section by a flexible bundle of wires designed for high mechanical loads. dynamic range of this sensor was found to be insufficient for very high amplitudes (the limit arises from exceeding the maximum allowed linear velocity), and so an additional strain gauge bridge was mounted on the bottom flat spring and calibrated to provide a redundant signal of the vertical deflection. The measuring section of the wind tunnel is equipped by a Prandtl tube, measuring the static pressure p 1 and total pressure p 0 . The Prandtl tube is connected to a pressure scanner wired to a PC, which registers the pressures and evaluates the Mach number using relation . 1 1 Here κ = c p /c v is the heat capacity ratio, for ideal diatomic gas equal to 7/5. The flow velocity can be then calculated as is the specific gas constant for dry air and T is the absolute local static temperature calculated from the ambient temperature T 0 as . 2 The signals from the airfoil pressure transducers, strain gauge bridge and optical length gauge sensor are digitized, monitored online and stored by a Dewetron software system.
In standard configuration, the flow field around the airfoil is measured by Mach-Zehnder interferometer and a high-speed camera. Due to technical issues, the interferometer was not mounted during this set of measurements. To monitor the airfoil motion (the vertical deflection is measured by the length gauge, but for the pitch angle there is currently no sensor installed on the model), a high-speed camera (NanoSense MKIII, maximum resolution 1280x1024 pixels @1000 Hz) and a special light source were mounted in the perpendicular direction, using the same optical access as the interferometric setup. The camera is triggered from the Dewesoft module, which ensures perfect synchronisation of the pressure and deflection measurements with the high-speed camera frames. The region of interest was set to a rectangle 756x700 pixel to spare disk storage and to save time during the transfer of data from the camera to the PC, which was the most time-consuming operation during the measurement.

Numerical model
The CFD model is described in detail in [18]; here it is sufficient to state that the model is based on cell-centred finite volume discretization of incompressible Navier-Stokes equations solved on a moving mesh. Two variants of the model are used: a fully 3D model, and a much computationally cheaper 2D model used for fast testing and preliminary analysis. The numerical code uses the OpenFOAM package, which is available as a free opensource software under the GNU general public licence. The motion of the airfoil is prescribed as a 2-DOF harmonic oscillation, with kinematic data identified from the experimental results.
Most of the CFD simulations were run in parallel on a SGI Altix UV 100 shared-memory supercomputer of the Supercomputing Centre of the Czech Technical University in Prague. This shared-memory machine is built on the cache coherent non uniform memory access (cc-NUMA) architecture and offers 12 6-core Intel Xeon Nehalem processors, with 8 GB per core RAM. The nodes are interconnected by the SGI NUMAlink 5 interconnects, providing low latency and 15 GB/s bandwidth through two 7.5 GB unidirectional links. The OpenFOAM library is compiled using standard gcc, and the parallel solvers run via SGI MPT implementation of the MPI standard.

Evaluation of the airfoil motion from the high-speed camera images
The measurement was performed for 32 cases: four configurations of additional masses and various values of inlet flow velocities. In all cases, first the airfoil was forced an initial deflection of y = 3 mm. Once the flow velocity in the wind tunnel was adjusted as needed, the high-speed camera started recording and shortly afterwards the airfoil was released.
During the 32 measurements, altogether more than 30 000 frames were recorded by the high-speed camera. The position of the airfoil in all these frames (vertical shift and pitch angle) were evaluated after the wind tunnel measurements using the positions of three key points (see figure 4), which were tracked automatically using Olympus i-Speed Advanced software suite. On a standard laptop, the software was able to process about 10 frames per second.  Figure 5 demonstrates a case where the flow velocity is subcritical. From the initial deflection (y = 3 mm, pitch = 2º), the vertical displacement oscillates with exponential attenuation to zero. As the two modes are coupled (the elastic axis is not located in the centre of gravity), a short transient in the pitching mode is excited, too. The offset of the vertical deflection waveform does not have physical meaning; it is caused by an arbitrary position of the coordinate system origin when processing the camera frames. The offset of the pitch waveform (dark red) of about 2º corresponds to a real initial inclination of the airfoil, caused by imperfect adjustment of the torsion bar.
When the flow velocity is increased (measurements 2792-13, 2792-26..31 -not shown here), the behaviour of the system remains qualitatively the same, but the damping factor decreases due to negative damping contribution from the aerodynamic damping matrix. A measurement at a supercritical flow velocity, when the flutter instability is present, is demonstrated in figure 6. In this case, the real part of the complex frequency is not negative any more, and within the validity of a linear approximation the amplitudes of the oscillations would exponentially grow until system destruction. Due to nonlinearities of the system (both mechanical and aerodynamic), the amplitudes eventually stabilize at a finite value.
To visualize the amplitudes and phase differences between the two modes, the airfoil motion can be also displayed in phase plane. Tab. 1 reports the phase plots of all the cases measured, with vertical displacement y on the bottom axis and pitch angle on the left axis. In the cases where the system reached the flutter instability, the frequency and amplitudes of both modes are listed in the table, too. For configurations #1, 2 and 3, no measurements were performed beyond the lowest flow velocity able to trigger the flutter instability in order not to destroy the model: it was not known a priori, whether the vibration amplitude stabilizes due to system nonlinearities earlier, than the torsion bar breaks. The last measurements in configuration #4, however, recorded the system behaviour at several supercritical flow velocities up to Ma = 0.2.
In general, from the data given in Tab. 1 it can be concluded that the current model behaves in a predictable way, consistent with the classical theory of 2-DOF aeroelastic systems. The flutter instability occurs at certain critical flow velocity, with a frequency in between the natural frequencies of the mechanical system (in our case somewhat closer to the frequency of the translational mode). When additional mass is applied to the frame, increasing the corresponding translational natural frequency and widening the gap between the two frequencies, the flutter instability is postponed to higher flow velocities. The effect is, however, not very distinct: the critical flow velocity for configuration #1 (m = 0 g) is Ma = 0.145, for configuration #4 (m = 350 g) flutter occurs at Ma = 0.160. When the flow velocity is increased beyond the critical value, the oscillation amplitudes grow significantly (from 7.1º / 0.81 mm at Ma = 0. 16

Pressure measured on the airfoil surface
In the current set of measurements, due to technical complications it was unfortunately not possible to install the interferometric setup. Thus, the only measurements of the flow field are the signals from the discrete pressure sensors mounted on the surface of the oscillating airfoil. The relatively cheap MPXH6115 miniature pressure sensors provided unexpectedly good and accurate signals: prior to wind tunnel measurements, they were verified against a high-end Kulite XCQ-080 pressure transducer, with no significant discrepancy in amplitude and phase up to at least 1 kHz. The time history of the pressure from these two sensors is shown in figure 7 for two cases, both in the flutter regime. The ambient atmospheric pressure was 98.32 kPa. As can be seen from comparison with

Numerical simulation results
The numerical simulations were run using the simplified 2D model, with the airfoil kinematic data and inlet Mach numbers matching values of measurements #2792-25 and #2792-32. The frequencies and amplitudes of the oscillation modes are listed in Tab. 1. For #2792-25, the phase lag of the pitching mode behind the plunging mode was identified to 1.182 π, in #2792-32 it was 1.134 π. Figure 8 shows the simulated pressure field for a time instant approximately in the middle of the 6 th period of vibration, where the vertical deflection of the airfoil passes through zero position (cf. figure 9). The pitching angle being 9.85º, the highest pressure (about 4 kPa) is located on the bottom surface near the leading edge. On the top surface, the airflow separates and forms a highly turbulent separation region, with vortex structures forming zones of very low pressure. As the simulation is performed in 2D, the vortex dynamics are completely different from the real 3D case. In particular, the vortex structures cannot interact and change orientation, and thus disperse with a much slower rate than they would in 3D. This is a known drawback of the 2D simulations; the predictions of the 2D model in the separation region must be interpreted cautiously.
The time history of the pressure in two probes (V1h, V2h), located at the places where the pressure transducers are mounted on the physical model, are shown in figure 10. These results can be directly compared to the values measured in the wind tunnel and displayed in figure 7. As expected, the pressure maxima correspond to the time instants of minimum pitching angle (cf. figure 9 for the simulation #2792-32).   In both simulations, the waveform of pressure probe V1h (in 18% of the chord length) compares very well to the measurement, only with the minimum values slightly lower and more influenced by pressure disturbances than the measured data. The pressure probe V2h is located at 47% chord. During almost half of the oscillation period (when pitching angles are positive), the probe already lies in the separation region. Thus, the pressure predictions of the numerical model are negatively influenced by the fact that the 2D simulation overestimates the turbulent pressure fluctuations, as discussed earlier. However, during negative pitching angles the V2h signal still compares well with the experiments.

Conclusions
Based on the previous experience and measurements of the authors, a physical model of a NACA0015 airfoil with two degrees of freedom has been modified and repaired, equipped with sensors and measured in a wind tunnel at various values of subcritical and supercritical flow velocities. The mechanical vibration and static pressure at two discrete locations of the upper airfoil surface were measured, with synchronous high-speed imaging of the airfoil motion. The stability boundaries for different values of additional masses (i.e., different eigenfrequencies of the plunging mode) have been identified, and the post-flutter behaviour of the system analysed. The results served well as a frame of reference for subsequent set of measurements, where an interferometric system was mounted. However, the interferometric data has not yet been processed and evaluated.
The measured signals from the pressure sensors have also been used as validation data for a concurrently developed CFD model of airflow past an oscillating airfoil. The results of the numerical simulations match well in the regions of attached flow. However, at high angles of attack when the airflow separates, the 2D CFD predictions strongly overestimate the turbulent pressure fluctuations. In this separated flow regime, both numerical simulations and interpretation of the interferometric data remain a challenging goal for the future.