Get 20M+ Full-Text Papers For Less Than $1.50/day. Start a 14-Day Trial for You or Your Team.

Learn More →

An effective computer modelling approach to the study of aeroelastic characteristics of an aircraft composite wing with high aspect ratio

An effective computer modelling approach to the study of aeroelastic characteristics of an... Mathematical and Computer Modelling of Dynamical Systems, 2015 Vol. 21, No. 1, 58–76, http://dx.doi.org/10.1080/13873954.2014.903283 An effective computer modelling approach to the study of aeroelastic characteristics of an aircraft composite wing with high aspect ratio a b a a a Fusheng Wang *, Shihui Huo , Shengjun Qiao , Junran Zhang and Zhufeng Yue Advanced Materials Test Center, School of Mechanics, Civil Engineering and Architecture, Northwestern Polytechnical University, 710129 Xi’an, China; Mechanics & Environment Research Center of Rocket Engine, Xi’an Aerospace Propulsion Institute, 710100 Xi’an, China (Received 6 October 2013; accepted 7 March 2014) Static aeroelastic and flutter characteristics of an aircraft composite wing with high aspect ratio were analysed by an effective Computational Fluid Dynamics and Computational Structure Dynamics coupled method. Effects of stiffness distribution on aeroelastic characteristics were considered. Honeycomb core sandwich composite was considered to be equivalent to an orthotropic material by stiffness and inertance equivalent method to allow highly efficient numerical simulation, which was used for analysis of bending and torsional stiffness distribution. The results showed that the redistributed aerodynamic load leads to a decrease of pressure difference between the upper and lower airfoils. The flutter speed of the composite wing is near 0.64 Ma. Both bending and torsional stiffness increases with a small increase of beam size. Stiffness of the wing root has a major influence generally on the static aeroelastic characteristics. Both the lift coefficient and the loss percent decrease with a small increase of beam size. Effects of stiffness distribution on frequency are not obvious. Flutter speed remains close to the initial value when the beam size is changed. Keywords: honeycomb core sandwich composite; high-aspect ratio wing; stiffness distribution; static aeroelasticity; flutter; CFD/CSD coupled method 1. Introduction Composite airplanes with high-aspect ratio wings have received widespread attention in recent years throughout the world due to their high flight altitude and long endurance. A high-aspect ratio wing usually has lighter weight and more flexibility, which results in many special and complex aeroelasticity problems. The stiffness distribution needs to be considered in the design of a high-aspect ratio wing, which is a beneficial reference for further structural and aeroelastic analysis. Structural and aerodynamic characteristics are both affected significantly by the stiffness distribution of a high-aspect ratio wing. Many scholars have conducted research in this field and several aeroelastic analytical methods and stiffness distribution discus- sions have been put forward. Haddadpour et al. [1], Demasi and Livne [2] and Dang et al. [3] investigated the aeroelastic stability of an aircraft wing, which was modelled as an anisotropic composite thin-walled beam. The coupled method of Computational Fluid Dynamics (CFD) and Computational Structure Dynamics (CSD) is an effective approach to aeroelastic analysis, in which the structural and aerodynamic characteristics are *Corresponding author. Email: fswang@nwpu.edu.cn © 2014 Taylor & Francis Mathematical and Computer Modelling of Dynamical Systems 59 analysed independently followed by information transfer of load and deformation between the two analysis tools. Slone et al. [4] considered the application of a full physics simulation in the context of the AGARD 445.6 wing geometry. Bartoli and Mannini [5] investigated flutter instability of flexible bridge decks in framework of a semi-empirical Scanlan’s approach based on flutter derivatives. Borglund and Kuttenkeuler [6] investi- gated aeroservoelastic behaviour of a thin rectangular wing with a controllable trailing edge flap. Guo [7] aimed at presenting an investigation into optimal design about the minimum weight and aeroelastic tailoring of an aerobatic aircraft wing. Information transfer between different disciplines and dynamic mesh methods are concerned in the coupled CFD and CSD methods. Many scholars have studied this field and several methods have been put forward such as spring analogy method [8–10], elastic solid method [11–13] and layered elastic solid method [14]. Since stiffness distribution influences aeroelastic characteristics of a high-aspect ratio wing greatly, appropriate stiffness distribution corresponding to ideal aeroelastic charac- teristic can be obtained. In the present study, analysis on stiffness distribution and aeroelastic characteristic of a high-aspect ratio wing was carried out by a highly efficient computer modelling approach. And then effects of stiffness distribution on aeroelastic characteristics of the high-aspect ratio wing were considered. 2. Stiffness-distribution analysis Stiffness-distribution analysis is carried out based on an aircraft composite wing with high aspect ratio. A finite-element (FE) model of the high aspect ratio wing is shown in Figure 1, in which typical cross sections are illustrated, representing different rib locations with red colour and used to study stiffness distribution along the spanwise wing. Here, x is the spanwise coordinate, y is the chordwise coordinate, direction of z coordinate is vertical Figure 1. Finite element model of high-aspect ratio wing. 60 F. Wang et al. to the plane made of x and y coordinates. The wing structure is mainly composed of skin, longitudinal and transverse components. Longitudinal components include beams and stringers, while transverse components include ribs and thin-walled jointed structures. The elastic constraint of the wing root is simulated through an extension used for reinforcement support of root cross section. The upper and lower wing skins are made of honeycomb core sandwich composite. When the FE software NASTRAN (MSC Inc., Elkhart, IN, USA) is adopted for structural analysis, no relative element types are selected for the honeycomb sandwich plate. Although both three-dimensional FE model and laminated elements can be used to analyse honeycomb core sandwich composite, a three dimensional FE model will be costly in terms of the computation time, and honeycomb core sandwich composite cannot be simulated well by laminated elements. An effective modelling approach will allow the input parameter values for the FE software NASTRAN to be obtained (such as equivalent density, Young’ modulus, the shear modulus and Poisson’sratio), through which honeycomb core sandwich composite is equivalent to an orthotropic material. In the present study, the honeycomb core sandwich composite was regarded as equivalent to an orthotropic material by the stiffness and inertance equivalent method [15], the validity of which has been established by comparing results obtained from theory, those calculated by the equivalent method. It is found that differences between results from theory and from the equivalent method are less than 3%. This facilitates the numerical simulation process. The initial honeycomb core sandwich composite is shown in Figure 2(a). The total thickness of honeycomb sandwich plate is 2h, in which thickness of honeycomb structure is 2h and that of upper or lower surface is h . The equivalent orthotropic material of the honeycomb core sandwich composite is shown in Figure 2(b). Material properties can be obtained by a stiffness and inertance equivalence method expressed as 0 1 0 1 E ðe e  e Þ=e 11 22 22 1 12 B C B E C ðe e  e Þ=e 11 22 11 2 12 B C B C B C B G C e B C ¼ (1) B C B C B G C e B C @ A @ A e =e 12 22 Figure 2. Honeycomb core sandwich composite. (a) Initial honeycomb core sandwich composite; (b) Equivalent orthotropic material. Mathematical and Computer Modelling of Dynamical Systems 61 where 0  . 1 3 3 3 f 3 c ðh þ h Þ  h e þ h e ðh þ h Þ c f c f c 11 c 11 0 1 B C 1 B C 3 3 3 f 3 c B C ðh þ h Þ  h e þ h e ðh þ h Þ c f c f B C c 12 c 12 e B C B C  . B C B C 3 3 3 f 3 c B C ðh þ h Þ  h e þ h e ðh þ h Þ B C c f c f c 22 c 22 ¼ (2) B C B C B C f c B C B ðh e þ h e Þ ðh þ h Þ C f c c f 44 44 @ A B C f c B C ðh e þ h e Þ ðh þ h Þ f c c f 55 55 @  . A 3 3 3 f 3 c ðh þ h Þ  h e þ h e ðh þ h Þ c f c f c 66 c 66 f c In Equation (2), e and e are the stiffness coefficients of the surface and honeycomb ij ij structure, respectively. The density of the equivalent orthotropic material can be expressed as h ρ þ h ρ f c f c ρ ¼ (3) h þ h f c where ρ and ρ are the densities of the surface and honeycomb structures, respectively. f c Based on Equations (1) and (3), the honeycomb core sandwich composite can be equivalent to an orthotropic material. For different parts of the high-aspect ratio wing, stiffness coefficients of the honeycomb core sandwich composite and properties of the orthotropic material on the upper skin are given in Tables 1 and 2, respectively. Also, the stiffness coefficients of the honeycomb core sandwich composite and properties of the orthotropic material on lower skin are given in Tables 3 and 4, respectively. In Tables 2 and 4, E and E are Young’ moduli, respectively. G , G and G are the shear 11 22 12 13 23 moduli, respectively. μ is the Poisson’s ratio. The stiffness distribution is analysed based on different cross sections of the high- aspect ratio wing, which is selected according to rib positions. In the present study, Table 1. Stiffness coefficients of the honeycomb core sandwich composite on upper skin (MPa). Stiffness Ribs Ribs Ribs Ribs Ribs Ribs Ribs coefficient (1–2) (2–4) (4–5) (5–6) (6–7) (7–10) (10–13) e 23,550 21,830 15,105 21,848 21,913 19,167 22,317 e 7500 6599 3285 2834 2143 5191 5658 e 11,212 11,040 11,609 10,671 10,925 10,409 10,578 e 28,322 27,108 20,179 29,807 30,131 25,000 28,164 e 28,322 27,108 20,179 29,807 30,131 25,000 28,164 e 26,588 25,393 18,869 27,568 27,815 23,332 26,303 e 13,116 12,720 11,174 12,724 12,740 12,108 12,832 e 9225 9218 8456 8353 8193 8894 9002 e 10,279 10,239 10,370 10,154 10,213 10,094 10,133 e 11,309 11,222 10,727 11,415 11,438 11,071 11,297 e 11,309 11,222 10,727 11,415 11,438 11,071 11,297 e 13,815 13,540 12,039 14,040 14,096 13,066 13,749 66 62 F. Wang et al. Table 2. Properties of the orthotropic material on upper skin. Ribs Ribs Ribs Ribs Ribs Ribs Ribs Properties (1–2) (2–4) (4–5) (5–6) (6–7) (7–10) (10–13) E (MPa) 58,600 61,100 62,400 58,600 69,100 65,000 53,700 E (MPa) 24,600 25,500 33,700 31,500 25,300 26,800 30,600 G = G = G (MPa) 16,200 15,300 12,600 14,500 12,900 13,800 16,200 12 13 23 μ 0.529 0.487 0.313 0.38 0.421 0.424 0.43 Table 3. Stiffness coefficients of the honeycomb core sandwich composite on lower skin (MPa). Stiffness coefficient Ribs (1–4) Ribs (4–5) Ribs (5–6) Ribs (6–7) Ribs (7–9) Ribs (9–13) e 23,289 26,040 21,913 19,167 22,860 22,317 e 6517 7934 2143 5191 3305 5658 e 11,216 11,062 10,925 10,409 10,615 10,578 e 28,491 30,214 30,131 25,000 30,120 28,164 e 28,491 30,214 30,131 25,000 30,120 28,164 e 26,682 28,117 27,816 23,332 27,906 26,304 e 13,056 13,689 12,739 12,108 12,957 12,833 e 9180 8974 8193 8894 8462 9002 e 10,280 10,244 10,213 10,094 10,142 10,133 e 11,321 11,444 11,434 11,071 11,437 11,297 e 11,321 11,444 11,434 11,071 11,437 11,297 e 13,836 14,166 14,096 13,066 14,118 13,749 Table 4. Properties of the orthotropic material on lower skin. Properties Ribs (1–4) Ribs (4–5) Ribs (5–6) Ribs (6–7) Ribs (7–9) Ribs (9–13) E (MPa) 57,100 53,900 69,100 65,000 59,900 53,700 E (MPa) 26,600 24,500 25,300 26,800 28,600 30,600 G = G =G (MPa) 16,200 17,600 12,900 13,800 14,800 16,200 12 13 23 μ 0.491 0.562 0.421 0.424 0.427 0.43 bending and torsional stiffness of each section are analysed when 13 different cross sections of wing structure are selected. Bending and torsional stiffness can be expressed as K ¼ EI (4) bending Ω Gd K ¼ (5) torsional where E is the Young’s modulus, I is the inertia moment relative to the main inertia along the chordwise wing, G is the shear modulus, Ω is two times of the cross-sectional area enclosed, d is the effective thickness of cross section and b is the width of cross section. Mathematical and Computer Modelling of Dynamical Systems 63 Bending and torsional stiffness of cross sections are shown in Figures 3 and 4, respectively. It can be seen that both bending and torsional stiffness decrease gradually from wing root to the tip. Sizes of skin, stringer and beam decrease from wing root to the tip, which lead to this trend in terms of stiffness change. The trends in term of bending and torsional stiffness are not exactly coherent, which is due to the different influencing factors of bending and torsional stiffness. Comparing bending stiffness with torsional stiffness, it can be seen that the value of bending stiffness is obviously larger than that of torsional stiffness at the same cross section. 3. Aeroelastic characteristics analysis Aeroelasticity is always a major concern in aircraft design, which includes two aspects such as static aeroelasticity and flutter. Here, an effective CFD/CSD coupled aeroelastic analysis method was used, in which the solving processes of aerodynamic equations and structural Figure 3. Bending stiffness of cross sections. Figure 4. Torsional stiffness of cross sections. 64 F. Wang et al. equations are independent, but additional information transfer must be considered between the analytical tools in terms of load and deformation. The CFD software FLUENT (ANSYS Inc., Canonsburg, PA, USA) is adopted for aerodynamic analysis and FE software NASTRAN is adopted for structural analysis. The Navier–Stokes equations are used in modelling the fluid dynamics aspect in the aerodynamic model. The computation zone is selected according to the wing dimension and generally 10–15 times of wing dimension. Here, the computation zone is selected as 80,000 × 40,000 × 40,000 mm. The whole computation zone and wing surface grids for aerodynamic analysis are shown in Figures 5 and 6, respectively. Grids near the wing model are dense while those far from it are sparse. There are 109, 255 nodes and 622, 915 elements within the aerodynamic model, in which unstructured grids and tetrahedron elements are adopted for the aerodynamic analysis. There are 6439 nodes and 9542 elements within the structural model, in which shell elements and uniform three-node or four-node grids are adopted for structural analysis. The pressure distribution on the coupled surface can be obtained through aerodynamic analysis and the load for structural analysis needs to be applied directly to structural nodes or elements. Locations of nodes and elements for both models are different and informa- tion transfer between two models is necessary. First, pressure distribution of the wing surface is extracted from the FLUENT computational results. And then forces of CFD nodes can be obtained through the following equation. 0 1 0 1 G n x x @ A @ A G ¼ S  P  n (6) y y G n z z Figure 5. The whole computation zone. Figure 6. Wing surface grids. Mathematical and Computer Modelling of Dynamical Systems 65 Figure 7. Load information transfer between CFD and CSD nodes. where G , G and G are forces in x, y and z directions, respectively. S is the applied area x y z of pressure P. The quantities n , n and n are normal vectors of the coupled surface in x, y x y z and z directions, respectively. Load information transfer between CFD and CSD nodes is shown in Figure 7. For each CFD node, its projection node on the triangle decided by three CSD nodes can be obtained. In general, the distance between the CFD projection node and each CSD node is selected to be as small as possible. Because each CSD node is associated with several CFD nodes, so the ultimate load information of this CSD node is the combination of information transfer results for all relative CFD nodes. It can be expressed as F ¼ F (7) i¼1 where F is the ultimate load of each CSD node, F is transfer load of the ith CFD node and can be obtained through the following equation. 0 1 0 10 1 F A =A 00 G x;i i x @ A @ A@ A ¼ 0 A =A 0 G (8) y;i i y 00 A =A G z;i i z where F , F and F are loads applied on wing structure in x, y and z directions, x;i y;i z;i respectively. A =A is the area coordinate of the ith CSD node. Flow field meshes for the next aerodynamic analysis are generated based on the deformed wing configuration and a mesh regeneration method is adopted. Model building and mesh generating are carried out through the preprocessing software GAMBIT (ANSYS Inc.) for the initial CFD model. Seven typical lines are selected based on the change of wing configuration. A number of the controlled point coordinates are obtained based on initial coordinates and deformation. Wing configuration is regenerated by fitting these typical lines. Aerodynamic mesh regeneration is also carried out through GAMBIT software. The main computation processes of static aeroelasticity and flutter are almost the same. The only difference is that structural and aerodynamic analysis is independent for static aeroelastic analysis while structural analysis is inserted into aerodynamic analysis for flutter analysis. 66 F. Wang et al. 3.1. Static aeroelastic analysis The analysis process of static aeroelasticity is shown in Figure 8. Iterated for about 4–20 times, both aerodynamics and deformation is generally convergent. At this stage, the aerodynamic distribution represents the true condition when the aircraft is in the cruise state and deformation is the real jig-shape of aircraft wing at this flight attitude. Static aeroelasticity of the high-aspect ratio wing is analysed at the 13,000 m flight attitude, 0.54 Ma and 1.81° attack angle. The ultimate deformation of high-aspect ratio wing is shown in Figure 9. It can be seen that the maximal deformation value is located on trailing edge point of wing tip and is 0.313 m. The deformation value gradually increases from wing root to the tip. The load distribution for the aircraft wing will be regenerated under the effect of static aeroelasticity. For both rigid and flexible wings, pressure coefficient (Cp) distributions for 60% and 90% spanwise locations are shown in Figures 10 and 11, respectively. It can be seen that the pressure difference of the rigid wing is larger than that of the flexible one. Figure 8. Analysis process of static aeroelasticity. Figure 9. Ultimate deformation of the high-aspect ratio wing. Mathematical and Computer Modelling of Dynamical Systems 67 Figure 10. Cp distribution of 60% spanwise location. Figure 11. Cp distribution of 90% spanwise location. 68 F. Wang et al. The load is redistributed under the effect of static aeroelasticity, which leads to a decrease of pressure difference between the upper and lower airfoils. 3.2. Flutter analysis Futter-analysis process is shown in Figure 12. The user-defined function of FLUENT software is used for load and deformation interchange between aerodynamic and structural models. First, modal analysis of the structural model is carried out through the solving module SOL 103 of the NASTRAN software. The first three modal shapes of the high- aspect ratio wing are shown in Figure 13. It can be seen that the first two modal shapes exhibit two kinds of bending deformation while the third modal shape exhibits a kind of oscillation deformation. Here, the fourth modal shape mainly involves the aileron and not considered for further dynamic analysis. In modal space, dynamical aeroelasticity can be expressed as the following state equation [16]. M €qðtÞþ Cq _ðtÞþ KqðtÞ¼ FðtÞ (9) where M, C, K, qðtÞ and FðtÞ are generalized mass, generalized damper, generalized stiffness, generalized displacement and generalized load matrices, respectively. t is the time. M and K can be extracted directly from the FE software NASTRAN. While C cannot be extracted directly and are approximately expressed as C ¼ λMω (10) where ω is the natural frequency matrix. λ is the damper ratio selected between 0.01 and 0.02. The above state equation can also be expressed as Figure 12. Flutter analysis process. Mathematical and Computer Modelling of Dynamical Systems 69 Figure 13. The first three modal shapes of the high-aspect ratio wing. (a) The first modal shape; (b) The second modal shape; (c) The third modal shape. Y ¼ f ðt; YÞ¼ AY þ BFðtÞ (11) q_ðtÞ Y ¼ (12) qðtÞ 1 1 M C M K A ¼ (13) I 0 B ¼ (14) A constant step-length program of the third-order Runge–Kutta method is adopted to resolve the state equation and analyse time history of the generalized displacement. The following iterative equation is used to calculate Equation (11) from the initial value (t , Y ). 0 0 hðf þ 4f þ f Þ 1 2 3 Y ¼ Y þ (15) kþ1 k where f ¼ f ðt ; Y Þ 1 k k h h f ¼ft þ ; Y þ f (16) 2 k k 1 2 2 f ¼ f ðt þ h; Y  hf þ 2hf Þ 3 k k 1 2 The real displacement or deformation of wing structure can be expressed as 70 F. Wang et al. X ðtÞ¼ f qðtÞ (17) where X ðtÞ is the real displacement, f is the first m order modal matrix extracted directly from the FE software NASTRAN. Information transfer including load and deformation between both of different disciplines is inserted into the whole computational process. Flutter analysis is carried out at the initial 13 km flight altitude and 0° attack angle with an impulse load applied to the tip of the wing. In order to meet the requirement in terms of computational precision, the time step length is generally expressed using the following equation which has been derived from experience. Δt ¼ð0:01, 0:1ÞT (18) where Δt is the time step length and T is the period corresponding to the maximum frequency of all selected modal components. Here, the time step length is 0.008 s. Two cases with different flight speeds are studied, in which flight speeds are selected as 0.54 and 0.64 Ma, respectively. Time histories of the first three generalized displacements for both cases are shown in Figures 14 and 15, respectively. It can be seen that the first three generalized displacements are convergent for the flight speed of 0.54 Ma, while all of them show continuous oscillations of equal amplitude for flight speed of 0.64 Ma. Hence, the preliminary flutter characteristic can be obtained and flutter speed is near 0.64 Ma. Figure 14. Time histories of the first three generalized displacements with flight speed of 0.54 Ma. (a) The first generalized displacement X (m); (b) The second generalized displacement X (m); (c) 1 2 The third generalized displacement X (m). 3 Mathematical and Computer Modelling of Dynamical Systems 71 Figure 15. Time histories of the first three generalized displacements with flight speed of 0.64 Ma. (a) The first generalized displacement X (m); (b) The second generalized displacement X (m); (c) 1 2 The third generalized displacement X (m). 4. Effects on the aeroelastic characteristics The stiffness distribution influences the structural characteristics directly, which can consequently change the aeroelastic characteristics. Here, the aeroelastic characteristics of wing structure with small and different stiffness distribution are discussed. Three kinds of models named as initial model, larger beam size model and smaller beam size model will be studied, in which web thicknesses are 3, 4 and 2 mm, respectively. The bending and torsional stiffness of cross sections are shown in Figures 16 and 17, respectively. It can be seen that both bending and torsional stiffness increases with a small increase of beam size. 4.1. Effects on static aeroelastic characteristics Static aeroelasticity of these three models with different beam sizes is analysed by the above CFD/CSD coupling method. Effects of stiffness distribution on static aeroelastic characteristics are discussed through analysing deformation and load redistribution of aircraft wing. The maximum deformation values and load redistribution results of the three models are given in Table 5. The stiffness distribution influences the static aero- elastic characteristics, especially for wing root stiffness. To weaken the effect of static aeroelasticity, an appropriate way is to reinforce the wing root stiffness. The lift coefficient 72 F. Wang et al. Figure 16. Bending stiffness of cross sections. Figure 17. Torsional stiffness of cross sections. Table 5. The maximum deformation and load redistribution results for the three models considered. The maximum C of rigid C of flexible Lift l l Model deformation (mm) wing wing loss (%) Initial model 313 0.273 0.259 5.1 Larger beam size model 308 0.273 0.261 4.4 Smaller beam size model 332 0.273 0.257 5.9 Mathematical and Computer Modelling of Dynamical Systems 73 C of the flexible wing decreases with a small increase of beam size. The lift loss percent also has the same trend, which aggravates the effect of static aeroelasticity. 4.2. Effects on flutter characteristics Modal analysis of these three models with different beam sizes is carried out. Frequency values of different models are given in Table 6. It can be seen that effects of stiffness distribution on the frequency are not obvious. The first-order frequency decreases when beam size is changed. The second-order frequency increases with a small increase of beam size while it has an opposite trend for the third-order frequency. Table 6. Frequency values for the three models considered. The first-order The second-order The third-order Model frequency (Hz) frequency (Hz) frequency (Hz) Initial model 1.051 2.943 3.672 Larger beam size model 1.037 2.957 3.657 Smaller beam size model 1.024 2.893 3.673 Figure 18. Time histories of the first three generalized displacements for the larger beam size model. (a) The first generalized displacement X (m); (b) The second generalized displacement X (m); (c) The third generalized displacement X (m). 2 3 74 F. Wang et al. Figure 19. Time histories of the first three generalized displacements for the smaller beam size model. (a) The first generalized displacement X (m); (b) The second generalized displacement X (m); (c) The third generalized displacement X (m). 2 3 The flutter speed of the initial model is about 0.64 Ma and it is selected as an initial setting speed to discuss effects of stiffness distribution on flutter characteristic. Time histories of the first three generalized displacements for the larger beam size model are shown in Figure 18 when flight speed is 0.64 Ma. Comparing Figure 18 with Figure 15,it can be seen that effects of beam size increase on flutter characteristic are not obvious. The wing structure of the second model is still in a condition of continuous oscillation with equal amplitude when flight speed is 0.64 Ma. Flutter speed remains almost at the initial value when the beam size is slightly increased. Time histories of the first three generalized displacements for the smaller beam-size model are shown in Figure 19 when flight speed is 0.64 Ma. Comparing Figure 19 with Figure 15, it can be seen that effects of the beam size decrease on flutter characteristic are also not obvious. The wing structure of the third model is still in a condition of continuous oscillation with equal amplitude when the flight speed is 0.64 Ma. Flutter speed remains almost at the initial value when the beam size is slightly decreased. 5. Conclusions The following useful conclusions can be drawn through combing an effective computer modelling of honeycomb core sandwich composite and a highly efficient CFD/CSD Mathematical and Computer Modelling of Dynamical Systems 75 coupled numerical simulation to study stiffness distribution and its effects on the aero- elastic characteristics of an aircraft composite wing with high aspect ratio, in which wing skin is made of honeycomb core sandwich composite. (1) The aerodynamic load is redistributed under the effect of static aeroelasticity, which leads to a decrease of pressure difference between the upper and lower airfoils. The preliminary flutter characteristic was obtained and flutter speed was near 0.64 Ma. (2) Both bending and torsional stiffness increases with a small increase of beam size. The stiffness distribution influences the static aeroelastic characteristics directly, especially for wing root stiffness. To weaken the effect of static aeroelasticity, an appropriate way is to reinforce wing root stiffness. Lift coefficient and lift loss percent of the flexible wing decreases with a small increase of beam size. (3) Effects of the stiffness distribution on frequency are not obvious. The first-order frequency decreases when beam size is changed. The second-order frequency increases with a small increase of beam size while it has an opposite trend for the third-order frequency. Flutter speed remains almost at the initial value when the beam size is changed. Funding This study is supported by the National Natural Science Foundation (No: 51210008), the 111 Project (No: B07050), Aviation Science Foundation (No: 2013ZF53068) and the Basic Research Foundation of Northwestern Polytechnical University (No: JC20110257). References [1] H. Haddadpour, M.A. Kouchakzadeh, and F. Shadmehri, Aeroelastic instability of aircraft composite wings in an incompressible flow, Compos. Struct. 83 (2008), pp. 93–99. doi:10.1016/j.compstruct.2007.04.012 [2] L. Demasi and E. Livne, Aeroelastic coupling of geometrically nonlinear structures and linear unsteady aerodynamics: two formulations, J. Fluid Struct. 25 (5) (2009), pp. 918–935. doi:10.1016/j.jfluidstructs.2009.03.001 [3] H.X. Dang, Z.C. Yang, and Y. Li, Accelerated loosely-coupled CFD/CSD method for non- linear static aeroelasticity analysis, Aerosp. Sci. Technol. 14 (4) (2010), pp. 250–258. doi:10.1016/j.ast.2010.01.004 [4] A.K. Slone, K. Pericleous, C. Bailey, M. Cross, and C. Bennett, A finite volume unstructured mesh approach to dynamic fluid-structure interaction: An assessment of the challenge of predicting the onset of flutter, Appl. Math. Modeling. 28 (2004), pp. 211–239. doi:10.1016/ S0307-904X(03)00142-2 [5] G. Bartoli and C. Mannini, A simplified approach to bridge deck flutter, J. Wind. Eng. Ind. Aerod. 96 (2) (2008), pp. 229–256. doi:10.1016/j.jweia.2007.06.001 [6] D. Borglund and J. Kuttenkeuler, Active wing flutter suppression using a trailing edge flap,J. Fluid Struct. 16 (3) (2002), pp. 271–294. doi:10.1006/jfls.2001.0426 [7] S.J. Guo, Aeroelastic optimization of an aerobatic aircraft wing structure, Aerosp. Sci. Technol. 11 (2007), pp. 396–404. doi:10.1016/j.ast.2007.01.003 [8] T.J. Blom, Considerations on the spring analogy, J. Aircraft. 32 (2000), pp. 647–668. [9] C. Farhat, A. Rajasekharan, and B. Koobus, A dynamic variational multiscale method for large eddy simulations on unstructured meshes, Comput. Method. Appl. M. 195 (2006), pp. 1667– 1691. doi:10.1016/j.cma.2005.05.045 [10] S.H. Huo, F.S. Wang, and Z.F. Yue, Spring analogy method for generating of 2D unstructured dynamic meshes, J. Vib. Shock. 30 (10) (2011), pp. 177–182. [11] T.J. Baker and P.A. Cavallo, Dynamic adaptation for deforming tetrahedral meshes,14th AIAA Computational Fluid Dynamics Conference, Norfolk, VA, 28 June–1 July 1999. 76 F. Wang et al. [12] M. Mitsuhiro, K. Nakahashi, and K. Matsushima, Unstructured dynamic mesh for largemove- ment and deformation, 40th AIAA Aerospace Science Meeting and Exhibt, Reno, NV, 14–17 January 2002. [13] K. Stein, T.E. Tezduyar, and R. Benney, Automatic mesh update with the solid extension mesh moving technique, Comput. Method. Appl. M. 193 (2004), pp. 2019–2032. doi:10.1016/j. cma.2003.12.046 [14] S.H. Huo, F.S. Wang, W.Z. Yan, and Z.F. Yue, Layered elastic solid method for the generation of unstructured dynamic mesh, Finite. Elem. Anal. Des. 46 (2010), pp. 949–955. doi:10.1016/ j.finel.2010.06.006 [15] S.J. Xu, X.R. Kong, B.L. Wang, X.R. Ma, and X.C. Zhang, Method of equivalent analysis for static and dynamics behavior of orthotropic honeycomb sandwich plates, Acta Mater. Compos. Sinica. 17 (3) (2000), pp. 92–95. [16] F.S. Wang, S.H. Huo, Z.F. Yue, and W. Qian, Further study on aeroelastic response of aircraft wing with freeplay, J. Vib. Control. 18 (11) (2012), pp. 1690–1697. doi:10.1177/ http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Mathematical and Computer Modelling of Dynamical Systems Taylor & Francis

