# Design of LPV-Based Sliding Mode Controller with Finite Time Convergence for a Morphing Aircraft

Authors: Nuan Wen; Zhenghua Liu; Yang Sun; Lingpu Zhu
Journal: International Journal of Aerospace Engineering
 (2017)

---

## Abstract

This paper proposes a finite time convergence sliding mode control (FSMC) strategy based on linear parameter-varying (LPV) methodology for the stability control of a morphing aircraft subject to parameter uncertainties and external disturbances. Based on the Kane method, a longitudinal dynamic model of the morphing aircraft is built. Furthermore, the linearized LPV model of the aircraft in the wing transition process is obtained, whose scheduling parameters are wing sweep angle and wingspan. The FSMC scheme is developed into LPV systems by applying the previous results for linear time-invariant (LTI) systems. The sufficient condition in form of linear matrix inequality (LMI) constraints is derived for the existence of a reduced-order sliding mode, in which the dynamics can be ensured to keep robust stability andL2 gain performance. The tensor-product (TP) model transformation approach can be directly applied to solve infinite LMIs belonging to the polynomial parameter-dependent LPV system. Then, by the parameter-dependent Lyapunov function stability analysis, the synthesized FSMC is proved to drive the LPV system trajectories toward the predefined switching surface with a finite time arrival. Comparative simulation results in the nonlinear model demonstrate the robustness and effectiveness of this approach.


---

## Body

## 1. Introduction

Morphing aircrafts, as a broad range of air vehicles and vehicle components that can make shape transition with internal morphing devices, can adapt to multimission requirements [1]. As a result, they will be a substitution, covering the roles of several different aircrafts that can allow optimized flight over a large flight envelope instead of merely one flight condition. However, this property will certainly bring great difficulties to the traditional modeling and control methods. Due to the massive change of aerodynamic configuration during flight, the morphing aircraft cannot be considered under full geometry variation and it must take the morphing structures as well as multirigid body variations into consideration. Therefore, unlike most conventional aircrafts which are only concentrated on a fixed-structure and treat themselves as a single rigid body, the consideration of morphing demands a combination with several research areas such as aerodynamic modeling [2], multirigid body dynamics [3], and flight control based on a large range of reference points [4]. Obviously, it is important to incorporate a wide range performance into modeling and controller for the dynamical systems that describe morphing aircrafts. LPV synthesis techniques naturally fit into this characteristic. By selecting appropriate operating conditions of the original nonlinear model, the complex dynamics in morphing process can be represented by LPV dynamic models [5]. LPV systems, being the specific instance of linear time-varying (LTV) systems, are the representation that the entries of the state-space matrices continuously depend on a time-varying parameter vector that belongs to a bounded compact set [6, 7]. Using LPV techniques, the dynamics of the original nonlinear systems can be reduced to the linear equations. Meanwhile, the controller outputs will be continuously “scheduled” according to the system operating conditions. Comparing with the classical gain schedule techniques [8], LPV control can display prominent advantages as it can theoretically guarantee global stability and robustness over the whole operating envelope.While its controller synthesis methods are theoretically well-founded based on LMI, LPV control is confronted with an inevitable challenge: direct generation of these methods to LPV systems will result in infinitely many LMI constraints that have to be satisfied [9]. There are three most common approaches to fixing this problem. An early possible way is to impose the constraints on a particular grid of the region of scheduled parameters [10]. But this discrete numerical method, which is similar to the classical gain schedule technique, will undoubtedly jeopardize robust stability and performance. Later developments in multiplier-based linear fractional transformation (LFT) synthesis make it possible to obtain finite-dimensional LMI constraints. In this case, the stability and performance are preserved, yet potential conservatism is introduced in the full-block S-procedure [11].The third approach to LPV controller synthesis is based on the properties of convex sets, in which it requires the LPV system to admit an affine or polytopic representation. Then, it is unnecessary to consider the full trajectories of scheduling parameter values, only the vertices of the convex hull need to satisfy the LMIs [12].However, a general LPV system may not be a desired polytopic or affine representation. Several kinds of mapping transformation are introduced to a model polytopic LPV system [13].The TP model transformation, as an equivalent transformation methodology in essence, is introduced in [14]. It can be directly employed to obtain a convex polytopic form with the known modeling accuracy, which greatly decreases the computing load of the controller design. Another key issue of the LPV system is to construct parameter-dependent Lyapunov functions to overcome potential conservatism brought by the slowly varying parameters. Aiming at providing the controller with less conservatism and better performance, Wu and Dong [11] incorporate quadratic LFT parameter-dependent Lyapunov functions and full-block multipliers into controller synthesis. Another way to construct Lyapunov functions with time-varying scheduled parameters is employing different Lyapunov functions at each of the vertices of the polytope LPV systems. By introducing slack variables, Peaucelle and Arzelier [15] separate the system matrices and the Lyapunov matrix to follow this method. Obviously, such LMI conditions will reduce conservatism at the expense of computational complexity.For a nonlinear system of aircrafts, the Jacobian linearization approach is the most widespread methodology to obtain LPV models [5]. Generally, the designers built the global LPV model by interpolation and fitting of the LTI models which is formulated at different frozen points of interest throughout the operational envelope. Unfortunately, such steps inherently introduce uncertainties and model mismatches [16]. Moreover, the morphing transition will bring challenges for exact modeling of flight dynamics for its own sake. Both of these must be considered in the control design phase, which means that robust LPV control is attractive. Most of the current researches for LPV systems mainly focus on H∞ control techniques; these methods can reduce the impact of disturbance to some extent. However, unexpected uncertainties and model mismatches will deteriorate the performance of the controller [17]. Robust methods designed for uncertain LPV models have appeared in the literature. Considering actuator fault tolerant, a robust control strategy specified for LPV systems with polytopic uncertainty is proposed in [18]. For uncertain scheduling parameters appearing in LPV systems, the work in [19] develops a gain-scheduled output-feedback controller. In [20], a more general polytopic model describing LPV systems with state-dependent uncertain parameters is considered, and then the output-feedback controllers are designed by introducing the Lagrange multiplier term and adding some weight matrix variables. However, none of them consider general mismatched uncertainties appearing in LPV systems.Sliding mode control is an attractive from ‘uncertainty’ standpoint, since the sliding surfaces are insensitive to matched or mismatched uncertainties [21]. To deal with uncertain systems with general mismatched uncertainties, a new reduced sliding mode is proposed in terms of LMI in [22]. Using the same method, the work in [23] ensures the robust stability and L2 performance in the presence of an external disturbance. In [24], a robust adaptive sliding strategy based on model reference is presented for the tracking control problem of a special aircraft. However, all the methods cited above are based on LTI system descriptions of the plant and are restricted to the single operational condition. There are very limited researches on the use of sliding mode techniques for LPV systems. Most of them are restricted to sliding mode observers for a specific class of LPV systems [25, 26]. By way of exception, the work in [27] proposes a time-varying sliding mode hyperplane for LPV systems. Also in [7], LPV-based integral sliding modes and control allocation are presented in a realistic operational environment. More recently, a standard second-order sliding mode control approach is directly explored for the control of LPV systems in [28]; however, the approach for turning the infinite LMI constraints into a finite set of LMIs is simply a gridding-based synthesis without any theoretical proof.Inspired by these works, the robust sliding mode control problem for uncertain LPV systems of a kind of large-scale morphing aircrafts is studied in this paper. First of all, according to multibody modeling and linearization, longitudinal nonlinear dynamic equations of the morphing aircraft in wing transition are simplified and transformed into an LPV model with mismatched uncertainties and external disturbance. Then, in the following LPV controller synthesis, we give a less conservative LMI existence condition of the reduced sliding dynamics. In order to avoid solving an infinite number of LMI constraints, the obtained LPV system is transformed into a convex polytopic by the TP model transformation. Finally, the synthesis of a parameter-dependent reaching law is investigated for driving the system trajectories toward the predefined switching surface with a finite time arrival.The exact contributions of this paper over previous works are concluded as follows:(i)
Compared with the works in [4, 17], we extend the morphing dimension into two degree of freedom: wingspan and sweep angle. Both detailed modeling process based on the Kane method, and explicit expressions of polynomial parameter-dependent LPV system are given.(ii)
In this paper, we explore and evaluate the potential that the sliding mode control may have in the management of LPV systems with mismatched uncertainties. Unlike the controller in [28], which is just a simple mixture of the existed second-order sliding mode with fixed-structure for uncertain LPV systems, we adopt the reduced-order sliding mode controller to deal with uncertainties in LPV systems. Meantime, by applying the slack matrix method to derive a sufficient condition in the form of finite set of LMIs, we guarantee that the reduced-order dynamics in the sliding mode are robust and stable with the bounded L2 gain performance. Then a parameter-dependent sliding mode control law can make the LPV system evolve into its sliding mode in a finite time.The remainder of the paper is organized as follows: some preliminaries are given in Section2. In Section 3, both a new morphing model and the corresponding LPV system are presented. In Section 4, the design process and stability analysis of FSMC law are demonstrated. In Section 5, numerical simulations are implemented to justify the performance of the proposed controller. The conclusion is shown in Section 6.The notation is standard. The superscript “T” denotes the transpose of a matrix. P>0P≥0 means that P is the real symmetric and positive-definite (semipositive definite); A<BA≤B stands for that the matrix A−B is the negative definite (semi-negative definite); I and 0 represent the identity matrix and a zero matrix, respectively; and • refers to the 2-norm for a vector or a matrix. In symmetric block matrices, we use ∗ to represent a term that is induced by symmetry. diag represents a diagonal matrix. Both a′ and ȧ denote the first order derivative. a″ denotes the second-order derivative. × denotes the outer product of vectors.


