2.1 Static stability

The total aircraft pitching moment is given by Ref. Reference Roskam20:

(1) \begin{equation}{{\bf{C}}_{\bf{m}}} = {{\bf{C}}_{{{\bf{m}}_{\bf 0}}}} + {\rm{\;}}{{\bf{C}}_{{{\bf{m}}_{\boldsymbol{\alpha }}}}} + {\rm{\;}}{{\bf{C}}_{{{\bf{m}}_{{{\bf{i}}_{\bf{h}}}}}}} + {{\bf{C}}_{{{\bf{m}}_{{{\boldsymbol{\delta }}_{\bf{e}}}}}}}\end{equation}

where ${C_{{M_0}}}$ is the value of C _M for α = i _h = δ _e = 0

(2) \begin{equation}{{\bf{C}}_{{{\bf{m}}_{\bf 0}}}} = {{\bf{C}}_{{{\bf{m}}_{{\bf{a}}{{\bf{c}}_{{\bf{wf}}}}}}}} + {\rm{\;}}{{\bf{C}}_{{{\bf{L}}_{{{\bf 0}_{{\bf{wf}}}}}}}}\!\left( {{\bf{h}} - {{\bf{h}}_{\bf 0}}} \right) + {{\bf{C}}_{{{\bf{L}}_{{{\boldsymbol{\alpha }}_{\bf{h}}}}}}}{\boldsymbol{\eta }}\frac{{{{\bf{S}}_{\bf{T}}}}}{{\bf{S}}}\!\left( {{{\bf{h}}_{\bf 0}} - {\bf{h}}} \right){{\boldsymbol{\varepsilon }}_{\bf 0}}\end{equation}

and ${C_{{M_\alpha }}}$ = , the change in total aircraft pitching moment coefficient with a change in angle-of-attack

(3) \begin{equation}{{\bf{C}}_{{{\bf{m}}_{\boldsymbol{\alpha }}}}} = {{\bf{C}}_{{{\bf{L}}_{{{\boldsymbol{\alpha }}_{{\bf{wf}}}}}}}}\!\left( {{\bf{h}} - {{\bf{h}}_{\bf 0}}} \right) - {{\bf{C}}_{{{\bf{L}}_{{{\boldsymbol{\alpha }}_{\bf{h}}}}}}}{\boldsymbol{\eta }}\frac{{{{\bf{S}}_{\bf{T}}}}}{{\bf{S}}}\!\left( {{{\bf{h}}_{\bf 0}} - {\bf{h}}} \right)\left( {{\bf 1} - \frac{{{\bf{d}{\boldsymbol\varepsilon }}}}{{{\bf{d}\boldsymbol{\alpha }}}}} \right)\end{equation}

and ${C_{{M_{{i_h}}}}}$ = , the change in total aircraft pitching moment coefficient with a change of horizontal tailplane angle, i_h for α = δ _e = 0

(4) \begin{equation}{{\bf{C}}_{{{\bf{m}}_{{{\bf{i}}_{\bf{h}}}}}}} = - {{\bf{C}}_{{{\bf{L}}_{{{\boldsymbol{\alpha }}_{\bf{h}}}}}}}{\boldsymbol{\eta }}\frac{{{{\bf{S}}_{\bf{T}}}}}{{\bf{S}}}\!\left( {{{\bf{h}}_{\bf 0}} - {\bf{h}}} \right) = - {{\bf{C}}_{{{\bf{L}}_{{{\boldsymbol{\alpha }}_{\bf{h}}}}}}}{\boldsymbol{\eta }}{\bar{\textbf{V}}_{\bf{h}}}\end{equation}

where ${{\bar{\textbf{V}}}_{\bf{h}}}$ , the tail volume ratio is given by:

(5) \begin{equation}{{\bar{\textbf{V}}}_{\bf{h}}} = \frac{{{{\bf{S}}_{\bf{T}}}}}{{\bf{S}}}\!\left( {{{\bf{h}}_{\bf 0}} - {\bf{h}}} \right)\end{equation}

