Modelling fraMework for flight dynaMics of flexible aircraft

. The flight dynamics and handling qualities of any flexible aircraft can be analysed within the Cranfield Air-craft Accelerated Loads Model (CA 2 LM) framework. The modelling techniques and methods used to develop the framework are presented. The aerodynamic surfaces were modelled using the Modified Strip Theory (MST) and a state-space representation to model unsteady aerodynamics. With a modal approach, the structural flexibility and each mode’s influence on the structure deflections are analysed. To supplement the general overview of the framework equations of motion, models of atmosphere, gravity, fuselage and engines are introduced. The AX-1 general transport aircraft model is analysed as an example of the CA 2 LM framework capabilities. The results showed that, according to the Gibson Dropback criterion, the aircraft with no control system lacks the stability and its longitudinal handling qualities are unsatisfactory. Finally, the steps for future developments of the CA 2 LM framework are listed within conclusions.


introduction
Airframe flexibility effects have always been of concern to aircraft designers.As a consequence, manufacturers have developed extensive loads and aeroelastic analysis processes aimed to minimise the airframe weight, develop technologies to achieve environmental targets (European Commission 2011;Tollefson 2016) and satisfy the safety requirements set by the regulatory authorities.However, for the design of traditional aircraft, these processes are usually decoupled from the flight dynamic analysis and assessments.This has been justified by the relatively small size and high stiffness of the traditional airframe.With the advent of modern large transport and high altitude long endurance (HALE) aircraft, where the extensive use of advanced materials has led to large and light weight flexible airframes, the interaction between the flight dynamics and aeroelasticity has become a more significant design driver.Flight dynamics analysis methods can no longer assume a rigid airframe and the aeroelasticity practices cannot ignore the rigid body flight dynamics.
Modelling frameworks of various complexity have been developed both in the industry and academia.Industrial frameworks are highly complex and aimed at supporting certification activities.These often couple computational fluid dynamics (CFD) with computational structural mechanics (CSM) resulting in processes that provide the desired insight, but are computationally very expensive (Cooper et al. 2016;Lindhorst et al. 2014;Wang et al. 2015).Reduced order models such as VARLOADS (Kier et al. 2005) have also been developed, but these have only seen limited research usage.In academia, Palacios et al. (Palacios et al. 2010;Palacios, Cesnik 2008;Simpson et al. 2015) have shown the capability to link aeroelasticity with flight control and develop novel approaches to the aeroservoelastic analysis of highly flexible configurations.Structural flexibility effects were modelled through the implementation of a nonlinear structural dynamics formulation.Aerodynamic contributions were captured through the implementation of an unsteady vortex lattice method code.Although the approach adopted by Palacios et al. is computationally cheaper than those used in the industry, a real time simulation is still not possible.
The Cranfield Aircraft Accelerated Loads Model (CA 2 LM) framework was initially developed for the evaluation of handling qualities of large flexible aircraft (Andrews 2011;Lone 2013).It also provides the capability for the flight control law design and a reduced order aeroservoelastic analysis of user-defined airframe configurations.This article provides a brief overview of this modelling framework and its components, along with examples demonstrating its use for the flight loads and handling qualities analysis.