## 2. Preliminary Results

Lemma 1. (see [23]).
The following block matrix inequality holds(1)A11A12A12TA22≤A~11A12A12TA22,
ifA11≤A~11 holds.Lemma 2 (see [24]).
X and Y are real matrices (or vectors) of appropriate dimensions; for any scalar ε>0, one has
(2)XTY+YTX≤εXTX+ε−1YTY.Lemma 3 (see [12]).
The following statements are equivalent:(3)INTMIN<0,M+NT−IGT+GN−I<0.Lemma 4 (see [29]).
Assume that a continuous, nonnegative functionVt satisfies the differential inequality V̇t≤−αVηt, where α>0 and 0<η<1 are constants. Then, for any given t0(4)V1−ηt≤V1−ηt0−α1−ηt−t0,t0≤t≤tr,and Vt=0 and t≥tr with tr given by
(5)tr=t0+V1−ηt0α1−η.


## 3. Model and Problem Formulation

Flight dynamics of a morphing aircraft is a key component in the design process and must therefore be modeled accurately as well as practically. In this section, the longitudinal nonlinear dynamic equations of the morphing aircraft are formulated first, and then the linear parameter-dependent model is developed for controller synthesis.


### 3.1. Longitudinal Nonlinear Dynamic Modeling for the Morphing Aircraft

While numerous morphing designs exist for advanced aircrafts, we consider a new type of large-scale morphing aircrafts [30] depicted in Figure 1(a), whose wings can extend, contract, and rotate to a certain extent.Figure 1
(a) The Morphing aircraft with variable sweep angle and variable span. (b) Depiction of the multibody representation.
(a)
(b)The Teledyne Ryan BQM-34 “Firebee”, which is designed as a high-speed target drone, is chosen as the baseline platform [31]. The aircraft is equipped with an easily removable wing that could be modified into a morphing wing. Also, many variants of this aircraft have been constructed to accommodate mission requirements such as subsonic or supersonic flight, reconnaissance, and combat. Next, we suppose that the wing of each side of the aircraft is equipped with an inner wing and an outer wing. In wing morphing process, the inner can rotate with Ob and the outer will keep level and move linearly along with the airfoil (see Figure 1(a)). Owing to this specific structure, the morphing aircraft can change the wing sweep angle and wingspan to accommodate mission requirements. The maximum wing sweep angle (Λmax) is 45° and the maximum extension (Λmax) is 2 m. The geometry of the other parts (fuselage, vertical tail, horizontal tail, etc.) remains unchanged in this paper.To simplify the modeling process, we make the following assumption:Assumption 1.
In the morphing process, the left wing and the right wing will realize the synchronous motion. The symmetry axis of the fuselage lies in the same plane as the wings on both sides.The traditional approach to modeling a conventional aircraft is to treat the whole plane as a single rigid body without regard for the complications brought by the aeroelastic problems. The morphing aircraft with large-scale varations mentioned above, however, cannot generally be modeled as a single rigid body due to the varying inertial properties caused by the changes of geometric parameters of the wings. To describe the aircraft dynamic behavior accurately, the whole structural changes of each part must be taken into account properly, which means that multibody modeling is an appropriate solution. The Kane method is a simple and efficient approach to derive the dynamic model of a multibody system [32]. Thereby, the longitudinal dynamics of the morphing motion is described by using Kane method in this paper.The aircraft under consideration consists of five separate rigid bodies (see Figure1(b)): fuselage (body 1), inner part of left wing (body 2), outer part of left wing (body 3), inner part of right wing (body 4), and outer part of right wing (body 5), and their masses are mb,m1,m2,m1,m2, respectively. The important parameters of the morphing aircraft are given in Table 1.Table 1
b
0.3 m
Λ
0~45 deg
Δ
0~2 m
Jy
3107.53 kg⋅m2
l
1
2.5 m
l
2
2.0 m
mb
771.79 kg
m
1
54 kg
m
2
13.605 kg
T
δ
T
129.27 N/%Firstly, as shown in Figure2, we set the origin of the whole aircraft body coordinate frame Obxbybzb to locate at the center the wing (the unit vectors corresponding to the axes xb,yb and zb are exb,eyb and ezb), and the ground coordinate system is described as Ogxgygzg. The coordinate frames of morphing aircraft’s fuselage and the other several moving parts of the wing are described as Oixiyizii=1,…,5, taking the center of mass of each body as the origin correspondingly. In what follows, Kane’s equations will be established in the reference frame Obxbybzb.Assumption 2.
Both the aerodynamic force and the thrust are merely acting on the fuselage (body 1). And the other parts suffer the effects of gravity only.Figure 2
Coordinate axis system.There are four steps to establish the Kane equations, as depicted in Figure3. Here, we only consider the longitudinal responses in morphing process. Thus, u,w,q are selected as independent generalized speed variables.Figure 3
Kane modeling procedure.In Figure3, FObu,FObw and FObq are generalized active forces corresponding to u,w,q. FObi and MObi are generalized active forces and active moments of the ith body projected into Obxbybzb. F~Obu,F~Obw, and F~Obq are generalized inertial forces relative to u,w,q. F~Obi and M~Obi stand for the generalized inertial force and inertial moment of the ith body. UObi~ and WObi~ are partial velocities of the ith body located into Obxbybzb. The specific definitions of these variables are given in Appendix A.According to the fourth step in Figure3, the detailed expressions of the Kane equations of longitudinal dynamics of the morphing aircraft can be written as follows:
(6)f1u̇,w,q,α,θ,Λ,Λ′,Λ″,Δ′,Δ″=0,f2ẇ,q̇,q,u,α,θ,Λ,Λ′,Δ′=0,f3ẇ,q̇,q,u,α,θ,Λ,Λ′,Δ′=0.The explicit expressions of (6) are given in Appendix B.By introducing (7), (6) can be transformed into a wind coordinate frame. Considering the kinematic equation (8), we finally get the full longitudinal dynamic equation (9).
(7)u=Vcosα,u̇=V̇cosα−Vsinαα̇,w=Vsinα,ẇ=V̇sinα+Vcosαα̇,(8)θ̇=q,ḣ=Vsinθ−α,(9)EV̇α̇q̇θ̇ḣ=F1F2F3qVsinθ−α+W1W2W300.The explicit expressions ofE,F1,F2,F3,W1,W2, and W3 are given in Appendix C.Remark 1.
By adopting the Kane modeling, only the position, velocity and acceleration ofΛ and Δ appear in the dynamic modeling equation (9). They can be computed and manipulated more practically than the position, velocity, and acceleration of the center of mass or instantaneous inertia tensor and its derivative [17]. This shows another major superiority of the Kane method. Notice that we will apply our designed LPV-based controller on the original nonlinear equation (9), and the exact values of Λ′,Λ″ and Δ′,Δ″ at every moment must be obtained so that we can substitute them into (9). Therefore, the motion trajectories of sweep angle and span need to be scheduled appropriately in advance. We will discuss it in the following section.Remark 2.
W
1
,
W
2, and W3 including the first and the second-order derivative of Λ and Δ are the force and moment variations affected by the morphing movement, which can be considered as the additional disturbances in the wing morphing dynamic response [17]. However, in the following section, these extra forces and moments are proved to have such slight influence to the entire dynamic equations that they can be ignored in the following linearization procedure.


### 3.2. Longitudinal LPV Model for the Morphing Aircraft

