Extreme aircraft maneuver under sudden lateral CG movement: Modeling and control

JID:AESCTE AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.1 (1-15)

Aerospace Science and Technology ••• (••••) •••–•••

1

Contents lists available at ScienceDirect

2

67 68

3

Aerospace Science and Technology

4 5

69 70 71

6

72

www.elsevier.com/locate/aescte

7

73

8

74

9

75

10

76

11 12 13

Extreme aircraft maneuver under sudden lateral CG movement: Modeling and control

14 15 16 17

81

Department of Aerospace Engineering, IIT, Kharagpur, India

83

82 84 85

a r t i c l e

i n f o

a b s t r a c t

86

21 23 24 25

87

Article history: Received 9 August 2016 Received in revised form 19 March 2017 Accepted 28 April 2017 Available online xxxx

26 27 28 29 30 31 32

79 80

19

22

78

Bijoy K. Mukherjee ∗ , Manoranjan Sinha

18 20

77

Keywords: 6-DOF dynamics Asymmetric dynamics Center-of-gravity Cobra maneuver Herbst maneuver Sliding mode control

Hitherto unaddressed issue of six degree-of-freedom transient dynamics during asymmetric ejection of stores with finite velocity, onboard a combat aircraft, is addressed and modeled from the first principle. Further, the effect of asymmetric center-of-gravity shift, post ejection of the store, on some complex high angle-of-attack maneuvers such as cobra and Herbst is also investigated. It is shown that the performance of the maneuvers drastically deteriorates when carried out with controller designed for the pre-ejection symmetric c . g. based dynamics. In order to improve the deteriorated performance, two new control schemes based on the standard sliding mode technique are proposed. The first sliding control is designed based on a simple ad-hoc model for the asymmetric dynamics, whereas the states are propagated using the exact model developed. It is shown that using this scheme the lost maneuver performance can be reasonably recovered. The second control scheme is formulated using an accurate asymmetric dynamics. This proposed control scheme almost completely recovers the original maneuver performance. © 2017 Elsevier Masson SAS. All rights reserved.

88 89 90 91 92 93 94 95 96 97 98

33

99

34

100

35

101

36 37

102

1. Introduction

38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59

Highly maneuverable and agile fighter aircraft are required to carry out complex maneuvers in extremely challenging combat scenarios. In such situations it is quite probable that the stores onboard the aircraft may not be fired in pairs due to malfunctioning of the ejection system. Quite often, a single store may be loaded along with a dummy on a wing of a combat aircraft putting unnecessary weight penalty leading to reduced performance of the aircraft. Sometimes paired firing of the stores may not be required at all. It is also possible that aircraft may sustain wing-damage under various circumstances. Asymmetric fuel consumption may further aggravate the situation. All these will contribute considerably to the lateral mass asymmetry thereby shifting the center-ofgravity (c . g.) of the aircraft out of its plane of symmetry. Although it is not conventional to carry out extreme maneuvers with such mass asymmetry, however it can not be downright ruled out that such a situation may not arise before the pilot where he needs to go for some extreme maneuvers with such mass asymmetry. This makes it an important research problem which needs to be investigated. Such ejection of stores with some finite velocity not only changes c . g. of the aircraft suddenly or continuously, but also pro-

60 61 62 63 64 65 66

*

Corresponding author. E-mail address: [email protected] (B.K. Mukherjee).

http://dx.doi.org/10.1016/j.ast.2017.04.030 1270-9638/© 2017 Elsevier Masson SAS. All rights reserved.

duces unwanted reaction forces and moments leading to transition dynamics which has not been addressed in literature yet. Modeling of post c . g. shift aircraft dynamics was addressed in the literature [1–3]. These works, however being the post ejection modeling, do not reflect the transition dynamics in the dynamics equations making them incomplete. Aircraft dynamics while refueling onboard as reported in [4,5] concentrates on the formation flying of the dispenser and the accepter aircraft. To the best knowledge of the authors, modeling of the transition dynamics during the sudden shift in the c . g. of a fighter aircraft due to asymmetric firing of a store with finite ejection velocity has not been addressed in literature. The transient/impulsive effects are usually neglected on the assumption that the store will be first released with zero relative velocity and then fired. This is indeed required in order to avoid the impulsive reactions. However, malfunctioning of the ejection system may result in such impulsive reaction which will unexpectedly change the intended motion of the aircraft. Therefore, modeling the transient dynamics is necessary for the real-time closed loop control purpose. The controller should be robust enough so that it can effectively reject the adverse impact of such a faulty ejection. Investigating the effects of such sudden considerable real-time lateral c . g. movement becomes all that more significant if some extreme maneuvers are carried out following the ejection. Work on the control design of aircraft for the asymmetric c . g. is scarce [2,6,7]. To handle the asymmetric dynamics in the post damage condition, a neural network and a classical adaptive con-

103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132

JID:AESCTE

AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.2 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

2

1 2

67

Nomenclature

68

3 4 5 6 7 8 9 10 11 12

C .G ., c . g . center-of-gravity (or center-of-mass in the present context) origin of the body frame (in the plane of symmetry) o 6-DOF six degree-of-freedom δe , δa , δr stabilator, aileron and rudder deflections . . . . . . . . . deg δ p , δ y , δT pitch and yaw thrust vector settings (deg), and throttle setting T r˜e  [ xe y e ze ] position of the store from o . . . . . . . . . . . . m

13

r˜cg  [ xcg

14

V˜ o  [ u

15 16 17 18 19 20

[ u

v

ω˜  [ p

y cg

zcg ]

T

position of the c . g. from o . . . . . . . . . m

T

v w ] velocity of the aircraft in the frame attached to o . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . m/s T w  ] velocity components in the frame attached to the shifted c.g. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . m/s T q r ] angular velocity vector with components along the body axes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . deg/s

V˜ R

˜ F˜ , M [I ] m, m0 S , b, c¯ q¯

α, β μ, γ , χ φ, θ, ψ sgn(.) sat (.) S MC

69

ejection velocity of the store relative to the aircraft . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . m/s total external force (N) and moment . . . . . . . . . . . . . . N m inertia matrix of the body and store in the frame attached to o . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . kg m2 mass of the aircraft and mass of the store . . . . . . . . . kg wing planform area (m2 ), wing span (b), and mean aerodynamic chord . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . m free stream dynamic pressure . . . . . . . . . . . . . . . . . . . N/m2 angle-of-attack (AOA) and sideslip angle . . . . . . . . . . deg velocity axis bank, flight path, and heading angles . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . deg roll, pitch, and yaw Euler angles (deg) (bank, elevation, and azimuth) signum function saturation function sliding mode control

21

24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87

22 23

70

88

troller were proposed in [2] and [6] respectively. However, both the controllers were designed for the normal flight objectives such as steady wings-level trim, thereby allowing them to use the linearized dynamic model for the purpose of control design and analysis. Nonlinear control design for carrying out aggressive high angle-of-attack (AOA) maneuvers with lateral mass asymmetry, as reported in the present paper, is a completely novel work. The proposed six degree-of-freedom (6-DOF) equation of motion is derived from the first principle for a variable mass aircraft assuming a continuous or discrete mass ejection, from a known location on aircraft and with finite ejection/relative velocity, as presented in section 2. The equation of motion is expressed in the body frame fixed at the initial pre-ejection c . g. location lying in the plane of symmetry. This is followed by control implementation in section 3, which is divided into four subsections. In the first two subsections, a brief review of the sliding mode control (SMC) implementation of the two popular post-stall maneuvers cobra and Herbst in the nominal situation, i.e. without any c . g. movement is presented (termed as nominal SMC). Underlying logic for SMC parameter tuning based on the time scale separation between the actuator dynamics, and the inner and outer loop dynamics along with stability margin computations are also presented therein. In the next two subsections, it is first shown that the nominal SMC scheme severely under performs when used for carrying out the maneuvers with sudden faulty ejection of a store. This is followed by the two proposed new closed loop control formulations using the same standard SMC technique to improve the maneuver performance. A novel ad-hoc model is proposed in section 3.3, which simply adds an extra moment due to gravity term to the standard symmetric 6-DOF equations of motion to approximately represent the asymmetric dynamics. Sliding control is formulated using this model in the ad-hoc SMC scheme. It is shown that this scheme reasonably improves the maneuver performance under asymmetric c . g. variation without conceding any extra computational complexity for control design and implementation. To further improve the maneuver performance, in section 3.4, control design is proposed based on the aircraft dynamics expressed in the body frame attached to the post-ejection shifted c . g. location rather than the pre-ejection c . g. location as considered in section 2. This is necessary in order to avoid complexity in control formulation and computation involved with the use of the model as derived in section 2. This strategy (termed accurate SMC) though almost completely recovers the lost maneuver performance, also requires the measurements to be transferred to the new c . g . location. In both