Overview of the CA 2 LM framework
The CA 2 LM framework provides an environment for the modelling and simulation of flexible aircraft (of various configurations) in MATLAB/Simulink.This not only allows the framework to be easily linked with in-house flight control toolboxes and open source codes such as SIDPAC (Morelli 2002) for system identification, but also allows the potential for connections with the flight simulation facilities available at Cranfield University.The framework was initially developed for modelling Fig. 1.AX-1 aeroplane model and its specifications the AX-1 configuration (Fig. 1).Since then, the framework has seen numerous upgrades and is now known as the CA 2 LM framework.This section discusses the high level structure of the framework and the techniques implemented to model aerodynamics, structural dynamics and the equations of motion.The AX-1 configuration will be used to demonstrate the capabilities of this simulation framework.
The overall structure of the CA 2 LM framework is shown in Figure 2. The user can provide the time domain signals representing inputs such as aileron, elevator, rudder and throttle variations.The outputs are aircraft rigid body states, such as aircraft position in the Earth's axis, its angular and translational velocities and the attitude in the body axis.Internal structural loads, such as bending moments and torsion, can also be the output.The core of the framework consists of the aerodynamic, structural, gravity and equations of motion blocks.These are discussed separately below.The gust/turbulence block provides an environment for modelling atmospheric disturbances and allows the implementation of a continuous turbulence and discrete gusts.The non-aerodynamic loading block allows the specification of specific mass properties.Fuel, cargo and passenger loadings can be specified in detail and this information is used to calculate aircraft mass, inertia tensor and the centre of gravity position.

Modelling of aerodynamic surfaces
The aerodynamic modelling process is further detailed in Figure 3.The wings, tailplanes and the fin are modelled in very similar ways.However, a block, modelling interference effects between the lifting surfaces and the fuselage, is added to the wing aerodynamics.Steady aerodynamic forces and moments are modelled using the Modified Strip Theory approach that relies on the input of appropriate aerofoil aerodynamic characteristics as a function of airspeed and angle of attack.This enables the calculation of the aerodynamic forces on a user defined wing planform and takes into account the compressibility effects via the Prandtl-Glauert correction factor (DeYoung, Harper 1948;Weissinger 1947).A Leishman-Beddoes unsteady aerodynamic model has been implemented in a state-space form.Therefore, the entire airframe is divided into strips and each strip has a focal point about which the forces and moments are calculated.These are referred to as the aerodynamic nodes.
The implementation of an unsteady aerodynamics model is considerably more involved than that of the steady model.Therefore, a brief summary of the modelling is provided here.The unsteady aerodynamics model is programmed in the following state-space form: (1) The state matrix A is a square matrix that may be represented as follows: 1,1 1,15 15,1 15,15 where the non-zero terms are defined as follows: ; ; where V is the airspeed, β is the Prandtl-Glauert compressibility correction factor, c aero is the chord of an aerofoil, b 1 …b 6 are the exponents of indicial functions (Leishman 1988).Within the state matrix and later in the output matrix C, the following non-circulatory time constants (Leishman 1988) are also used: where A 1 …A 4 are the coefficients of the various indicial functions (Leishman 1988), a is the speed of sound, M is the Mach number, F 1 …F 11 are wing control surface geometric properties (Theodorsen 1949), x e represents the hinge location of control surface as a percentage of the chord.The input matrix B takes the following form: The non-zero terms of the matrix B are as follows: The output matrix C is represented in the following form: The non-zero terms of the matrix C are as follows: ( ) 1 2 1,10 2,1 2 2 ; ; 2 2 ; ; 0;  ; ; (8) The feedthrough matrix D takes the following form: And its non-zero terms are as follows: ( ) ( )( ) The state vector x and input vector u are as follows: where n is the number of states.
Each aerodynamic node has a 15 element state vector x associated with it, together with an input vector u consisting of the angle of attack, the angle of the control surface and their rates of change (Andrews 2011;Lone 2013).For the AX-1 model, the surfaces generating lift are modelled using 58 aerodynamic nodes that result in 870 unsteady aerodynamic states.
The steady aerodynamic coefficients for each section of the lifting surfaces are found from pre-programmed look-up tables (LUTs).Therefore, parameters such as viscous drag, zero lift drag, aerofoil profile drag and zero lift pitching moment coefficients and profile drag increase due to flaps are obtained through the simple interpolation for a specified Mach number and Reynolds number.To take into account the 3D effects, an indicial angle of attack (α ind ) is added to the steady state angle of attack (α) and the effective angle of attack (α eff ) is calculated: (13) The Modified Strip Theory is then applied to provide the forces acting on aerodynamic surfaces in the wind axes system.These are transferred into the body axes system via the application of the following direction cosine matrix (DCM) that considers local deformation along with the relative changes in the orientation of the two axes systems:  The aerodynamic model also estimates the wing downwash effect on the tailplanes.Downwash circulation strength (Γ) to estimate this influence is calculated as follows: where s is the span of a wing, A is the coefficients matrix of the Modified Strip Theory, is the zero lift angle of attack.Circulation Γ is then evaluated through a reduced order state-space model to get the indicial angle of attack for the tailplane.This implementation of the Modified Strip Theory and unsteady aerodynamic modelling has been found to provide a satisfactory balance between precision and computational cost (Andrews 2011).