In this paper, the Jacobian linearization approach is adopted to transform the nonlinear model of the morphing aircraft into an LPV model. The flight condition of interest is selected as the altitudeh0=9144m and the Mach number Ma=0.5. By defining the scheduling variables λ=Λ/45° and ξ=Δ/l1, we get λ∈0,1 and ξ∈0,0.8.Assumption 3 (see [17]).
In this paper, unsteady aerodynamic effects are neglected. We only focus on the aerodynamic forces and moments introduced by the quasi-steady assumption; that is, the aerodynamic forces and moments of the morphing process in different wing configurations are nearly identical to those of the corresponding static configuration.The lift, drag, pitching moments, and thrust in (9) are given as
(10)D=12ρV2SwCDαα,L=12ρV2SwCL=12ρV2SwCLα=0+CLαα+CLδeδe,M=12ρV2Swc¯CM=12ρV2Swc¯CMα=0+CMαα+CMδeδe+c¯2VCmqq,T=Tδtδt.Choosingλ=0,0.2,0.4,…,1.0 and ξ=0,0.2,0.4,…,0.8 as the reference points, then the corresponding aerodynamic parameters for these configurations can be approximately calculated through computational fluid dynamics (CFD). Firstly, the outline dimensions of the aircraft are obtained by observation. Then the three-dimensional model of the morphing aircraft under different configurations is established by Catia. After gridding model and importing it into the fluent software of ANSYS, the aircraft aerodynamic force and moment can be calculated.Finally, we integrate, analyze, and calculate the obtained data of those static configurations and the expressions of aerodynamic derivatives of the morphing aircraft during the whole morphing process can be derived by least square fitting of MATLAB. Therefore, we get the continuous aerodynamic parameter expressions in AppendixC.Next, we take the following morphing process (shown in Figure4) as an example. Configuration I stands for loiter mode with λ=0 and ξ=0.8. With this configuration, the aircraft is provided with a larger lift-to-drag ratio that can reduce energy consumption during loiter. Configuration II represents dash mode with λ=1 and ξ=0, and the aircraft can obtain more speed with smaller drag force.Figure 4
Morphing process.In Figure5, where the values of aerodynamic coefficients are varying along with an angle of attack from configuration I to configuration II, we can conclude that the aerodynamic coefficients of the morphing aircraft change greatly during the morphing process. For example, the lift coefficient decreases by about 53% with α=5∘ when the aircraft changes from configuration I to configuration II. The drag coefficient has little changes in the range of the small angle of attack. Furthermore, the lift-to-drag ratio varies greatly under different aerodynamic configurations, which provides the possibility for the morphing aircraft to complete different tasks. The pitching moment coefficient has a dramatic change during the morphing process, which is caused by the difference of the distance between the aerodynamic center and the mass center.Figure 5
Aerodynamic coefficients.After linearizing the morphing model equation (9) at each equilibrium reference point, the longitudinal-directional dynamics can be expressed as
(11)ẋ=Aλ,ξx+B2λ,ξu,where x=ΔV,Δα,Δq,Δθ,ΔhT and u=Δδe,ΔδtT.It is easy to verify that the elements in the state-space matrices depend on the scheduling variables, and hence we get(12)Aλ,ξ=A11λ,ξA12λ,ξ0−g00A22λ,ξ1000A32λ,ξ000001000−V00V00,B2λ,ξ=B11λ,ξ0.1425B21λ,ξ0B31λ,ξ00000.The explicit expressions of the elements in the above matrixes are given in the AppendixD.Notice that the LPV model equation (11) is constructed on a family of LTI plants linearized with respect to a set of equilibrium points of a nonlinear system that represents the flight envelope of interest. Therefore, the open-loop dynamic response of morphing process produced by these two models should be similar, which guarantees that LPV model can be a reliable alternative to describe the original nonlinear system (9).However, considering the existence ofΛ′,Λ″,Δ′, and Δ″ in a nonlinear system (6), the motion trajectories of the sweep angle and span need to be scheduled. According to the work in [2], we know that the variation rate of the morphing movement will not change the trend of dynamic responses, but it will have a marked impact on the pace of change of the dynamic responses. Taking the long-period mode V and the short-period mode q as examples, Figure 6 shows that the dynamic responses vary along with span and sweep angle in different morphing velocities. Also, we must notice that too much fast variation rate of morphing movement will definitely bring a nonignorable unsteady aerodynamic effect.Figure 6
Response variation in different morphing velocities.In what follows, the morphing time is assumed to start at 0 s and end at 10 s. Notice that the morphing bodies are driven by servo actuators, then a real-time energy-optimal trajectory generation for a servomotor system is adopted to arrange morphing mechanism motions. The servomotor system of the morphing mechanism is subject to acceleration and speed constraints. The trajectory generation is formulated as a linear-constrained optimal control problem and the Pontryagin’s maximum principle can be applied to derive necessary optimality conditions. Instead of solving multipoint boundary value problems directly, we can solve this constrained optimal control by transforming it into an optimal time interval of the speed-constrained arc and a specific acceleration-constrained optimal control problem, which can be rapidly solved using the numerical method and analytic method. Here, the energy-optimal morphing trajectory is presented without concrete derivation and the details are given in [33]. The proposed trajectory not only can guarantee the smooth morphing movement but also will save the additional energy caused by morphing mechanisms.Here we define(13)Λmax′=7deg/s,Λmax″=2deg/s2,Λmin″=−2deg/s2,Δmax′=0.31m/s,Δmax″=0.09m/s2,Δmin″=−0.09m/s2.Then, the energy-optimal morphing trajectories from configuration I to configuration II are designed and shown in Figure7.Figure 7
Energy-optimal trajectories.The distinct responses of the two kinds of models are shown in Figure8. We can see that the altitude gradually decreases and at the same time the speed increases accordingly. This is mainly due to the reduction of the wing area brought by the contraction and sweep motion of the wing. The value of lift coefficient is becoming smaller along with a reduction of the wing area. Therefore, the lift force belonging to the vertical direction is not large enough to balance the gravity, which provides the aircraft with the downward acceleration. Meantime, with the decrease of pitch angle, the component of the engine thrust has changed from the upward direction into the downward direction. And this will certainly result in the downward acceleration. During the overall morphing process, the angle of attack increases as the value of the vertical velocity component becomes larger. At the end of the morphing process, the aircraft is going to seek a new balance state where the speed and the angle of attack are larger to counteract the smaller value of lift force brought by the reduction of the wing area. Furthermore, the discrepancies of dynamic responses between the LPV model and the nonlinear model are small. The LPV model is able to capture the dynamic behavior and match the nonlinear model to some extent. Therefore, the validated LPV model can become an object for the controller synthesis.Figure 8
Nonlinear model versus the LPV model.Figure9 also shows the proportion of the additional disturbances W1, W2, and W3 in the whole dynamic model. It is easy to see that W1, W2, and W3 are not in the same order magnitude comparing with F1, F2, and F3. Therefore, we can simplify the linearization procedure by omitting these extra forces and moments. However, this will certainly propose a higher requirement for the robustness of the controller.Figure 9
Conventional dynamic versus additional dynamic.Moreover, the aerodynamics in wing morphing obtained by CFD is regarded as being quasi steady. Unfortunately, the unsteady aerodynamic characteristic caused by dynamic hysteresis effect does exist in the practical morphing process [34]. This will bring modeling uncertainties of the aerodynamics in the nonlinear model to a certain extent. Meanwhile, external disturbances such as gust wind or turbulent interference exert considerable influence during flight. Therefore, now we consider the more challenging case in which the LPV system is affected by parametric perturbation in the system matrices and external disturbance, that is,
(14)ẋ=Aλ,ξ+ΔAλ,ξx+B1w+B2λ,ξ+ΔB2λ,ξu.In view of this, robust control techniques for uncertain LPV system are indispensable to maintain nominal performance and robustness in the face of massive uncertain factors.


## 3.1. Longitudinal Nonlinear Dynamic Modeling for the Morphing Aircraft