and ${C_{{M_{{\delta _e}}}}}$ = , the change in total aircraft pitching moment coefficient with a change of elevator angle, δ _e for α = i _h = 0

(6) \begin{equation}{{\bf{C}}_{{{\bf{m}}_{\boldsymbol{\delta }}}_{\bf{e}}}} = - {{\bf{C}}_{{{\bf{L}}_{{{\boldsymbol{\alpha }}_{\bf{h}}}}}}}{\boldsymbol{\eta }}{{\bar{\textbf{V}}}_{\bf{h}}}{{\boldsymbol{\tau }}_{\bf{e}}}\end{equation}

and horizontal tailplane efficiency, η ranges from 0 to 1.

Thus, the contributions to the pitching moment fall into two major categories [Reference Gudmundsson21], contributions that vary with angle-of-attack and this that do not vary with angle-of-attack. The tailplane efficiency (η) can be used to model the degradation of lift due as a result of tailplane icing.

The total aircraft pitching moment about the aircraft’s centre of gravity [Reference Nelson22] consists of contributions from wing, tail and fuselage (Fig. 3). Each contribution generates moments that vary with angle-of-attack and contributions that are independent of angle-of-attack (constant). A negative total pitching moment slope represents positive static stability – the aircraft returning to the trimmed flight condition following a disturbance (e.g., sudden change of tail lift).

Figure 3. Contributions to total aircraft pitching moment about the CG.

As ice builds up on the tail, it becomes less aerodynamically efficient, increased negative tail lift is needed to keep maintain the trimmed flight condition and prevent the aircraft tending to nose down and overspeed with flaps extended. This is achieved by increasing elevator trailing edge up to compensate (Fig. 4).

Figure 4. Forces and moments on an aircraft and effect of elevator deflection (adapted from Ref. 23).

When flaps are retracted (Fig. 5), the lift generated by the wing decreases and the point of lift moves forward. The aircraft tends to nose UP as the nose down pitching moment due to wing lift decreases. Increased trailing edge down (TED) elevator is required to compensate. At the tailplane, downwash angle decreases hence negative tail lift decreases. The aircraft tends to nose DOWN as the pitching moment due to negative tail lift decreases. This time, increased trailing edge UP (TEU) elevator is required to compensate. The resultant pitching moment about the CG is the net effect of these contributions. The effects due to tailplane are usually dominant (Fig. 3).

Figure 5. Forces and moments on an aircraft and effect of flap retraction (adapted from Ref. 23).

2.2 Dynamic stability

The longitudinal dynamic stability characteristics of the aircraft are dependent upon the airspeed, pitching moment variation with angle-of-attack, moment of inertia, lift to drag ratio and aerodynamic tail damping. The system may be exhibit positive, neutral or negative dynamic stability [Reference Hurt24].

Given the trimmed flight condition and static stability of the generic business jet, dynamic analysis can be undertaken to consider the effects of small disturbances (perturbations) such as turbulence (external) or control inputs (internal). Using defined aircraft notation and axes state variable, control inputs and matrix/vector notations are defined. Systems with more than one input and more than one output are known as multi-input multi-output systems (MIMO). Systems that have only a single-input and a single-output are defined (SISO). The aircraft in longitudinal and pitching motion maybe modelled using a state-space model. The longitudinal aircraft dynamics are linearised about the setpoint (airspeed, pitch angle, pitch rate, elevator deflection and throttle setting) and can be written in state space form [Reference Nelson22]:

(7) \begin{equation}{\dot{\textbf{X}}} = {{\bf{A}}_{{\bf{long}}}}{\bf{X}} + {{\bf{B}}_{{\bf{long}}}}{\bf{u}}\end{equation}

The elements of matrix A are stability derivatives describing the effect of state variables on forces and moments. The elements of matrix B are control derivatives representing the effects of elevator and throttle commands on the body referenced forces and moments.

2.3 Transfer functions, elevator to pitch

Using this matrix method, transfer functions (input-output relationships) are derived using desktop computing mathematical modelling tool Matlab to determine the relationship between Input: Elevator Deflection ( ${{\rm{\delta }}_e}$ ) and Output: Pitch Angle ( $\theta$ ).

