#### Abstract

A folding wing morphing aircraft should complete the folding and unfolding process of its wings while in flight. Calculating the hinge moments during the morphing process is a critical aspect of a folding wing design. Most previous studies on this problem have adopted steady-state or quasi-steady-state methods, which do not simulate the free-flying morphing process. In this study, we construct an aeroelastic flight simulation platform based on the secondary development of ADAMS software to simulate the flight-folding process of a folding wing aircraft. A flexible multibody dynamic model of the folding wing structure is established in ADAMS using modal neutral files, and the doublet lattice method is developed to generate aerodynamic influence coefficient matrices that are suitable for the flight-folding process. The user subroutine is utilized, aerodynamic loading is realized in ADAMS, and an aeroelastic flight simulation platform of a folding wing aircraft is built. On the basis of this platform, the flight-folding process of the aircraft is simulated, the hinge moments of the folding wings are calculated, and the influences of the folding rate and the aircraft’s center of gravity (c.g.) position on the results are investigated. Results show that the steady-state method is applicable to the slow folding process. For the fast folding process, the steady-state simulation errors of the hinge moments are substantially large, and a transient method is required to simulate the flight-folding process. In addition, the c.g. position considerably affects the hinge moments during the folding process. Given that the c.g. position moves aft, the maximum hinge moments of the inner and outer wings constantly increase.

#### 1. Introduction

As a typical morphing aircraft, the folding wing aircraft can change its configuration autonomously by folding its wings during flight to fit variable flying environments and satisfy multimission demands. Its operational performance is much higher than that of a conventional fixed configuration aircraft [1, 2] and thus has elicited widespread attention.

Many studies on a folding wing aircraft have been conducted in recent years, and the focus of research includes two aspects. One is the aeroelastic characteristic in fixed configuration. Related works have undergone a phase from linear to nonlinear and from simulation to experiment [3–5]. This aspect has been the subject of considerable research. The other is the dynamic response during the folding process. Related studies are based on simulation and require further development. Zhao and Hu [6] established the flexible multibody dynamic model of a folding wing aircraft by combining the Craig–Bampton synthesis technique with the flexible multibody dynamic approach, and the accuracy of this method was verified by simulating the folding process. Jung and Kim [7] investigated the aerodynamic characteristics of a wing during the folding process using the unsteady vortex lattice method and discussed the influences of the folding angle and angular velocity in detail. Hu et al. [8] developed aerodynamic influence coefficient (AIC) matrices in the time domain based on the doublet lattice method and coupled this approach with the structural model to simulate the aeroelastic response of a body-fixed folding wing aircraft during the folding process. Lee and Weisshaar [9] regarded the folding process of the folding wing aircraft as a steady-state process and used the steady-state approach of static trimming to calculate the hinge moments at different folding angles on the basis of the ZAERO software. Reich et al. and Bowman et al. [10, 11] developed an integrated aeroelastic multibody morphing simulation tool by using flexible multibody dynamics (ADAMS) and vortex lattice model aerodynamic loading (an in-house code). On the basis of this tool, Scarlett et al. [12] conducted a series of wind tunnel simulations and studied the change in load paths and hinge moments under specific motions.

Given the aforementioned research status, the calculation of the hinge moment, as a representative issue of the folding wing aircraft, has attracted extensive attention. Considering the complexity of flight-folding simulation, the majority of current studies on this problem adopt steady-state or quasi-steady-state simplified methods, which are suitable for a slow folding process. However, the applicability of such methods to a fast folding process has yet to be verified.

In the current work, we construct an aeroelastic flight simulation platform to simulate the flight-folding process of a folding wing aircraft. A flexible multibody dynamic model of the folding wing structure is established using ADAMS software, and the aerodynamic force is calculated via the AIC matrices and loaded in ADAMS through the user subroutine. The dynamic response during the flight-folding process of the aircraft is simulated on the basis of this platform, the hinge moments of the folding wings are calculated, and the influences of the folding rate and the center of gravity (c.g.) position on the results are investigated.

#### 2. Simulation Platform

##### 2.1. Model of Folding Wing Aircraft