Fuselage and engines modelling
The fuselage and engines make a significant contribution of forces and moments acting on aircraft.A sketch of the fuselage and engines with the corresponding sources of modelling methods is shown in Figure 5.
For aircraft such as the Airbus A340 or the Boeing 777 the fuselage flexibility effects on the flight dynamics and handling qualities cannot be ignored.Within this framework, the fuselage flexibility is taken into account through the definition of elastic angles of attack and sideslip.The fuselage is divided into three parts -the nose, the tail and the central section which consists of the wing fuselage junction.Flexibility is considered via the changes in the angles of attack and sideslip for the nose and tail parts due to their deflection as shown in Figure 5.The forebody of the fuselage is modelled as an axisymmetric slender body (ESDU 89008 and ESDU 89014) and the aftbody is modelled as an axisymmetric conical body (ESDU 87033).
Engine dynamics is also modelled in this framework.The forces and moments from each engine are split into two parts -the nacelle aerodynamics and the thrust producing unit.The nacelles are modelled as annular aerofoils (ESDU 77012).The flexibility is taken into account, as in the case of the fuselage, through additional terms for the angles of attack and sideslip.Within the AX-1 implementation, the thrust producing unit is modelled as a turbofan engine.The forces and moments of each engine are calculated in the engine axes system.However, these are transferred to the body axes system to be included in the total forces and moments.

Differences between aerodynamic and structural frames
The forces and moments from aerodynamic surfaces, fuselage and engines are calculated at the aerodynamic nodes and the structural forces and moments are evaluated at the structural nodes.Therefore, an aeroelastic simulation requires a transformation of the aerodynamic forces and moments in aerodynamic nodes to the structural nodes and vice versa.However, a typical implementation, such as the AX-1 model, has the aerodynamic contributions being calculated at a higher resolution than the structural dynamic contributions.Thus, the number of aerodynamic nodes often exceeds the number of the structural nodes and these nodes are not collocated in space.Loads therefore need to be transformed from the aerodynamic frame to the structural frame and it is very important to analyse the difference between those two frames.
Figure 6 shows the scheme as applied to the AX-1 model.Along the wing of the AX-1 there are 35 aerodynamic nodes, which correspond to 35 aerodynamic strips.It is shown (Fig. 7) that this is the optimal number of strips providing the desired balance between the model accuracy and the computational cost.The same analysis was done for the tailplane and fin, resulting in the selection of 15 and 8 aerodynamic strips respectively.On the other hand, the structural layout is modelled with 21, 7 and 4 nodes for the wing, fuselage and tailplane, respectively.10 nodes are used for the fuselage modelling, 2 of which coincide with the central nodes of the wing and tailplane.Hence, additional operations are done converting aerodynamic loads to structural loads and then structural frame deflections to aerodynamic frame deflections.