2.4 Effects of flap retraction

The effects of flap retraction are simulated using a Simulink switching model. This type of model enables the dynamic analysis to account for changes in stability and control derivatives as a result of flap configuration changes. Given stability and control derivatives for the selected aircraft, transfer functions are estimated for different flap configurations and tailplane efficiencies.

3.1 Static stability

For the generic business jet, using the given flight condition of V= 204 KTAS, H = 4,000ft pressure height, $\textrm{T} = 29^{\circ}\textrm{F}$ , aft CG/low MTOM with flap = 7 degrees, the variation of pitching moment with angle-of-attack was obtained for a range of horizontal tailplane efficiencies from 1.0 (100%) to 0.5 (50%) to simulate the effects of icing on the aircraft tailplane (Fig. 6).

Figure 6. Pitch stability – pitching moment vs angle-of-attack @204kts.

The results show that the pitching moment versus angle-of-attack gradient decreases as the horizontal tailplane efficiency decreases, hence aircraft static stability also decreases as horizontal tailplane efficiency decreases.

The results show that the pitching moment versus angle-of-attack gradient decreases as the horizontal tailplane efficiency decreases, hence aircraft static stability also decreases as horizontal tailplane efficiency decreases. This trend is also apparent by examining Equation (3) showing that the tail is associated with a negative (hence stabilising) contribution to the pitching moment coefficient, C _Mα.

The results also show that the angle-of-attack for trimmed flight conditions increases as tailplane efficiency decreases. Analysis of horizontal tailplane efficiency versus stick-fixed static margin suggests that static margin decreases as the horizontal tailplane efficiency decreases and that the aircraft is neutrally (statically) stable when the horizontal tailplane efficiency is approximately 0.68 for the generic business jet model in given flight conditions (Fig. 7).

Figure 7. Tailplane efficiency factor vs static margin.

Further analysis of elevator deflection versus horizontal tailplane efficiency (Fig. 8) suggests that increasing UP elevator (–ve elevator deflection according to convention) is required to maintain trimmed flight as the horizontal tailplane efficiency decreases. The range of elevator deflection for the generic business jet model was 20 degrees UP and 15 degrees DOWN. The results suggest that as horizontal tailplane efficiency decreases below approximately 0.2 (20%), there is insufficient UP elevator to maintain trimmed flight.

Figure 8. Elevator deflection vs tailplane efficiency.

3.2 Dynamic stability

Having established the trimmed flight condition and static stability of the generic business jet the dynamic analysis was undertaken to consider the effects of small disturbances (perturbations) such as turbulence (external) or control inputs (internal). Before conducting the dynamic analysis, the aircraft, the system of notation and axes were defined using the right-hand rule [Reference Hopkin25]. State variable, control inputs and matrix/vector notations were defined as follows:

(8) \begin{equation}\boldsymbol{m}{\dot{\boldsymbol{V}}_{\boldsymbol{I}}} = \boldsymbol{m}\!\left( {{{\dot{\boldsymbol{V}}}_{\boldsymbol{B}}} + {\boldsymbol\omega _{\boldsymbol{B}/\boldsymbol{I}}} \times {\boldsymbol{V}_{\boldsymbol{B}}}} \right) = {\boldsymbol{F}_{\boldsymbol{Aero}}} + {{\boldsymbol{F}}_{\boldsymbol{gravity}}} + {\boldsymbol{F}_{\boldsymbol{T}}}\end{equation}

(9) \begin{equation}\boldsymbol{I}{\dot{\boldsymbol\omega}_{\boldsymbol{I}}} = \boldsymbol{I}{\dot{\boldsymbol\omega}_{\boldsymbol{B}}} + {\boldsymbol\omega _{\boldsymbol{B}/\boldsymbol{I}}} \times \boldsymbol{I}{\boldsymbol\omega _{\boldsymbol{B}}} = {\boldsymbol{M}_{\boldsymbol{Aero}}} + {\boldsymbol{M}_{\boldsymbol{gravity}}}\end{equation}