A typical half model of a folding wing aircraft is shown in Figure 1. This half model can be treated as a flexible multibody structure, as each flexible body is a substructure (i.e., fuselage (I), inner wing (II), and outer wing (III)). In addition, an aileron (IV) is arranged at the trailing edge of the outer wing to maintain stability during the folding process. The folding angle is defined as the angle between the fuselage and the inner wing. To guarantee the lift performance of the folded wing, the outer wing is designed to be parallel to the fuselage.

##### 2.2. Framework of the Platform

The aeroelastic flight simulation platform constructed in this paper is shown in Figure 2. To effectively complete the dynamic simulation of the folding wing aircraft during the flight-folding process, the simulation platform has the following functions: flexible multibody dynamic modeling, unsteady aerodynamic load modeling, and flight control modeling.

The platform integrates the capabilities of ADAMS and NASTRAN in structural modeling, combining advanced finite element models with multibody dynamic models to simulate complex flexible multibody systems. ADAMS currently allows incorporating flexible bodies from finite element analyses by the Craig–Bampton component mode synthesis method [13]. By using the NASTRAN software, a modal analysis of each aircraft component is conducted with the classical Craig–Bampton modal synthesis method to generate modal neutral files suitable for importing into the ADAMS modeling environment. Each component must be created separately in order for ADAMS to treat them as separate parts. Subsequently, the components introduced into the system are assembled by revolute joints in ADAMS, and rotational motion is applied to each revolute joint. The entire model (assembled machine) for simulation is shown in Figure 3.

Thus far, we have completed the flexible multibody dynamic modeling of the folding wing aircraft, and it needs to integrate the unsteady aerodynamic force to achieve dynamic simulation. However, ADAMS itself does not have an aerodynamic calculation module and cannot directly calculate complex unsteady aerodynamic loads. Fortunately, ADAMS provides a friendly second-development interface, and we can use the user subroutine to realize aerodynamic updating and loading during the folding process. Thus, we set GFORCE-generalized forces at the structural nodes and write the SUBROUTINE GFOSUB user subroutine accordingly. The subroutine flow is shown in Figure 4. The subroutine uses SYSARY and SYSFNC macros to read the displacement, velocity, and acceleration information of the structural nodes from ADAMS solver and passes them to the aerodynamic calculation program. By using the AIC method, the aerodynamic load is calculated and then fed back to ADAMS solver.

In addition, mass distribution and aerodynamic load distribution change considerably during the folding process. To maintain flight stability during folding, we establish a longitudinal stabilization control system in the ADAMS built-in control module. The stability control system of the aircraft uses the angle of attack, pitch rate, altitude, and plunging velocity as the feedback input signals and then generates a control law by using the PID control method.

##### 2.3. Aerodynamic Model

The modeling of the unsteady aerodynamic load is a critical work in the dynamic simulation of the folding process. Based on the idea of rational approximation, we develop the doublet lattice method and propose a calculation method for the AIC matrices in the time domain suitable for the calculation of unsteady aerodynamic loads during the flight-folding process.

The conventional unsteady AIC matrix based on the doublet lattice method can be obtained in simple harmonic conditions. The matrix is a function of the reduced frequency and related to the Mach number Ma [14]. The aerodynamic model can be expressed in the following form: where is the dynamic pressure, is the AIC matrix for the aerodynamic nodes, is the Mach number, is the reduced frequency, and and are the aerodynamic force vector and the displacement vector, respectively, for aerodynamic nodes (for the folding wing aircraft, containing components in both directions of - and -axes).

###### 2.3.1. Interpolation Technique

In aeroelastic analysis, to adapt to the specific theory of structural dynamics and aerodynamics, the division of the aerodynamic mesh and the structural mesh are conducted independently. To complete the analysis of aeroelasticity, two transformations are needed: one transforms the displacement of the structural nodes into the displacement of the aerodynamic nodes, while the other transforms the aerodynamic force of the aerodynamic nodes into an equivalent force that acts on the structural nodes. This transformation can be achieved by applying spline matrix [15, 16], which can be expressed as where and are the aerodynamic force vector and the displacement vector, respectively, for the structural nodes. By substituting equation (2) into equation (1), we can obtain the relationship between aerodynamic force and displacement for the structural nodes. where is the AIC matrix for the structural nodes.

###### 2.3.2. Rational Approximation