While numerous morphing designs exist for advanced aircrafts, we consider a new type of large-scale morphing aircrafts [30] depicted in Figure 1(a), whose wings can extend, contract, and rotate to a certain extent.Figure 1
(a) The Morphing aircraft with variable sweep angle and variable span. (b) Depiction of the multibody representation.
(a)
(b)The Teledyne Ryan BQM-34 “Firebee”, which is designed as a high-speed target drone, is chosen as the baseline platform [31]. The aircraft is equipped with an easily removable wing that could be modified into a morphing wing. Also, many variants of this aircraft have been constructed to accommodate mission requirements such as subsonic or supersonic flight, reconnaissance, and combat. Next, we suppose that the wing of each side of the aircraft is equipped with an inner wing and an outer wing. In wing morphing process, the inner can rotate with Ob and the outer will keep level and move linearly along with the airfoil (see Figure 1(a)). Owing to this specific structure, the morphing aircraft can change the wing sweep angle and wingspan to accommodate mission requirements. The maximum wing sweep angle (Λmax) is 45° and the maximum extension (Λmax) is 2 m. The geometry of the other parts (fuselage, vertical tail, horizontal tail, etc.) remains unchanged in this paper.To simplify the modeling process, we make the following assumption:Assumption 1.
In the morphing process, the left wing and the right wing will realize the synchronous motion. The symmetry axis of the fuselage lies in the same plane as the wings on both sides.The traditional approach to modeling a conventional aircraft is to treat the whole plane as a single rigid body without regard for the complications brought by the aeroelastic problems. The morphing aircraft with large-scale varations mentioned above, however, cannot generally be modeled as a single rigid body due to the varying inertial properties caused by the changes of geometric parameters of the wings. To describe the aircraft dynamic behavior accurately, the whole structural changes of each part must be taken into account properly, which means that multibody modeling is an appropriate solution. The Kane method is a simple and efficient approach to derive the dynamic model of a multibody system [32]. Thereby, the longitudinal dynamics of the morphing motion is described by using Kane method in this paper.The aircraft under consideration consists of five separate rigid bodies (see Figure1(b)): fuselage (body 1), inner part of left wing (body 2), outer part of left wing (body 3), inner part of right wing (body 4), and outer part of right wing (body 5), and their masses are mb,m1,m2,m1,m2, respectively. The important parameters of the morphing aircraft are given in Table 1.Table 1
b
0.3 m
Λ
0~45 deg
Δ
0~2 m
Jy
3107.53 kg⋅m2
l
1
2.5 m
l
2
2.0 m
mb
771.79 kg
m
1
54 kg
m
2
13.605 kg
T
δ
T
129.27 N/%Firstly, as shown in Figure2, we set the origin of the whole aircraft body coordinate frame Obxbybzb to locate at the center the wing (the unit vectors corresponding to the axes xb,yb and zb are exb,eyb and ezb), and the ground coordinate system is described as Ogxgygzg. The coordinate frames of morphing aircraft’s fuselage and the other several moving parts of the wing are described as Oixiyizii=1,…,5, taking the center of mass of each body as the origin correspondingly. In what follows, Kane’s equations will be established in the reference frame Obxbybzb.Assumption 2.
Both the aerodynamic force and the thrust are merely acting on the fuselage (body 1). And the other parts suffer the effects of gravity only.Figure 2
Coordinate axis system.There are four steps to establish the Kane equations, as depicted in Figure3. Here, we only consider the longitudinal responses in morphing process. Thus, u,w,q are selected as independent generalized speed variables.Figure 3
Kane modeling procedure.In Figure3, FObu,FObw and FObq are generalized active forces corresponding to u,w,q. FObi and MObi are generalized active forces and active moments of the ith body projected into Obxbybzb. F~Obu,F~Obw, and F~Obq are generalized inertial forces relative to u,w,q. F~Obi and M~Obi stand for the generalized inertial force and inertial moment of the ith body. UObi~ and WObi~ are partial velocities of the ith body located into Obxbybzb. The specific definitions of these variables are given in Appendix A.According to the fourth step in Figure3, the detailed expressions of the Kane equations of longitudinal dynamics of the morphing aircraft can be written as follows:
(6)f1u̇,w,q,α,θ,Λ,Λ′,Λ″,Δ′,Δ″=0,f2ẇ,q̇,q,u,α,θ,Λ,Λ′,Δ′=0,f3ẇ,q̇,q,u,α,θ,Λ,Λ′,Δ′=0.The explicit expressions of (6) are given in Appendix B.By introducing (7), (6) can be transformed into a wind coordinate frame. Considering the kinematic equation (8), we finally get the full longitudinal dynamic equation (9).
(7)u=Vcosα,u̇=V̇cosα−Vsinαα̇,w=Vsinα,ẇ=V̇sinα+Vcosαα̇,(8)θ̇=q,ḣ=Vsinθ−α,(9)EV̇α̇q̇θ̇ḣ=F1F2F3qVsinθ−α+W1W2W300.The explicit expressions ofE,F1,F2,F3,W1,W2, and W3 are given in Appendix C.Remark 1.
By adopting the Kane modeling, only the position, velocity and acceleration ofΛ and Δ appear in the dynamic modeling equation (9). They can be computed and manipulated more practically than the position, velocity, and acceleration of the center of mass or instantaneous inertia tensor and its derivative [17]. This shows another major superiority of the Kane method. Notice that we will apply our designed LPV-based controller on the original nonlinear equation (9), and the exact values of Λ′,Λ″ and Δ′,Δ″ at every moment must be obtained so that we can substitute them into (9). Therefore, the motion trajectories of sweep angle and span need to be scheduled appropriately in advance. We will discuss it in the following section.Remark 2.
W
1
,
W
2, and W3 including the first and the second-order derivative of Λ and Δ are the force and moment variations affected by the morphing movement, which can be considered as the additional disturbances in the wing morphing dynamic response [17]. However, in the following section, these extra forces and moments are proved to have such slight influence to the entire dynamic equations that they can be ignored in the following linearization procedure.


## 3.2. Longitudinal LPV Model for the Morphing Aircraft

In this paper, the Jacobian linearization approach is adopted to transform the nonlinear model of the morphing aircraft into an LPV model. The flight condition of interest is selected as the altitudeh0=9144m and the Mach number Ma=0.5. By defining the scheduling variables λ=Λ/45° and ξ=Δ/l1, we get λ∈0,1 and ξ∈0,0.8.Assumption 3 (see [17]).
In this paper, unsteady aerodynamic effects are neglected. We only focus on the aerodynamic forces and moments introduced by the quasi-steady assumption; that is, the aerodynamic forces and moments of the morphing process in different wing configurations are nearly identical to those of the corresponding static configuration.The lift, drag, pitching moments, and thrust in (9) are given as
(10)D=12ρV2SwCDαα,L=12ρV2SwCL=12ρV2SwCLα=0+CLαα+CLδeδe,M=12ρV2Swc¯CM=12ρV2Swc¯CMα=0+CMαα+CMδeδe+c¯2VCmqq,T=Tδtδt.Choosingλ=0,0.2,0.4,…,1.0 and ξ=0,0.2,0.4,…,0.8 as the reference points, then the corresponding aerodynamic parameters for these configurations can be approximately calculated through computational fluid dynamics (CFD). Firstly, the outline dimensions of the aircraft are obtained by observation. Then the three-dimensional model of the morphing aircraft under different configurations is established by Catia. After gridding model and importing it into the fluent software of ANSYS, the aircraft aerodynamic force and moment can be calculated.Finally, we integrate, analyze, and calculate the obtained data of those static configurations and the expressions of aerodynamic derivatives of the morphing aircraft during the whole morphing process can be derived by least square fitting of MATLAB. Therefore, we get the continuous aerodynamic parameter expressions in AppendixC.Next, we take the following morphing process (shown in Figure4) as an example. Configuration I stands for loiter mode with λ=0 and ξ=0.8. With this configuration, the aircraft is provided with a larger lift-to-drag ratio that can reduce energy consumption during loiter. Configuration II represents dash mode with λ=1 and ξ=0, and the aircraft can obtain more speed with smaller drag force.Figure 4
Morphing process.In Figure5, where the values of aerodynamic coefficients are varying along with an angle of attack from configuration I to configuration II, we can conclude that the aerodynamic coefficients of the morphing aircraft change greatly during the morphing process. For example, the lift coefficient decreases by about 53% with α=5∘ when the aircraft changes from configuration I to configuration II. The drag coefficient has little changes in the range of the small angle of attack. Furthermore, the lift-to-drag ratio varies greatly under different aerodynamic configurations, which provides the possibility for the morphing aircraft to complete different tasks. The pitching moment coefficient has a dramatic change during the morphing process, which is caused by the difference of the distance between the aerodynamic center and the mass center.Figure 5
Aerodynamic coefficients.After linearizing the morphing model equation (9) at each equilibrium reference point, the longitudinal-directional dynamics can be expressed as
(11)ẋ=Aλ,ξx+B2λ,ξu,where x=ΔV,Δα,Δq,Δθ,ΔhT and u=Δδe,ΔδtT.It is easy to verify that the elements in the state-space matrices depend on the scheduling variables, and hence we get(12)Aλ,ξ=A11λ,ξA12λ,ξ0−g00A22λ,ξ1000A32λ,ξ000001000−V00V00,B2λ,ξ=B11λ,ξ0.1425B21λ,ξ0B31λ,ξ00000.The explicit expressions of the elements in the above matrixes are given in the AppendixD.Notice that the LPV model equation (11) is constructed on a family of LTI plants linearized with respect to a set of equilibrium points of a nonlinear system that represents the flight envelope of interest. Therefore, the open-loop dynamic response of morphing process produced by these two models should be similar, which guarantees that LPV model can be a reliable alternative to describe the original nonlinear system (9).However, considering the existence ofΛ′,Λ″,Δ′, and Δ″ in a nonlinear system (6), the motion trajectories of the sweep angle and span need to be scheduled. According to the work in [2], we know that the variation rate of the morphing movement will not change the trend of dynamic responses, but it will have a marked impact on the pace of change of the dynamic responses. Taking the long-period mode V and the short-period mode q as examples, Figure 6 shows that the dynamic responses vary along with span and sweep angle in different morphing velocities. Also, we must notice that too much fast variation rate of morphing movement will definitely bring a nonignorable unsteady aerodynamic effect.Figure 6
Response variation in different morphing velocities.In what follows, the morphing time is assumed to start at 0 s and end at 10 s. Notice that the morphing bodies are driven by servo actuators, then a real-time energy-optimal trajectory generation for a servomotor system is adopted to arrange morphing mechanism motions. The servomotor system of the morphing mechanism is subject to acceleration and speed constraints. The trajectory generation is formulated as a linear-constrained optimal control problem and the Pontryagin’s maximum principle can be applied to derive necessary optimality conditions. Instead of solving multipoint boundary value problems directly, we can solve this constrained optimal control by transforming it into an optimal time interval of the speed-constrained arc and a specific acceleration-constrained optimal control problem, which can be rapidly solved using the numerical method and analytic method. Here, the energy-optimal morphing trajectory is presented without concrete derivation and the details are given in [33]. The proposed trajectory not only can guarantee the smooth morphing movement but also will save the additional energy caused by morphing mechanisms.Here we define(13)Λmax′=7deg/s,Λmax″=2deg/s2,Λmin″=−2deg/s2,Δmax′=0.31m/s,Δmax″=0.09m/s2,Δmin″=−0.09m/s2.Then, the energy-optimal morphing trajectories from configuration I to configuration II are designed and shown in Figure7.Figure 7
Energy-optimal trajectories.The distinct responses of the two kinds of models are shown in Figure8. We can see that the altitude gradually decreases and at the same time the speed increases accordingly. This is mainly due to the reduction of the wing area brought by the contraction and sweep motion of the wing. The value of lift coefficient is becoming smaller along with a reduction of the wing area. Therefore, the lift force belonging to the vertical direction is not large enough to balance the gravity, which provides the aircraft with the downward acceleration. Meantime, with the decrease of pitch angle, the component of the engine thrust has changed from the upward direction into the downward direction. And this will certainly result in the downward acceleration. During the overall morphing process, the angle of attack increases as the value of the vertical velocity component becomes larger. At the end of the morphing process, the aircraft is going to seek a new balance state where the speed and the angle of attack are larger to counteract the smaller value of lift force brought by the reduction of the wing area. Furthermore, the discrepancies of dynamic responses between the LPV model and the nonlinear model are small. The LPV model is able to capture the dynamic behavior and match the nonlinear model to some extent. Therefore, the validated LPV model can become an object for the controller synthesis.Figure 8
Nonlinear model versus the LPV model.Figure9 also shows the proportion of the additional disturbances W1, W2, and W3 in the whole dynamic model. It is easy to see that W1, W2, and W3 are not in the same order magnitude comparing with F1, F2, and F3. Therefore, we can simplify the linearization procedure by omitting these extra forces and moments. However, this will certainly propose a higher requirement for the robustness of the controller.Figure 9
Conventional dynamic versus additional dynamic.Moreover, the aerodynamics in wing morphing obtained by CFD is regarded as being quasi steady. Unfortunately, the unsteady aerodynamic characteristic caused by dynamic hysteresis effect does exist in the practical morphing process [34]. This will bring modeling uncertainties of the aerodynamics in the nonlinear model to a certain extent. Meanwhile, external disturbances such as gust wind or turbulent interference exert considerable influence during flight. Therefore, now we consider the more challenging case in which the LPV system is affected by parametric perturbation in the system matrices and external disturbance, that is,
(14)ẋ=Aλ,ξ+ΔAλ,ξx+B1w+B2λ,ξ+ΔB2λ,ξu.In view of this, robust control techniques for uncertain LPV system are indispensable to maintain nominal performance and robustness in the face of massive uncertain factors.