where F _T is the throttle force that only acts along the X axis (other components on thrust due to misalignment were verified to be negligible). The throttle acceleration will be denoted $\tau = \dfrac{{{F_T}}}{m}$ longitudinal motion only affects the X and Z axes and the pitch rotation (Fig. 9), with v=r=p=Y=L=N=0:

(10) \begin{equation}\dot{\boldsymbol{U}} = - {\boldsymbol{g}}{\bf{sin\boldsymbol\Theta }} - {\boldsymbol{QW}} + \frac{{{\boldsymbol{X}} + {{\boldsymbol{F}}_{\boldsymbol{T}}}}}{\boldsymbol{m}}\end{equation}

(11) \begin{equation}\dot{\boldsymbol{W}} = - {\boldsymbol{g}}{\bf{cos\boldsymbol\Theta }} + {\boldsymbol{QU}} + \frac{{\boldsymbol{Z}}}{{\boldsymbol{m}}}\end{equation}

(12) \begin{equation}\dot{\boldsymbol{q}} = \frac{{\boldsymbol{M}}}{{{{\boldsymbol{I}}_{{\boldsymbol{yy}}}}}}\end{equation}

Figure 9. Aircraft trim condition at the time of the event (Adapted from Ref. 23).

Equations are linearised around a constant pitch $\theta$ ₀ and constant longitudinal velocity u ₀, with q ₀=w ₀=0 representing a zero rate of change in pitch and altitude.

The trim condition is given by:

(13) \begin{equation}{\boldsymbol{X}_{\textbf{0}}} - \boldsymbol{mg}{\bf{sin}}\!\left( {{\boldsymbol\theta _{\bf{0}}}} \right) = \bf{0}\end{equation}

(14) \begin{equation}{\boldsymbol{Z}_{\bf{0}}} + \boldsymbol{mg}{\bf{cos}}\!\left( {{\boldsymbol\theta _{\bf 0}}} \right) = {\bf 0}\end{equation}

(15) \begin{equation}{\boldsymbol{M}_{\bf 0}} = {\bf 0}\end{equation}

More precisely, U,W,Q and Θ are linearised around u ₀, w ₀, q ₀, $\theta_0$ such that:

\begin{equation*}\boldsymbol{U} = {\boldsymbol{u}_{\bf 0}} + {\bf{u}}\end{equation*}

\begin{equation*}{\bf{W}} = {\boldsymbol{w}_{\bf 0}} + {\bf{w}}\end{equation*}

\begin{equation*}{\bf{Q}} = {\boldsymbol{q}_{\bf 0}} + {\bf{q}}\end{equation*}

\begin{equation*}{\boldsymbol{\Theta }} = {\boldsymbol\theta _{\bf 0}} + \theta \end{equation*}

Θ is replaced by $\theta$ ₀+ $\theta$

By taking the differences and substituting Θ by $\theta$ ₀+ $\theta$ , we have:

(16) \begin{align}\dot{\boldsymbol{u}} & = \frac{\bf 1}{\boldsymbol{m}}\!\left( {\left( {\frac{{\boldsymbol\partial \boldsymbol{X} + \boldsymbol\partial {\boldsymbol{X}_{\boldsymbol{T}}}}}{{\boldsymbol\partial \boldsymbol{u}}}} \right)\boldsymbol{u} + \frac{{\boldsymbol\partial \boldsymbol{X}}}{{\boldsymbol\partial \boldsymbol{w}}}\boldsymbol{w} - \boldsymbol{mg}\boldsymbol\theta {\bf{cos}}{\boldsymbol\theta _{\bf 0}}} \right)\nonumber\\[4pt]\dot{\boldsymbol{w}} & = \frac{\bf 1}{\boldsymbol{m}}\!\left( {\frac{{\boldsymbol\partial \boldsymbol{Z}}}{{\boldsymbol\partial \boldsymbol{u}}}\boldsymbol{u} + \frac{{\boldsymbol\partial \boldsymbol{Z}}}{{\boldsymbol\partial \boldsymbol{w}}}\boldsymbol{w} + \frac{{\boldsymbol\partial \boldsymbol{Z}}}{{\boldsymbol\partial \boldsymbol{q}}}\boldsymbol{q} - \boldsymbol{mg}\boldsymbol\theta {\bf{sin}}{\boldsymbol\theta _{\bf 0}} + {\boldsymbol{u}_{\bf 0}}\boldsymbol{q}} \right)\nonumber\\[4pt]\dot{\boldsymbol{q}} & = \frac{\bf 1}{\boldsymbol{I}_{\boldsymbol{yy}}}\!\left( \frac{\boldsymbol\partial \boldsymbol{M}}{\boldsymbol\partial \boldsymbol{u}}\boldsymbol{u} + \frac{\boldsymbol\partial \boldsymbol{M}}{\boldsymbol\partial \boldsymbol{w}}\boldsymbol{w} + \frac{\boldsymbol\partial \boldsymbol{M}}{\boldsymbol\partial \dot{\boldsymbol{w}}} \dot{\boldsymbol{w}} + \frac{\boldsymbol\partial \boldsymbol{M}}{\boldsymbol\partial \boldsymbol{q}}\boldsymbol{q} \right)\\[4pt]\dot{\boldsymbol\theta} & = \boldsymbol{q}\nonumber\end{align}