An effective computer modelling approach to the study of aeroelastic characteristics of an aircraft composite wing with high aspect ratio

Download PDF
19 pages

Loading next page...
 
/lp/taylor-francis/an-effective-computer-modelling-approach-to-the-study-of-aeroelastic-zTZbNFrtda

References (17)

Publisher
Taylor & Francis
Copyright
© 2014 Taylor & Francis
ISSN
1744-5051
eISSN
1387-3954
DOI
10.1080/13873954.2014.903283
Publisher site
See Article on Publisher Site

Abstract

Mathematical and Computer Modelling of Dynamical Systems, 2015 Vol. 21, No. 1, 58–76, http://dx.doi.org/10.1080/13873954.2014.903283 An effective computer modelling approach to the study of aeroelastic characteristics of an aircraft composite wing with high aspect ratio a b a a a Fusheng Wang *, Shihui Huo , Shengjun Qiao , Junran Zhang and Zhufeng Yue Advanced Materials Test Center, School of Mechanics, Civil Engineering and Architecture, Northwestern Polytechnical University, 710129 Xi’an, China; Mechanics & Environment Research Center of Rocket Engine, Xi’an Aerospace Propulsion Institute, 710100 Xi’an, China (Received 6 October 2013; accepted 7 March 2014) Static aeroelastic and flutter characteristics of an aircraft composite wing with high aspect ratio were analysed by an effective Computational Fluid Dynamics and Computational Structure Dynamics coupled method. Effects of stiffness distribution on aeroelastic characteristics were considered. Honeycomb core sandwich composite was considered to be equivalent to an orthotropic material by stiffness and inertance equivalent method to allow highly efficient numerical simulation, which was used for analysis of bending and torsional stiffness distribution. The results showed that the redistributed aerodynamic load leads to a decrease of pressure difference between the upper and lower airfoils. The flutter speed of the composite wing is near 0.64 Ma. Both bending and torsional stiffness increases with a small increase of beam size. Stiffness of the wing root has a major influence generally on the static aeroelastic characteristics. Both the lift coefficient and the loss percent decrease with a small increase of beam size. Effects of stiffness distribution on frequency are not obvious. Flutter speed remains close to the initial value when the beam size is changed. Keywords: honeycomb core sandwich composite; high-aspect ratio wing; stiffness distribution; static aeroelasticity; flutter; CFD/CSD coupled method 1. Introduction Composite airplanes with high-aspect ratio wings have received widespread attention in recent years throughout the world due to their high flight altitude and long endurance. A high-aspect ratio wing usually has lighter weight and more flexibility, which results in many special and complex aeroelasticity problems. The stiffness distribution needs to be considered in the design of a high-aspect ratio wing, which is a beneficial reference for further structural and aeroelastic analysis. Structural and aerodynamic characteristics are both affected significantly by the stiffness distribution of a high-aspect ratio wing. Many scholars have conducted research in this field and several aeroelastic analytical methods and stiffness distribution discus- sions have been put forward. Haddadpour et al. [1], Demasi and Livne [2] and Dang et al. [3] investigated the aeroelastic stability of an aircraft wing, which was modelled as an anisotropic composite thin-walled beam. The coupled method of Computational Fluid Dynamics (CFD) and Computational Structure Dynamics (CSD) is an effective approach to aeroelastic analysis, in which the structural and aerodynamic characteristics are *Corresponding author. Email: fswang@nwpu.edu.cn © 2014 Taylor & Francis Mathematical and Computer Modelling of Dynamical Systems 59 analysed independently followed by information transfer of load and deformation between the two analysis tools. Slone et al. [4] considered the application of a full physics simulation in the context of the AGARD 445.6 wing geometry. Bartoli and Mannini [5] investigated flutter instability of flexible bridge decks in framework of a semi-empirical Scanlan’s approach based on flutter derivatives. Borglund and Kuttenkeuler [6] investi- gated aeroservoelastic behaviour of a thin rectangular wing with a controllable trailing edge flap. Guo [7] aimed at presenting an investigation into optimal design about the minimum weight and aeroelastic tailoring of an aerobatic aircraft wing. Information transfer between different disciplines and dynamic mesh methods are concerned in the coupled CFD and CSD methods. Many scholars have studied this field and several methods have been put forward such as spring analogy method [8–10], elastic solid method [11–13] and layered elastic solid method [14]. Since stiffness distribution influences aeroelastic characteristics of a high-aspect ratio wing greatly, appropriate stiffness distribution corresponding to ideal aeroelastic charac- teristic can be obtained. In the present study, analysis on stiffness distribution and aeroelastic characteristic of a high-aspect ratio wing was carried out by a highly efficient computer modelling approach. And then effects of stiffness distribution on aeroelastic characteristics of the high-aspect ratio wing were considered. 2. Stiffness-distribution analysis Stiffness-distribution analysis is carried out based on an aircraft composite wing with high aspect ratio. A finite-element (FE) model of the high aspect ratio wing is shown in Figure 1, in which typical cross sections are illustrated, representing different rib locations with red colour and used to study stiffness distribution along the spanwise wing. Here, x is the spanwise coordinate, y is the chordwise coordinate, direction of z coordinate is vertical Figure 1. Finite element model of high-aspect ratio wing. 60 F. Wang et al. to the plane made of x and y coordinates. The wing structure is mainly composed of skin, longitudinal and transverse components. Longitudinal components include beams and stringers, while transverse components include ribs and thin-walled jointed structures. The elastic constraint of the wing root is simulated through an extension used for reinforcement support of root cross section. The upper and lower wing skins are made of honeycomb core sandwich composite. When the FE software NASTRAN (MSC Inc., Elkhart, IN, USA) is adopted for structural analysis, no relative element types are selected for the honeycomb sandwich plate. Although both three-dimensional FE model and laminated elements can be used to analyse honeycomb core sandwich composite, a three dimensional FE model will be costly in terms of the computation time, and honeycomb core sandwich composite cannot be simulated well by laminated elements. An effective modelling approach will allow the input parameter values for the FE software NASTRAN to be obtained (such as equivalent density, Young’ modulus, the shear modulus and Poisson’sratio), through which honeycomb core sandwich composite is equivalent to an orthotropic material. In the present study, the honeycomb core sandwich composite was regarded as equivalent to an orthotropic material by the stiffness and inertance equivalent method [15], the validity of which has been established by comparing results obtained from theory, those calculated by the equivalent method. It is found that differences between results from theory and from the equivalent method are less than 3%. This facilitates the numerical simulation process. The initial honeycomb core sandwich composite is shown in Figure 2(a). The total thickness of honeycomb sandwich plate is 2h, in which thickness of honeycomb structure is 2h and that of upper or lower surface is h . The equivalent orthotropic material of the honeycomb core sandwich composite is shown in Figure 2(b). Material properties can be obtained by a stiffness and inertance equivalence method expressed as 0 1 0 1 E ðe e  e Þ=e 11 22 22 1 12 B C B E C ðe e  e Þ=e 11 22 11 2 12 B C B C B C B G C e B C ¼ (1) B C B C B G C e B C @ A @ A e =e 12 22 Figure 2. Honeycomb core sandwich composite. (a) Initial honeycomb core sandwich composite; (b) Equivalent orthotropic material. Mathematical and Computer Modelling of Dynamical Systems 61 where 0  . 1 3 3 3 f 3 c ðh þ h Þ  h e þ h e ðh þ h Þ c f c f c 11 c 11 0 1 B C 1 B C 3 3 3 f 3 c B C ðh þ h Þ  h e þ h e ðh þ h Þ c f c f B C c 12 c 12 e B C B C  . B C B C 3 3 3 f 3 c B C ðh þ h Þ  h e þ h e ðh þ h Þ B C c f c f c 22 c 22 ¼ (2) B C B C B C f c B C B ðh e þ h e Þ ðh þ h Þ C f c c f 44 44 @ A B C f c B C ðh e þ h e Þ ðh þ h Þ f c c f 55 55 @  . A 3 3 3 f 3 c ðh þ h Þ  h e þ h e ðh þ h Þ c f c f c 66 c 66 f c In Equation (2), e and e are the stiffness coefficients of the surface and honeycomb ij ij structure, respectively. The density of the equivalent orthotropic material can be expressed as h ρ þ h ρ f c f c ρ ¼ (3) h þ h f c where ρ and ρ are the densities of the surface and honeycomb structures, respectively. f c Based on Equations (1) and (3), the honeycomb core sandwich composite can be equivalent to an orthotropic material. For different parts of the high-aspect ratio wing, stiffness coefficients of the honeycomb core sandwich composite and properties of the orthotropic material on the upper skin are given in Tables 1 and 2, respectively. Also, the stiffness coefficients of the honeycomb core sandwich composite and properties of the orthotropic material on lower skin are given in Tables 3 and 4, respectively. In Tables 2 and 4, E and E are Young’ moduli, respectively. G , G and G are the shear 11 22 12 13 23 moduli, respectively. μ is the Poisson’s ratio. The stiffness distribution is analysed based on different cross sections of the high- aspect ratio wing, which is selected according to rib positions. In the present study, Table 1. Stiffness coefficients of the honeycomb core sandwich composite on upper skin (MPa). Stiffness Ribs Ribs Ribs Ribs Ribs Ribs Ribs coefficient (1–2) (2–4) (4–5) (5–6) (6–7) (7–10) (10–13) e 23,550 21,830 15,105 21,848 21,913 19,167 22,317 e 7500 6599 3285 2834 2143 5191 5658 e 11,212 11,040 11,609 10,671 10,925 10,409 10,578 e 28,322 27,108 20,179 29,807 30,131 25,000 28,164 e 28,322 27,108 20,179 29,807 30,131 25,000 28,164 e 26,588 25,393 18,869 27,568 27,815 23,332 26,303 e 13,116 12,720 11,174 12,724 12,740 12,108 12,832 e 9225 9218 8456 8353 8193 8894 9002 e 10,279 10,239 10,370 10,154 10,213 10,094 10,133 e 11,309 11,222 10,727 11,415 11,438 11,071 11,297 e 11,309 11,222 10,727 11,415 11,438 11,071 11,297 e 13,815 13,540 12,039 14,040 14,096 13,066 13,749 66 62 F. Wang et al. Table 2. Properties of the orthotropic material on upper skin. Ribs Ribs Ribs Ribs Ribs Ribs Ribs Properties (1–2) (2–4) (4–5) (5–6) (6–7) (7–10) (10–13) E (MPa) 58,600 61,100 62,400 58,600 69,100 65,000 53,700 E (MPa) 24,600 25,500 33,700 31,500 25,300 26,800 30,600 G = G = G (MPa) 16,200 15,300 12,600 14,500 12,900 13,800 16,200 12 13 23 μ 0.529 0.487 0.313 0.38 0.421 0.424 0.43 Table 3. Stiffness coefficients of the honeycomb core sandwich composite on lower skin (MPa). Stiffness coefficient Ribs (1–4) Ribs (4–5) Ribs (5–6) Ribs (6–7) Ribs (7–9) Ribs (9–13) e 23,289 26,040 21,913 19,167 22,860 22,317 e 6517 7934 2143 5191 3305 5658 e 11,216 11,062 10,925 10,409 10,615 10,578 e 28,491 30,214 30,131 25,000 30,120 28,164 e 28,491 30,214 30,131 25,000 30,120 28,164 e 26,682 28,117 27,816 23,332 27,906 26,304 e 13,056 13,689 12,739 12,108 12,957 12,833 e 9180 8974 8193 8894 8462 9002 e 10,280 10,244 10,213 10,094 10,142 10,133 e 11,321 11,444 11,434 11,071 11,437 11,297 e 11,321 11,444 11,434 11,071 11,437 11,297 e 13,836 14,166 14,096 13,066 14,118 13,749 Table 4. Properties of the orthotropic material on lower skin. Properties Ribs (1–4) Ribs (4–5) Ribs (5–6) Ribs (6–7) Ribs (7–9) Ribs (9–13) E (MPa) 57,100 53,900 69,100 65,000 59,900 53,700 E (MPa) 26,600 24,500 25,300 26,800 28,600 30,600 G = G =G (MPa) 16,200 17,600 12,900 13,800 14,800 16,200 12 13 23 μ 0.491 0.562 0.421 0.424 0.427 0.43 bending and torsional stiffness of each section are analysed when 13 different cross sections of wing structure are selected. Bending and torsional stiffness can be expressed as K ¼ EI (4) bending Ω Gd K ¼ (5) torsional where E is the Young’s modulus, I is the inertia moment relative to the main inertia along the chordwise wing, G is the shear modulus, Ω is two times of the cross-sectional area enclosed, d is the effective thickness of cross section and b is the width of cross section. Mathematical and Computer Modelling of Dynamical Systems 63 Bending and torsional stiffness of cross sections are shown in Figures 3 and 4, respectively. It can be seen that both bending and torsional stiffness decrease gradually from wing root to the tip. Sizes of skin, stringer and beam decrease from wing root to the tip, which lead to this trend in terms of stiffness change. The trends in term of bending and torsional stiffness are not exactly coherent, which is due to the different influencing factors of bending and torsional stiffness. Comparing bending stiffness with torsional stiffness, it can be seen that the value of bending stiffness is obviously larger than that of torsional stiffness at the same cross section. 3. Aeroelastic characteristics analysis Aeroelasticity is always a major concern in aircraft design, which includes two aspects such as static aeroelasticity and flutter. Here, an effective CFD/CSD coupled aeroelastic analysis method was used, in which the solving processes of aerodynamic equations and structural Figure 3. Bending stiffness of cross sections. Figure 4. Torsional stiffness of cross sections. 64 F. Wang et al. equations are independent, but additional information transfer must be considered between the analytical tools in terms of load and deformation. The CFD software FLUENT (ANSYS Inc., Canonsburg, PA, USA) is adopted for aerodynamic analysis and FE software NASTRAN is adopted for structural analysis. The Navier–Stokes equations are used in modelling the fluid dynamics aspect in the aerodynamic model. The computation zone is selected according to the wing dimension and generally 10–15 times of wing dimension. Here, the computation zone is selected as 80,000 × 40,000 × 40,000 mm. The whole computation zone and wing surface grids for aerodynamic analysis are shown in Figures 5 and 6, respectively. Grids near the wing model are dense while those far from it are sparse. There are 109, 255 nodes and 622, 915 elements within the aerodynamic model, in which unstructured grids and tetrahedron elements are adopted for the aerodynamic analysis. There are 6439 nodes and 9542 elements within the structural model, in which shell elements and uniform three-node or four-node grids are adopted for structural analysis. The pressure distribution on the coupled surface can be obtained through aerodynamic analysis and the load for structural analysis needs to be applied directly to structural nodes or elements. Locations of nodes and elements for both models are different and informa- tion transfer between two models is necessary. First, pressure distribution of the wing surface is extracted from the FLUENT computational results. And then forces of CFD nodes can be obtained through the following equation. 0 1 0 1 G n x x @ A @ A G ¼ S  P  n (6) y y G n z z Figure 5. The whole computation zone. Figure 6. Wing surface grids. Mathematical and Computer Modelling of Dynamical Systems 65 Figure 7. Load information transfer between CFD and CSD nodes. where G , G and G are forces in x, y and z directions, respectively. S is the applied area x y z of pressure P. The quantities n , n and n are normal vectors of the coupled surface in x, y x y z and z directions, respectively. Load information transfer between CFD and CSD nodes is shown in Figure 7. For each CFD node, its projection node on the triangle decided by three CSD nodes can be obtained. In general, the distance between the CFD projection node and each CSD node is selected to be as small as possible. Because each CSD node is associated with several CFD nodes, so the ultimate load information of this CSD node is the combination of information transfer results for all relative CFD nodes. It can be expressed as F ¼ F (7) i¼1 where F is the ultimate load of each CSD node, F is transfer load of the ith CFD node and can be obtained through the following equation. 0 1 0 10 1 F A =A 00 G x;i i x @ A @ A@ A ¼ 0 A =A 0 G (8) y;i i y 00 A =A G z;i i z where F , F and F are loads applied on wing structure in x, y and z directions, x;i y;i z;i respectively. A =A is the area coordinate of the ith CSD node. Flow field meshes for the next aerodynamic analysis are generated based on the deformed wing configuration and a mesh regeneration method is adopted. Model building and mesh generating are carried out through the preprocessing software GAMBIT (ANSYS Inc.) for the initial CFD model. Seven typical lines are selected based on the change of wing configuration. A number of the controlled point coordinates are obtained based on initial coordinates and deformation. Wing configuration is regenerated by fitting these typical lines. Aerodynamic mesh regeneration is also carried out through GAMBIT software. The main computation processes of static aeroelasticity and flutter are almost the same. The only difference is that structural and aerodynamic analysis is independent for static aeroelastic analysis while structural analysis is inserted into aerodynamic analysis for flutter analysis. 66 F. Wang et al. 3.1. Static aeroelastic analysis The analysis process of static aeroelasticity is shown in Figure 8. Iterated for about 4–20 times, both aerodynamics and deformation is generally convergent. At this stage, the aerodynamic distribution represents the true condition when the aircraft is in the cruise state and deformation is the real jig-shape of aircraft wing at this flight attitude. Static aeroelasticity of the high-aspect ratio wing is analysed at the 13,000 m flight attitude, 0.54 Ma and 1.81° attack angle. The ultimate deformation of high-aspect ratio wing is shown in Figure 9. It can be seen that the maximal deformation value is located on trailing edge point of wing tip and is 0.313 m. The deformation value gradually increases from wing root to the tip. The load distribution for the aircraft wing will be regenerated under the effect of static aeroelasticity. For both rigid and flexible wings, pressure coefficient (Cp) distributions for 60% and 90% spanwise locations are shown in Figures 10 and 11, respectively. It can be seen that the pressure difference of the rigid wing is larger than that of the flexible one. Figure 8. Analysis process of static aeroelasticity. Figure 9. Ultimate deformation of the high-aspect ratio wing. Mathematical and Computer Modelling of Dynamical Systems 67 Figure 10. Cp distribution of 60% spanwise location. Figure 11. Cp distribution of 90% spanwise location. 68 F. Wang et al. The load is redistributed under the effect of static aeroelasticity, which leads to a decrease of pressure difference between the upper and lower airfoils. 3.2. Flutter analysis Futter-analysis process is shown in Figure 12. The user-defined function of FLUENT software is used for load and deformation interchange between aerodynamic and structural models. First, modal analysis of the structural model is carried out through the solving module SOL 103 of the NASTRAN software. The first three modal shapes of the high- aspect ratio wing are shown in Figure 13. It can be seen that the first two modal shapes exhibit two kinds of bending deformation while the third modal shape exhibits a kind of oscillation deformation. Here, the fourth modal shape mainly involves the aileron and not considered for further dynamic analysis. In modal space, dynamical aeroelasticity can be expressed as the following state equation [16]. M €qðtÞþ Cq _ðtÞþ KqðtÞ¼ FðtÞ (9) where M, C, K, qðtÞ and FðtÞ are generalized mass, generalized damper, generalized stiffness, generalized displacement and generalized load matrices, respectively. t is the time. M and K can be extracted directly from the FE software NASTRAN. While C cannot be extracted directly and are approximately expressed as C ¼ λMω (10) where ω is the natural frequency matrix. λ is the damper ratio selected between 0.01 and 0.02. The above state equation can also be expressed as Figure 12. Flutter analysis process. Mathematical and Computer Modelling of Dynamical Systems 69 Figure 13. The first three modal shapes of the high-aspect ratio wing. (a) The first modal shape; (b) The second modal shape; (c) The third modal shape. Y ¼ f ðt; YÞ¼ AY þ BFðtÞ (11) q_ðtÞ Y ¼ (12) qðtÞ 1 1 M C M K A ¼ (13) I 0 B ¼ (14) A constant step-length program of the third-order Runge–Kutta method is adopted to resolve the state equation and analyse time history of the generalized displacement. The following iterative equation is used to calculate Equation (11) from the initial value (t , Y ). 0 0 hðf þ 4f þ f Þ 1 2 3 Y ¼ Y þ (15) kþ1 k where f ¼ f ðt ; Y Þ 1 k k h h f ¼ft þ ; Y þ f (16) 2 k k 1 2 2 f ¼ f ðt þ h; Y  hf þ 2hf Þ 3 k k 1 2 The real displacement or deformation of wing structure can be expressed as 70 F. Wang et al. X ðtÞ¼ f qðtÞ (17) where X ðtÞ is the real displacement, f is the first m order modal matrix extracted directly from the FE software NASTRAN. Information transfer including load and deformation between both of different disciplines is inserted into the whole computational process. Flutter analysis is carried out at the initial 13 km flight altitude and 0° attack angle with an impulse load applied to the tip of the wing. In order to meet the requirement in terms of computational precision, the time step length is generally expressed using the following equation which has been derived from experience. Δt ¼ð0:01, 0:1ÞT (18) where Δt is the time step length and T is the period corresponding to the maximum frequency of all selected modal components. Here, the time step length is 0.008 s. Two cases with different flight speeds are studied, in which flight speeds are selected as 0.54 and 0.64 Ma, respectively. Time histories of the first three generalized displacements for both cases are shown in Figures 14 and 15, respectively. It can be seen that the first three generalized displacements are convergent for the flight speed of 0.54 Ma, while all of them show continuous oscillations of equal amplitude for flight speed of 0.64 Ma. Hence, the preliminary flutter characteristic can be obtained and flutter speed is near 0.64 Ma. Figure 14. Time histories of the first three generalized displacements with flight speed of 0.54 Ma. (a) The first generalized displacement X (m); (b) The second generalized displacement X (m); (c) 1 2 The third generalized displacement X (m). 3 Mathematical and Computer Modelling of Dynamical Systems 71 Figure 15. Time histories of the first three generalized displacements with flight speed of 0.64 Ma. (a) The first generalized displacement X (m); (b) The second generalized displacement X (m); (c) 1 2 The third generalized displacement X (m). 4. Effects on the aeroelastic characteristics The stiffness distribution influences the structural characteristics directly, which can consequently change the aeroelastic characteristics. Here, the aeroelastic characteristics of wing structure with small and different stiffness distribution are discussed. Three kinds of models named as initial model, larger beam size model and smaller beam size model will be studied, in which web thicknesses are 3, 4 and 2 mm, respectively. The bending and torsional stiffness of cross sections are shown in Figures 16 and 17, respectively. It can be seen that both bending and torsional stiffness increases with a small increase of beam size. 4.1. Effects on static aeroelastic characteristics Static aeroelasticity of these three models with different beam sizes is analysed by the above CFD/CSD coupling method. Effects of stiffness distribution on static aeroelastic characteristics are discussed through analysing deformation and load redistribution of aircraft wing. The maximum deformation values and load redistribution results of the three models are given in Table 5. The stiffness distribution influences the static aero- elastic characteristics, especially for wing root stiffness. To weaken the effect of static aeroelasticity, an appropriate way is to reinforce the wing root stiffness. The lift coefficient 72 F. Wang et al. Figure 16. Bending stiffness of cross sections. Figure 17. Torsional stiffness of cross sections. Table 5. The maximum deformation and load redistribution results for the three models considered. The maximum C of rigid C of flexible Lift l l Model deformation (mm) wing wing loss (%) Initial model 313 0.273 0.259 5.1 Larger beam size model 308 0.273 0.261 4.4 Smaller beam size model 332 0.273 0.257 5.9 Mathematical and Computer Modelling of Dynamical Systems 73 C of the flexible wing decreases with a small increase of beam size. The lift loss percent also has the same trend, which aggravates the effect of static aeroelasticity. 4.2. Effects on flutter characteristics Modal analysis of these three models with different beam sizes is carried out. Frequency values of different models are given in Table 6. It can be seen that effects of stiffness distribution on the frequency are not obvious. The first-order frequency decreases when beam size is changed. The second-order frequency increases with a small increase of beam size while it has an opposite trend for the third-order frequency. Table 6. Frequency values for the three models considered. The first-order The second-order The third-order Model frequency (Hz) frequency (Hz) frequency (Hz) Initial model 1.051 2.943 3.672 Larger beam size model 1.037 2.957 3.657 Smaller beam size model 1.024 2.893 3.673 Figure 18. Time histories of the first three generalized displacements for the larger beam size model. (a) The first generalized displacement X (m); (b) The second generalized displacement X (m); (c) The third generalized displacement X (m). 2 3 74 F. Wang et al. Figure 19. Time histories of the first three generalized displacements for the smaller beam size model. (a) The first generalized displacement X (m); (b) The second generalized displacement X (m); (c) The third generalized displacement X (m). 2 3 The flutter speed of the initial model is about 0.64 Ma and it is selected as an initial setting speed to discuss effects of stiffness distribution on flutter characteristic. Time histories of the first three generalized displacements for the larger beam size model are shown in Figure 18 when flight speed is 0.64 Ma. Comparing Figure 18 with Figure 15,it can be seen that effects of beam size increase on flutter characteristic are not obvious. The wing structure of the second model is still in a condition of continuous oscillation with equal amplitude when flight speed is 0.64 Ma. Flutter speed remains almost at the initial value when the beam size is slightly increased. Time histories of the first three generalized displacements for the smaller beam-size model are shown in Figure 19 when flight speed is 0.64 Ma. Comparing Figure 19 with Figure 15, it can be seen that effects of the beam size decrease on flutter characteristic are also not obvious. The wing structure of the third model is still in a condition of continuous oscillation with equal amplitude when the flight speed is 0.64 Ma. Flutter speed remains almost at the initial value when the beam size is slightly decreased. 5. Conclusions The following useful conclusions can be drawn through combing an effective computer modelling of honeycomb core sandwich composite and a highly efficient CFD/CSD Mathematical and Computer Modelling of Dynamical Systems 75 coupled numerical simulation to study stiffness distribution and its effects on the aero- elastic characteristics of an aircraft composite wing with high aspect ratio, in which wing skin is made of honeycomb core sandwich composite. (1) The aerodynamic load is redistributed under the effect of static aeroelasticity, which leads to a decrease of pressure difference between the upper and lower airfoils. The preliminary flutter characteristic was obtained and flutter speed was near 0.64 Ma. (2) Both bending and torsional stiffness increases with a small increase of beam size. The stiffness distribution influences the static aeroelastic characteristics directly, especially for wing root stiffness. To weaken the effect of static aeroelasticity, an appropriate way is to reinforce wing root stiffness. Lift coefficient and lift loss percent of the flexible wing decreases with a small increase of beam size. (3) Effects of the stiffness distribution on frequency are not obvious. The first-order frequency decreases when beam size is changed. The second-order frequency increases with a small increase of beam size while it has an opposite trend for the third-order frequency. Flutter speed remains almost at the initial value when the beam size is changed. Funding This study is supported by the National Natural Science Foundation (No: 51210008), the 111 Project (No: B07050), Aviation Science Foundation (No: 2013ZF53068) and the Basic Research Foundation of Northwestern Polytechnical University (No: JC20110257). References [1] H. Haddadpour, M.A. Kouchakzadeh, and F. Shadmehri, Aeroelastic instability of aircraft composite wings in an incompressible flow, Compos. Struct. 83 (2008), pp. 93–99. doi:10.1016/j.compstruct.2007.04.012 [2] L. Demasi and E. Livne, Aeroelastic coupling of geometrically nonlinear structures and linear unsteady aerodynamics: two formulations, J. Fluid Struct. 25 (5) (2009), pp. 918–935. doi:10.1016/j.jfluidstructs.2009.03.001 [3] H.X. Dang, Z.C. Yang, and Y. Li, Accelerated loosely-coupled CFD/CSD method for non- linear static aeroelasticity analysis, Aerosp. Sci. Technol. 14 (4) (2010), pp. 250–258. doi:10.1016/j.ast.2010.01.004 [4] A.K. Slone, K. Pericleous, C. Bailey, M. Cross, and C. Bennett, A finite volume unstructured mesh approach to dynamic fluid-structure interaction: An assessment of the challenge of predicting the onset of flutter, Appl. Math. Modeling. 28 (2004), pp. 211–239. doi:10.1016/ S0307-904X(03)00142-2 [5] G. Bartoli and C. Mannini, A simplified approach to bridge deck flutter, J. Wind. Eng. Ind. Aerod. 96 (2) (2008), pp. 229–256. doi:10.1016/j.jweia.2007.06.001 [6] D. Borglund and J. Kuttenkeuler, Active wing flutter suppression using a trailing edge flap,J. Fluid Struct. 16 (3) (2002), pp. 271–294. doi:10.1006/jfls.2001.0426 [7] S.J. Guo, Aeroelastic optimization of an aerobatic aircraft wing structure, Aerosp. Sci. Technol. 11 (2007), pp. 396–404. doi:10.1016/j.ast.2007.01.003 [8] T.J. Blom, Considerations on the spring analogy, J. Aircraft. 32 (2000), pp. 647–668. [9] C. Farhat, A. Rajasekharan, and B. Koobus, A dynamic variational multiscale method for large eddy simulations on unstructured meshes, Comput. Method. Appl. M. 195 (2006), pp. 1667– 1691. doi:10.1016/j.cma.2005.05.045 [10] S.H. Huo, F.S. Wang, and Z.F. Yue, Spring analogy method for generating of 2D unstructured dynamic meshes, J. Vib. Shock. 30 (10) (2011), pp. 177–182. [11] T.J. Baker and P.A. Cavallo, Dynamic adaptation for deforming tetrahedral meshes,14th AIAA Computational Fluid Dynamics Conference, Norfolk, VA, 28 June–1 July 1999. 76 F. Wang et al. [12] M. Mitsuhiro, K. Nakahashi, and K. Matsushima, Unstructured dynamic mesh for largemove- ment and deformation, 40th AIAA Aerospace Science Meeting and Exhibt, Reno, NV, 14–17 January 2002. [13] K. Stein, T.E. Tezduyar, and R. Benney, Automatic mesh update with the solid extension mesh moving technique, Comput. Method. Appl. M. 193 (2004), pp. 2019–2032. doi:10.1016/j. cma.2003.12.046 [14] S.H. Huo, F.S. Wang, W.Z. Yan, and Z.F. Yue, Layered elastic solid method for the generation of unstructured dynamic mesh, Finite. Elem. Anal. Des. 46 (2010), pp. 949–955. doi:10.1016/ j.finel.2010.06.006 [15] S.J. Xu, X.R. Kong, B.L. Wang, X.R. Ma, and X.C. Zhang, Method of equivalent analysis for static and dynamics behavior of orthotropic honeycomb sandwich plates, Acta Mater. Compos. Sinica. 17 (3) (2000), pp. 92–95. [16] F.S. Wang, S.H. Huo, Z.F. Yue, and W. Qian, Further study on aeroelastic response of aircraft wing with freeplay, J. Vib. Control. 18 (11) (2012), pp. 1690–1697. doi:10.1177/

Journal

Mathematical and Computer Modelling of Dynamical Systems – Taylor & Francis

Published: Jan 2, 2015

Keywords: honeycomb core sandwich composite; high-aspect ratio wing; stiffness distribution; static aeroelasticity; flutter; CFD/CSD coupled method

There are no references for this article.