## 4. Sliding Mode LPV Control

In order to guarantee the stability of the morphing aircraft, especially the stability during the morphing process, an appropriate controller has to be designed to achieve the following objectives:(i)
Not only at any static morphing completion state but also during the dynamic morphing process, the closed loop of the flight system should stable globally.(ii)
Keep the altitude and speed constant during the transition from configuration I to configuration II.To accomplish a smooth transition process against a great deal of uncertainty, we design the FSMC for the uncertain LPV system (14).


### 4.1. Sliding Mode Dynamics Analysis

Consider the system with the following state-space representation for the LPV system(15)ẋt=AΘ+ΔAΘxt+B1Θwt+B2ut,zt=C1Θxt+D1Θwt,where x∈Rnx is the state, w∈Rnw is the external disturbance, u∈Rnu is the control input, and z∈Rnz is the performance of the system. We denote continuous, measurable quantities that range in some admissible set Θ as scheduling variables, Θt=θ1t,θ2t,…,θnθtT. ΔAΘ represents the system matrix uncertainty.The following assumptions are made.Assumption 4.
The pairAΘ,B2 is stabilizable for all Θ, and the state x is available.Assumption 5.
The input matrixB2 has full column rank and B⊥ is any basis of the null space of B2T.Assumption 6.
There exists a known constantεΔ such that for all Θ, ΔAΘ≤εΔ. As for wt, we set εw=wtmax is the upper bound of wt which is assumed to be obtained in advance.Remark 3.
The matrixB2 is assumed to be constant. When this is not the case, low-pass filters with sufficiently large bandwidth can be used to filter the system’s inputs and hence to move all the time-varying parameters into the state matrix [12]. Therefore, even in the case when B2 is the functions of Θ, the augmented model of the system can be converted into the form (15). In consideration of the system (14), we define a new control input u~ by
(16)ẋu=Auxu+Buu~,u=Cuxu,where Au is stable. Then the augmented LPV plant is described by
(17)ẋẋu=Aλ,ξ+ΔAλ,ξB2λ,ξ+ΔB2λ,ξCu0Auxxu+B10w+0Buu~.
For now, note that the control matrixB2 is now parameter-free. The bandwidth of filter of (16) must be chosen larger than the desired system bandwidth. With this constraint, the proposed prefiltering will not significantly alter the original system (14). Therefore, the specific assumption of B2 in the system (15) is withheld.Remark 4.
As described in Assumption5, the uncertainty ΔAΘ is not necessary to satisfy the matched conditions.
To stabilize the LPV system, the FSMC technique is utilized here. First, a sliding mode surface function is designed as(18)st=B2TP−1xt,where P∈Rnx⋅nx is a positive-definite matrix which needs to be designed later.
Define a transformation matrix (M) and the associated vector (v) as
(19)M=M1M2=B⊥TPB⊥−1B⊥TB2TP−1B2−1B2TP−1,v=v1v2=Mx,where v1∈Rnx−nu and v2∈Rnu.
The “regular form” of the original system (15) is then represented as follows
(20)v̇1tv̇2t=M1AΘ+ΔAΘPB⊥M1AΘ+ΔAΘB2M2AΘ+ΔAΘPB⊥M2AΘ+ΔAΘB2v1tv2t+M1B1ΘM2B1Θwt+0M2B2ut.
Corresponding to the regular form, the reduced-order sliding mode dynamics on the sliding surfacest=0 with dimension of nx−nu can be obtained as [22, 23]
(21)v̇1t=M1AΘ+ΔAΘPB⊥⏟Acv1t+M1B1Θ⏟Bcwt,zt=C1ΘPB⊥⏟Ccv1t+D1Θwt.
In the following, a robust stability criterion will be proposed to guarantee the existence of the sliding surface (18) under the existing of ΔAΘ and wt.Theorem 1.
Given a scalarγ>0, the switching surface st exists and the sliding mode dynamics in (21) is admissible with bounded L2 gain performance γ, that is, the system in the sliding mode is robust stable and has an L2 gain less than γ for all admissible uncertainties ΔAΘ and for all exogenous input wt∈L20,∞, if there exist matrices P>0, G11,G21,…,G51,G61, G12,G22,…,G52,G62, and a scalar σ>0 such that for all Θ(22)Ω11ΘΩ12ΘΩ13ΘΩ14ΘΩ15ΘΩ16Θ∗−γI+G21C1TΘ+C1ΘG21TD1Θ+C1ΘG31TG21+C1ΘG41T−G21+C1ΘG51T−G22+C1ΘG61T∗∗−γIG31−G31−G32∗∗∗−σI+G41+G41T−G41+G51T−G42+G61T∗∗∗∗−G51−G51T−G52−G61T∗∗∗∗∗−G62−G62T<0,where
(23)Ω11Θ=G11AΘTB⊥+B⊥TAΘG11T+εΔ2G12B⊥+εΔ2B⊥TG12T,Ω12Θ=B⊥TAΘG21T+εΔ2B⊥TG22T+G11C1ΘT,Ω13Θ=B⊥TB1Θ+B⊥TAΘG31T+εΔ2B⊥TG32T,Ω15Θ=B⊥TP+εΔ2B⊥TG52T+B⊥TAΘG51T−G11,Ω14Θ=B⊥TAΘG41T+εΔ2B⊥TG42T+G11,Ω16Θ=12σB⊥T+B⊥TAΘG61T+εΔ2B⊥TG62T−G12.Proof.
By using the bounded real lemma [12], the system (21) satisfies the robust L2 performance requirement if and only if there exists a positive-define matrix Xc, such that
(24)AcTXc+XcAcXcBcCcT∗−γID1ΘT∗∗−γI<0.
By selectingXc=B⊥TPB⊥, substituting (21) into inequality (24), and applying Lemmas 1 and 2, inequality (24) holds if
(25)B⊥TAΘP+PAΘTB⊥+σεΔ2B⊥TB⊥+1σB⊥TPPB⊥B⊥TB1ΘB⊥TPC1ΘT∗−γID1ΘT∗∗−γI<0.
According to Schur complement formula and elementary transformation, inequality (25) is equivalent to
(26)B⊥TAΘP+PAΘTB⊥+σεΔ2B⊥TB⊥B⊥TPC1ΘTB⊥TB1ΘB⊥TP∗−γID1Θ0∗∗−γI0∗∗∗−σI<0.
Inequality (26) can be further expressed as
(27)INT00B⊥TB1Θ0B⊥TP12σB⊥T∗−γID1Θ000∗∗−γI000∗∗∗−σI00∗∗∗∗00∗∗∗∗∗0IN<0,where
(28)N=AΘTB⊥C1ΘT0IεΔ2B⊥000.
By using Lemma3, we have
(29)00B⊥TB1Θ0B⊥TP12σB⊥T∗−γID1Θ000∗∗−γI000∗∗∗−σI00∗∗∗∗00∗∗∗∗∗0+N−IGT+GN−I<0where
(30)G=G11TG21TG31TG41TG51TG61TG12TG22TG32TG42TG52TG62TT.
After matrix decomposition and combination, we can obtain the final inequality (22). This completes the proof.Remark 5.
In inequality (22), G is considered as a full-block matrix to get the feasible solution of LMIs. We may also simplify G by setting Gij=0 or Gij=I partly to reduce computational complexity. But this will bring conservatism definitely.
We notice that the inequality (22) has the infinite number of LMI constraints. Therefore, the main additional challenge in the LPV system sliding surface existence is hence the reduction of these to finite-dimensional LMIs that can be solved via semidefinite programming (SDP). Here, a possible approach is to transform the nominal LPV system (15) into a convex polytopic one by TP model transformation. Then based on the properties of a convex set, it can be concluded that these infinite constraints will hold for all parameter trajectories if and only if it holds for the finite vertex systems.
In TP model transformation, the step of making high-order singular value decomposition (HOSVD) to a tensor is regarded as the crucial point [35]. When constructing a polytopic LPV model by HOSVD, its computational complexity is directly proportional to the number of the reserved singular values. In the work of [36], we can find the details of HOSVD about computational complexity relaxation and the trade-off between approximation accuracy and complexity.
Now, we consider the following system matrix(31)SΘ=AΘB1ΘC1ΘD1Θ.
Using TP model transformation and discarding the smaller singular values as well as their corresponding singular vectors, one obtains a higher order form as(32)SΘ≈∑i1=1I1∑i1=1I2⋯∑iNIN∏n=1nθwn,inθnSi1,i2,…,iN,with ∑in=1Inwn,inθn=1,∀n,θn, and wn,inθn≥0,∀n,θn,in. In (32), I1,I2,…,IN stand for the number of the reserved singular values along each dimension.
The tensorSΘ is then approximately constructed from finite LTI vertex systems Si1,i2,…,iN and the weight function wn,inθn which is obtained through reserved singular vectors. Here, to guarantee that the resulting LTI vertex systems form a convex hull of the LPV system, we use the canonical form to transform the weight functions.
Consequently, we only need to consider finite number of LMI constraints corresponding to the obtained vertex systemSi1,i2,…,iN for the controller design, which means that inequality (22) in Theorem 1 can be rewritten as
(33)Ω11iΩ12iΩ13iΩ14iΩ15iΩ16i∗−γI+G21C1iT+C1iG21TD1i+C1iG31TG21+C1iG41T−G21+C1iG51T−G22+C1iG61T∗∗−γIG31−G31−G32∗∗∗−σI+G41+G41T−G41+G51T−G42+G61T∗∗∗∗−G51−G51T−G52−G61T∗∗∗∗∗−G62−G62T<0,where
(34)Ω11i=G11AiTB⊥+B⊥TAiG11T+εΔ2G12B⊥+εΔ2B⊥TG12T,Ω12i=B⊥TAiG21T+εΔ2B⊥TG22T+G11C1iT,Ω13i=B⊥TB1i+B⊥TAiG31T+εΔ2B⊥TG32T,Ω14i=B⊥TAiG41T+εΔ2B⊥TG42T+G11,Ω15i=B⊥TPi+εΔ2B⊥TG52T+B⊥TAiG51T−G11,Ω16i=12σB⊥T+B⊥TAiG61T+εΔ2B⊥TG62T−G12.Remark 6.
We apply Lemma3 to eliminate the coupling between the Lyapunov function matrix and the system matrix. Then for the purpose of reducing conservatism, different Lyapunov functions at each of the vertices of the polytope system should be employed. The details and verifications can be found in [35, 37]. Furthermore, if the solution Pi1,i2,…,iN exists at each of the vertices, then the parameter-dependent Lyapunov function can be obtained as
(35)PΘ=∑i1=1I1∑i1=1I2⋯∑iNIN∏n=1Nwn,inθnPi1,i2,…,iN.