3.3 MIMO (multiple input multiple output) aircraft longitudinal state space model

Systems with more than one input and more than one output are known as multi-input multi-output systems (MIMO). Systems that have only a single input and a single output are defined (SISO). The aircraft in longitudinal and pitching motion maybe modelled using a state-space model. The longitudinal aircraft dynamics, linearised about the setpoint (V= 204 KTAS, $\theta$ = 24 degrees pitch), can be written as a state space model:

(17) \begin{equation}{\dot{\textbf{X}}} = {{\bf{A}}_{{\bf{long}}}}{\bf{X}} + {{\bf{B}}_{{\bf{long}}}}{\bf{u}}\end{equation}

Matrix/vector notations in bold

A: State matrix, depends on stability derivatives and setpoint condition

B: Input matrix, depends on control derivatives and actuator layout

X: State vector, components are u, w, q, $\theta$

u: Control input vector, components are η, τ

where:

\begin{equation*}{\bf{X}} = \left[ {\begin{array}{*{20}{c}}u\\[2pt]w\\[2pt]q\\[2pt]\theta \end{array}} \right],{\bf{u}} = \left[ {\begin{array}{*{20}{c}}{{\delta _e}}\\[2pt]\tau \end{array}} \right]\end{equation*}

u: X axis speed, w : Z axis speed.

q: pitch rate, $\theta$ : pitch, ${{\rm{\delta }}_e}$ : elevator deflection, τ : throttle command.

Throttle history was recorded within flight data but elevator history was not.

Using matrix notation, the state transition is given by:

(18) \begin{equation}{\bf{A}} = \left[ {\begin{array}{c@{\quad}c@{\quad}c@{\quad}c}{{\boldsymbol{X}_{\boldsymbol{Tu}}} + {\boldsymbol{X}_{\boldsymbol{u}}}} & {{\boldsymbol{X}_{\boldsymbol{w}}}} & {{{\bf{X}}_{\boldsymbol{q}}} - \boldsymbol{m}{\boldsymbol{w}_{\bf 0}}}& { - \boldsymbol{g}{\bf{cos}}\!\left( {{\boldsymbol\theta _{\bf 0}}} \right)}\\[5pt]{\dfrac{{{\boldsymbol{Z}_{\boldsymbol{u}}}}}{\boldsymbol{f}}}& {\dfrac{{{{\bf{Z}}_{\boldsymbol{w}}}}}{\boldsymbol{f}}} &{{\boldsymbol{Z}_{\boldsymbol{q}}} + \dfrac{{{\boldsymbol{u}_{\bf 0}}}}{\boldsymbol{f}}}& { - \dfrac{{\boldsymbol{g}{\bf{sin}}\!\left( {{\boldsymbol\theta _{\bf 0}}} \right)}}{\boldsymbol{f}}}\\[13pt]{{{\bf{M}}_{{\bf{Tu}}}} + {{\bf{M}}_{\boldsymbol{u}}} + {{\bf{M}}_{{\dot{\textbf{w}}}}}\dfrac{{{{\bf{Z}}_{\boldsymbol{u}}}}}{\boldsymbol{f}}{\rm{\;}}}& {{{\bf{M}}_{\boldsymbol{w}}} + \dfrac{{{{\bf{M}}_{{\bf{T\alpha }}}}}}{{{\boldsymbol{u}_{\bf 0}}}} + {{\bf{M}}_{\dot{\boldsymbol{w}}}}\dfrac{{{{\bf{Z}}_{\boldsymbol{w}}}}}{\boldsymbol{f}}}& {{{\bf{M}}_{\dot{\boldsymbol{w}}}}\dfrac{{{{\bf{Z}}_{\boldsymbol{w}}}}}{{\bf{f}}} {{\bf{M}}_{\boldsymbol{q}}} + {{\bf{M}}_{\dot{\boldsymbol{w}}}}\dfrac{{{{\bf{Z}}_{\boldsymbol{q}}} + {{\bf{u}}_{\bf 0}}}}{\boldsymbol{f}}}& { - {{\bf{g}}{\dot{\textbf{w}}}}\dfrac{{\bf{sin}\! \left( {{\boldsymbol\theta _{\bf 0}}} \right)}}{\boldsymbol{f}}}\\[13pt]{\bf 0} &{\bf 0} &{\bf 1} &{\bf 0}\end{array}} \right]\end{equation}

where ${\boldsymbol{f}} = {\bf 1} - {{\boldsymbol{Z}}_{\dot{\boldsymbol{w}}}}$

All stability derivatives are in dimensional form and represent the partial derivatives of a force or moment with respect to a state variable, divided by mass in the case of forces or inertia in the case of moments, for example:

\begin{equation*}{{\boldsymbol{X}}_{\boldsymbol{u}}} = \frac{\bf 1}{\boldsymbol{m}}\frac{{\boldsymbol\partial \boldsymbol{X}}}{{\boldsymbol\partial \boldsymbol{u}}},{\rm{\;}}{{\boldsymbol{M}}_\boldsymbol{q}} = \frac{\bf 1}{{{\boldsymbol{I}_{\boldsymbol{yy}}}}}\frac{{\boldsymbol\partial \boldsymbol{M}}}{{\boldsymbol\partial \boldsymbol{q}}}, \textrm{etc}.\end{equation*}

The dimensional stability derivatives were derived from the non-dimensional stability derivatives of the aircraft using the equations defined by Roskam [Reference Roskam26].

The non-dimensional derivatives were lift, drag, pitching moment coefficients ${C_{L1}},{C_{D1}},{C_{M1}}$ respectively, defined at speed setpoint u₀ and pitch setpoint $\theta$ ₀, and ${C_{L0}},{C_{D0}},{C_{M0}}$ that are constants depending on aircraft geometry.

1. The throttle coefficients: ${C_{TX0}},{C_{MT0}}$

2. Speed dependent coefficients: $\;{C_{Lu}},{C_{Du}},{C_{Mu}},{C_{TXu}},{C_{MTu}}$

3. Force coefficients defined with respect to angle-of-attack and pitch rate: $\;{C_{D\alpha }},{C_{L\alpha }},{C_{L\dot \alpha }},{C_{Lq}}$

4. Moment coefficients defined with respect to angle-of-attack and pitch rate: $\;{C_{M\alpha }},{C_{M\dot \alpha }},{C_{Mq}},{C_{MT\alpha }}$

5. These were derived using aircraft design software for stability and control analysis [27].

The control matrix is given by:

(19) \begin{equation}{\bf{B}} = \left[ {\begin{array}{c@{\quad}c}{{{\boldsymbol{X}}_{{\boldsymbol\delta _{\boldsymbol{e}}}}}} & {{{\boldsymbol{X}}_{\boldsymbol\tau} }}\\[9pt]{\dfrac{{{{\boldsymbol{Z}}_{{\boldsymbol\delta _{\boldsymbol{e}}}}}}}{{\boldsymbol{f}}}}& {\bf 0}\\[13pt]{{{\bf{M}}_{\boldsymbol\eta} } + {{\bf{M}}_{\dot{\boldsymbol{w}}}}\dfrac{{{{\boldsymbol{Z}}_{{\boldsymbol\delta _{\boldsymbol{e}}}}}}}{{\boldsymbol{f}}}} & {\bf 0}\\[13pt]{\bf 0} & {\bf 0}\end{array}} \right]\end{equation}