Structural modelling
Structural modelling is done in a structural dynamics block and the process is shown in Figure 8.This block converts aerodynamic and gravitational loads to structural loads.The structural dynamics for the AX-1 implementation is done in the modal domain, thus stiffness and mass matrices are generated to obtain structural mode shapes.The first 4 structural modes of the AX-1 model are shown in Figure 9.
The following structural equations of motion are solved to acquire modal accelerations ( x  ), velocities ( x  ) and displacements ( x ): where F i is modal force, m i is modal mass, ω n is modal natural frequency, ζ is damping ratio, i is the number of mode.
Finally, the transition from modal to nodal displacements, velocities and accelerations is done.The results of each mode influence are acquired and summed with other modal influences to give the resultant displacements, velocities and accelerations.The AX-1 implementation only considers the first 12 modes because the model aims to analyse the flight dynamics phenomena that are typical at low frequencies.
It is important to note that only small deflections (less than 10% of a wing semi-span) are modelled within the CA 2 LM framework, as it is assumed that the properties of each beam vary linearly.However, recent developments in highly flexible aircraft (Patil, Hodges 2006) have introduced wing deflections of more than 25% of the wing semi-span.As a result, a non-linear approach to model structural dynamics is currently under investigation.

Equations of motion
The forces and moments acting on the aircraft are concentrated at the centre of gravity (CG), about which the accelerations, attitude, position and velocities are calculated.However, airframe flexibility is taken into account through the recalculation of the moments for a constantly changing CG position.This method has been considered as an appropriate way of taking flexibility effects into account.The equations solved for the body forces (F b ) and the moments (M b ) are in vector form as follows (Stengel 2004;Cook 2007): T is linear velocities in x-, y-and z-axes, ω = [p q r] T is angular velocities around x-, y-and z-axes, Yet, it should be noted that significant changes in the CG position are expected because of high structural deformations.A constantly changing CG position will result in a time-varying inertia tensor I. Hence, a contribution of each node should be taken into account in the equations of motion and a new approach is currently under development.

Gravity and atmosphere modelling
Gravity is modelled according to the WGS-84 reference (WGS-84 1991).Gravitational constant (γ h ) is calculated using the following equation: , (19)   where γ e is theoretical gravity at the equator, k is theoretical gravity formula constant, e is the first ellipsoidal eccentricity, φ is geodetic latitude, a is the semi-major axis, f is ellipsoidal flattening, ω is the angular velocity of the Earth, b is the semi-minor axis, GM is the Earth's gravitational constant, h is height.The gravitational constant is then applied at the CG position for solving equations of motion.Additionally, it is applied to each structural node to solve the structural equation of motion.Atmospheric properties such as air density and temperature are modelled as the International Standard Atmosphere (ISA) according to ESDU 77021.

case studies utilising the ca 2 lM framework
This section briefly presents two case studies demonstrating the capabilities of the CA 2 LM framework.The first case study focuses on the handling qualities analysis and the second demonstrates the capability of performing failure case assessments.Both case studies are based on the AX-1 model, which is representative of a large transport aircraft.

Time domain handling qualities analysis
The Gibson Dropback Criterion (Gibson 1982) is a well-known approach developed to predict longitudinal handling qualities and assist in the design of command and stability augmentation systems.The key advantage of this approach is that it is based in the time domain, so the effects of nonlinear dynamics arising due to nonlinear flight control can be considered in the handling qualities analysis.Such effects cannot be captured through approaches based on low order equivalent systems (LOES).The key parameters for evaluating the Dropback criterion are: 1. Pitch rate overshoot ratio, which is expressed as a ratio between the maximum pitch rate (q max ) and the steady state pitch rate (q ss ). 2. Attitude dropback (DB) to the steady state pitch rate (q ss ) ratio.These parameters are illustrated graphically in Figure 10.The criterion is based on these ratios and the extensive pilot opinion gathered to outline the regions of satisfactory and undesirable response characteristics, as shown in Figure 11.The boundaries shown in Figure 11 are based on the research conducted by Mooij (Mooij 1985), which focused on large transport aircraft.
In this case study, the AX-1 model was trimmed at an altitude of 10000 ft and the Dropback criterion was evaluated at several airspeeds.This was carried by specifying an elevator pulse input of ±5°. Figure 11 shows the variation of longitudinal handling qualities with varying airspeed.It should be noted that no stability augmentation system has been implemented, and, consequently,