The conventional AIC matrix based on the doublet lattice method can be obtained in simple harmonic conditions. However, the output can only be used to calculate frequency domain, and it cannot be directly applied to the aerodynamic calculations in the time domain. To calculate the unsteady aerodynamic force in arbitrary motion, the Roger method [17] is used in this study to derive a rational approximation of the unsteady AIC matrices in the Laplace domain. The fitting formula is where is the Laplace variable, is the reference length, is the air speed, and are the undetermined matrices, and are the undetermined coefficients. The values in the reduced frequency range of interest are taken as . For the general fitting accuracy requirement, is 4.

By transforming equation (4) for application in the time domain environment, we express the aerodynamic force for the structural nodes as

For in equation (5), the following relationship exists:

By iterating the above equation into time steps, the equation can be transformed into

Thus far, we have obtained the AIC matrices in the time domain for a specific configuration. On the basis of the small-perturbation assumption of the doublet lattice method, the obtained AIC matrices in the time domain have become suitable for the calculation of aerodynamic force in the small perturbation range near the given configuration. Here, the small perturbation range includes the elastic deformation of the wings and a small angle of folding and pitching motion. In addition, considering that the aerodynamic force is calculated relative to the reference plane, the reference plane can be dynamically selected as the plane at the aircraft’s altitude. Thus, the AIC matrices can also be used to calculate aerodynamic force in heavy plunging motion.

###### 2.3.3. Interpolation of the Folding Angle

The folding process of the folding wing aircraft is a large-scale morphing process. The AIC matrices in the time domain in single configuration cannot fully describe the aerodynamic load during the entire folding process. Thus, we calculate a set of AIC matrices every 15° of the folding angle. The AIC matrices and of each folding angle are obtained by linearly interpolating two sets of AIC matrices (i.e., before and after ). The interpolation algorithm is defined as follows:

Figure 5 shows the flowchart of building the aerodynamic model of the folding wing aircraft. By using the method, we can obtain the AIC matrices in the time domain that continuously changes with the folding angle. The aerodynamic force at any folding angle can be expressed as where can be obtained in time steps with equation (7).

#### 3. Numerical Examples

##### 3.1. Geometrical Analysis

The geometries and dimensions of the folding wing in the present study are illustrated in Figure 6. The fuselage and wings are composed of a skin and inner skeleton and assume an NACA0006 airfoil shape that defines the upper and lower surfaces. The aileron is set as the control surface with 20% of the outer-wing chord and is located between 11% and 83% of the outer-wing span [9]. These substructures are composed of aluminum and connected by rotating hinges, as indicated by the circles in Figure 6(c). The drivers are placed at the hinges of the folding wings. Counterweights are added to the fuselage and wings to represent the weight of the engine, actuators, and transmissions, as indicated by the triangles in the figure. The intersections of the spars and ribs are selected as the structural nodes where aerodynamic forces are applied and are identified in Figure 6(c). The aerodynamic panel model is presented in Figure 6(b). The aircraft mass is 823 kg, and other main parameters are listed in Table 1.

##### 3.2. Dynamic Response during the Flight-Folding Process

Based on the aeroelastic flight simulation platform built in this study, a simulation of dynamic response of the folding wing aircraft during the flight-folding process is conducted. We set the flight altitude to 2 km, maintain the airspeed to 100 m/s, and add a driver to the aileron for the control system. The aircraft is designed to maintain an unfolded wing configuration at the initial 10 s for the preliminary trimming calculations, and then, the wings begin to fold at 10 s. When the wings are folded in place, the aircraft maintains a leveled flight in the folded wing configuration.

In the simulation, we first investigate the general variation of the dynamic parameters of the aircraft during the folding process, and then, we study the influence of the folding rate and the c.g. position on the results. Thus, we first simulate the folding process of the aircraft with a folding time of 30 s. To reduce the impact at the start and end of the folding process, a cosine driving law was applied to the folding. The main dynamic parameters of the aircraft are depicted in Figure 7.

**(a) Folding angle**

**(b) Angle of attack**

**(c) Altitude**

**(d) Deflection angle of aileron**

**(e) Hinge moment of the inner wing**

**(f) Hinge moment of the outer wing**