where the dimensional control derivatives were obtained from the non-dimensional control derivatives. ${C_{D{\delta _e}}},{C_{L{\delta _e}}},{C_{M{\delta _e}}}$ of the aircraft at pitch and speed setpoint.

The elements of matrix A are stability derivatives describing the effect of state variables on forces and moments. The elements of matrix B are control derivatives representing the effects of elevator and throttle commands on the body referenced forces and moments.

Note: The B1 type model requires computation of flap aerodynamic parameters, which are not available. The simulations were generated by computing the A and B for different tailplane stall efficiency settings.

3.4 Transfer functions, elevator to pitch

Aircraft parameters (Appendix A, Table A1) were used define elements in respective matrices. Using this matrix method, transfer functions were derived using Matlab [28] to determine the relationship between Input: Elevator Deflection ( ${{\rm{\delta }}_e}$ ) and Output: Pitch Angle ( $\theta$ ). Transfer function coefficients were derived from expressions for the aircraft stability derivatives [Reference Roskam26].

For the given flight condition (Appendix A, Table A1):

The elevator to pitch transfer function, denoted G(s), with 100% horizontal tailplane efficiency was found to be:

(20) \begin{equation}{\boldsymbol{G}_{\boldsymbol{elevator} \to \boldsymbol{pitch}{\rm{\;}}\left( {\bf{100}{\rm{\% \;}}{\bf{tail}}{\rm{\;}}{\bf{efficiency}}} \right)}}\!\left( \boldsymbol{s} \right) = \frac{{ - \bf{17.39}{\boldsymbol{s}^{\bf 2}} - {\bf{52.13}}\boldsymbol{s} - \bf{1.128}}}{{{\boldsymbol{s}^{\bf 4}} + \bf{5.521}{\boldsymbol{s}^{\bf 3}} + \bf{12.5}{\boldsymbol{s}^{\bf 2}} + \bf{0.3387}\boldsymbol{s} + \bf{0.2689}}}\end{equation}

where the negative numerator coefficients are consistent with a pitch down response to elevator down deflections and the positive coefficients are consistent with a stable response.

Decreasing the horizontal tailplane efficiency factor was found to significantly change the intrinsic characteristics of the transfer function. For example, a 50% tailplane efficiency factor results in the following transfer function:

(21) \begin{equation}{\boldsymbol{G}_{\boldsymbol{elevator} \to \boldsymbol{pitch}{\rm{\;}}\left( {50{\rm{\% \;}}{\bf{tail}}{\rm{\;}}{\bf{efficiency}}} \right)}}\!\left( \boldsymbol{s} \right) = \frac{{\bf{8.736}{\boldsymbol{s}^{\bf 2}} - \bf{26.17}\boldsymbol{s} - \bf{0.6315}}}{{{\boldsymbol{s}^{\bf 4}} + \bf{4.226}{\boldsymbol{s}^{\bf 3}} - \bf{0.03886}{\boldsymbol{s}^{\bf 2}} + \bf{0.0396}\boldsymbol{s} + \bf{0.1132}}}\end{equation}

The negative coefficient (–0.03886) in the denominator for Equation (15) is linked to an unstable pole, and the pitch response to elevator commands is therefore unstable when horizontal tailplane efficiency is reduced to 50%. Lower values of the tailplane efficiency factor were found to further increase instability.

Figure 10. Short period oscillation (SPO) with $100\%$ horizontal tailplane efficiency.

Figure 11. Long period oscillation (LPO) with $100\%$ horizontal tailplane efficiency.

3.5 Short period oscillation (SPO) and long period oscillation (LPO)