Aileron soft failure simulation
A control surface failure scenario is one of many extreme cases that need to be considered for the flight loads evaluation.Here a soft aileron failure is simulated, where the port aileron undergoes an actuation failure whilst the starboard aileron remains in the original trim setting.The main results obtained from the simulation of the AX-1 model are shown in Figure 12.The port aileron is forced to effectively undergo a limit cycle oscillation at a constant frequency of 1.16 Hz, which corresponds to The frequency content of the roll rate (p) and yaw rate (r) signals shows that the failure has excited a low frequency lateral-directional mode corresponding to the periods of T p = 10.24 sec and T r = 10.92 sec in roll and yaw, respectively.These correspond to the usual frequencies of the aircraft's Dutch roll mode.The highest peaks, just above 1 Hz, are the direct result of the simulated aileron forcing function.The load factor (n) only exhibits large transients when the aileron failure is initiated.
Figure 13 shows the frequency content of the wing root bending moment M root at different aileron excitation frequencies.At a frequency of 1.245 Hz, slightly higher than the frequency of the first structural mode of the wing (1.1634Hz), the first aeroelastic mode appears and a resulting resonance is observed.Upon magnification (bottom right subfigure), other two peaks can be observed at 2.5 Hz and 3 Hz.These correspond to the aeroelastic modes associated with the 5 th and 11 th structural wing bending modes.At the frequency of 0.9 Hz, the M root is higher than at the frequency of 1.1 Hz, which can be explained by the fact that the forcing function frequency is getting closer to the rigid body frequencies.

conclusions
A brief overview of the CA 2 LM framework designed to model flexible aircraft has been presented in this paper.Structural deformations are obtained through a linear modal formulation of the aircraft structure.An assumption of linearity limits the model to small deformations that are less than 10% of the wing semi-span.The aerodynamics is modelled by coupling the steady Modified Strip Theory with the Leishman-Beddoes unsteady model in the state-space form.The CA 2 LM framework effectively combines these methods in a MATLAB/ Simulink environment.The capabilities of such an environment are demonstrated through two case studies.These cases have focused on the AX-1 model, which represents a generic large transport aircraft.The first case study focuses on the handling qualities analysis based on the Dropback criterion.It demonstrates that the AX-1 model's response to a longitudinal control input is unsatisfactory without a stability augmentation system.The second case study simulates a port aileron failure case and its impact on structural loads.It shows that the coupling between aeroelastic modes and rigid body flight dynamic modes appears when the aileron undergoes a limit cycle oscillation at a slightly higher frequency than the first wing bending mode.
Recent developments in highly flexible aircraft have introduced wing deflections of more than 25% of wing semi-span.Thus, a new approach to structural modelling is currently being developed.Moreover, such a flexible aircraft cannot be assumed as a rigid body when solving the flight dynamic equations of motion.Hence, a new approach including additional terms due to the flexibility into the equations of motion is being investigated.

acknowledgements
The author would like to acknowledge Stuart P. Andrews for developing the AX-1 aircraft model and tools for its analysis during his PhD degree studies at Cranfield University.These tools were the basis for developing the current CA 2 LM framework.

Fig. 2 .
Fig. 2. Structure of the CA 2 LM framework where θ i is local twist angle, λ i is local sweep angle, γ i is local dihedral angle.The various axes systems used in the model are shown in Figure4.

Fig. 4 .
Fig. 4. Different axes systems used in the model

Fig. 10 .
Fig.10.Visualisation of q max , q ss and DB terms used in the Gibson Dropback criterion

Fig. 12 .
Fig. 12. Ailerons deflection δ A , angular rate, load factor n, wing root bending moment M root and wing root torsion T root and the roll, pitch and yaw rates at the aileron excitation frequency f = 1.1634Hz