Figures 7(a)–7(d) show that the aircraft is in equilibrium at approximately 6 s. The initial trimmed angle of attack is 1.87°, and the deflection angle of the ailerons is −5.15°. At 10 s, the wings begin to fold. As the folding angle increases, the effective wing area decreases and the aircraft drops (Figure 7(c)). To recover the altitude, the ailerons are designed to deflect upward movement (Figure 7(d)), which then generates a head-up torque that can increase the angle of attack (Figure 7(b)) and maintain balance of lift. The wings are folded in place at 40 s, followed by a short stabilization process. After approximately 2 s, the aircraft stabilizes to reach a new equilibrium state; that is, the trimmed condition of the folded wing configuration is achieved. Relative to the trimmed condition of the initial unfolded wing configuration, the angle of attack increases to 6.83°, and the ailerons deflect upward to −8.78°.

In addition, we focus on the change in hinge moments during the folding process (see Figures 7(e) and 7(f)). In the initial unfolded state, the driver of the inner wing must simultaneously maneuver the inner and outer wings. The hinge moment in this state is relatively large at 1537 N m, whereas the hinge moment of the outer wing is relatively small at 351 N m. With the folding of the wings, two important changes occur in the aerodynamic force that acts on the wings. First, the arm of the aerodynamic force that acts on the outer wing to the inner wing hinges is reduced, thereby also reducing the hinge moment of the inner wing. Second, the aerodynamic force that acts on the outer wing increases with the angle of attack, thus also increasing the hinge moment of the inner wing. Under the common variation of the two changes at 30.75 s at an approximately 94° folding angle, the hinge moment of the inner wing reaches a maximum of 2637 N m. The hinge moment of the outer wing increases with the aerodynamic load in the outer wing and reaches a maximum value of 1762 N m when the wings are folded in place.

Considering the complexity of the dynamic behavior during the flight-folding process, we verify the above conclusions only for the equilibrium state of the unfolded and folded wing configurations. The statically trimmed results calculated by the aeroelastic flight simulation platform (ADAMS platform) and by the NASTRAN platform are listed in Table 2.

The comparison result shows that the statically trimmed calculation based on the NASTRAN platform supports the conclusion of the ADAMS platform, which also demonstrates the reliability of the dynamic simulation results derived from the simulation platform proposed in this paper.

##### 3.3. Effect of the Folding Rate

To complete the folding process at a faster speed and within a shorter time, achieving rapid attack-defense conversion is undoubtedly one of the objectives of a folding wing design. However, most previous studies have regarded the folding process as slow and have adopted steady-state or quasi-steady-state simplified methods. When folding speed is increased, the applicability of such methods has yet to be verified. Here, we simulated the folding process with four different folding times (5, 10, 20, and 30 s) to study the influence of the folding rate on the simulation results and the applicability of the steady-state method. The dynamic response of the flight-folding process is shown in Figure 8. The response indicates that as the folding speed increases, the time required for an aircraft to complete the folding process is shortened, while the transient effect of the folding process becomes increasingly evident and a longer stabilization process is required after folding in place, which also makes the steady-state method face a test.

**(a) Folding angle**

**(b) Angle of attack**

**(c) Altitude**

**(d) Deflection angle of aileron**

**(e) Hinge moment of the inner wing**

**(f) Hinge moment of the outer wing**

Figure 9 shows the specific time required to complete the folding and reach an equilibrium state. The blue bar in the figure indicates the time required for folding, and the red one indicates the time spent in the stabilization process. The figure shows that the time required for the stabilization process increases as the folding speed increases, whereas the total time decreases. The total time required to complete the folding and reach an equilibrium state is 32 s when the folding time is 30 s. The total time required is decreased to 14.8 s when the folding time is 5 s.

**(a) Description of the time spent**

**(b) Time spent**

The hinge moments of the intermediate configurations at different folding speeds are compared to illustrate the transient effect during the folding process and the applicability of the steady-state method. The result is presented in Figure 10, where the hinge moment at a folding angle of 120° is the maximum hinge moment after folding in place. When the folding process is relatively slow (folding time: 30 s), the transient effect during the process is weak, and the steady-state and transient simulation results of the hinge moments are similar. The steady-state simulation error (the difference between the hinge moments calculated by the steady-state and transient method) in this case is less than 1%. When the folding process is fast (folding time: 5 s), the transient effect during the process is evident, and the steady-state method results in considerable simulation errors. The hinge moments of the inner and outer wings are 16.4% and 12.1% larger than the transient simulation results, respectively, at a folding angle of 60°, and 12.4% and 12.2% smaller than the transient simulation results, respectively, at a folding angle of 120°.