For model validation purposes, the longitudinal dynamics were analysed by excitation of the short period (Fig. 10) and long period (Fig. 11) modes using the eigenvectors of A associated with each mode to specify initial conditions as follows (i is the mode number):

(22) \begin{equation}{{\bf{x}}_{\bf{0}}} = {{\boldsymbol{E}}_{\boldsymbol{i}}}{{\boldsymbol{T}}_{\boldsymbol{i}}} + {\boldsymbol{E}}_{\boldsymbol{i}}^{\rm{*}}{\boldsymbol{T}}_{\boldsymbol{i}}^{\rm{*}}\end{equation}

where E _i and T _i are respectively the eigenvalue and eigenvector associated with the mode i, with i=1 for the phugoid mode and i=2 for the short period mode.

The variables simulated were relative changes (variations with respect to a trim condition). Zero value indicates the initial trimmed flight condition and all angles are in radians.

The results show that for a horizontal tailplane efficiency of 100%, the aircraft is statically and dynamically stable with heavy and moderate damping for the SPO and LPO, respectively. In addition, the model was independently verified by comparing the modal characteristics of the full state space model to the reduced order models of the two modes.

When the efficiency is reduced to 80%, the response for the LPO is also stable although with less damping of oscillations than the 100% efficiency case (Fig. 12). The results of the dynamic stability are in agreement with those of the static stability analysis.

Figure 12. Long period oscillation (LPO) with 80% horizontal tailplane efficiency.

Instability was found to start from approximately 50% tailplane efficiency factor, which was in line with the static stability analysis.

3.6 The effects of flap retraction

The effects of flap retraction were simulated using a Simulink switching model. This model enables the dynamic analysis to account for changes in stability and control derivatives as a result of flap configuration changes. Stability and control derivatives were determined for the generic business jet model using the commercial aircraft design software as before. Transfer functions were estimated for four different conditions using the switching model (Table 1).

Table 1. Switching conditions

Figure 13. Pitch response with tailplane stall then flap retraction.

For flap retraction following a tailplane stall, initially, the tailplane efficiency is 100% and the flap angle is 7 degrees (Fig. 13). At t=4 seconds, a tailplane stall is applied by a step reduction of 50% to horizontal tailplane efficiency, which destabilises the system response. The aircraft pitches down a further -30 degrees within 2 seconds. At t= 6 seconds, the retraction of flaps from 7 degrees to 0 degrees helps initially, but insufficiently to stabilise the response with only 50% tailplane efficiency. Note that the pitch response in Fig. 13 effectively represents an asymptotic pitch plot. The true pitch is limited by the upper and lower bounds of this asymptotic plot at the model transition times (t=4 seconds and t=6 seconds) and approaches this asymptotic plot at other times. This is a limitation of using a switched linear system model with a Matlab/Simulink implementation in this instance. The model switches between the discrete stability and control derivatives for each condition (tailplane efficiency and flap angle) resulting in discontinuities. This leads to discrete change in the elevator to pitch transfer function parameters. The asymptotic plot however still provides a good understanding of the pitch variations around the transition times, with an increasingly accurate description of the pitch variations and dynamics at other times.

3.7 Tailplane efficiency factor effects

At 100% horizontal tailplane efficiency, a large 20deg elevator down during 1 second is needed to initiate a large magnitude but stable phugoid response, which is only shown during the first 20 seconds for comparison purposes (Fig. 14). At 50% efficiency, a 1 second 2-degree elevator down input destabilises the system with growing oscillations and at 20% efficiency, a 1-degree elevator down is sufficient to produce very fast divergence without oscillations.

Figure 14. Pitch response with varying tailplane efficiency and elevator inputs.

Load factor changes with 1 second elevator down commands for varying horizontal tailplane efficiencies were also presented (Fig. 15). With 100% efficiency, load factor remains very close to 1 with a small elevator command, as expected. With 50% tailplane efficiency, load factor is reduced but changes are small. Major changes to load factor are however obtained with 20% tailplane efficiency and negative G is quickly reached in this case. This is believed to be closer to the conditions during the flight incident.

Figure 15. Load factor with varying tailplane efficiency and elevator inputs.