ad-hoc and accurate SMCs, the states are, however, propagated in the simulation using the equation of motion as derived in section 2. Finally, conclusions are presented in section 4.

89 90 91 92

2. Modeling of aircraft transient dynamics under lateral c . g. variation

93 94 95

Conventionally, aircraft dynamics modeling is carried out assuming c . g. to be located in the plane of symmetry, which makes the cross products of inertia vanish except that related to the x and z body axes, resulting in simplified aircraft dynamics. However, when the c . g. lies off the plane of symmetry, the dynamics of an aircraft can be modeled either in the body reference frame centered about the actual c . g. or fixed at the reference position in the plane of symmetry. If the c . g. is located off the plane of symmetry, then the cross products of inertia about the body axes fixed in the plane of symmetry do not vanish and lead to complicated aircraft dynamics equations. On the other hand, fixing the body reference frame at the actual c . g. when it lies off the plane of symmetry is not a good idea (though conventionally followed) as its position will change continuously or discretely with time due to mass ejection. This also requires the state variables to be defined about this new location, whereas the aerodynamic quantities such as AOA, sideslip-angle, and aerodynamic forces/moments are defined with respect to the body/wind axes located in plane of symmetry (comes from wind tunnel measurements). Transforming them to the c . g. location off the plane of symmetry is not convenient for physical comprehension. Therefore, in the present section the body axis is assumed to be fixed and remain at the reference/original c . g. position in the plane of symmetry. The standard 6-DOF aircraft equations of motion are usually derived assuming the c . g. to lie on the body centerline [8,9]. This assumption is relaxed here and a general situation is considered where the c . g. is off the plane of symmetry whereas the body axes (xb y b zb ) are fixed in the plane of symmetry at the reference point o , and then a mass is ejected from it. Fig. 1 shows this situation schematically. Let (−m) the mass to be ejected be located at re from the origin of the body frame o  . Let the c . g. of the whole body be located at rcg from o (not shown in Fig. 1 for the sake of clarity) and the mass of the whole body be m together with the mass (−m) to be ejected. Let X Y Z denote the Earth fixed  be the angular velocity of the aircraft with inertial frame and ω components p, q, and r along the body axes. All other usual assumptions such as flat non-rotating Earth and rigid aircraft body

96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132

JID:AESCTE AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.3 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

3

1

67

2

68

3

69

4

70

5

71

6

72

7

73

Fig. 2. c . g. shift in the body frame due to the ejection of mass (−m) at re from o .

8 9

75

the body frame, when a mass (−m) is ejected from the body at a distance re from o , can be written as

10 11

   (m + m) rcg b = (−m) rcg − re    m rcg  ≈ (m) re − rcg

12 13 14 15 16

Fig. 1. Schematic of the primary and the ejected bodies in the inertial and body frames with the body frame fixed at o .

19 20 21 22 23 24

are, however, retained. Let dm be an elemental mass at a distance  from o and r from the origin (o) of the inertial frame. R Before ejection, linear momentum of the element dm with respect to the inertial frame is given by

 i = r˙ dm dP

    dro  dR  = dm +  dm dt  I dt 

25 26

 × R dm = V o dm + ω

28 30

33 34 35

(3)

 is fixed with respect to the body frame. Integrating over the as R entire body yields

31 32

(2)

I

27 29

(1)

i = P



× V o dm + ω



 dm R

(4)

 × rcg = m V o + m ω

(5)

36

Let the velocity of o and the angular velocity of the body frame 

37

fixed in the parent body be

38

tively after the ejection. Let the changedposition vector of the c . g.  in the body frame after the ejection be rcg + rcg then the total linear momentum can be written as

39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56



65 66



(6)

(7)

Dividing throughout by t and letting t → 0

F ext

(11)

F ext

I

˙ω ˙ V e ˙ × rcg + m  × rcg − m +mω

 drcg   ˙ V o + m ω  × V o + m × = m V˙ o b + m ω dt b   ˙ω ˙ V e ˙ × rcg + m  × rcg + m ω × ω  × rcg − m +mω

Therefore, in the limit t → 0

(8)

(9)

˙ where V o b is the acceleration of point o in the body frame. Referring to Fig. 2, the change in the position vector of the c . g. in

77 79 80 81 82

   drcg   =m ˙ re − rcg m  dt

83

(12)

84 85

b

86

˙ < 0. Substitution of Eq. (12) into Eq. (9) yields where m

   ˙ ˙ V o + ω ˙ re − rcg × m  × V o + m F ext = m V o b + m ω   ˙ω ˙ V e ˙ × rcg + m  × rcg + m ω × ω  × rcg − m +mω ˙ r ˙ω ˙ V o + m ˙ω  × V o + m  × re −  × F ext = m V o b + m ω m  cg    ˙ rcg − m ˙ω ˙ V e  × rcg + m ω × ω ×  × rcg +  +mω m    ˙  × rcg + m ω  × V o + m ω × ω ˙ × rcg F ext = m V o b + m ω 

 ˙ V e − V o + ω  × re −m   ˙  × rcg  × V o + m ω × ω F ext = m V o b + m ω ˙ V R ˙ × rcg − m +mω

87 88 89

(13)

90 91 92

(14)

93 94 95

(15)

96 97 98 99

(16)

100

where V R is the ejection velocity relative to the parent body. The equation of motion for the rotational dynamics can be derived as follows. The equation for the rate of change of angular momentum can be written as

101

     dh r × d F =  dt 

102 103 104 105 106

(17)

107 108

I

  is the angular momentum of the body about the oriwhere h gin (o) of the inertial frame. The term r × d F denotes the total external moment acting on the body about the origin (o) of the inertial frame. The angular momentum of the elemental mass dm about (o) (Fig. 1a) before ejection is given by

109

  = r × V dm dh i    × R dm = r × V o + ω

(18)

115

(19)

117

110 111 112 113 114 116 118





  r × ω  × R dm (20)         ro + R × V o dm + ro + R × ω  × R dm = (21)  



  dm × V o + ro × ω  dm × = m ro × V o + R R    × ω  × R dm R (22) +

 = h i

     drcg  d Vo  dP   ˙ V o + m ω × =  =m  +m dt  dt  dt  I I

60

64



ω + ω respec-

 × rcg  × rcg + m ω  P = m  V o + m V o + m ω  × rcg + (−m) V e + m ω

59

63



where V e is the velocity of the mass (−m) with respect to the inertial frame after the ejection and m < 0. Subtracting Eq. (5) from Eq. (6) and neglecting the second and higher order terms yields

58

62

and

   f = (m + m) V o +  V o + (m + m) ω  + ω  P   × rcg + rcg + (−m) V e

57

61

V o +  V o

76 78

(10)

b

17 18

74

r × V o dm +





 

  × rcg + I ω = m ro × V o + m rcg × V o + m ro × ω 

(23)

where I  is the inertia matrix of the whole body with respect to the body frame fixed in the plane of symmetry. Post ejection

119 120 121 122 123 124 125 126 127 128 129 130 131 132

JID:AESCTE

AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.4 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

4

1 2

angular momentum of the whole system, i.e. the parent body and the ejected body combined is given by

3 4 5

 = h f