**(a) Hinge moment of the inner wing**

**(b) Hinge moment of the outer wing**

The preceding analysis shows that increasing folding speed can achieve fast attack-defense conversion. However, the flight-folding process is complex and dynamic and involves the coupling of aerodynamics, inertia, and control. The steady-state method only considers the balance of the constant aerodynamic loads and gravity. When the folding speed is increased to a certain level, the transient effect during the folding process cannot be ignored. At this point, the steady-state method will generate a large simulation error, and a transient method is required to simulate the flight-folding process.

##### 3.4. Effect of the c.g. Position

We know that the c.g. position is an important factor that affects the distribution of aerodynamic forces on the wings. The distributed aerodynamic forces and corresponding moments also have an effect on the hinge moments during the folding process. Therefore, studying the effect of the c.g. position on the abovementioned simulation results is important. Considering the maneuverability of the folded wing configuration, we examine a wide range of the c.g. position. We adjust the fore and aft locations of the fuselage counterweights so that the aircraft c.g. position varies from 2.12 m to 2.80 m (the aerodynamic center of the extended and fully folded wing configurations is 2.47 and 2.65 m, respectively). Figure 11 depicts the time history of the main dynamic parameters of the aircraft during the folding process at three c.g. positions.

**(a) Angle of attack**

**(b) Deflection angle of aileron**

**(c) Hinge moment of the inner wing**

**(d) Hinge moment of the outer wing**

Figure 11 demonstrates that the c.g. position exhibits great influence on the simulation. Given that the c.g. position moves aft, the angle of attack decreases, and the deflection angle of aileron increases. This result is mainly due to the following reasons. First, when the c.g. position moves aft, it requires substantial aileron deflection to balance the moment in the pitch direction. Second, the downward aileron deflection increases the lift coefficient of the outer wing, and a small angle of attack can satisfy the balance of the normal load. Finally, the reduced angle of attack reduces the aerodynamic load of the fuselage and inner wing, at which point a considerable aerodynamic force will act on the outer wing, thereby resulting in an increase in the hinge moments of the inner and outer wings. Interestingly, we found that for the hinge moment of the inner wing, unlike the three other sets of measurements, the changing trend with the folding angle changes dramatically when the c.g. position moves aft and fore. The results show an increasing trend for a fore location of the c.g. but a decreasing trend for an aft c.g. location. Here, we provide an explanation. As mentioned previously, with the folding of the wings, two important changes occur in the aerodynamic force that acts on the wings. When the c.g. position is located fore, the aerodynamic load on the aircraft mainly acts on the fuselage and the inner wing at the initial unfolded state. As the wings are folded, the aerodynamic load of the outer wing increases sharply. At this time, the second change dominates, thus eventually increasing the hinge moment of the inner wing. When the c.g. position is located aft, the aerodynamic load of the aircraft mainly acts on the outer wing at the initial unfolding state. The aerodynamic load on the outer wing increases limitedly with the folding of the wings, but the arm to the hinges of the inner wing decreases. At this point, the first change dominates, thereby causing the hinge moment of the inner wing to continue to drop.

Figure 12 exhibits the effect of the c.g. position on the maximum hinge moments. At the c.g. , the maximum hinge moments of the inner and outer wings are the smallest at 2637 and 1762 N m, respectively. When the c.g. position moves aft, the maximum hinge moments of the inner and outer wings increase continuously. When the c.g. , the maximum hinge moments of the inner and outer wings are the largest at 3706 and 2563 N m, correspondingly.

#### 4. Conclusion

To simulate the flight-folding process and calculate the hinge moments of a folding wing aircraft, this study develops a calculation method of the AIC matrices in the time domain and constructs an aeroelastic flight simulation platform based on ADAMS. The effectiveness of the platform is verified by numerical examples, and the influences of the folding rate and the c.g. position on the results are investigated. Results show that, for the fast folding process, the steady-state simulation errors of the hinge moments are substantially large, and a transient method is required to simulate the flight-folding process. In addition, the c.g. position considerably affects the hinge moments during the folding process. Given that the c.g. position moves aft, the maximum hinge moments of the inner and outer wings constantly increase.

#### Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

#### Conflicts of Interest

The authors declare no conflict of interest.

#### Acknowledgments

This work was supported by the National Natural Science Foundation of China (Grant No. 11472133).