Robust Autoland Design by Multi-Model H ∞ Synthesis with a Focus on the Flare Phase

Recent advances in the resolution of multi-model and multi-objective control problems via non-smooth optimization are exploited to provide a novel methodology in the challenging context of autoland design. Based on the structured H∞ control framework, this paper focuses on the demanding flare phase under strong wind conditions and parametric uncertainties. More precisely, the objective is to control the vertical speed of the aircraft before touchdown while minimizing the impact of windshear, ground effects, and airspeed variations. The latter is indeed no longer controlled accurately during flare and strongly affected by wind. In addition, parametric uncertainties are to be considered when designing the control laws. To this purpose, extending previous results published by the authors in a conference paper, a specific multi-model strategy taking into account variations of mass and center-of-gravity location is considered. The methodology is illustrated on a realistic aircraft benchmark proposed by the authors, which is fully described in this paper and freely available from the SMAC (Systems Modeling Analysis & Control) toolbox website (http://w3.onera.fr/smac).


Introduction
The steady growth of air traffic in recent years has led to drastic safety standards with the goal of limiting the number of accidents.Since approach and landing remain the most critical flight phases (almost 50% of fatal accidents and 75% of non-fatal hull losses between 1997 and 2016 [1]), particular attention has recently focused on improving autoland systems in adverse conditions.With the help of CAT III instrument landing systems (ILS), which are now available in a rapidly growing list of airports, automatic landing control laws have helped to secure these two phases notably in degraded weather conditions (such as fog and crosswinds).However, despite numerous methodological works [2][3][4][5] over the past two decades, the design, tuning, and validation process of final approach and flare control systems remains a challenging and time-consuming task.As is observed in [4], where a complete design framework together with a dedicated software is proposed, the tuning phase requires rather tricky multi-objective optimization.Note, however, that optimizing the parameters of the flare control system (such as initial altitude and vertical velocity profile) becomes much easier when the internal loops are correctly designed and tuned.In our context, a natural choice is to track the vertical velocity.The main difficulty is to obtain fast and accurate responses despite a significant airspeed decrease, possible windshear, and ground effects.Recall indeed that the flare segment generally lasts less than 7 s and that, to ensure good robust performance properties, the desired vertical speed should be reached at least 1 or 2 s before touchdown.Unlike the approaches detailed in [4,6,7], either using robust nonlinear dynamic inversion or adaptive control schemes to design the aforementioned inner loops, linear-oriented techniques are often preferred by a majority of contributions.It can be observed indeed that, during approach and landing, airspeed and altitude variations remain rather small so that the aircraft behavior is almost linear.In the field of flight control design, the most popular methods are still based today on eigenstructure assignment [8] or LQR (Linear Quadratic Regulator) control techniques [5].Both approaches have been successfully used in Airbus and Boeing design offices and have contributed for nearly 30 years to considerable improvements in the flight control design process.In the meantime, H ∞ control techniques have been progressively developed and evaluated on various flight control problems-see, for example, [2,3], where the flare phase receives particular attention.Yet, despite promising results even in flight tests [9], this third approach has not become as popular as the other two in the industry.Things are, however, likely to change in the near future with the emergence of new tools based on non-smooth optimization techniques [10,11].With these approaches, it becomes possible to impose constraints on the structure and the order of the controller.Although convexity is unfortunately lost in that case, the aforementioned algorithms converge to local solutions, which are (in a large majority of standard applications) not so far from the global (non-structured) optimum.Another interesting feature of these new tools is their capacity to handle multiple models and multiple separate channels [12].This last feature offers new ways to define H ∞ design models, which will be used in this paper to solve the flare control problem despite modeling errors and parametric uncertainties, thus extending the results of [13].
The paper is organized as follows.ILS-based automatic landing issues are briefly reviewed in Section 2. The flare control problem is then detailed and solved by a multi-model and multi-channel H ∞ design approach in Section 3. Specific attention is devoted in this section to robustness.Implementation issues and nonlinear simulation results are then presented in Section 4. Finally, Section 5 concludes the paper and details a few perspectives.For the sake of completeness, a thorough description of the aircraft model used in the nonlinear simulations is presented in Appendix A. This aircraft model together with a complete description of the design problem is freely available from the SMAC toolbox (dedicated to Systems Modeling, Analysis and Control) website (http://w3.onera.fr/smac).

ILS-Based Automatic Landing
As illustrated in Figure 1, automatic landing in the vertical plane can be divided into two main phases: the final approach during which the aircraft must follow a descent path (glide) and the flare segment, which is activated when the landing gear height H LG falls below a threshold value H FLA .The latter is most often fixed around H FLA ≈ 50 ft but might be slightly updated as a function of ground speed.Similarly, in the horizontal plane, the aircraft trajectory must coincide with the runway axis (localizer phase) as long as H LG ≥ H DEC = 30 ft.The alignment phase (or decrab mode) is then activated in order to minimize the lateral efforts on the landing gears at touchdown.Automatic landing control systems are generally designed on simplified models, using decoupling hypotheses between the longitudinal and the lateral axes.These assumptions, particularly during the landing phase, are usually satisfied without any severe restriction for a large majority of civil and military aircraft.

Final Approach
During the final approach, as already clarified, the trajectory of the aircraft must be kept as close as possible to the ILS beam.This is achieved by simultaneously minimizing the norm of the longitudinal and lateral errors ∆ z = f ( gld ) and ∆ y = g( loc ), where gld is the angular error between the nominal glide path and that of the aircraft (see Figure 1), and loc denotes the lateral angular error in the horizontal plane.In both cases, f (.) and g(.) are elementary altitude-dependent arctangent-based nonlinear functions that convert angular errors into metric deviations.
Moreover, during this phase, the calibrated airspeed V c is to be kept constant, and the aerodynamic sideslip angle β must remain zero.Based on the decoupling assumption and the above constraints, most autoland systems are usually based on two inner control loops.In this paper, the following structures are used for the longitudinal and the lateral axis, respectively: where v z , n z , n y , p, q, r, φ denote the vertical speed, the longitudinal and lateral load factors, the roll, the pitch and yaw rates, and the bank angle, respectively.δ th , δ e , δ a , δ r are the control inputs (thrust and elevator/aileron/rudder deflections), while the subscript c stands for a commanded value.The static gains K LON and K LAT are tuned on linearized models by a standard eigenstructure assigment method.The eigenvectors are constrained so as to ensure decoupling between v c and v z along the longitudinal axis and between n y and φ along the lateral axis.Note the use of small letters in Equation (1) to denote the longitudinal and the vertical speeds, which correspond here to variations about their nominal or initial values.This notation will further be used in the remaining of the paper.Thus, the outer control loops are easily tuned.Observing that ∆z ≈ v z , a fairly standard proportional controller k p z can be used to ensure a nominal trajectory in the vertical plane, while the commanded variations v c c = ∆V c c on the calibrated airspeed are simply set to 0: Next, using a bank-to-turn strategy, the lateral deviation ∆ y is controlled via a proportional-derivative action on the bank angle.Finally, the sideslip angle, which is not directly measured, is controlled via the lateral load factor: Remark 1.The derivative term ∆y is not directly accessible to the control law and a pseudo-derivation of ∆ y would introduce a very high noise level since loc is highly affected by noise as well.Then, the following approximation can be used: where χ and V gr denote the route angle and the ground speed, respectively.

Flare and Align Phases
The nominal slope during approach is generally fixed to γ gld = −3 deg, and the airspeed is a function of mass and mean wind.For the considered aircraft, whose mass belongs to the interval [120 t , 180 t], the nominal approach airspeed is around 70 ms −1 .The vertical speed V z then remains above 3.5 ms −1 during this phase but should fall below 1 ms −1 before touchdown.The flare phase, which starts at approximately 50 ft above ground and generally lasts less than 7 s, is then essential for reducing the vertical speed as fast as possible from 3.5 to less than 1 ms −1 .To achieve this goal, the simplest strategy is to update v z c in the longitudinal outer control loop (Equation ( 3)).Unfortunately, the inner control loop (Equation ( 1)) is generally much too slow.Moreover, the airspeed can no longer be controlled accurately during this phase and the throttle must be reduced to idle.A completely new design, to be further detailed in Section 3, is then required.
Initiated 30 ft above ground, the objective of the align phase is to make sure that the longitudinal axis of the plane coincides with that of the runway at touchdown.This is simply realized by relaxing the constraint on the aerodynamic sideslip angle.Thus, the outer-loop n y c = 0 is simply replaced by where ∆ψ denotes the heading angle error to be reduced to 0. The parameter k align is generally tuned empirically until the decrab objectives are met.

Overview of the Autoland Architecture
The above four elements consisting of glide, localizer, flare, and alignment control subsystems are now assembled as shown in Figure 2 to define the global autoland system.These four subsystems are coordinated by two decoupled longitudinal and lateral outer-loops, which, respectively, deliver the desired signals to be followed on each axis ((v cc , v zc ) for the longitudinal axis and (n yc , φ c ) for the lateral axis).Note that these signals are adapted as a function of the landing gear height.Along the longitudinal axis, the desired vertical speed profile v zc is modified when H LG ≤ H FLA to ensure a soft landing while controlling the impact point.Simultaneously, the control loop gains and architecture are modified by switching from the glide to the flare control system.The same architecture is used along the lateral axis.When H LG ≤ H DEC , the desired lateral acceleration n yc is adapted to minimize the landing gear sideslip angle at touchdown despite crosswind.Simultaneously, the localizer tracking system is deactivated in favor of the alignment control device.

Overall Controller Structure
As briefly discussed above, the objective of the flare phase is to provide enhanced vertical speed control capacities so that the aircraft hits the runway about 400 m after threshold with a nominal vertical speed of 2.5 ft/s.During this short maneuver, rarely exceeding 7 s, the engines are set to the idle position and the airspeed V c is no longer accurately controlled (provided, of course that any stall risk is avoided).Along the longitudinal axis, the only remaining active control input is then the elevator deflection.As shown in Figure 3, the general structure of the flare system consists of two nested loops.A reference vertical speed profile is generated by the outer-loop.It depends on both the initialization altitude (generally close to 50 ft) and the ground speed.Despite perturbations and uncertainties, the generated profile should then be tracked accurately thanks to an appropriate tuning of the vertical speed inner-loop controller, whose design is described in the following section.

Vertical Speed Control via Multi-Channel Structured H ∞ Optimization
The central objectives and constraints of the flare phase are fundamentally different from those encountered during the final approach.From a control viewpoint, the design problem appears to be simpler since only one variable (the vertical speed) is being tracked and only one control input (the elevator deflection) is being managed.However, the fluctuations of the (no longer controlled) longitudinal airspeed V c and the engine deceleration generate perturbations, which, in addition to the ground and wind effects, ultimately aggravate the control problem.It appears that modal control techniques (which have been successfully used to tune the glide controller) are no longer the best option to address this specific problem.The latter now falls within the scope of H ∞ design methods.
The key issue is then to define an appropriate H ∞ design model, which begins with the choice of relevant linear approximations of the aircraft longitudinal dynamics near the ground.Such models are easily obtained with trimming and linearization tools (available with the benchmark package as mentioned in Section A.6) for various points in the flare two-dimensional trajectory (height above runway and kinematic slope).After standard matrix manipulations of the state-space data, the following parameterized set of third-order models is considered: where θ is the pitch angle, and H = H LG denotes the height of the landing gear above ground (to alleviate notation, especially in Figure 4).Note that, in this model, v c and δ th are no longer viewed as a state and a control input signal, respectively, but as external perturbations.However, unlike wind inputs w d , these are measured and can then be used for feedback.This is illustrated in the design diagram depicted in Figure 4.Note that the controller to be optimized reduces to a static gain K ∈ R 7 .
The proposed vertical speed control law structure then reads where, as shown in the diagram, the signal v z r coincides with the output of a linear reference model R(s) to be tracked.The main interest of the proposed structure is to remain close to the one used during the glide segment (see Equation ( 1)), which will not only simplify the switching strategy but also contribute to an improved safety level.Let us now focus on the optimization criteria.Beyond the standard closed-loop stability requirements, the objective is to determine the best controller K such that the control activity is such that magnitude and rate saturations are avoided, and • the perturbations w pert have a negligible effect on the vertical speed variations.
In the H ∞ framework [14], the first two objectives are realized through the norm minimization of the following transfer function (using the notation of Figure 4): where P 1 (s, H) denotes the parameterized (height-dependent) weighted standard interconnection from w 1 to [z 11 z 12 ] , and F l (., .)stands for the lower linear fractional transformation (LFT).The weighting functions W (s) and W u (s) are chosen following standard rules to be further discussed in Section 3.4.The third objective is then achieved by the H ∞ norm minimization of an independent transfer from w 2 (with where the three-dimensional filter W pert (s) is used to shape the input perturbations.However, as this transfer is considered independently of Equation ( 9), it is not restrictive to consider a static weighting function W pert (s) = diag(w p vc , w p δth , w p wind ), where each static gain is used to scale the relative weights of each input perturbation.The three transfers then (Equation ( 10)) are weighted by a single scalar dynamic filter W v z (s).In order to remove any static error (mostly induced by the speed variations) on v z , which could be critical near ground, a high-gain low-pass filter must be used.Summarizing the above discussion, the controller optimization may now be stated as the following multi-model & multi-channel structured H ∞ control design problem:

Robustness Against Parametric Uncertainties and Modeling Errors
By simultaneously considering several design models for different values of the height of the landing gear above ground, interesting robustness properties against ground effects modeling errors are enforced in the proposed design process.This is, however, not sufficient in practice, where parametric robustness against mass and center-of-gravity location is also required.To this end, the above multi-objective H ∞ design framework is easily generalized to take these additional parametric uncertainties explicitly into account.Thus, the design problem (11) becomes where the new integer j ∈ {1, . . ., N} indexes a selection of different models in the (mass × centering) hypercube.A common choice is to consider the vertices of the hypercube.Such a design may then be followed by an analysis phase in order to detect any potential worst-case that could appear inside the hypercube.An iterative design process described in [15] can therefore be considered by incorporating the worst-case scenarios into the initial selection of design models.

Weight Selection & Resolution Aspects
The above problems (Equations ( 11) and ( 12)) are strongly non-convex and cannot be solved with the standard H ∞ control optimization tools, which provide unstructured full-order controllers and are unable to optimize separate channels on multiple models.However, using specialized non-smooth optimization techniques, some efficient algorithms have been proposed in [10] and then improved when they were implemented in the HINFSTRUCT and SYSTUNE routines of the MATLAB R Robust Control Toolbox [16].
In this application, the HINFSTRUCT routine was used on an augmented design model with a block-diagonal structure including 16 (2 × 4 × 2) models associated with

•
two values of the landing gear height: H LG ∈ {15 m, 5 m}, • four models (indexed by j) at the vertices of the (mass × centering) hypercube, and • two separate optimization channels (indexed by i).
Despite potential numerical difficulties, for a given set of weighting functions, a solution is obtained in less than 2 min on a standard computer thanks to the efficiency of the specialized non-smooth optimization algorithm implemented in the HINFSTRUCT routine.Such a limited computational time allows for different trials in the selection of the weighting filters, which are detailed in the next section.

Reference Model on the Vertical Speed
To ensure good time response while remaining compatible with the aircraft dynamics, the reference model is chosen to be where τ = 0.4 s, ω = 1 rad/s, and ξ = 0.75.

Performance & Robustness Oriented Weighting Filters
A low-pass filter (enforcing a high performance tracking of slowly varying inputs) is selected for W (s) while a high-pass filter (here approximated by a constant) is used for W u (s).After a short trial-and-error process, the weights are tuned as follows: Note that a constant weighting function is sufficient in our context since the optimized controller is a static gain K.This may not be the case with dynamic controllers, for which roll-off properties generally have to be enforced via high-pass dynamic weighting functions.
Finally, as is already clarified in Section 3.2, static input weighting functions and a high gain low-pass output filter are tuned for the second channel, whose objective is to minimize the effects of perturbations (induced by airspeed, thrust, and wind variations) on the vertical speed:

Optimization Results
With the above filters, a worst-case H ∞ norm is obtained around 1, and the following gains are computed: In the following, two flare control laws (denoted Law I and Law II) are compared.In the first case (Law I), the airspeed feedback (V c ) is not used (K 7 = 0), while in the second case (Law II), the gains are rigorously implemented as they appear in Equations ( 8) and ( 16).

Nonlinear Implementation & Simulation Results
The structured flare control law resulting from the above optimization process was implemented in a nonlinear closed-loop simulation diagram, as shown in Figure 5.The flare system, highlighted in green, is activated as soon as the height H LG of the landing gear is lower than a prescribed value H FLA , which has been classically set here to 50 ft.To illustrate the benefits of using the longitudinal speed variations during the flare phase, two different controllers are compared with deterministic windsteps first and then with turbulences and ILS noises.Next, specific simulations also illustrate the robustness of the proposed controller by comparing the results of different values of the mass and the center-of-gravity location.A comparison of the upper-right plots clearly shows that the second controller (using airspeed) outperforms the first one.The second controller is thus used for the remainder of the simulations.In this last case indeed, the vertical speed profile is very close to the reference represented by the dash-dotted red line.

Wind & Turbulence Effects
Deterministic windsteps and turbulence effects are then analyzed in the following simulations displayed in Figures 8-10.In the first two figures, the global trajectory (final approach and flare) is visualized.In the first case, a windstep is applied at 30 s (see the upper-right subplot of Figure 8).In the second case (Figure 9), strong turbulence conditions are added.In both situations, the proposed guidance and control system behaves well.The impact point is slightly 400 m after the hreshold meets the requirements.
Moreover, as is displayed in Figure 10, a zoom on the flare trajectory reveals that the vertical speed remains close to the reference trajectory (upper-right subplot).It is also noted that, during the flare maneuver, the pitch angle is barely decreasing.This property is of high interest to make the autoland system compatible with standard pilots habits.

Illustration of Robustness
To conclude this section, the robustness of the proposed autoland control system is now evaluated through various simulations corresponding to the four extreme configurations of the mass (varying between 120 t and 180 t) and the center-of-gravity location (x cg ∈ [15%, 40%]).The results (showing the two most critical variables in the longitudinal axis, namely, the height of the landing gear and the vertical speed) are presented in Figure 11.A reasonable dispersion of the impact point (lying in the interval [400 m , 600 m]) is observed.It can also be pointed out that the vertical speed is well tracked whatever the flight condition.

Conclusions
A new multi-model H ∞ -based design approach has been proposed and evaluated in this paper to improve and shorten the tuning process of flare control systems.Exploiting new capabilities of H ∞ optimization tools, the proposed methodology allows to design a fast and accurate vertical controller despite measured (longitudinal speed variations) and external (wind) perturbations.Thanks to this new controller, the flare outer-loop design becomes much easier and the safety of the whole maneuver is improved.Moreover, thanks to a flexible framework enabled by a multi-model strategy, robustness characteristics against parametric uncertainties and modeling errors are easily taken into account.This is achieved here to enforce robustness against inaccurate modeling of the ground effects and significant variations of the mass and the center-of-gravity location.
To comply with Equation (A1), these must also be expressed in the body-axis frame.Assuming that the thrust is aligned with the longitudinal axis, one obtains The gravity forces in the body-axis frame are given by Finally, the aerodynamic forces, initially expressed in a stability-axis frame, are obtained as where α, q d , and S denote the angle-of-attack, the dynamic pressure (see Equation (A15)), and the reference surface, respectively.The drag (C X = −C D ), lateral (C Y ), and lift (C Z = −C L ) coefficients are detailed in Section A. 4.
The moment M about the center of gravity G of the aircraft results from the engines thrust and aerodynamic actions: Assuming that both engines deliver the same thrust, its location can be reduced to a single point E whose coordinate along the y-axis is zero.The moment resulting from the thrust is Note that the engines are located below the center of gravity, such that z eng > 0 and one observes a pitching moment.The aerodynamic moment exhibits two terms: The first one contains the main contribution and is directly proportional to the moment coefficients C l , C m , and C n (about roll, pitch, and yaw axes, respectively) and to the aerodynamic mean chord L. The second term is linked to the fact that the point of application A of the aerodynamic forces F a possibly differs from the center of gravity G (which may change during flight with fuel consumption).

Appendix A.3. Wind and Atmosphere Effects
Aerodynamic moments and forces strongly depend on the airspeed V a = V a , but also on the air density ρ, via the dynamic pressure: Remark A2.Along the longitudinal axis, the aircraft is trimmed by the horizontal stabilizers whose dynamics are much slower than those of the elevators.It is assumed that the stabilizers are not used during the final approach, so they are not represented here.For control and simulation purposes, the model is delivered with a MATLAB R routine ACStrim.m to enable fast trimming and linearization of the nonlinear model.Various linear representations can thus be easily obtained according to the selected approach speed, initial slope, runway altitude, temperature, mass, and center-of-gravity location.This model and the associated tools can be downloaded from the benchmark section of the SMAC project (http://w3.onera.fr/smac),where a technical note co-signed by Onera and Airbus is also available [17].The latter includes a brief description of the aircraft model but is mainly focused on autoland design guidelines and requirements.It was written to offer control researchers a realistic and challenging autoland design problem and to promote competition in this field.It is quite interesting to note that a number of high-quality contributions have already been reported [18,19].

Figure 1 .
Figure 1.Basic segments of the approach and landing maneuver in a vertical plane.

4. 1 .
On the Interest of Using the Airspeed During Flare: Comparison of Laws I & II In a nominal mass and centering configuration (mass = 150 t and x cg = 21%), the last 15 s of simulation before touchdown (thus focusing on the flare trajectory) are displayed in Figures 6 and 7 .

Figure 6 .
Figure 6.Flare trajectory with Control Law I (without speed feedback).

Figure 7 .
Figure 7. Flare trajectory with Control Law II (with speed feedback).

Figure 8 .
Figure 8. Final approach and flare trajectory with deterministic windsteps.

Figure 9 .
Figure 9. Final approach and flare trajectory with deterministic windsteps and turbulences.

Figure 10 .
Figure 10.Zoom on the flare trajectory with deterministic windsteps and turbulences.

Figure 11 .
Figure 11.Robustness to mass and center-of-gravity variations.

Appendix A. 6 .
SIMULINK R Implementation The above equations are implemented in an open-access SIMULINK R diagram ACS.slx (Figure A1) whose main block "FALSIM" (for Final Approach & Landing SImulation Model) exhibits 9 inputs (4 control inputs, 3 wind inputs, and 2 ILS noises inputs) and 40 outputs.The first 16 ones correspond to the state vector (including actuators dynamics).The next 5 correspond to simulation oriented signals.These are not available for feedback.Finally, the last 19 outputs are measured signals.Note that the last two (∆ y and ∆ z ) correspond, respectively, to the lateral and vertical deviations from the nominal trajectory.These are delivered by the ILS device.

Figure A1 .
Figure A1.A general SIMULINK R view of the aircraft simulation model.