3.1. 1st Test Case - Forward swept wing
The first test case presented is the aeroelastic analysis of a forward-swept wing. This type of wing, unlike conventional wings, presents a criticality due to the phenomenon of static aeroelasticity of divergence. In other words, the negative sweep angle causes the divergence speed to be lower than the equivalent straight or sweep back wing [
19,
20].
Figure 2 shows the aeroelastic model: aerodynamic paneling or Double Lattice panel model connected to the structural beam-like model. In a stick model, the wing is represented by beams capable of bending, shear, torsional and axial deformations.
In this example
Figure 3 and
Figure 4 represent the wingspan distribution of the stiffness and mass respectively. Each beam is divided into several elements with a flexural stiffness
EI and torsional stiffness
GJ estimated from the member section properties by classical structural analysis methods.
A considerable degree of detail is known for the mass distribution on the wing structure, but only part of the mass is structural, and therefore associated with the FE structural model, whereas a significant amount of mass is linked with non-structural elements such as fuel [
21].
To represent the mass distribution, the wing is divided into several strips, centered on the nodes of the beam-like model. For each section, the mass is lumped at the reference positions and attached to the beam axis node by a rigid link element that allows the lumped mass to be represented by section mass, a moment of inertia, and mass moments [
22,
23]. The aerodynamic panels distribution for the wing is carried out by taking into account the criteria suggested by [
24]. The first ten structural elastic modes are used for flutter analysis.
Figure 5 shows the mode shapes and natural frequencies only for the first 4 modes for the sake of brevity.
Figure 6 show the V-g and V-f diagrams of the P-K-method and the V-g and V-f diagrams of the G method for the wing with 10 modes at Mach number
M = 0.
Figure 6 shows the NASTRAN™ P-K solutions of the wing previously described at different given velocities from 50 to 500 m/s.
From the damping plot versus airspeed, it is possible to observe that the aerodynamic damping of mode 5, corresponding to the wing torsional mode, is nil at a speed of about 430 m/s.
In the frequency versus airspeed plot, at the same velocity, a coalescence of the frequencies of mode 3 and mode 5 is observed: there is a modal coupling between the aforementioned wing torsional mode and the 2nd wing bending mode. A flutter condition is attributed to this speed.
A significant aeroelastic observation is needed where divergence speed instability appears at 300 m/s i.e. a speed lower than the flutter velocity. This divergence speed instability is evident by its associated zero frequency. Comparing the upper plots with the lower ones in
Figure 6, it can be seen that the agreement between the damping computed by the P-K method and the g-method is expected. Good agreement in terms of the overall V-g and V-f comparisons between these methods is obtained except the g-method predicts one extra aerodynamic lag root (represented by the crosses in
Figure 6).
This aerodynamic lag root appears at divergence speed (about V=300 m/s) with stable damping. This is an interesting phenomenon because it indicates that the divergence speed could be a bifurcation point. This result suggests that the divergence speed is caused by the coupling of a structural mode and an aerodynamic lag root and should be considered as a special case of flutter instability, the so-called dynamic divergence. This is supported by the g-method results shown in the lower plots of the figures where the frequency coalescence of the first bending mode and the aerodynamic lag root is seen. On the other hand, the plots of the P-K method do not show this result due to its incapacity of generating the non-zero-frequency aerodynamic lag root.
Table 1 shows the flutter analysis results obtained by means ZAERO, in terms of the mode associated with flutter, natural frequency, and flutter speed.
For this first comparison, two types of instability are predicted by both the P-K method and the g-method: a flutter speed at about V = 430 m/s and a divergence speed at about V = 300 m/s. This agreement is expected since at damping g = 0 the flutter equation of both methods reduces to the same form. The damping curves of the structural modes computed by both methods are in excellent agreement. The frequency curves of the two structural modes computed by both methods are also in good agreement except for the absence of the aerodynamic lag roots of the P-K method.
For a question of investigation completeness, a second step involved an analysis of the same wing by adding the mass due to the fuel. The V-g and V-f diagrams relating to the case of full fuel tanks and a table with the flutter data are reported.
The V-g diagram shows that the damping of mode 3 crosses the zero damping axis, indicating a fluctuation limit of the wing structure, at higher values than in the previous case.
From the V-g and V-f diagrams, calculated for the different load conditions we can draw the following conclusions:
in all the calculations the divergence condition for the wing is present, the latter being independent of the mass contributions from the model;
in both cases, with and without fuel, the divergence rate does not change;
the most critical condition for flutter occurs when the wing is unloaded (the addition of the masses of the fuel involves the removal of the bending frequencies from the torsional ones);
when the wing tanks are full there is the annulment of the aerodynamic damping of the first torsional frequency at a speed of about 44% higher than in the case of empty tanks. From the V-f diagram, we can ob-serve a coalescence of the aforesaid frequency with the second flexural frequency.
Upon completion of the study of the first test case, a dynamic flutter analysis was performed as the stiffness distribution varies along the wingspan. In particular, a wing section of about 1 meter was chosen, a representative for example of an inspection window for maintenance or the extraction of the landing gear. It has been hypothesized that, in the notch area, with the same mass distribution, the flexural and torsional stiffnesses decreased by 30% of the initial value, i.e. wing without the door.
To study how the flutter speed varies as a function of the location of the trap door, for the two stiffness variations the first area was considered respectively between 20 and 30% of the wingspan
Figure 8a and the second between 30 and 40%
Figure 8b. For the first variation of stiffness distribution, a flutter speed of 410 m/s is obtained (
Figure 8c) or a decrease of about 4.7% compared to the starting wing. For the second variation of stiffness distribution, on the other hand, the flutter speed is about 418 m/s (
Figure 8d) and therefore is decreased by 2.8% compared to the reference wing.
Figure 7.
V-g and V-f curves P-k method (Nastran) and G method (ZAERO) - wing with full fuel tank.
Figure 7.
V-g and V-f curves P-k method (Nastran) and G method (ZAERO) - wing with full fuel tank.
Table 2.
ZAERO Aeroelastic Analysis results - wing with full fuel tanks.
Table 2.
ZAERO Aeroelastic Analysis results - wing with full fuel tanks.
|
K method |
|
|
G method |
|
mode |
vf(m/s) |
Freq. (Hz) |
mode |
vf(m/s) |
Freq. (Hz) |
2 |
326.1 |
6.53 |
2 |
324.5 |
6.53 |
3 |
671.0 |
5.60 |
4 |
665.8 |
5.77 |
5 |
591.9 |
19.59 |
5 |
581.8 |
19.59 |
3.2. 2nd Test Case - Forward swept wing
The second test case has been selected due to the availability of the experimental modal analysis results and offered the opportunity to study a realistic case. The aircraft under study is an ultralight twin-seat aircraft whose fuselage has been changed from metallic to a carbon-fiber material (max cruise speed = 59 m/s and max take-off weight = 600 kg). The main dimensions are shown in
Figure 9. A complex FEM model of the aircraft under examination made it possible to perform preliminary dynamic flutter analyzes. This model and the mode shapes deriving from the dynamic analysis, is shown in the
Figure 10.
For the study of the comparison between the two different codes, it was instead decided to adopt a simplified model of the horizontal tailplane, which would allow obtaining adequate results with the advantage of saving time.
Table 3.
Results of the aeroelastic analysis as the stiffness varies along the wingspan.
Table 3.
Results of the aeroelastic analysis as the stiffness varies along the wingspan.
Reference |
Stiffness var. |
Stiffness var. |
stiffness distribution |
@ 20-30% |
@ 30-40% |
vf(m/s) |
vf(∆%) |
vf(∆%) |
430 |
−4.7 |
−2.8 |
In particular, the flutter analysis refers only to the horizontal tailplane. As regards the first, simplified stabilator FEM model, the structure of the fuselage is simplified and idealized through six elastic elements that constitute the degrees of freedom. In more detail, the central node of the stabilator is connected to a fixed node of the fuselage utilizing six elastic elements, which represent the stiffness of the fuselage. Rigid elements have also been introduced in the model to better highlight the modes, facilitating their visualization in the post-processing phase and therefore the distinction between bending and torsional mode shapes. The aeroelastic FEM model is shown in
Figure 11, in which it is possible to appreciate the constraint, the concentrated masses, the rigid elements, the aerodynamic panels, the set of nodes used for the interpolation, the spline itself.
For the analysis of flutter of the stabilator, various study cases are taken into account and they differ in the constraint condition (fixed) or in the unconstrained degree of freedom (rotation around the longitudinal axis x, rotation around the transverse axis y, and rotation around the vertical axis z) and for the value assumed for the stiffness of the elastic element. From the calculations performed by NASTRANTM, the Flutter condition is not found in any of the previously mentioned cases. A comparison is then made with the results obtained using the ZAERO code. First, however, it was deemed necessary to carry out a sort of sensitivity analysis to better understand how the aeroelastic responses of the structure are influenced by changes to the variables of the ZAERO code such as the number of aerodynamic panels along chord and along wingspan, nodes used for structural and aerodynamic interpolation, type of spline: linear or surface spline.
The purpose is to identify the parameters that make the analysis more conservative concerning the identification of the Flutter condition.
By neglecting to report the results of the numerous sensitivity analysis, from an evaluation of the values and the trend of the V-g and V-f graphs obtained from the analysis, it is possible to conclude that the best way to proceed is to use a surface-type spline that allows the end nodes of the rigid elements to be included for interpolation. In any case, the numerical flutter modes obtained through this model have too high and unrealistic natural frequencies.
When ensuring the safety of aircraft against aeroelasticity phenomena, one of the important steps of investigation is Ground Vibration Test.
GVT can be performed to obtain the vibration amplitude distributions and the frequencies of the structure. The GVT results can be used for checking FEM results or used directly in the prediction of the flutter problems [
25].
The next phase concerns the verification and tuning of the analytical results to the experimental ones. At this point, aerodynamic loads are modeled and measured using a wind tunnel test and combined with the dynamic properties in an analytical flutter prediction.
Finally, test flights are carried out for validation of the reached results and final certification, according to the requirement of the airworthiness standard [
26].
The vibrations test has been performed by measuring the frequency response functions (FRF) between several points of the stabilator structure, each concerning a fixed position of the excitation point.
The excitation point is located in the middle of the rear spar of the horizontal tail. A maximum of 12 points for each side of the stabilator was measured. The excitation signal used during these tests has been a Sine Sweep Excitation, obtained through an electro-dynamic exciter Modal Shop Model 2100 E11 while the structural response has been measured with PCB 333B32 ICP piezoelectric accelerometers.
The evaluation of the FRF was carried out by using a standard 8 channels FFT analyzer, the LMS Scadas 05 Mobile. The LMS TestLab software provides the subsequent steps for the analysis and the extraction of the modal parameters by means LMS PolyMAX algorithm. This modal parameters identification method yields extremely clear stabilization charts. Following the vibration tests, based on the experimental modal analysis results (4), the idea of updating the model of the horizontal tailplane was evaluated, thus making it better representative of the test article on which the tests were carried out.
The best modeling consists of discretizing the tailplane in as many points as are those considered in the test setup, in which the accelerometers have been placed (
Figure 12). As regards the masses, it has been assumed that the stabilator weighs a total of 0.8 kg. This weight is distributed equally over the 24 nodes in which the stabilator is discretized while for the fuselage section weight of about 2 kg is assumed, equally distributed over the 5 nodes of the fuselage.
The aerodynamic model consists of 10 boxes along the chord and 15 along the wingspan. The new aeroelastic model is shown in
Figure 13 in which the distribution of concentrated masses and the nodes used for the aerodynamic-structural spline are highlighted.
Preliminarily, a modal analysis and a Dynamic Flutter analysis are carried out of the new model using the P-K method. The results of these analyzes are reported for the first 20 modes.
Table 4 indicates the natural frequencies obtained while
Figure 13 indicates the mode shapes.
As regards the dynamic flutter analysis, as can be seen from
Figure 14, employing the P-K method of Nastran, no Flutter instability occurs for zero structural damping up to
Vmax = 200 m/s. The analysis results obtained using the ZAERO software, on the other hand, highlights a flutter instability associated with the 9
th mode of a torsional nature, as shown in
Figure 15 and
Figure 16. The flutter speed calculated with both the k and g methods is shown in
Table 5. In the same table it is also observed how by increasing the value of the structural damping, the flutter speed increases.