### 4.2. SMC Law Synthesis

In this section, we shall synthesize an FSMC law, by which the trajectories of the LPV system in (15) can be driven onto the prespecified switching surface st=0 in finite time and then are maintained there for all subsequent time.Theorem 2.
Suppose that the switching surface function (18) exists, namely, the LMIs in (33) is feasible for Pi1,i2,…,iN. Then the trajectories of system (15) can be driven onto the switching surface st=0 by the following sliding mode control law:
(36)ut=−M2ΘAΘ+P′ΘPΘ−1B2M2Θ−Ixt−ηΘsignM2Θxt,where M2Θ is described in (19) with P−1 turning into PΘ−1. P′Θ has the form as P′Θ=∑i1=1I1∑i2=1I2⋯∑iN=1IN∏n=1Nwn,inθn′Pi1,i2,…,iN and ηΘ is designed as
(37)ηΘ=μ+εΔM2Θx+εwM2ΘB1Θ,
whereμ is a positive scalar.Proof.
As the surface is defined asst=0, we choose the following parameter-dependent Lyapunov function
(38)Vt=12stTB2TPΘ−1B2−1TB2TPΘ−1B2−1st=12xTM2ΘTM2Θx.
Taking the derivative ofVt, considering the SMC law in (36) and noting M2ΘB2=I and M2′Θ=M2ΘP′ΘPΘ−1B2M2Θ−I, we have
(39)V̇t=xTM2ΘTM2ΘΔAΘx+M2ΘB1Θw−ηΘsignM2Θx.
In addition, we notice thatxTM2ΘTsignM2Θx≥M2Θx and ηΘ>0. Thus,
(40)V̇t≤xTM2ΘTεΔM2Θx+εwM2ΘB1Θ−ηΘ.
Substituting (37), inequality (40) can be further reduced to
(41)V̇t≤−μxTM2ΘT=−μV1/2t.
By using Lemma4, it can be seen from (41) that there exists a time t∗∈0,V1/20/μ such that Vt=0 and consequently st=0 when t≥t∗. This means that the system trajectories can be driven onto the predefined switching surface in a finite time, thereby completing the proof.Remark 7.
To deal with the problems introduced by chattering phenomenon, we replace the switching termsign• with •/•+β, where β is an adjustable scalar. Then, one can tune the gain μ and β to get a better robust performance and reduce the chattering simultaneously.


## 4.1. Sliding Mode Dynamics Analysis