8 9

rne w



  × V ne w dm + ro + re × −m V e

 (24)

p .b.

6 7



where p .b. denotes the parent body (excluding the ejected mass). The new position and velocity of the element with respect to the inertial frame rne w and V ne w are given by

10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30

rne w V ne w

= ro + ro + R      × R  + ω = V o +  V o + ω

 = h f

        ro + ro + R × V o +  V o + ω  × R dm  + ω

p .b.

  + (−m) ro + re × V e       = ro + ro + R × V o +  V o dm h f       ro + ro + R × ω  × R dm  + ω + p .b.

  + (−m) ro + re × V e       = ro + ro × V o +  V o h dm f ⎛ ⎜ +⎝

35



38

⎡

 ⎢    + ω  × + ro + ro × ⎣ ω

39



41

+

× R







44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66



= ro × F ext +

 p .b.

ω + ω × R dm

p .b.

  + (−m) ro + re × V e (29)       = (m + m) ro + ro × V o +  V o h f     + (m + m) rcg + rcg × V o +  V o        + ω  × rcg + rcg (m + m) + ro + ro × ω     + ω  + I  + I  ω   + (−m) ro + re × V e (30) Neglecting the second and higher order terms, the above equation gets reduced to





  = m ro × V o + m ro × V o + m ro ×  V o + m ro × V o h f  + m rcg × V o + m rcg × V o + m rcg ×  V o    × rcg  × rcg + m ω + m rcg × V o + ro × m ω   × rcg + m ω  × rcg + mω         + I  ω  + I  ω   × rcg + I ω + m ro × ω   + (−m) ro + re × V e



 × d F R

(35) (36)

       drcg  d Vo  dh  ˙ V o + m ω ×  = ro × m  +m dt  dt  dt  I I I  ˙  ext ˙ω ˙ V e + M  × rcg − m  × rcg + m + mω

73 74 75 76 77 79 80 81 82 83 84 85 86 87 88 89 90 91 92

96 98 99 100 101 102 103 104 105 106 107 108

(37)

109

Comparing Eqs. (33) and (37) and canceling their common terms results in

111

 ext M

    drcg  d Vo   × V o + m rcg ×  × rcg =m  + m V o × ω  dt I dt  I     d I  ˙ rcg × V o − m ˙ re × V e ˙ + m  + I ω + (38)  ω dt  I

Expressing the above equation in body frame

 ext M

˙ re × V e −m

110 112 113 114 115 116 117 118 119 120

   drcg   × V o + m ω  × rcg × V o + m rcg × V˙ o b =m dt b      × V o + m V o × ω  × rcg + m rcg × ω     d I   ˙ rcg × V o  + I ω ˙ + m  × I ω  +ω +  ω dt 

121 122 123 124 125 126 127 128

b

(31)

72

97

 ext is the resultant external moment acting on the body where M and computed about o . Substituting the expression for F ext from Eq. (8) into Eq. (36)

⎥

71

95

(34)

 ext = ro × F ext + M

 dm⎦ R

70

94

r × d F

I

69

93

Now from Eq. (17)



68

78

      drcg  dh d Vo   × V o + m rcg    = m V o × V o + m ro ×  +m dt  dt  dt  I I I    d Vo  ˙ × rcg ×  + m ro × ω dt  I 

  drcg   + m V o × ω ×  × rcg + m ro × ω dt  I    d I   ˙ ro × V o + m ˙ rcg × V o ˙ + m  + I ω +  ω  dt I   ˙ ro × V e − m ˙ ro × ω ˙ re × V e  × rcg − m +m (33)

 = dt 

42 43

Dividing throughout by t and letting t → 0 yields

   dh

⎤

67

h  = m ro × V o + m ro ×  V o + m rcg × V o    × rcg + m rcg ×  V o + m ro × ω    × rcg + m ro × ω      + I  ω  + m ro × V o  × rcg + I  ω + m ro × ω    × rcg + m rcg × V o + m ro × ω   + (−m) ro + re × V e (32)

   dm⎟ R ⎠ × V o +  V o

p .b.

40

(28)

p .b.

⎞

36 37

(27)

p .b.

32 34

(26)

Hence

31 33

(25)

Subtracting Eq. (23) from Eq. (31) yields

(39)

