Abstract
This paper describes a framework for an active control technique applied to gust load alleviation (GLA) of a flexible wing, including geometric nonlinearities. Nonlinear structure reduced order model (ROM) and nonplanar double-lattice method (DLM) are used for structural and aerodynamic modeling. The structural modeling method presented herein describes stiffness nonlinearities in polynomial formulation. Nonlinear stiffness can be derived by stepwise regression. Inertia terms are constant with linear approximation. Boundary conditions and kernel functions in the nonplanar DLM are determined by structural deformation to reflect a nonlinear effect. However, the governing equation is still linear. A state-space equation is established in a dynamic linearized system around the prescribed static equilibrium state after nonlinear static aeroelastic analysis. Gust response analysis can be conducted subsequently. For GLA analysis, a classic proportional-integral-derivative (PID) controller treats a servo as an actuator and acceleration as the feedback signal. Moreover, a wind tunnel test has been completed and the effectiveness of the control technology is validated. A remote-controlled (RC) model servo is chosen in the wind tunnel test. Numerical simulation results of gust response analysis reach agreement with test results. Furthermore, the control system gives GLA efficacy of vertical acceleration and root bending moment with the reduction rate being over 20%. The method described in this paper is suitable for gust response analysis and control strategy design for large flexible wings.
1. Introduction
With the development of high-altitude long-endurance (HALE) aircraft, the geometric nonlinear effects in aeroelastic analysis are becoming more and more significant, because they severely influence aerodynamic configuration and structural dynamic characteristics. Light weight and large flexibility have been necessary for the improvement of flight performance which causes a large deformation of the structure. The characteristics of both structural dynamics and aerodynamics need to be solved based on the deformed configuration as the geometric nonlinear aeroelasticity of a slender wing is analyzed. On the one hand, the lifting surface of aerodynamics is a curved surface. The geometrical relationship and boundary conditions in aerodynamic calculations are different from the linear method, which we call the nonplanar effect. On the other hand, the stiffness of the structure varies with a large deformation which is neglected in a linear condition. According to Patil et al., the concept of geometric nonlinear aeroelasticity was proposed in 1999 [1, 2]. The aeroelastic behavior turns out to be very different after considering nonlinear flexibility. Scores of research illustrate that geometric nonlinearities notably influence static aeroelastic configuration, flutter speed, and dynamic response. Patil et al. applied the exact intrinsic beam model in nonlinear aeroelastic analysis with linear state-space aerodynamics [3]. Tang and Dowell used the Hodge-Dowell equations in response to the problem of a high-aspect-ratio wing with the ONERA model of dynamic stall [4]. Frulla et al. investigated the dynamic characteristics of a flexible wing with a nonlinear beam model and an indicial function [5, 6]. The procedure took into account the large static deflection of the wing under aerodynamic loads. They point out that chordwise stiffness and elastic axis eccentricity parameters play an important role in nonlinear aeroelastic behavior. Abbas et al. performed the response analysis of HALE aircraft considering the imbalanced effects originating from the noncoincidence of the center of mass with the elastic axis [7]. Jian and Jinwu have studied the LCO, velocity, and drag effect with structural geometric nonlinearities on the static and dynamic aeroelastic characteristics [8]. Xie and Yang used nonlinear finite element method (FEM) and strip theory to simulate aeroelastic behavior of a metal wing under large deformation [9].
For large flexible aircraft, more attention needs to be paid on gust load due to the variation of stiffness with large deformation. Su and Cesnik studied the dynamic response of a highly flexible wing under a spatially-distributed discrete gust model in the time domain [10]. Guo et al. studied the nonlinear gust response of free flexible aircraft using a CFD/CSD method [11]. Wang et al. used the exact intrinsic beam model for discrete gust response simulation of very flexible aircraft with unsteady vortex lattice method [12]. Liu et al. established a theoretical geometrically nonlinear aeroelastic analysis framework and validated results with wind tunnel tests [13]. On the part of gust load alleviation (GLA), Dillsaver et al. performed linear-quadratic-Gaussian (LQG) GLA control of the rigid-elastic coupling in a flexible aircraft by reducing wing deflection by 47% [14]. Cook et al. conducted robust GLA control and stability analysis of a flexible wing. A comparison of closed-loop responses with the open-loop dynamics for both linearized and nonlinear systems was made to discrete gust distributions [15]. Wang et al. employed static output feedback to reduce the wing tip deflection of a flexible unmanned aerial vehicle by 33% [16]. Ricci et al. summarized the GLAMOUR project, which was a project that worked on GLA techniques for advanced regional aircraft. The wind tunnel model design and control strategy were presented as well [17].
Because of aerodynamic configuration and structural stiffness variation with deformation, large flexible structure modeling is vital in geometric aeroelastic analysis. Nonlinear FEM as a traditional modeling method has been taken widely in geometric nonlinear analysis with an iterative process [3, 18–20]. The intrinsic beam model is an exact model that describes the beam dynamics evaluated by Hodges [21]. It can be coupled with the finite-state inflow theory in a state-space form, which is convenient for the nonlinear analysis of dynamic behavior [22]. The multibody dynamic simulation is another tool for aeroelastic system analysis with arbitrary types of nonlinearities. Kruger et al. and Zhao and Ren researched the stability and control of the HALE configuration aircraft [23, 24]. The advantage of multibody simulation is the convenience for considering the interactions of rigid body motion and aeroelastic deformation [25]. While beam modeling may not be sufficient to capture complex structural details present in aircraft wings, nonlinear FEM may be very costly for dynamic simulations.
Modal analysis is the standard method in linear aeroelastic analysis. However, a large deformation brings the changing of stiffness and the deformation cannot be represented by the linear modal approach, especially the foreshortening effect. Geometric nonlinearities cannot be captured like the linear way. The reduced order model (ROM) based on the modal form was developed for resolving that problem. Mignolet and Soize showed that the nonlinear stiffness can be described in an equation with the quadratic and cubic terms of the basis function [26]. Based on that formulation, McEwan et al. performed the modal/FE method (MFE) by static analysis with numbers of specified static load cases. It is useful in the dynamic simulation of shells under intense acoustic excitation [27]. Hollkmap and Gordon improved that method to recover the in-plane deformation with the outer one [28]. Harmin and Cooper implied it for modeling the geometric nonlinearity of a large-aspect-ratio wing model [29]. Although the vertical deflection of the wing can be solved exactly, the foreshortening effect cannot be solved. Cestino et al. reduced the partial differential equations of the nonlinear beam model to a dimensionless form in terms of three ordinary differential equations with Galerkin’s method. Nonlinear flutter analysis is presented based on a nonlinear equilibrium condition [30]. Kim et al. used Proper Orthogonal Decomposition (POD) to capture the foreshortening effect. Spanwise deflection can be solved, but the foreshortening effect cannot be reflected in a nonlinear stiffness matrix directly [31].
The purpose of the present paper is to give an analysis framework of gust response analysis and GLA based on the dynamic linearization of a structural ROM. Nonlinear stiffness coefficients are described in a specified modal formulation. Except for linear modes, additional orthogonal spanwise modes are used to describe the foreshortening effects of the wing. Displacements under the follower forces as load cases are solved by nonlinear FEM. Unknown nonlinear stiffness coefficients are obtained from a curve-fitting procedure, and useless terms are rejected by stepwise regression. The method presented here achieves good agreement with those from the nonlinear FEM.
After the structural ROM was obtained, the dynamic linearization of ROM was implemented for stiffness description under a static equilibrium state. The unsteady nonplanar aerodynamics of the linearized system is coupled with a linearized ROM. A state-space model for a deformed wing has been established for gust response and GLA analysis. Feedback control is established for GLA, including the acceleration of wing tip and root bending moments. Significantly, a wind tunnel test is conducted with a remote-controlled (RC) model servo. The wing model and the wind tunnel test have also been introduced in a previous work [32]. Gust response results show that the method presented in this paper is feasible and the GLA strategy is effective in a numerical simulation and in a wind tunnel test.
2. Formulation
2.1. Structure of ROM
2.1.1. Nonlinear Structure Equation
The development of the structural ROM method is based on equations derived from the Galerkin approach to solve the geometric nonlinear dynamics in a weak form [26]. The equation of motion of a structure may often be given in dynamic equations as follows: where the tensor is the Second Piola-Kirchoff stress tensor. The tensor is the deformation gradient tensor, is the mass density of the structure, and is the volumetric force. is the position vector of the structure in the reference configuration, and is the deformed one so that the displacement vector is . The deformed gradient tensor can be defined as follows: where is the Kronecker symbol.
In order to approximate the solution with the Galerkin approach, a set of basis functions that describes the displacement can be written as follows: where satisfies all the boundary conditions and represents the component of the th basis function, and are the associated generalized coordinates.
Applying equation (3) to equation (1) which should be expressed in a weak form and proceeding with the Galerkin approach as shown by Mignolet and Soize [26], a third degree polynomial describes the nonlinear relationship of the dynamic equation corresponding to mode in terms of the generalized coordinates: where are the terms of the reduced mass matrix, is the modal component of the external force, and , , and are the components of tensors of reduced stiffness. The Einstein summation convention is applied in formulation. The analytical expression of the mass matrix and stiffness matrix can be achieved being directly related to the basis functions and relations of materials. But it is not practical, especially the FEM model which is complex. It should be noted that the geometric nonlinear effects are assumed not to impact the structural inertia, and therefore, the reduced mass matrix.
When a truncated basis of the linear modes is chosen as the basis function, and can be expressed as follows: where is the modal mass term of the th basis function and is the modal stiffness term of the th basis function. The formulation of nonlinear dynamic equations corresponding to mode can be
The stiffness is related to the deformation of structure. The tangent stiffness equation for the incremental forces and displacement can be written as follows: where is the vector of modal force and is the vector of modal coordinates. Considering equation (6), the terms of tangent stiffness matrix under an equilibrium deformation can be described as:
The mass matrix is assembled as equation (6) including the modal mass term:
It is generally a good approximation even for a very large structural deformation.
As a consequence, introducing the harmonic oscillating assumption, the eigenvalue problem of a structure under a stable deformation can be expressed as follows: where is the vibration circular frequency under the equilibrium deformation. The eigenvectors, also called quasimodes, can also be solved. They are useful in the stability analysis of geometric nonlinear aeroelasticity.
2.1.2. Foreshortening Effect under Large Deformation
Foreshortening effects under large deformations have a great influence in geometric nonlinear problems, especially aeroelastic problems, because the shortening effect of the spanwise projection of a wing will change the aerodynamic distribution. A truncated basis of the linear modes has been chosen as the basis function in our ROM research. A few linear modes that are truncated will not show the shortening effect of the wing in ROM analysis. Two orthogonal spanwise modes are taken into the establishment of ROM to describe the foreshortening effects of the wing. The combination of truncated linear modes and orthogonal spanwise modes is generated as a basis function in the nonlinear structure of ROM described in Section 2.1.1.
2.1.3. Regression Analysis for Solving the Nonlinear Stiffness Coefficient
As described before, the explicit calculation of nonlinear stiffness is not practical. Regression analysis is used to identify the nonlinear stiffness coefficients and . The static formulation of equation (6) is as follows:
Evidently, if there is a set of specified static loads and corresponding structural deformations, the unknown nonlinear stiffness terms related to the applied loads and the structural displacement resultant can be solved by using regression analysis. The set of specified loads and corresponding structural deformations can be denoted as the static test load case and calculated by a commercial FEM software package.
Consider that there aresets of static test load cases. We can get sets of corresponding displacement after completing static FEM analysis of the sets of test loads on the model in a commercial finite element software. Loads and deformations are translated into modal space. For the sake of finding the unknown nonlinear modal stiffness terms, the left side of equation (11) can be curve fitted in a regression progress. A regression problem can be presented as follows:
The Einstein summation convention is applied in equation (11) and equation (12). It should be noted that the superscripts without brackets denote the serial number of the test cases from 1 to instead of the power value. Simplify equation (12) and the regression problem can be presented in matrix short notation as follows:
Here, is a vector containing the modal restoring forces corresponding to mode in each of the test sets, and is a vector containing curve-fitted values of modal restoring forces with the caret indicating the fitted rather than the exact values. is a vector including all unknown stiffness terms corresponding to mode . is a matrix corresponding to mode denoted as the design matrix in the regression problem.
In our research, two bending modes, two twist modes, one chordwise bending mode, and two spanwise modes are used to establish nonlinear ROM. The regression problem should be solved seven times with those seven modes.
2.1.4. Strategy for Generating Test Load Cases
Through the abovementioned analysis, regression analysis is performed using the commercial software package on the actual deformation and load testing after FEM analysis; thus, the accuracy of the nonlinear stiffness coefficients directly depends on the rationality of the selected test load case, which is related to the success of the recovery of the nonlinear structure equation. The selection of test load cases must emphasize that the aerodynamic force on the wing should be a follower force, which more closely resembles the actual characteristics of the aerodynamic force. That is, taking an oriented load as the load test case cannot satisfy the requirements. In this paper, the aerodynamic force under the deformation-combined bending and torsion modes is chosen as the test load case. The formulation of the wing deformation, which generates aerodynamic forces, should be as follows:
Here, and represent the investigated bending and torsion modes, respectively, and represents the scalar mode weight factors, through which the selected test cases contain the nonlinear characteristics of the structure investigated in our research.
2.1.5. Selection of Nonlinear Stiffness Coefficients
The selection of nonlinear stiffness coefficients is the main problem in structural ROM establishment. In general, for the accuracy of ROM to be high, the nonlinear stiffness coefficients should have enough flexibility to approximate the training points. However, if the number of coefficients is large, an overfitting problem may obviously occur. The number of nonlinear stiffness coefficients utilizing modes is as follows: where With 7 modes used in our research, the number of nonlinear stiffness coefficients is 112. The overfitting problem cannot be neglected in the nonlinear structural ROM method.
The stepwise regression procedure is a statistical procedure to select variables in a regression space. It is computationally efficient because it gives intermediate statistical information at each stage of the calculation which is used to select the most appropriate coefficients to be added into the model. The procedure of stepwise regression can be summarized as follows: (1)The regression analysis problem corresponding to equation (12) has the variables . At first, all variables are not introduced into the regression equation(2)Chose the term in unintroduced variables which has the maximum value of the partial regression quadratic sum as the introduced term for the regression equation, and do the significance test. If it passes the test, the term can be preserved in the equation. Otherwise, it will not be considered(3)Choose the term in the regression equation which has the minimal value of the partial regression quadratic sum, and do the significance test. If it passes the test, the term can be preserved in the equation. Otherwise, it will be removed(4)Repeat (2) and (3) above, until no terms can be introduced or removed
By the stepwise regression procedure, the best nonlinear ROM model will be determined. The joint hypotheses test (-test) is used for the specified criterion in the significance test and is expressed as follows [33]: where is the fitted value of the regression, is the observed value of the regression, and is the mean value of the observed values.
2.2. Unsteady Aerodynamics
In this paper, unsteady aerodynamics in gust response analysis and GLA analysis are calculated by the nonplanar double-lattice method (DLM) [32, 34]. Aerodynamics mesh and double-lattice element location should be determined on the deformed surface and updated with structure deflection as shown in Figure 1(a).
(a) Nonplanar lattice on a deformed lifting surface
(b) Typical doublet-lattice line and its normal vector
In traditional DLM, unsteady aerodynamic pressure at the lattice element is as follows: where is the density of airflow, is the air speed, is the aerodynamic coefficient matrix, and is the boundary condition expressed as follows: where is the mode shape matrix at the control points, is the mode shape slope along the inflow direction at the control point, is the general coordinate of modes, is the reference chord length, and is the reduced frequency.
Another critical problem of nonplanar DLM is the implementation of exact geometric boundary conditions. The local normal wash velocity is shown in Figure 1(b). Due to large deformations, the lifting surface cannot be treated as oscillating about the -plane. The boundary condition related to the local geometric nonlinearity can be written based on modes as follows: where is the nondimensional induced velocity in the normal direction of the local lattice and is the modal vibration shape function of the structure in a local normal direction and it should be different under various equilibrium states. According to the deformed configuration, is updated related to and described in equations (7)–(9) before.
2.3. Gust Model
The discrete gust model is often used in aircraft design, flight quality assessment, and control system design. It is a simple and effective gust modeling form. The “1-cos” gust model is chosen as the gust model in this paper. The indication of the “1-cos” gust model is shown in Figure 2. The gust section can be expressed as follows: where is the gust amplitude, is the location of aircraft, and is the gust scale.
2.4. Synthesized Modeling
Unsteady aerodynamic forces excited by structure vibration, deflection of the control surface, and gust in the frequency domain are calculated by nonplanar DLM, and we have where are the vibration-induced equivalent nodal loads; are the nodal loads caused by control surfaces; are the gust-induced equivalent nodal loads; and , , and are the aerodynamic influence coefficient matrices with respect to the structure of the quasimodes of the static equilibrium configuration, control surface modes, and gust modes, respectively, which are complex functions of reduced frequency. denote the modal coordinate vector, denote deflection vector of control surface, and denote the gust velocity vector.
In order to express the motion equation in a state-space form, the frequency domain unsteady aerodynamics should be described in the time domain. The time domain aerodynamics can be obtained with Karpel’s minimum-state rational function approximation [32]: where , is the Laplace variable, and is the reference chord length. , , and are the polynomial fitting coefficient matrices solved by a least square fitting, and is a diagonal matrix with aerodynamic lag roots.
By combining quasimodes around the equilibrium state with aerodynamic modeling and servo modeling, the aeroelastic equations of motion in generalized coordinates can be expressed as follows: where , , and are generalized unsteady aerodynamic influence coefficient matrices corresponding to the static equilibrium configuration modes, control surfaces, and gust in equation (20). is the mass matrix expressed with in equation (9), is the mass matrix coupled with structural modes and control surface modes, is the damping matrix, and is the tangent stiffness matrix expressed with in equation (7).
By applying the Laplace transformation to equation (22), the open-loop aeroelastic system with gust excitation can be expressed in the following state-space form: where , , , and . is the aerodynamic lag state. To eliminate the coefficient associated with the second-order derivative of the gust velocity, we use
The output can be expressed as follows:
The procedure of synthesized modeling can be seen in Figure 3. Nonlinear static aeroelastic analysis is implemented first. Then, the static equilibrium state would be obtained using the structural ROM and nonplanar vortex lattice method (VLM) [35]. Under this deformed configuration, quasimodes around the equilibrium state are solved by the linearization of the structural ROM. The state-space equation for gust analysis can be generated by coupling with quasimodes and unsteady aerodynamic influence coefficients.
2.5. Active Control Strategy
Based on the state-space model established in Section 2.4, the control law to alleviate gust loads would be designed in the time domain. The diagram of the control system of the GLA of a flexible wing is shown in Figure 4.
An aileron with an RC model servo is settled. A servo can be activated by a pwm signal generated from a microcontrol unit (MCU) through A/pwm translation from the control computer. The vertical acceleration measured is fed back to the servo as a signal with a proportional-integral-derivative (PID) control algorithm. On the feedback route, the vertical acceleration passes through a first-order low-pass filter which improves the quality of the sensor signal and then is sent to a PID controller to obtain an external signal used to alleviate gust response loads. One degree deflection of the servo corresponds to one volt of control voltage.
The reduction rate is given as follows: where is the open-loop response peak value and is the close-loop response peak value.
3. Wing Model and Wind Tunnel Test System
Wind tunnel tests are carried out to study the gust response of a flexible wing and to verify the active gust load alleviation techniques. All tests of the flexible wing model are performed in an FD-09 low-speed wind tunnel which is 12 m long, 3 m high, and 3 m wide in a close test environment. Velocity from 10 m/s to 100 m/s can be performed in the wind tunnel, and the test velocity needed is from 25 m/s to 45 m/s. The maximum test velocity of 45 m/s corresponds to a Reynolds number per meters of about 3.1 million. The turbulence coefficient of the wind tunnel is lower than 0.1%.
3.1. Wing Model
The flexible wing model designed in the wind tunnel test is shown in Figure 5 and Figure 6. Its design parameters are shown in Table 1.
(a) 3-D CATIA model
(b) FEM model
A gradually varied cross section has been chosen for the beam simulation of the wing. The beam is located at the 40% chord line of the wing, and the dimensions of the cross section decreases from the root to the tip of the wing. The dimensions of the root section and tip section are given in Table 2. The density of the beam materials is kg/m3, and the module is 70 Gpa. The wing frame includes 11 light wood boxes and the gap between the boxes is about 2 mm. The box is connected to the main beam by a single point, and the skin of the model is tissue paper just as a shape contour. As a result, the influence of the stiffness coming from the wing box has been reduced to the lowest.
The aileron prepared for the analysis and wind tunnel test is located 804.5 mm from the root of the wing. The area of the control surface is 16,800 mm2, and the average chord is 40 mm.
The main linear modes of the wing model before ROM is presented in Table 3. It can be seen that the frequencies of the first two modes are low, which means that the flexibility of the model is high.
3.2. Design of Wind Tunnel Test Subsystem
3.2.1. Support System
The wing model is vertically fixed in the wind tunnel with a lamped root. The angle of attack (AOA) which is related to the undeformed wing model can be changed from 0 deg to 5 deg. At the joint between the support system and the wing model, a force balance measures the force and moments in six directions.
3.2.2. Gust Generator
A gust generator is applied to create the expected gust disturbance during the tunnel test. Two rectangular blades deflect sinusoidally and synchronously. The span length of the blades is 2000 mm and the chord length of the blades is 300 mm. The distance between two blades is 600 mm, and the distance from the model and the generator is 2000 mm. The airfoil of the blades is NACA0015. The deflecting angle range of the blades is -7 deg to 7 deg, and the deflecting frequency is 1-7 Hz. The gust generator in the wind tunnel is shown in Figure 7.
When the blades are deflected sinusoidally at a certain frequency, an approximately sinusoidal lateral gust can be generated in the test field, and the gust velocity can be written as follows: where is the amplitude of the blade deflecting angle, and is the gust disturbance coefficient which is relevant to the velocity and the blade deflecting frequency .
3.2.3. Measure-Control System
A measure-control system is required to contain GLA in a wind tunnel test, as shown in Figure 8. During the test, the GLA control loop can be switched on and off in order to validate the effect of the GLA system. Due to the particularity of using the RC model servo, the MCU needs to translate the analog signal and pwm signal. Figure 9 shows the parts of the hardware devices used in the wind tunnel test. The accelerometer we chose is an ADXL345 located at about 33% of the span from the wing tip.
(a) Accelerometer and servo
(b) MCU for servo
3.3. Frequency Sweep Test of Servo
Because the servo we chose cannot output its deflection directly, we use a laser rangefinder to measure the vertical displacement of a point on the rudder, then we calculate the angle of deflection of the servo from that vertical displacement. Figure 10 shows the frequency sweep test of the RC model servo [32].
The frequency range is 0-6 Hz in a sweep test, and a deflection of 14.4 deg is given to the RC model servo with a pwm signal. The time domain result is given in Figure 11. The black curve is a deflection signal input, and the red curve is a deflection signal output. It can be seen that the dynamic response characteristics of the RC model servo became worse with the increase of frequency. The reason is that the RC model servo is usually used for the static control of an RC model airplane. The frequency-domain result is given in Figure 12. The black curve describes the derivation of the deflection input and deflection output in gain and phase. The red curve with a circle is the fit curve of the transfer function. It can be seen that the lagging derivation of the phase can be nearly 150 deg when the frequency is close to 6 Hz.
(a) Gain
(b) Phase
The transfer function of the servo in an aileron is as follows: where is the Laplace variable.
4. Theoretical Analysis and Wind Tunnel Test Results
4.1. Validation of Nonlinear ROM
A nonlinear ROM has been obtained by regression analysis from test load cases chosen before and from corresponding deformations. The validation of the nonlinear ROM calculation from a structural equation must be established. Then, ROM can be applied to static and response analysis reasonably. The numbers of nonlinear stiffness coefficients are less than 10 by stepwise regression procedure. The regression procedure of the equation corresponding to the 1st bending mode as an example is presented in Table 4. Modal coordinates from to stand for the 1st bending mode, 2nd bending mode, 1st chordwise bending mode, 1st torsion mode, 2nd torsion mode, 1st spanwise mode, and 2nd spanwise mode in that order. There are 11 steps included in the operation, where 9 terms were introduced and 2 terms were removed with a decreasing root mean squared error (RMSE). The final equation has 9 nonlinear stiffness coefficients.
The validation procedure is illustrated in Figure 13. The nonlinear FEM analysis is completed in MSC.NASTRAN. The module SOL 106 is used for static nonlinear FEM analysis, and the module SOL 400 is used for dynamic nonlinear FEM analysis. The FEM model has 339 CBEAM elements, and it is presented in Figure 5.
4.1.1. Static Validation with Nonlinear FEM
In static validation of nonlinear ROM with the data from nonlinear FEM, we take 50 sets of a growing validation load as examples which all have a distribution of the actual aerodynamic load. The displacement of a wing tip under these sets of validation loads is presented in Figure 14(a). The choice of validation loads makes the vertical displacement of a wing tip within the range of 10% to 25%, which is in the demand of nonlinear analysis. Spanwise displacement of a wing tip is also considered. It is important in geometric nonlinear aeroelastic analysis, since it will change aerodynamic distribution.
(a) Displacement of a wing tip
(b) Relative deviation between FEM and ROM
The calculation results of those 50 cases are shown in Figure 14(b) along with a comparison of ROM solutions and FEM solutions. Here, the black line with a circle represents the relative deviation of a wing tip vertical displacement and the red line with a rectangle represents the relative deviation of a wing tip spanwise displacement. It can be seen that the relative deviation between ROM solutions and nonlinear FEM solutions under these validation loads are no greater than 1%. ROM solutions can reach good agreement with nonlinear FEM solutions.
4.1.2. Dynamic Response Validation with Nonlinear FEM
In the dynamic response verification of nonlinear ROM with the data from nonlinear FEM, two load cases are set for comparison.
Case 1. A 0.5 N follower force is applied at the tip of the beam in the vertical direction as a step input at the initial time.
Case 2. A set of follower forces has the distribution of an actual aerodynamic load as a step input at the initial time; the load can make a vertical deflection nearly 20% the length of the span.
A comparison of the FEM solution results and the ROM solution results is shown in Figure 15 and Figure 16. Nonlinear ROM solution results in both the vertical deflection and spanwise deflection under Case 1 and Case 2 coincide with nonlinear FEM solution results. The structural ROM is reliable.
(a) Vertical deflection
(b) Spanwise deflection
(a) Vertical deflection
(b) Spanwise deflection
4.2. Static Deformation and Flutter Analysis
Before the gust response and GLA analysis, nonlinear aeroelastic stable deformation and flutter speed are calculated. The nonlinear analysis used nonplanar VLM and nonplanar DLM coupled with structural ROM. As described before, the stiffness of the structure has a relationship with large deformations which will have a big influence on gust analysis. All the methods implemented are self-programmed.
Figure 17 shows the nonlinear static aeroelastic analysis results of a wing model under 4 deg AOA. The vertical deflection of a wing tip grows with the increase of velocity with a nonlinear trend. The vertical deflection of a wing tip reaches 615 mm, nearly 40% of the wingspan, which is a large deformation. Meanwhile, the spanwise deflections are comparable with the vertical deflections, which are not reflected in the linear analysis.
(a) Vertical deflection
(b) Spanwise deflection
The distributions of lift, drag, and side force under the static aeroelastic equilibrium state with 3 deg AOA and 4 deg AOA in 37 m/s velocity are presented in Figure 18. It is obvious that there are differences between considering and not considering spanwise displacement. The total lift has 10% and 7.5% overall errors for the 4 deg AOA condition and the 3 deg AOA condition, respectively, which is computed as follows: where is the total lift of the condition with spanwise displacement and is the total lift of the condition without spanwise displacement. The lifting surfaces are shown in Figure 19. The red one is the lifting surface with spanwise displacement and the blue one is the lifting surface without spanwise displacement.
(a) Lift
(b) Drag
(c) Side force
Quasimode frequencies change a lot with growing deformation. Figure 20 shows the changing of the first six order mode frequencies with velocity variations from 25 m/s to 43 m/s under 2 deg AOA. It can be seen that the first five mode frequencies decrease as the velocity increases; meanwhile, the deformations of the wing also increase. The 2nd order mode frequency decreases nearly 12%, which is the 1st chordwise bending mode in the linear analysis. In nonlinear aeroelastic analysis, the chordwise bending mode plays an important role caused by the geometric nonlinear effect as shown by Frulla et al. [5, 6].
In linear flutter analysis of this wing model, it is a typical bend twist coupling flutter with a flutter speed of 70.2 m/s. It is not interrelated with structure deformation. The nonlinear flutter speed has a correlation with the structure deformation described above. The flutter speed is 57.5 m/s in a nonlinear flutter analysis under 2 deg AOA in tune with the quasimode frequency analysis condition. Coupling modes have changed to the bending mode and the chordwise bending mode. The chordwise bending mode may be significant in nonlinear flutter analysis.
4.3. Gust Response Analysis and Test Results
Both the gust response analysis simulation results and the wind tunnel test results are presented and compared to illustrate the consistency. Acceleration responses in the test were obtained from an ADXL345 accelerometer. Root bending moments in the test were obtained from the force balance measures fixed on the support system. Figure 21 and Figure 22 show the comparison between the simulation and test results corresponding to the vertical acceleration and root bending moment. The results presented are peak values of frequency domain analysis. Most of the test results reach agreement with the simulation results. The instability of the sensors may bring out bad data in the wind tunnel test. However, the tendencies of the simulation results are obvious. Because the frequencies of the gusts in the test did not reach the first mode frequency of the wing model, the acceleration and root bending moment grew with the increase of velocity and frequency. It can be said that the theoretical method is valid.
4.4. Gust Load Alleviation Analysis and Test Results
Two sets of conditions have been given as examples to present GLA results under 2 deg AOA with velocities of 28 m/s and 31 m/s. Figure 23 shows the simulation results in the frequency domain under velocity conditions of 28 m/s and 31 m/s, including both open-loop and close-loop results. It can be seen that acceleration alleviation and root bending moment alleviation are valid under the low-frequency conditions with reduction rates of over 35%. As the frequency grows, reduction rates evidently decrease. When the frequency reaches 4.4 Hz, the control system is almost ineffective. The main reason is that the RC model servo we used has a large lagging of phase and decaying of amplitude. The same situation comes to the wind tunnel test results. Figure 24 shows the wind tunnel test results in the frequency domain under the velocity conditions of 28 m/s and 31 m/s. The control system is valid when the frequency is kept under low frequencies. But if the frequency increases, the effectiveness of control system is quickly lost. This will be more obvious in test results. Moreover, the signal of the small ADXL345 accelerometer would have heavy noises in the 2.8 Hz test condition which leads to the small rate of reduction. This situation would not occur in theoretical simulation calculations. Table 5 presents the alleviation efficacy with the reduction rate.
(a) Acceleration response
(b) Bending moment response
(a) Acceleration response
(b) Bending moment response
5. Conclusions
In this paper, an effective method for gust response analysis and GLA analysis of a flexible wing including geometric nonlinearities was performed. A nonlinear structural ROM and nonplanar aerodynamic calculation with the linearized system under a static equilibrium state is applied for gust response and alleviation solving.
The structural modeling method describes the geometric nonlinear stiffness coefficient in an equation with quadratic and cubic terms of a linear modal function. Two orthogonal spanwise modes are used to capture the foreshortening effects, and the spanwise deflection can be solved directly. The influence of the foreshortening effects in an aerodynamic calculation has been discussed. Stepwise regression ensures the accuracy of nonlinear stiffness identification, and follower force effect is also considered.
The linearized formulation of structural ROM is given, and a framework of gust response analysis and GLA is derived based on linearized structural ROM and nonplanar DLM. A wing model for the validation of this method was designed, and a wind tunnel test was carried out successfully. An acceleration feedback control strategy with PID controllers is used to realize alleviation. We used an RC model servo in our first attempt in the wind tunnel test. This would serve as a guide for applications of this kind or those using a servo in a wind tunnel test or a flight test.
Gust response analysis shows that most of the test results reach agreement with the simulation results. The acceleration and root bending moment grew with the increase of velocities and frequencies. The method described is valid. In GLA analysis and in the wind tunnel test, the reduction rate of vertical acceleration and root bending moment could be more than 20% with lower frequencies. When the frequency reaches a higher value, the alleviation efficacy will seriously decrease due to a defect of the dynamic response characteristics of the RC model servo. This is important to consider in similar tests or flight tests. Using this actuator gives the control and test foundation of our upcoming geometric nonlinear aeroelastic flight test with a lightweight airplane model which cannot afford a fin actuator. A small ADXL345 accelerometer is chosen for the same reason. The alleviation technology described in this paper is validated.
Data Availability
The simulation results and wind tunnel test results data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
The study was supported by the National Key Research and Development Program (2016YFB 0200703).