Consider the system with the following state-space representation for the LPV system(15)ẋt=AΘ+ΔAΘxt+B1Θwt+B2ut,zt=C1Θxt+D1Θwt,where x∈Rnx is the state, w∈Rnw is the external disturbance, u∈Rnu is the control input, and z∈Rnz is the performance of the system. We denote continuous, measurable quantities that range in some admissible set Θ as scheduling variables, Θt=θ1t,θ2t,…,θnθtT. ΔAΘ represents the system matrix uncertainty.The following assumptions are made.Assumption 4.
The pairAΘ,B2 is stabilizable for all Θ, and the state x is available.Assumption 5.
The input matrixB2 has full column rank and B⊥ is any basis of the null space of B2T.Assumption 6.
There exists a known constantεΔ such that for all Θ, ΔAΘ≤εΔ. As for wt, we set εw=wtmax is the upper bound of wt which is assumed to be obtained in advance.Remark 3.
The matrixB2 is assumed to be constant. When this is not the case, low-pass filters with sufficiently large bandwidth can be used to filter the system’s inputs and hence to move all the time-varying parameters into the state matrix [12]. Therefore, even in the case when B2 is the functions of Θ, the augmented model of the system can be converted into the form (15). In consideration of the system (14), we define a new control input u~ by
(16)ẋu=Auxu+Buu~,u=Cuxu,where Au is stable. Then the augmented LPV plant is described by
(17)ẋẋu=Aλ,ξ+ΔAλ,ξB2λ,ξ+ΔB2λ,ξCu0Auxxu+B10w+0Buu~.
For now, note that the control matrixB2 is now parameter-free. The bandwidth of filter of (16) must be chosen larger than the desired system bandwidth. With this constraint, the proposed prefiltering will not significantly alter the original system (14). Therefore, the specific assumption of B2 in the system (15) is withheld.Remark 4.
As described in Assumption5, the uncertainty ΔAΘ is not necessary to satisfy the matched conditions.
To stabilize the LPV system, the FSMC technique is utilized here. First, a sliding mode surface function is designed as(18)st=B2TP−1xt,where P∈Rnx⋅nx is a positive-definite matrix which needs to be designed later.
Define a transformation matrix (M) and the associated vector (v) as
(19)M=M1M2=B⊥TPB⊥−1B⊥TB2TP−1B2−1B2TP−1,v=v1v2=Mx,where v1∈Rnx−nu and v2∈Rnu.
The “regular form” of the original system (15) is then represented as follows
(20)v̇1tv̇2t=M1AΘ+ΔAΘPB⊥M1AΘ+ΔAΘB2M2AΘ+ΔAΘPB⊥M2AΘ+ΔAΘB2v1tv2t+M1B1ΘM2B1Θwt+0M2B2ut.
Corresponding to the regular form, the reduced-order sliding mode dynamics on the sliding surfacest=0 with dimension of nx−nu can be obtained as [22, 23]
(21)v̇1t=M1AΘ+ΔAΘPB⊥⏟Acv1t+M1B1Θ⏟Bcwt,zt=C1ΘPB⊥⏟Ccv1t+D1Θwt.
In the following, a robust stability criterion will be proposed to guarantee the existence of the sliding surface (18) under the existing of ΔAΘ and wt.Theorem 1.
Given a scalarγ>0, the switching surface st exists and the sliding mode dynamics in (21) is admissible with bounded L2 gain performance γ, that is, the system in the sliding mode is robust stable and has an L2 gain less than γ for all admissible uncertainties ΔAΘ and for all exogenous input wt∈L20,∞, if there exist matrices P>0, G11,G21,…,G51,G61, G12,G22,…,G52,G62, and a scalar σ>0 such that for all Θ(22)Ω11ΘΩ12ΘΩ13ΘΩ14ΘΩ15ΘΩ16Θ∗−γI+G21C1TΘ+C1ΘG21TD1Θ+C1ΘG31TG21+C1ΘG41T−G21+C1ΘG51T−G22+C1ΘG61T∗∗−γIG31−G31−G32∗∗∗−σI+G41+G41T−G41+G51T−G42+G61T∗∗∗∗−G51−G51T−G52−G61T∗∗∗∗∗−G62−G62T<0,where
(23)Ω11Θ=G11AΘTB⊥+B⊥TAΘG11T+εΔ2G12B⊥+εΔ2B⊥TG12T,Ω12Θ=B⊥TAΘG21T+εΔ2B⊥TG22T+G11C1ΘT,Ω13Θ=B⊥TB1Θ+B⊥TAΘG31T+εΔ2B⊥TG32T,Ω15Θ=B⊥TP+εΔ2B⊥TG52T+B⊥TAΘG51T−G11,Ω14Θ=B⊥TAΘG41T+εΔ2B⊥TG42T+G11,Ω16Θ=12σB⊥T+B⊥TAΘG61T+εΔ2B⊥TG62T−G12.Proof.
By using the bounded real lemma [12], the system (21) satisfies the robust L2 performance requirement if and only if there exists a positive-define matrix Xc, such that
(24)AcTXc+XcAcXcBcCcT∗−γID1ΘT∗∗−γI<0.
By selectingXc=B⊥TPB⊥, substituting (21) into inequality (24), and applying Lemmas 1 and 2, inequality (24) holds if
(25)B⊥TAΘP+PAΘTB⊥+σεΔ2B⊥TB⊥+1σB⊥TPPB⊥B⊥TB1ΘB⊥TPC1ΘT∗−γID1ΘT∗∗−γI<0.
According to Schur complement formula and elementary transformation, inequality (25) is equivalent to
(26)B⊥TAΘP+PAΘTB⊥+σεΔ2B⊥TB⊥B⊥TPC1ΘTB⊥TB1ΘB⊥TP∗−γID1Θ0∗∗−γI0∗∗∗−σI<0.
Inequality (26) can be further expressed as
(27)INT00B⊥TB1Θ0B⊥TP12σB⊥T∗−γID1Θ000∗∗−γI000∗∗∗−σI00∗∗∗∗00∗∗∗∗∗0IN<0,where
(28)N=AΘTB⊥C1ΘT0IεΔ2B⊥000.
By using Lemma3, we have
(29)00B⊥TB1Θ0B⊥TP12σB⊥T∗−γID1Θ000∗∗−γI000∗∗∗−σI00∗∗∗∗00∗∗∗∗∗0+N−IGT+GN−I<0where
(30)G=G11TG21TG31TG41TG51TG61TG12TG22TG32TG42TG52TG62TT.
After matrix decomposition and combination, we can obtain the final inequality (22). This completes the proof.Remark 5.
In inequality (22), G is considered as a full-block matrix to get the feasible solution of LMIs. We may also simplify G by setting Gij=0 or Gij=I partly to reduce computational complexity. But this will bring conservatism definitely.
We notice that the inequality (22) has the infinite number of LMI constraints. Therefore, the main additional challenge in the LPV system sliding surface existence is hence the reduction of these to finite-dimensional LMIs that can be solved via semidefinite programming (SDP). Here, a possible approach is to transform the nominal LPV system (15) into a convex polytopic one by TP model transformation. Then based on the properties of a convex set, it can be concluded that these infinite constraints will hold for all parameter trajectories if and only if it holds for the finite vertex systems.
In TP model transformation, the step of making high-order singular value decomposition (HOSVD) to a tensor is regarded as the crucial point [35]. When constructing a polytopic LPV model by HOSVD, its computational complexity is directly proportional to the number of the reserved singular values. In the work of [36], we can find the details of HOSVD about computational complexity relaxation and the trade-off between approximation accuracy and complexity.
Now, we consider the following system matrix(31)SΘ=AΘB1ΘC1ΘD1Θ.
Using TP model transformation and discarding the smaller singular values as well as their corresponding singular vectors, one obtains a higher order form as(32)SΘ≈∑i1=1I1∑i1=1I2⋯∑iNIN∏n=1nθwn,inθnSi1,i2,…,iN,with ∑in=1Inwn,inθn=1,∀n,θn, and wn,inθn≥0,∀n,θn,in. In (32), I1,I2,…,IN stand for the number of the reserved singular values along each dimension.
The tensorSΘ is then approximately constructed from finite LTI vertex systems Si1,i2,…,iN and the weight function wn,inθn which is obtained through reserved singular vectors. Here, to guarantee that the resulting LTI vertex systems form a convex hull of the LPV system, we use the canonical form to transform the weight functions.
Consequently, we only need to consider finite number of LMI constraints corresponding to the obtained vertex systemSi1,i2,…,iN for the controller design, which means that inequality (22) in Theorem 1 can be rewritten as
(33)Ω11iΩ12iΩ13iΩ14iΩ15iΩ16i∗−γI+G21C1iT+C1iG21TD1i+C1iG31TG21+C1iG41T−G21+C1iG51T−G22+C1iG61T∗∗−γIG31−G31−G32∗∗∗−σI+G41+G41T−G41+G51T−G42+G61T∗∗∗∗−G51−G51T−G52−G61T∗∗∗∗∗−G62−G62T<0,where
(34)Ω11i=G11AiTB⊥+B⊥TAiG11T+εΔ2G12B⊥+εΔ2B⊥TG12T,Ω12i=B⊥TAiG21T+εΔ2B⊥TG22T+G11C1iT,Ω13i=B⊥TB1i+B⊥TAiG31T+εΔ2B⊥TG32T,Ω14i=B⊥TAiG41T+εΔ2B⊥TG42T+G11,Ω15i=B⊥TPi+εΔ2B⊥TG52T+B⊥TAiG51T−G11,Ω16i=12σB⊥T+B⊥TAiG61T+εΔ2B⊥TG62T−G12.Remark 6.
We apply Lemma3 to eliminate the coupling between the Lyapunov function matrix and the system matrix. Then for the purpose of reducing conservatism, different Lyapunov functions at each of the vertices of the polytope system should be employed. The details and verifications can be found in [35, 37]. Furthermore, if the solution Pi1,i2,…,iN exists at each of the vertices, then the parameter-dependent Lyapunov function can be obtained as
(35)PΘ=∑i1=1I1∑i1=1I2⋯∑iNIN∏n=1Nwn,inθnPi1,i2,…,iN.