Ejecting a mass (−m) from the body at a distance r e from the  (in the origin (o ) of the body frame will change the term I  ω

129 130 131 132

JID:AESCTE AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.5 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

1 2 3 4 5 6 7 8 9 10 11 12 13 14



    × re  = − (−m) re × ω  I b ω     d I   ˙ re × ω  × re  =m  ω dt 

(41)

    (( ( ˙ rcg(× ˙ re − rcg × V o + ( =m m( ω( ×( V o + m rcg × V o b    (((  o(  × V o + m V  × rcg + m rcg × ω ×(ω (

(





 

18 19

 ext M

22 23 24

 ext M

26

˙ re × V e −m      × V o + ω  × I ω = m rcg × V˙ o b + m rcg × ω

   ˙ re × V e − V o + ω ˙ − m  × re + I ω    ˙ + m rcg × ω  × V o = m rcg × V˙ o b + I  ω

30 31 32 33

(44)

where V R is the ejection velocity relative to the parent body as defined earlier. It may be noted that in this formulation the parameters m, r cg , and [ I  ] are all varying with time. They are obtained at each instant by solving the following three additional ordinary ˙ profile say, f (t ) differential equations for any given m

34 35 36 37 38 39 40 41 42 43 44 45

˙ = f (t ) m

(45)

˙ (r e − r cg ) m r˙cg ,b = m

(46)

˙I  = m ˙ S (˜re ) S T (˜re ) b

(47)

where S (.) denotes the skew-symmetric matrix representing the vector cross product operator as given below.

⎡

0   S r˜e = ⎣ − ze ye

46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

ze 0 −xe

− ye

⎤

xe ⎦ 0

(48)

Once the above equations are solved for m, r cg , and [ I  ] they are used in solving Eqs. (16) and (44). Combining all the equations yields the 6-DOF dynamics described by fifteen coupled nonlinear differential equations rather than the standard twelve ones. The force and moment balance equations can be combined and expressed in a compact matrix form as

⎤−1

⎡

⎡ ˙ ⎤ 1 0 0 0 zcg − y cg u 0 1 0 − zcg 0 xcg ⎥ ⎢ v˙ ⎥ ⎢ ⎥ ⎢ ⎢ ⎥ ⎢ 0 0 1 y − x 0 ⎥ ˙ w ⎢ ⎥ ⎢ cg cg ⎢ ˙ ⎥=⎢ 0    ⎥ ⎥ × −mzcg my cg I xx − I xy − I xz ⎢ p⎥ ⎢ ⎥ ⎣ q˙ ⎦ ⎣ mz    0 −mxcg − I xy I y y − I yz ⎦ cg  r˙ −my cg mxcg 0 − I xz − I yz I zz ⎛⎡ ⎤ 2 2 ⎜⎢ ⎜⎢ ⎜⎢ ⎜⎢ ⎝⎣

 pq + ( I  I xz yy

−qw + r v + (q + r )xcg − pqy cg − rpzcg −ru + p w − pqxcg + (r 2 + p 2 ) y cg − qrzcg − pv + qu − rpxcg − qr y cg + ( p 2 + q2 ) zcg  rp + I  (q2 − r 2 ) + m(qu − pv ) y + m(ru − p w ) z − I zz )qr − I xy cg cg yz

 qr + ( I  − I  )rp + I  (r 2 − p 2 ) + m( pv − qu )x + m(r v − qw ) z − I yz pq + I xy cg cg zz xx xz

(I 

xx

−

I

y y ) pq

−

 qr I xz

+

I yz rp

+

I

xy ( p

2

− q ) + m( p w − ru )xcg + m(qw − r v ) y cg 2

1 V m er x 1 V m er y 1 V m er z

⎤⎞

⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥+m ⎢ ˙ ⎥ ⎢ y V −z V ⎥ ⎢ e er z e er y ⎥ ⎢ ⎦ ⎣ ze V erx − xe V er z xe V er y − y e V er x

67

⎥⎟ ⎥⎟ ⎥⎟ ⎥⎟ ⎥⎟ ⎥⎟ ⎥⎟ ⎥⎟ ⎦⎠

68 69

(49)

⎡

F˜ ext

⎤

70 71 72 73 74

where the external forces and moments both having all the three components namely gravity, aerodynamic, and propulsive are given in matrix form by the following expressions

⎡

⎤

⎤

⎡

q¯ SC X TX −m g sinθ = F˜ G + F˜ A + F˜ T = ⎣ m g cosθ sinφ ⎦ + ⎣ q¯ SC Y ⎦ + ⎣ T Y ⎦ m g cosθ cosφ q¯ SC Z TZ (50)

75 76 77 78 79 80 81 82

⎥ ⎥ ⎥ ⎥ ⎦

84

˜ ext = M ˜G+M ˜ A+M ˜T M ⎡

⎤ −zcg m g cosθ sinφ + y cg m g cosθ cosφ = ⎣ −zcg m g sinθ − xcg m g cosθ cosφ ⎦ y cg m g sinθ + xcg m g cosθ sinφ ⎤ ⎡ ⎤ ⎡ q¯ SbC l MT X + ⎣ q¯ S c¯ C m ⎦ + ⎣ M T Y ⎦ q¯ SbC n MT Z

(43)

28 29

⎡

83

(42)

 ˙ re × V R  −m  × I ω +ω

27

⎤

 

˙ rcg × V o ˙ re × ω × I ω + I ω ˙ + m  × re + ω +m

17

25

M ext ,z

On substitution of Eqs. (12) and (41) into Eq. (39) yields

 ext M

1 F m ext ,x 1 F m ext , y 1 F m ext , z

⎢ ⎢ ⎢ ⎢ +⎢ ⎢ M ⎢ ext ,x ⎢ ⎣ M ext , y

(40)

b

16

21

⎡





15

20





 × re assuming that ω  − (−m) re × ω  rebody frame) to I  ω of mains unchanged before and after deletion of the mass (effect ˙ term  is already accommodated separately in the I  ω change in ω in the above equation). This can be written as

5

85 86 87 88 89

(51)

90 91 92

The kinematic equations (both translational and rotational), however, remain unchanged. They are derived in any standard textbook on the subject such as [8,9] and are given in the appendices. It is further observed that the impulsive terms do not appear in the equations of motion if the store is released without any velocity relative to the aircraft. From Eq. (49) it is clear that when a single store of mass −m0 > 0 is fired with some relative ejection velocity, the integration of the first two terms inside the parenthesis on the right hand side over an infinitesimal time interval will be zero, but the contribution due to the ejection of the store will not be zero ˙ term is present in the third term in the parenthesis. During as m the infinitesimal time interval of ejection, the change in the state variables due to the impulse will, therefore, is given by

93

⎡

106

⎤ u ⎡ m ⎢ v ⎥ ⎢ 0 ⎢ ⎥ ⎢ 0 ⎢ w ⎥ m = ⎢ ⎥ 0⎢ ⎢ p ⎥ ⎣ 0 m zcg ⎣ q ⎦ −m y cg r ⎡

0 m 0 −m zcg 0 m xcg

0 0 m m y cg −m xcg 0

V er x ⎢ V er y ⎢ ⎢ V er z ×⎢ ⎢ y e V er z − ze V er y ⎢ ⎣ ze V erx − xe V er z xe V er y − y e V er x

0 −m zcg m y cg  I xx  − I xy  − I xz

m zcg 0 −m xcg  − I xy I y y − I yz

−m y cg m xcg 0  − I xz − I yz I zz

⎤−1

95 96 97 98 99 100 101 102 103 104 105 107 108

⎥ ⎥ ⎥ ⎦

109 110 111 112

⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

94

113 114 115

(52)

116 117 118 119 120

3. Control design for high-alpha maneuvers under lateral c . g . variation

121 122 123

In this section, the effect of laterally asymmetric c . g. on the closed loop dynamics of the aircraft, especially for performing some aggressive maneuvers, is investigated. For this purpose, two well known post-stall maneuvers the cobra and the Herbst are considered. In the cobra maneuver, the aircraft is made to suddenly pitch up to 90◦ or even more and then return quickly to the initial pitch, within a few seconds. This helps the aircraft to suddenly brake, and therefore gain some advantage in a typical dog fight situation. In the Herbst maneuver, the aircraft first pitches to

124 125 126 127 128 129 130 131 132

JID:AESCTE

AID:4014 /FLA

6

[m5G; v1.215; Prn:9/05/2017; 17:04] P.6 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

1

67

2

68

3

69

4

70

5

71

6

72

7

73

8

74

9

75

10

76

11

77

12

78

13

79

14

80

15 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37

a sufficiently high AOA (60◦ –70◦ ), and thereafter is commanded to bank to a high angle at a high rate to initiate a quick turn. This is followed by returning the pitch and bank angles to their initial values after a few seconds. This helps the aircraft to rapidly reverse its flight direction (mostly in the vertical plane), and therefore quickly return to the base after performing some combat mission [10,11]. While the cobra maneuver takes place predominantly in the longitudinal plane, the Herbst maneuver involves both longitudinal as well as lateral–directional channels. As both are high AOA maneuvers, the aircraft dynamics becomes highly nonlinear because of aerodynamic and trigonometric nonlinearities, kinematic and inertial coupling etc. Therefore, nonlinear controllers need to be designed for carrying out such maneuvers automatically. Nonlinear control methods have traditionally been applied for aircraft flight control by various researchers over the past two decades [12–16]. Standard SMC design formulation is briefly reviewed in the next subsection. Afterwards, the effects of lateral c . g. movement on the cobra and the Herbst maneuver performance are demonstrated and two novel means to reformulate the control design scheme so as to nullify such adverse effects are proposed and validated.

38 39

3.1. Sliding mode control formulation

40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

81

Fig. 3. Block diagram of the closed loop control scheme.

16

SMC is a robust nonlinear control design method which is particularly suitable for dealing with structured or parametric uncertainties [17]. The basic idea is to design a linear asymptotically stable error surface in the error phase space so that once the system trajectory is on the surface, error will slide along the surface to zero. A discontinuous control is applied across this surface to make sure that the trajectories starting from any arbitrary initial condition reaches the surface in finite time. This discontinuous term also ensures that the trajectory returns to the sliding surface whenever external disturbances and system uncertainties try to take it away from the error surface. Therefore, the system trajectory experiences a rough slide along the sliding surface instead of a smooth one. This phenomenon, called chattering, poses the main challenge in this method as it causes unnecessary high control activity. Chattering may also excite the high frequency dynamics of the system. There are several ways to mitigate chattering as reported in the literature [17]. In the present work, the hard discontinuous signum function is replaced by a softer saturation function to tackle this problem. Use of saturation function, however, may give rise to some small steady state errors. It is shown in [12] that the rigid aircraft dynamics can be effectively split into fast and slow variables comprising of the body angular rates (p, q, r) and the angular variables such as (α , β , μ) respectively. The slow variables are controlled by an outer loop controller which generates the required reference signals (i.e. the desired body rates) for the faster inner loop controller in a cascade

82

architecture as shown in Fig. 3 [12]. Two sliding mode controllers are designed for the slow and fast channels on the two time scales for the proposed flight dynamics. It can be shown that such a system is controllable [14]. It is assumed that the control deflections do not cause any change in the external forces (they only produce moments) thereby removing the non-minimum phase characteristics of the system, and therefore possibility of unstable internal dynamics [15]. Suitable desired profiles for α , β , and μ are fed to the outer loop controller for carrying out the specific maneuvers.



For the outer loop controller, sliding surface S˜ out  S α S β S μ = 0 is chosen element-wise as

T



S a = (ad − a) + ka

(ad − a)dt , f or a = α , β, μ

(53)

83 84 85 86 87 88 89 90 91 92 93 94 95 96 97

where ka are strictly positive constants and the subscript d denotes the desired profile. For the reachability condition to be satisfied

S˙ a = −ka sgn ( S a ) 

(54)

where ka are positive gains. Hence,

98 99 100 101 102 103

(˙ad − a˙ ) + ka (ad − a) = −ka sgn ( S a )

(55)

104

˙ μ ˙ , β, ˙ ) from the system On substitution of expressions for a˙ (i .e . α dynamics as given in Appendix B into Eq. (55), and after rearrangement, the desired body rates are obtained as

105

⎡

⎤

⎤ −1

⎡

pd − cos α tan β ⎣ qd ⎦ = ⎣ sin α rd cos α sec β

⎤

⎛⎡

⎡

S p Sq Sr

T

111 112

⎤⎞

113

(56)

116

(˙cd − c˙ ) + kc (cd − c ) = −kc sgn ( S c )

118 120

(57)

where kc are strictly positive constants and Eq. (57) is differentiated with respect to time and the reachability condition is satisfied by equating 

117 119



(cd − c )dt , f or c = p , q, r

114 115

= 0 are chosen as

S c = (cd − c ) + kc

108 110

where F α , F β , and F μ are functions of state vector and desired profiles. Similarly, for the inner loop error surfaces S˜ in 



107 109

1 − sin α tan β 0 − cos α ⎦ 0 sin α sec β

kα sgn( S α ) Fα  ⎦ ⎣ ⎝ ⎣ F β + kβ sgn( S β ) ⎦⎠ ×  Fμ kμ sgn( S μ )

106

(58)

where kc are positive constants. After substitution of the expressions for c˙ , i.e. ( p˙ , q˙ , r˙ ) from the system dynamics (refer Appendix A) into Eq. (58) and after simplifications, the final control surface deflections are computed as

121 122 123 124 125 126 127 128 129 130 131 132

JID:AESCTE AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.7 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

7

1

67

2

68

3

69

4

70

5

71

6

72

7

73

8

74

9

75

10

76

11

77

12

78

13

79

14

80

15

81

16

82

17

83

18

84

19

85

20

86

21

87

22

88

23

89

24

90

25

91

26

92

27

93

28

94

29

95

30

96

31

97

32

98

33

99

34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50

Fig. 4. Time response and control signals for the cobra maneuver when the c . g. is located in the plane of symmetry.

⎡

⎤ δe ⎡ g g 12 g 13 ⎢ δa ⎥ ⎢ ⎥ ⎣ 11 ⎢ δr ⎥ = g 21 g 22 g 23 ⎣δ ⎦ g 31 g 32 g 33 p δy ⎤ ⎡ ⎛⎡

g 14 g 24 g 34

⎤#

g 15 g 25 ⎦ g 35

⎤⎞

kp sgn( S p ) Fp × ⎝⎣ F q ⎦ + ⎣ kq sgn( S q ) ⎦⎠ Fr kr sgn( S r )

(59)

where the symbol ‘#’ denotes matrix pseudo inverse and F p , F q , and F r are functions of the state vector and desired profiles. The terms g 11 , g 12 etc. also depend on the states and can be derived in line with [12]. Pseudo-inverse method solves the optimal control allocation problem by minimizing the total control energy.

51 52

3.2. Performance of the nominal SMC

53 54 55 56 57 58 59 60 61 62 63 64 65 66

In this subsection, the cobra and the Herbst maneuvers are implemented using SMC technique as formulated in the previous subsection when the c . g. of the aircraft is at the reference point in the plane of symmetry. Let this control scheme be called the Nominal SMC scheme. The corresponding time evolutions of states and controls are shown in Figs. 4 and 5. For the cobra maneuver, a bell-shaped curve for the desired AOA profile is considered keeping sideslip and bank angles zero. These time profiles are fed to the outer loop controller as shown in Fig. 3. Maneuver duration is considered to be 6 s starting at t = 5 s and the initial trim is considered at 0.6 Mach at an altitude of 3000 m. Similarly, for the Herbst maneuver, suitable bell-shaped desired profiles are considered for the AOA and bank angle for a total maneuver duration of

18 s starting at t = 5 s with the same initial trim as considered for the cobra maneuver. The sideslip angle is commanded to be zero throughout. The aircraft is assumed to possess thrust vector control in both pitch and yaw planes apart from the four conventional controls, namely the stabilator, rudder, aileron, and throttle. However, throttle is controlled separately in an open loop fashion while the rest are controlled automatically in a closed loop manner. All the actuators are assumed to have a first order dynamics with a bandwidth as shown in Table 1. Moreover, they are limited by position and rate saturations as shown in the same table. Throttle dynamics is assumed to be ten times slower than the control surface dynamics. Further, availability of full state feedback is assumed. Figs. 4 and 5 show very good tracking of the commanded signals with negligible sideslip buildup. Herbst maneuver shows some small chattering in the aileron command at the initial stage. A ±10% random uncertainty has been assumed in the aerodynamic data of the aircraft. A larger uncertainty level may lead to more pronounced chattering. Various SMC parameters, which are tuned maintaining the time scale separation of the inner and the outer loops and the bandwidth availability of the control surfaces are listed in Table 2. To illustrate how the controller parameters are designed let us consider the Herbst maneuver problem. To remove chattering, signum functions are replaced by saturation functions whose linear region is set at ±0.05 by trial and error. As the reference profiles considered are sufficiently smooth functions, error variables more often stay within this linear band. Therefore, for the α -channel (or the longitudinal channel), the error dynamics can be expressed from Eq. (55) as

100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132

JID:AESCTE

AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.8 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

8

1

67

2

68

3

69

4

70

5

71

6

72

7

73

8

74

9

75

10

76

11

77

12

78

13

79

14

80

15

81

16

82

17

83

18

84

19

85

20

86

21

87

22

88

23

89

24

90

25

91

26

92

27

93

28

94

29

95

30

96

31

97

32

98

33

99

Fig. 5. Time response and control signals for the Herbst maneuver when the c . g. is located in the plane of symmetry.

34 35

101

Table 1 Saturation levels and bandwidths of the control effectors.

36 37 38

Position limit Rate limit Bandwidth

39 40 41

102 103

Stabilator

Aileron

Rudder

Pitch nozzle

Yaw nozzle

104

±25 deg ±60 deg/s

±25 deg ±90 deg/s

±30 deg ±90 deg/s

±20 deg ±80 deg/s

±20 deg ±80 deg/s

105

10 Hz

10 Hz

10 Hz

10 Hz

10 Hz

108

43

109

Table 2 SMC parameters.

44

110

45

ka

46

[0.01 0.01 0.01] [0.05 0.01 0.01]T

Cobra maneuver Herbst maneuver

47

k

kc

a

T

[0.1 0.1 0.1] [0.1 0.1 0.1]T T

111

k

c

[0.05 0.01 0.05] [0.05 0.01 0.05]T T

[1.5 1.0 1.5] [1.5 1.0 1.5]T T

48 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

106 107

42

49

100

112 113 114

(α˙ d − α˙ ) + kα (αd − α ) + kα sat ( S α ) = 0

115

(60)

Substituting S α from Eq. (53) in Eq. (60) results in the following equations in the linear region of the saturation function  1

e˙ α + kα e α + kα b

e α + kα

1 e¨ α + (kα + kα )˙e α + b





kα kα b

e α dt eα = 0

=0

(61) (62)

where the linear region of the saturation function is ±b and e α denotes error in α . Substituting the gains from Table 2, the time constant (τ = 1/ζ ωn ) of the second order error dynamics (Eq. (62)) is found to be 0.97 s. In the same manner, the time constant for the β –μ channel (or the lateral–directional channel) is computed to be about 0.99 s. Similarly, the time constants for the inner longitudinal loop is about 0.1 s and for the inner lateral–directional loop is about 0.07 s. Therefore, the outer loop time constants are about

ten times more than the inner loop time constants. Moreover, the inner loop time constants are nearly four to five times more than the available actuator time constant which is about 0.02 s for the given bandwidth limit as shown in Table 1. It can be shown that similar time scale separation among the actuator, inner loop, and outer loop dynamics exists for the cobra maneuver also. It may be noted from Table 1 that approximately time scale separation between the inner and outer loop gains is approximately maintained on the sliding surface also. The gains associated with the signum or the saturation functions are chosen to be much higher in the inner loop as compared to the outer loop to increase the robustness of the inner loop controller. This is because the aerodynamic uncertainties affect primarily the inner loop; robustness of the outer loop is less significant as the outer loop is predominantly kinematics. The gain and phase margins of the controllers are computed to check if the closed loop dynamics has sufficient robustness property. To experimentally compute the gain and phase margins of the

116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132

JID:AESCTE AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.9 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

9

1

67

2

68

3

69

4

70

5

71

6

72

7

73

8

74

9

75

10

76

11

77

12

78

13

79

14

80

15

81

16

82

17

Fig. 6. Oscillations in AOA after insertion of (a) gain; (b) delay (dotted – commanded values, solid – actual values).

83

18

84

19

85

20

86

21

87

22

88

23

89

24

90

25

91

26

92

27

93

28

94

29

95

30

96

31

97

32

98

33

99

34

100

35 36

101

Fig. 7. Oscillations in bank angle after insertion of (a) gain; (b) delay (dotted – commanded values, solid – actual values).

37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

102 103

nonlinear controllers, a pure gain and a pure time delay are respectively inserted one at a time in the forward path and their gain and delay values are gradually increased until sustained oscillations are observed at the output [18,19]. For the Herbst maneuver problem in the longitudinal channel, a pure gain of 3.5 yielded oscillation in AOA as shown in Fig. 6a. Therefore, the longitudinal channel has a gain margin (GM) of around 11 dB. A pure time delay of 70 ms initiated oscillation of a time period of 1.2 s in the AOA as shown in Fig. 6b. This corresponds to a phase margin (PM) of at least 200 . Similarly, for the lateral–directional channel it is found that a pure gain of 2.5 and a pure delay of 50 ms kicked off sustained oscillations in bank angle as shown in Fig. 7. Fig. 7b shows oscillations of a time period of about 1.2 s. Hence, the lateral channel has a minimum GM and PM of 8 dB and 15◦ respectively. Thus, both the channels have sufficiently high gain margins and moderate phase margins. However, higher margins were not attempted because that would definitely reduce maneuver performance level or, in other words, increase maneuver duration. Fig. 4 reveals that excellent tracking of the commanded profile is achieved with the controls remaining mostly within the saturation limits. High pitch rate induced lateral–directional motions are observed to be negligible and hence not shown. However, the lateral–directional controls are still moderate as they are less effective in the high alpha regions forcing considerable deflections even for a very small amount of moment. From Fig. 5 it is observed that the closed loop performance for the Herbst maneuver is also very good with the flight direction changing exactly by 180◦ within about 18 s. In either of the maneuvers, throttle is approx-

imately doubled by the pilot during the maneuver time to boost the thrust vector control power.

104 105 106

3.3. Performance of the ad-hoc SMC

107 108

To investigate the effects of lateral movement of c . g. on the closed loop control performance during the cobra and Herbst maneuvers, the aircraft is assumed to initially carry two identical stores weighing 900 kg each located at a lateral distance of 176 cm (which is about one third of the semispan) on both sides from its plane of symmetry and a vertical location of 45 cm downward from the nominal c . g. location. One of the stores is assumed to be fired from the port end immediately prior to the initiation of the maneuvers with a forward relative ejection velocity of about 100 m/s (which is about half the forward velocity of the aircraft) and a small downward relative ejection velocity of about 10 m/s (about 10% of the forward ejection velocity). Usually, the stores are first ejected and then fired to avoid the impulsive reaction. However, in the present study the stores are assumed to be fired with some relative velocity due to some fault in the ejection mechanism or other. The plots with dotted lines in Figs. 8 and 9 show the departure in maneuver performance when the nominal SMC as designed in the previous subsection is retained. From these figures it is evident that the adverse effects of asymmetric c . g. variation on the maneuvers is significant. In the cobra maneuver, pitch angle remains less than 70◦ (instead of going upto 90◦ ) and there arises a large lateral deviation along with a large sideslip buildup. In the Herbst maneuver, tracking error in the bank angle is excessive. Moreover, the ground track of the trajectory initially deviates

109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132

JID:AESCTE

10

AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.10 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

1

67

2

68

3

69

4

70

5

71

6

72

7

73

8

74

9

75

10

76

11

77

12

78

13

79

14

80

15

81

16

82

17

83

18

84

19

85

20

86

21

87

22

88

23

89

24

90

25

91

26

92

27

93

28

94

29

95

30

96

31

97

32

98

33

99

34

100

35

101

36

102

37

103

38

104

39

105

40

106

41

107

42

108

43

109

44

110

45

111

46

112

47

113

48

114

49

115

50

Fig. 8. Departure in the cobra maneuver performance due to the lateral shift in c . g. (- - - nominal SMC, — ad-hoc SMC).

51

117

52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

116 118

in the opposite direction before finally reversing much more than the desired 180◦ defeating the very purpose of the maneuver. As the c . g offset produces mainly rolling and partly yawing moments about the body axes, the lateral asymmetric mass distribution scenario can be modeled in an ad-hoc (or approximate) manner by simply adding to the external moments, the corresponding moment due to the gravity (along the three body axes) in the standard symmetric 6-DOF equations of motion as given in Appendix A. This additional external moment term is given by the first term on the right hand side in Eq. (51). Because of the presence of the extra lump mass of the remaining store, the inertia matrix is also modified accordingly. Unlike the exact model, as derived in Section 2, this ad-hoc model does not add any further complexity or nonlinearity to the existing 6-DOF equations of mo-

tion (as the inertial effects of the lump mass is neglected). The sliding mode controllers are designed using this ad-hoc mathematical model (denoted as the Ad-hoc SMC) and the corresponding time simulations are shown with the solid lines in Figs. 8 and 9. The states are, however, propagated considering the exact equation of motion as derived in Section 2. Control parameters and uncertainty in the aero database are kept the same as in the previous subsection. From these figures it is readily observed that considerable improvement in maneuver performance can be achieved by the ad-hoc SMC without complicating control design and implementation. In the cobra maneuver, pitch angle now goes up to about 85◦ and the heading deviation is nearly one third the Nominal SMC. In the Herbst maneuver, though the trajectory still deviates initially in the opposite direction, it finally reverses by

119 120 121 122 123 124 125 126 127 128 129 130 131 132

JID:AESCTE AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.11 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

11

1

67

2

68

3

69

4

70

5

71

6

72

7

73

8

74

9

75

10

76

11

77

12

78

13

79

14

80

15

81

16

82

17

83

18

84

19

85

20

86

21

87

22

88

23

89

24

90

25

91

26

92

27

93

28

94

29

95

30

96

31

97

32

98

33

99

34

100

35

101

36

102

37

103

38

104

39

105

40

106

41

107

42

108

43

109

44

110

45

111

46

112

47

113

48

114

49

115

50

116

51

Fig. 9. Departure in the Herbst maneuver performance due to the lateral shift in c . g. (- - - nominal SMC, — ad-hoc SMC).

52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

117 118 119

nearly 180◦ . The tracking error in bank is also drastically reduced. However, sideslip buildup is still significant.

120 121 122

3.4. Performance of the accurate SMC As stated earlier, the asymmetric 6-DOF dynamics can be alternately represented in the body frame fixed at the new shifted c . g. location (Fig. 10). This is used to design the controllers although the states are propagated using exact equations modeled about the original c . g. position in the plane of symmetry. This helps avoid complexity in the control formulation and computation involved with the use of the asymmetric dynamics as derived in Section 2. In the body frame attached at the shifted c . g., the translational motion can be described by the following equation [3]:

123 124 125 126 127 128 129 130

Fig. 10. Schematic of a rigid body with the body axes fixed at the laterally shifted c . g. location.

131 132

JID:AESCTE

AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.12 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

12

1

67

2

68

3

69

4

70

5

71

6

72

7

73

8

74

9

75

10

76

11

77

12

78

13

79

14

80

15

81

16

82

17

83

18

84

19

85

20

86

21

87

22

88

23

89

24

90

25

91

26

92

27

93

28

94

29

95

30

96

31

97

32

98

33

99

34

100

35

101

36

102

37

103

38

104

39

105

40

106

41

107

42

108

43

109

44

110

45

111

46

112

47

113

48

114

49

115

50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

Fig. 11. Departure in the cobra maneuver performance due to the lateral shift in c . g. (- - - nominal SMC, — accurate SMC).

⎡

⎤

⎤

⎡

⎡

⎤

u˙  FX −qw  + r v   ⎣ v˙ ⎦ = ⎣ −ru  + p w  ⎦ + 1 ⎣ F Y ⎦ m ˙ w − pv  + qu  FZ where



(63)

T u v  w  denotes the velocity vector in the new body T

frame and they are related with their counterparts [u v w] in the original body frame fixed in the plane of symmetry as

⎡

⎤

⎡

⎤

⎡

u u 0 ⎣ v  ⎦ = ⎣ v ⎦ + ⎣ − zcg w w y cg

zcg 0 −xcg

− y cg xcg 0

⎤⎡

⎤

p ⎦⎣ q ⎦ r

(64)

It may be noted that the aerodynamic forces and moments should be computed from the look-up tables based onα , β and not α  , β  . Using Eq. (64), linear velocity measurements which are normally

available in the frame fixed at the initial c . g. position, can be converted to the frame attached to the shifted c . g. position. With



 T

u v  w available at every time step from the above computation, Eq. (63) which has the same form as the standard symmetric dynamics (given in Appendix A), can be used for the outer loop control design. Since the angular velocity remains the same in this new frame, the rotational dynamic equations needed for the inner loop control computations do not require any modification apart from the fact that the inertia matrix is required to be computed about this new frame using the parallel axis theorem. Let the controller designed using this formulation be denoted as the Accurate SMC. The time response for the cobra and the Herbst maneuvers using this Accurate SMC formulation are shown and compared with the performance of the nominal SMC in Figs. 11 and 12. From

116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132

JID:AESCTE AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.13 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

13

1

67

2

68

3

69

4

70

5

71

6

72

7

73

8

74

9

75

10

76

11

77

12

78

13

79

14

80

15

81

16

82

17

83

18

84

19

85

20

86

21

87

22

88

23

89

24

90

25

91

26

92

27

93

28

94

29

95

30

96

31

97

32

98

33

99

34

100

35

101

36

102

37

103

38

104

39

105

40

106

41

107

42

108

43

109

44

110

45

111

46

112

47

113

48

114

49

115

50

Fig. 12. Departure in the Herbst maneuver performance due to the lateral shift in c . g. (- - - nominal SMC, — accurate SMC).

116

51

117

52

118

53 54 55 56 57 58 59 60 61

these figures it is readily observed that the closed loop system is almost completely insensitive to the lateral c . g. movement. For the cobra maneuver, there is negligible lateral deviation or sideslip buildup, and for the Herbst maneuver, the trajectory neither deviates initially in the opposite direction nor ultimately reverses by more or less than the desired 180◦ . The controller parameters and aerodynamic uncertainties considered are again kept the same as in the previous two subsections.

62 63

4. Conclusion

64 65 66

The problem of significant practical importance of aircraft flight dynamics with laterally varying c . g. was addressed in the present

work. A detailed investigation of the consequences of such a situation, which may arise because of asymmetric firing of stores and/or asymmetric structural damage, required a mathematical model of the dynamics to be developed. Such a mathematical model was developed from the variable mass perspective assuming the body reference frame to be attached at the initial c . g. location in the plane of symmetry. This generalized model not only describes the pre-ejection and post-ejection dynamics but also captures the impulsive effects during the sudden transition. Moreover, it is valid for both continuous as well as discrete expulsion of masses. The effect of such laterally asymmetric c . g. on the aircraft flight was found to be significant on comparing the closed loop performance for executing two standard post-stall maneuvers such as the co-

119 120 121 122 123 124 125 126 127 128 129 130 131 132

JID:AESCTE

AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.14 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

14

⎡ 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19

bra and the Herbst when the c . g. lies in the plane of symmetry or shifts laterally. To make the closed loop dynamics insensitive to the lateral c . g. variation, a simple ad-hoc equation of motion, which adds a moment due to gravity term to the standard symmetric 6-DOF equation of motion to describe the asymmetric aircraft dynamics, was proposed. It was shown with the help of simulations that if the closed loop control is designed considering this ad-hoc model, the lost maneuver performance can be reasonably recovered without complicating control design and computations at all. It was further shown that alternately, for control design, the dynamics can be expressed about the shifted c . g. location. The resulting control scheme, although requires transformation of the measurements to the new shifted c . g. location, does not require the complex equation of motion as derived in this paper for the control computations thereby not putting any extra computational burden which may be a significant factor so far as the real time implementation of the closed loop control scheme is concerned. This control was found to recover back the nominal maneuver performance almost completely.

⎤

μ˙ ⎣ γ˙ ⎦ χ˙ ⎡

= ⎢⎣

67 68 69

p cos α sec β + r sin α sec β −

g V

− Vg cos γ 1 mV

where

⎡ ⎣ ⎡

F Xw F Yw F Zw

1 tan cos mV 1 1 − mV sin F Yw − mV cos 1 sec F Yw − mV sin sec

tan β cos μ cos γ +

μ F Yw −

γ

μ

cos μ

γ

1 mV

(tan β + sin μ tan γ ) F Zw

μ F Zw

⎤

70

⎥ ⎦

71

γ F Zw

μ

72 73

⎤



⎦=

cos α cos β − cos α sin β − sin α

⎤ ⎡ ⎤ X˙ V cos γ cos χ ⎣ Y˙ ⎦ = ⎣ V cos γ sin χ ⎦ V sin γ h˙

sin β cos β 0

sin α cos β − sin α sin β cos α



q¯ SC X + T X q¯ SC Y + T Y q¯ SC Z + T Z



23

None declared. Appendix A. Standard equations of motion in body axes

26 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52

Following are the standard symmetric 6-DOF equations of motion of a rigid aircraft expressed in body axes [8,9]. The notations have their usual meanings.

⎡

⎤

⎤

⎡

⎡

⎤

u˙ −mg sin θ + q¯ SC X + T X −qw + r v ⎣ v˙ ⎦ = ⎣ −ru + p w ⎦ + 1 ⎣ mg cos θ sin φ + q¯ SC Y + T Y ⎦ m ˙ w − pv + qu mg cos θ cos φ + q¯ SC Z + T Z

⎡

I xx ⎣ 0 − I xz

⎡

0

− I xz

I yy 0

0 I zz

⎤⎡

⎤

p˙ ⎦ ⎣ q˙ ⎦ r˙

80 81

denotes the force components along the wind axes. Equations for angular accelerations remain the same as given in Appendix A.

⎤

⎤

⎡

55

Appendix B. Standard equations of motion in wind-body axes

57 59

6-DOF equations of motion in combined wind and body axes as given in [9,20] can be alternately expressed as:

60

⎡

61 62 63 64 65 66

83 85

Appendix C. Aircraft model

V˙

⎤

⎣ α˙ ⎦ β˙  =

87

Loaded mass: Length: Wingspan: Wing area: Mean chord: Roll-axis inertia: Pitch-axis inertia: Yaw-axis inertia:

− g sin γ + m1 F Xw g q − p cos α tan β − r sin α tan β + V sec β cos μ cos γ + g 1 p sin α − r cos α + V sin μ cos γ + mV F Yw

 1 mV

sec β

F zw

88

24000 kg 17.20 m 11.80 m 39.90 m2 4.4 m 3.512 × 104 kg m2 2.643 × 105 kg m2 2.911 × 105 kg m2

89 90 91 92 93 94 95 96 97

Aerodynamic model [21]:

98

C X = C X ,α (α ) + C X ,β (α , β) + C X ,δe (α ) δe

+ C Y ,δr (α ) δr

u ×⎣ v ⎦ w

54

58

82 84

C Y = C Y ,β (α , β) + C Y , p (α )

( I y y − I zz )qr + I xz pq q¯ SbC l + M T X = ⎣ ( I zz − I xx )rp + I xz (r 2 − p 2 ) ⎦ + ⎣ q¯ S c¯ C m + M T Y ⎦ q¯ SbC n + M T Z ( I xx − I y y ) pq − I xz qr ⎡ ⎤ ⎡ ⎤⎡ ⎤ φ˙ 1 sin φ tan θ cos φ tan θ p ⎣ θ˙ ⎦ = ⎣ 0 cos φ − sin φ ⎦ ⎣ q ⎦ r 0 sin φ sec θ cos φ sec θ ψ˙ ⎡ ⎤ ˙X ⎣ Y˙ ⎦ h˙   cos ψ cos θ − sin ψ cos φ + cos ψ sin θ sin φ sin ψ sin φ + cos ψ sin θ cos φ = sin ψ cos θ cos ψ cos φ + sin ψ sin θ sin φ − cos ψ sin φ + sin ψ sin θ cos φ sin θ − cos θ sin φ − cos θ cos φ ⎡ ⎤

53

56

77 79

Geometric model:

24

27

76

86

Conflict of interest statement

22

25

75

78

20 21

74

pb 2V

+ C Y ,r (α )

C l = C l,β (α , β) + C l, p (α )

C m = C m,α (α ) + C m,q (α )

pb

q c¯

q c¯

2V

100

+ C Y ,δa (α ) δa

101 102

+ C m,α˙ (α )

2V pb 2V

+ C Z ,δe (α ) δe

2V rb

+ C l,r (α )

2V

C n = C n,β (α , β) + C n, p (α )

+ C n,δr (α ) δr

rb

103

C Z = C Z ,α (α ) + C Z ,β (α , β) + C Z ,q (α )

+ C l,δr (α ) δr

99

2V

105 106

+ C l,δa (α ) δa

α˙ c¯ 2V rb

+ C n,r (α )

104

2V

107 108 109

+ C m,δe (α ) δe + C n,δa (α ) δa

References [1] B.J. Bacon, I.M. Gregory, General equations of motion for a damaged asymmetric aircraft, in: AIAA Atmospheric Flight Mechanics Conference and Exhibit, South Carolina, 20–23 August 2007. [2] N. Nguyen, K. Krishnakumar, J. Kaneshige, P. Nespeca, Flight dynamics and hybrid adaptive control of damaged aircraft, J. Guid. Control Dyn. 31 (3) (2008) 751–764, http://dx.doi.org/10.2514/1.28142. [3] B.K. Mukherjee, M. Sinha, Modeling and bifurcation analysis of combat aircraft dynamics under lateral C.M. Shift, in: AIAA Atmospheric Flight Mechanics Conference, San Diego, 04–08 January 2016, http://dx.doi.org/10.2514/6.2016-1037. [4] J. Waishek, A. Dogan, W. Blake, Derivation of the dynamics equations of receiver aircraft in aerial refueling, J. Guid. Control Dyn. 32 (2) (2009) 585–597, http://dx.doi.org/10.2514/1.35892. [5] A. Dogan, E. Kim, W. Blake, Control and simulation of relative motion for aerial refueling in racetrack maneuvers, J. Guid. Control Dyn. 30 (5) (2007) 1551–1557, http://dx.doi.org/10.2514/1.29487. [6] Y. Liu, G. Tao, S.M. Joshi, Modeling and model reference adaptive control of aircraft with asymmetric damage, J. Guid. Control Dyn. 33 (5) (2010) 1500–1517, http://dx.doi.org/10.2514/1.47996.

110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132

JID:AESCTE AID:4014 /FLA

[m5G; v1.215; Prn:9/05/2017; 17:04] P.15 (1-15)

B.K. Mukherjee, M. Sinha / Aerospace Science and Technology ••• (••••) •••–•••

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18

[7] J. Guo, G. Tao, Y. Liu, Multivariable adaptive control of NASA generic transport aircraft model with damage, J. Guid. Control Dyn. 34 (5) (2011) 1495–1506, http://dx.doi.org/10.2514/1.53258. [8] B. Etkin, Dynamics of Atmospheric Flight, Dover Publications Inc., New York, 2005. [9] B.L. Stevens, F.L. Lewis, Aircraft Control and Simulation, 2nd ed., Wiley India Pvt. Ltd., 2003. [10] S.H. Pourtakdoust, J. Karimi, S. Shajiee, Design of a tracking control system for an optimal post-stall maneuver using dynamic inversion approach, in: 25th International Congress of Aeronautical Sciences, Hamburg, Germany, 3–8 Sept. 2006. [11] R.Y. Chiang, M.G. Safonov, K. Haiges, K. Madden, J. Tekawt, A Fixed H ∞ controller for a super maneuverable fighter performing the herbst maneuver, Automatica 29 (1) (1993) 111–127, http://dx.doi.org/10.1016/0005-1098(93) 90176-t. [12] S.A. Snell, D.F. Enns, W.L. Garrard Jr., Nonlinear inversion flight control for a supermaneuverable aircraft, J. Guid. Control Dyn. 15 (4) (1992) 976–984, http:// dx.doi.org/10.2514/3.20932. [13] A.K. Khatri, J. Singh, N.K. Sinha, Accessible regions for controlled aircraft maneuvering, J. Guid. Control Dyn. 36 (6) (2013) 1829–1834, http://dx.doi.org/ 10.2514/1.59592. [14] S.N. Singh, Asymptotically decoupled discontinuous control of systems and

[15]

[16]

[17] [18]

[19] [20]

[21]

15

nonlinear aircraft maneuver, IEEE Trans. Aerosp. Electron. Syst. 25 (3) (1989) 380–391, http://dx.doi.org/10.1109/7.30793. S.H. Lane, R.F. Stengel, Flight control design using nonlinear inverse dynamics, Automatica 24 (4) (1988) 471–483, http://dx.doi.org/10.1016/0005-1098 (88)90092-1. D.J. Bugajski, D.F. Enns, Nonlinear control law with application to high angleof-attack flight, J. Guid. Control Dyn. 15 (3) (1992) 761–767, http://dx.doi.org/ 10.2514/3.20902. J.J.E. Slotine, W. Li, Applied Nonlinear Control, Prentice-Hall, Englewood Cliffs, NJ, 1991. S.H. Kim, M.J. Tahk, Missile acceleration controller design using proportionalintegral and non-linear dynamic control design method, Proc. Inst. Mech. Eng., G J. Aerosp. Eng. 226 (2012) 882–897, http://dx.doi.org/10.1177/ 0954410011417510. P. Zarchan, Tactical and Strategic Missile Guidance, 4th ed., American Institute of Aeronautics and Astronautics, USA, 2002. S.H. Sadati, M.S. Parvar, N.B. Menhaj, M. Bahrami, Backstepping controller design using neural networks for a fighter aircraft, Eur. J. Control 13 (5) (2007) 516–526, http://dx.doi.org/10.3166/EJC.13.516-526. Technical Report No. ADA/AMC/100, Aeronautical Development Agency, Bangalore, India, 2014.

67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84

19

85

20

86

21

87

22

88

23

89

24

90

25

91

26

92

27

93

28

94

29

95

30

96

31

97

32

98

33

99

34

100

35

101

36

102

37

103

38

104

39

105

40

106

41

107

42

108

43

109

44

110

45

111

46

112

47

113

48

114

49

115

50

116

51

117

52

118

53

119

54

120

55

121

56

122

57

123

58

124

59

125

60

126

61

127

62

128

63

129

64

130

65

131

66

132