## 4.2. SMC Law Synthesis

In this section, we shall synthesize an FSMC law, by which the trajectories of the LPV system in (15) can be driven onto the prespecified switching surface st=0 in finite time and then are maintained there for all subsequent time.Theorem 2.
Suppose that the switching surface function (18) exists, namely, the LMIs in (33) is feasible for Pi1,i2,…,iN. Then the trajectories of system (15) can be driven onto the switching surface st=0 by the following sliding mode control law:
(36)ut=−M2ΘAΘ+P′ΘPΘ−1B2M2Θ−Ixt−ηΘsignM2Θxt,where M2Θ is described in (19) with P−1 turning into PΘ−1. P′Θ has the form as P′Θ=∑i1=1I1∑i2=1I2⋯∑iN=1IN∏n=1Nwn,inθn′Pi1,i2,…,iN and ηΘ is designed as
(37)ηΘ=μ+εΔM2Θx+εwM2ΘB1Θ,
whereμ is a positive scalar.Proof.
As the surface is defined asst=0, we choose the following parameter-dependent Lyapunov function
(38)Vt=12stTB2TPΘ−1B2−1TB2TPΘ−1B2−1st=12xTM2ΘTM2Θx.
Taking the derivative ofVt, considering the SMC law in (36) and noting M2ΘB2=I and M2′Θ=M2ΘP′ΘPΘ−1B2M2Θ−I, we have
(39)V̇t=xTM2ΘTM2ΘΔAΘx+M2ΘB1Θw−ηΘsignM2Θx.
In addition, we notice thatxTM2ΘTsignM2Θx≥M2Θx and ηΘ>0. Thus,
(40)V̇t≤xTM2ΘTεΔM2Θx+εwM2ΘB1Θ−ηΘ.
Substituting (37), inequality (40) can be further reduced to
(41)V̇t≤−μxTM2ΘT=−μV1/2t.
By using Lemma4, it can be seen from (41) that there exists a time t∗∈0,V1/20/μ such that Vt=0 and consequently st=0 when t≥t∗. This means that the system trajectories can be driven onto the predefined switching surface in a finite time, thereby completing the proof.Remark 7.
To deal with the problems introduced by chattering phenomenon, we replace the switching termsign• with •/•+β, where β is an adjustable scalar. Then, one can tune the gain μ and β to get a better robust performance and reduce the chattering simultaneously.


## 5. Simulation

Having formulated the approach to modeling of morphing dynamic equations and FSMC-LPV controller synthesis, we are in a position to verify the effectiveness and performance of the presented algorithm when applied in the longitudinal nonlinear model. To this end, this section is intended to illustrate that the controller can make the altitude and speed remain constant during morphing process when subjected to parameter uncertainty and external disturbance.To obtain a standard realization of the LPV plant from system (14) to system (15), we consider an L2 gain performance minimization problem with a performance output zt=IVIhT, where IV=∫eVdt and Ih=∫ehdt represent the error integrals of speed and altitude, respectively. Meantime, for turning B2 into a constant matrix, two additional state variables are introduced to filter inputs u=Δδe,ΔδtT. Here, we select the matrices of prefilter (16) as
(42)Au=−0.510−0.5,Bu=11T,Cu=I2.The rest of coefficient matrixes in (15) are set as B1Θ=I3×30T, C1Θ=0I2×2T, and D1Θ=0.In what follows, we calculate the parameter-dependent Lyapunov functionPΘ in Theorem 1. Considering the computational load and approximation accuracy, four singular values of λ dimension and three singular values of ξ dimension as well as their corresponding singular vectors are remained to obtain the minimal vertex systems. Therefore, it is concluded that the polynomial parameter-dependent LPV system (11) can be approximately presented in the convex polytopic form with minimum 4 × 3 = 12 LTI systems, that is,
(43)Sλ,ξ=Aλ,ξB2λ,ξ≈∑i4∑j3w1,iλw2,jξAi,jB2i,j,where the normalized weighting functions are given in Figure 10.Figure 10
Weight functions of the convex model.By using the YALMIP interface with the SDP solver SeduMi toolbox in MATLAB, we can calculate the minimal performanceγmin=4.6567 in Theorem 1. However, optimality is not desirable in this case since it is achieved at the expense of large entry values in the concerned matrixes. To prevent this phenomenon, we use a slightly higher value of γ. In this simulation, the desired performance is chosen as γsubopt=5.0. Then, the Lyapunov matrixes Pi,j can be calculated at each LTI vertex systems Sλ,ξi,j.For the aircraft during both morphing state and fixed state, the concerned model uncertain parameters and external disturbance are given in (44). Based on these uncertainties, the controller parameters are designed in (45).
(44)ΔCL=ΔCD=ΔCM=±30%,Δmt=ΔJy=±5%,Δ12ρV02=±20%,dFV=0.01sin2πt,dFα=0.02sin2πt,dMq=0.05sin2πt,(45)μ=5.6,εw=0.2445,εΔ=0.2,β=0.1.Here, we consider the following case: morphing movement starts from configuration I to configuration II (see Figure4) taking 10 seconds to complete. The energy-optimal motion trajectories of the sweep angle and span are just designed as Figure 7. The resulting FSMC with L2 gain performance is applied at the morphing aircraft’s wing transition phase. Comparative simulations are also developed to verify the robustness of the proposed controllers. Gain self-scheduled H∞ controllers (GS-H∞) can be designed for the same vertex LTI systems (43) based on the method which is derived from [12, 17], and, correspondingly, the suboptimal performance is chosen as γsubopt=6.0.In what follows, we perform 50 separate simulation experiments with different perturbation parameters for the two controllers. At each group of two comparative experiments; however, the values of perturbation parameters are identical. The response of the aircraft in a closed-loop simulation is shown in Figure11. In addition, Figure 12 gives the changes in elevator deflection and throttle.Figure 11
Closed-loop response during morphing process.Figure 12
Closed-loop controller outputs.It is clear that under parameter perturbation and external disturbances, the FSMC-LPV controller can guarantee the variation of speed being less than 0.1 m/s and the altitude increases or decreases 0.2 m to the maximum during the wing morphing process. They can both converge to the previous value in a short time. Due to the uncertainties, the elevator deflection and the rate of pitch may change drastically at the beginning, but they are both within acceptable ranges. Since the wing area is reducing, the angle of attack will increase smoothly to achieve a new balance. It is obvious that the FSMC-LPV controller based on the morphing aircraft’s linear model can eliminate the adverse effect of uncertainties caused by wing morphing and guarantee constant speed and altitude flight during the wing transition process.On the other hand, provided that there exist no parameters perturbation, conventional GS-H∞ LPV controller can perform as well as the FSMC-LPV controller. However, H∞ controller, although being a robust strategy especially for external disturbance, cannot restrain parameter perturbation effectively. This will cause the speed and altitude to fail and to remain in previous levels. Although the closed-loop system is stable, the changes in elevator deflection and throttle significantly increase during the morphing process, which may generate failures in the actuators. Overall, Figures 11 and 12 elaborate the control performances between the two controllers and we can conclude that the FSMC-LPV controller is more robust under parameter perturbation.


## 6. Conclusion

In this paper, a systematic method of LPV-based sliding mode controller synthesis is explored for a morphing aircraft with a variable sweep angle and a variable span. At first, the longitudinal dynamic model of the morphing aircraft is derived by using Kane method. Then the LPV model with mismatched uncertainties of the morphing aircraft is developed and it can describe the morphing phase’s complex behavior. Therefore, the following controller synthesis can be based on the obtained LPV model. Secondly, we extend the sliding mode control of existing works for LTI systems via an LPV methodology. Sufficient condition in the form of LMIs has been established so as to ensure that the dynamics in the sliding mode are robustly admissible with a boundedL2 gain performance. In order to avoid solving an infinite number of LMI constraints during the controller synthesis, we give a less conservative LMI existence condition by applying slack variables decoupling method and convex polytopic transformation. Furthermore, according to the parameter-dependent Lyapunov function stability analysis, it is verified that the system trajectories can be driven onto the predefined switching surface in a finite time by the synthesized FSMC. Finally, the sliding mode approach has been verified to maintain nominal performance and robustness in the face of parameter perturbation and external disturbance by the numerical simulations. The proposed sliding mode control synthesis for uncertain LPV system can be extended to cover a wide range of general LPV systems.



---

*Source: 8426348-2017-12-31.xml*