2. PRELIMINARIES
2.4. Mathematical Models of the Missile
2.4.1. Coordinate Systems
I) Definitions of Coordinate Frames
Before proceeding with the derivation, it is necessary to assume that the earth is an inertial reference, and unless otherwise stated the atmosphere is fixed with respect to the earth [46]. In addition, the coordinate systems adopted in the present discussion are right-handed axis systems.
i) Earth (Inertial) coordinate frame (Oxgygzg)
The origin Og is at the ground tracker. Oxg-axis is taken as north. The positive Oyg-axis points upward in the vertical plane including Oxg-axis. The positive Ozg-axis is the right direction or completes the right-handed coordinate system.
ii) Body coordinate frame (Oxbybzb)
The origin Ob is at the center of gravity of the missile. The positive Oxb-axis coincides with the center line (or longitudinal axis) of the missile or forward direction. The positive Oyb-axis points upward in the vertical plane including Oxb-axis. The Ozb-axis completes the right-handed coordinate system.
iii) Ballistic coordinate frame (Oxtytzt)
The origin Ot is at the center of gravity of the missile. The positive Oxt-axis coincides with the velocity of the missile. The positive Oyt-axis points upward in the vertical plane including Oxt-axis and perpendicular to the horizontal plane of the earth. The Ozt-axis completes the right-handed coordinate system.
iv) Wind coordinate frame (Oxvyvzv)
The origin Ov is at the center of gravity of the missile. The positive Oxv -axis coincides with the velocity of the missile. The positive Oyv-axis points upward in the vertical plane including Oxv-axis. The Ozv-axis completes the right-handed coordinate system. Note that if the wind coordinate frame is nonrotating with respect to Oxv-axis and the initial definition of the Oyv-axis is including the vertical plane of the inertial coordinate frame (or the direction
of Oyv-axis is always the same as Oyt-axis), it is equal to ballistic coordinate frame.
II) Definitions of Angles
Herein, in order to clearly understand the definitions of the angles, the plus or minus sign of each angle is according to a rule which the rotation axis directs to the reader.
i) Angles between wind frame and body frame
α: It is between Oxv-axis and the plane composed of Oxb-axis and Ozb-axis and defined the sign is positive when Oxv-axis is under that plane.
β: It is between Oxv-axis and the plane composed of Oxb-axis and Oyb-axis and defined the sign is positive when Oxv-axis is on the right of that plane.
Fig. 2.6. Definitions of α and β
ii) Angles between ballistic frame and wind frame
γv: It is between Ozv-axis and the plane composed of Oxt-axis and Oyt-axis and defined the sign is positive when Ozv-axis is on the left of that plane.
iii) Angles between inertial frame and ballistic frame
θ: It is between Oxt-axis and the plane composed of Oxg-axis and Ozg-axis and defined the sign is positive when Oxt-axis is on the above of that plane.
ψv: It is between Oxt-axis and the plane composed of Oxg-axis and Oyg-axis
Fig. 2.7. Definition of γv
and defined the sign is positive when Oxt-axis is on the left of that plane.
Fig. 2.8. Definitions of θ and ψv
iv) Angles between inertial frame and body frame
ϑ: It is between Oxb-axis and the plane composed of Oxg-axis and Oyg-axis and defined the sign is positive when Oxb-axis is on the above of that plane.
ψ: It is between Oxb-axis and the plane composed of Oxg-axis and Ozg-axis and defined the sign is positive when Oxb-axis is on the left of that plane.
γ: It is between Oyb-axis and the plane composed of Ox′
g-axis and Oyg-axis and defined the sign is positive when Oyb-axis is on the left of that plane.
III) Coordinate Transformation
Define Mk(ϕ) to be the rotation by k axis with angle ϕ which any one frame is rotated counterclockwise away from other one with. and it is called direct cosine
Fig. 2.9. Definitions of ϑ, ψ, and γ
matrix (DCM). So far as the surface-to-air missile is concerned, each one coordinate follows three rotated steps to other one: 1) yaw, 2) pitch, and 3) roll, and the derivation of the transformation only is established in the direction presented in Fig. 2.10. The DCM rotated by the axes y, z, and x will be:
My(ϕ) =
cϕ 0 −sϕ
0 1 0
sϕ 0 cϕ
Mz(ϕ) =
cϕ sϕ 0
−sϕ cϕ 0
0 0 1
(2.68)
Mx(ϕ) =
1 0 0
0 cϕ sϕ 0 −sϕ cϕ
where cϕ and sϕ denote cos ϕ and sin ϕ, respectively.
Fig. 2.10. The relationships among coordinates
And the each coordinate transformation is described separately as follows:
i) Wind frame transforms to body frame
ii) Ballistic frame transforms to wind frame
Ttv(γv) =
iii) Inertial frame transforms to ballistic frame
Tgt(θ, ψv) =
iv) Inertial frame transforms to body frame
Tgb(γ, ϑ, ψ) = 2.4.2 Rigid-Body Equations of Motion
In this section we will consider a typical missile and derive the equations of motion according to Newton’s laws. In deriving the rigid-body equations of motion, the following assumptions will be made [46], [47], [48]:
I) The missile is a rigid body, that is, the missile does not undergo any change in size and shape.
II) The missile is approximate a cylinder, that is, it is an axisymmetric or rotational symmetry missile.
III) The mass of the missile remains constant during any particular dynamic analysis.
IV) The aerodynamic forces and moments acting on the missile are invariant with the roll position of the missile relative to the free-stream velocity vector.
In addition, we note that in general, a vector Q can be transformed from a fixed frame OXY Z to a rotating coordinate system oxyz by the relation [49]
Q˙OXY Z = ˙QOxyz+ ωQ× Q (2.73)
where ωQ is the angular velocity of a rotating frame relatively to a fixed frame. Further-more, if the rotating frame stops rotating, the two frames will has the same time rate of the change of state variables. Herein, the equations of motion are derived by the kine-matics and dynamics. They will be presented in four formats: I) kinekine-matics equations of translation about mass center, II) kinematics equations of rotation about mass center, III) dynamics equations of translation about mass center, and IV) dynamics equations of rotation about mass center, respectively [46], [47], [48].
I) Kinematics of Translation about Mass Center
In engineering practice, it is the simplest criterion for describing the missile trans-lation in the ballistic frame. Denoting the angular velocity of ballistic frame rela-tively to inertial frame by Ω and the missile velocity V, the missile velocity expressed in the ballistic frame can be written in the form
mdV
dt = m (δV
δt + Ω× V )
(2.74) Let us first resolve the vector Ω and V into components Ωtx, Ωty, Ωtz and Vtx, Vty, Vtz, respectively, along the axes of the ballistic frame. Denoting by itx、jty、ktz
the corresponding to unit vectors of the ballistic frame (Oxtytzt), we write
Ω = Ωtxitx+ Ωtyjty+ Ωtzktz (2.75)
V = Vtxitx+ Vtyjty+ Vtzktz (2.76)
δV
δt = dVtx
dt itx+dVty
dt jty+dVtz
dt ktz (2.77)
where
Substituting from Eq. (2.78)into Eq. (2.76)， δV And, we have known that
Ω = ˙ψv + ˙θ (2.81)
Replacing Eq. (2.80) with Eq. (2.82), we have
Ω× V =
Hence, substituting from Eqs. (2.77) and (2.79) into Eq. (2.74), we obtain
where Ftx, Fty, Ftz are the components of net external force, which is formed by thrust, aerodynamic force, gravity, and lateral force, etc., with respect to the ballistic frame. Herein, we only analyze the first four sources of external force and the others are regarded as disturbances.
i) Thrust vector control (TVC)
The positive force of TVC Fpb is fixed in the direction of Oxb-axis; that is, Fpb =
Using Eqs. (2.69) and (2.70), the force can be projected onto the ballistic frame and denoted as Fpt.
The components of the force are defined as the positive drag force X along negative Oxv-axis, the lift force Y positive to the Oyb-axis, and the side force positive to the Ozb-axis in the wind frame. Using Eq. (2.70), we can project it onto the ballistic frame.
The components of the force are in the directions of Oyb-axis and Ozb-axis in the body frame, respectively. In the same manner as thrust vector force, we
have [
The negative force is in the direction of Oyg-axis. Using Eq. (2.71), the force can be projected onto the ballistic frame.
Substituting Eqs. (2.86)-(2.89) into Eq. (2.84), the kinematics equations of trans-lation about mass center are
II) Kinematics of Rotation about Mass Center
In engineering practice, it is the simplest criterion for describing the missile ro-tation in body frame. Denoting the angular velocity of body frame corresponding to inertial frame by ω and the angular momentum H, the kinematics equations of rotation about mass center has the form [49]
dH
dt = δH
δt + ω× H (2.91)
The vectors ω and H are divided into components ωbx, ωby, ωbz and Hbx, Hby, Hbz,respectively, along the axes of the body frame. Denoting by ibx, jby, kbz the corresponding to unit vectors of the body frame, we write
ω = ωbxibx + ωbyjby + ωbzkbz (2.92)
H = Hbxibx+ Hbyjby + Hbzkbz (2.93) The first term of the right side in Eq. (2.91) will be
δH
δt = dHbx
dt ibx+dHby
dt jby+ dHbz
dt kbz (2.94)
Besides, we have known that
H = I· ω (2.95)
where I is inertia tensor, including moments and products of inertia. According to the assumption 2, the products of inertia are zero in the body frame and Eq. (2.95) can be simplified as follows
Replacing the second term of the right side in Eq. (2.91) with Eq. (2.92) and Eq.
Substituting Eqs. (2.94) and (2.97) into Eq. (2.91), we write
where Mbx, Mby, Mbz, which are the rolling, yawing, and pitching moments, respec-tively, are the components of net external moment produced mainly by aerodynamic and lateral moments with respect to the body frame. Finally, we adjust Eq. (2.98), the kinematics equations of rotation about mass center are
III) Dynamics of Translation about Mass Center
These equations are defined in the inertial frame such that we will understand the trajectory of the missile clearly. Furthermore, we must consider the altitude of the missile when calculating air density, dynamic pressure, and aerodynamic force.
Hence, it is necessary to build these equations. In order to get these vectors, we must use Eqs. (2.69) and (2.72) to let Eq. (2.78) project into the inertial frame.
The procedure is
Expanding the above equation, we have dynamics equations of translation about
IV) Dynamics of Rotation about Mass Center
For the purpose of describing the attitude of the missile in the inertial frame, it is indispensable to construct these equations. According to the relationship between body frame and inertial frame, we have known
ω = ˙ψ + ˙ϑ + ˙γ (2.102)
where ˙ψ, ˙ϑ, and ˙γ are in the direction of Oyg-axis, Ozg′-axis, and Oxb-axis. We can modify Eq. (2.102)on the basis of the regulation of coordinate transformation in this thesis.
Inversing the above matrix and expanding Eq. (2.103), the dynamics equations of rotation about mass center will be
2.4.3 Equations of Attitude Dynamics
In guidance law design, inputs are the overloadings of the missle, while the autopilot must provide these overloadings to guidance law to successfully hit-to-kill the target.
Besides, the overloading is produced by angle of attack or sideslip angle. Herein, for convenience, we will build the equations of attitude of these two angles in the wind frame.
The angular velocity in the wind frame can be separated into two parts: 1) one is yielded by wind frame relatively to body frame; 2) the other is yielded by body frame relatively to inertial frame and projected into the wind frame. That is,
ωv = (ωvb)v+ (ωbg)v (2.105)
Then, we import Eq. (2.73) to express the acceleration of wind coordinate system relative to inertial coordinate system:
Substituting Eqs. (2.109) and (2.110) into Eq. (2.108) and multiplying m in each term
where Fvx, Fvy, Fvz are the components of net external force with respect to the wind frame, and we have presented the net external force is yielded by four parts. Now, we will analyze there force one by one as follows:
I) TVC
Using Eq. (2.69) to project into the wind frame, we obtain
The components of the aerodynamic force are below:
Using Eq. (2.69) to project into wind frame, we have
Using Eqs. (2.72) and (2.69), we get
where
mgbx mgby mgbz
=
−mgsϑ
−mgcϑcγ mgcϑsγ
(2.116)
Substituting Eqs. (2.112)-(2.115) into Eq. (2.111), the first equation ˙V is the same as Eq. (2.90). Thus, the attitude equations are
˙
α = ωbz − sβ
cβ(ωbxcα− ωbysα)− 1
mV cβ (Fpvy+ Favy+ Ftvy + mgvy) β = ω˙ bxsα+ ωbycα+ 1
mV (Fpvz+ Favz + Ftvz+ mgvz)
(2.117)
2.4.4 Model of Tail Fins
The X-tail is located in the bottom of the missile, which is based on the Patriot Advanced Capability-3 (PAC-3) or MIM-104F missile [13]. From the end to the head of the missile, the signs of the first fin deflection δ1 is on the top-left corner, and the others abide by the direction of clockwise are δ2, δ3, δ4, respectively. Besides, the relationship between fin deflections and total equivalent fin deflections of aileron deflection angle, rudder deflection angle, and elevator deflection angle is analyzed below, and Fig. 2.11 shows the X-tail physical characteristics [50].
Fig. 2.11. Force analysis of X-tail
δx = δ1+ δ2+ δ3+ δ4 (2.118)
2.4.5 Reaction-Jet Control System (RCS)
This technology has been successfully implemented in PAC-3 since the Iraq War in 2003. This system, installed in front of the center of gravity of the missile or between the center of gravity and the top of the missile, yields lateral thrust changing the missile’s attitude immediately for additional auxiliary thrust mounted [51]. It is contented 180 impulse attitude control motors (IACMs), arraying in 10 circles (each one composed of 18 IACMs), staggered distributing along the Oxb-axis equably. Note that the IACM is disposable. Define ith circle (i = 1, 2,· · · , 10) for each circle from top to the center of gravity and jth IACM (j = 1, 2,· · · , 18) for the number in each circle, the odd and even number circles are shown in Fig. 2.12.
In Fig. 2.12, the layout of the odd number circles presents that the first IACM is opposite direction to the Oyb-axis, and the number of the others follows the direction of counterclockwise, respectively. In the similar manner, the layout of the even number circles presents that the first IACM is on the left 20 degrees of the opposite direction to the Oyb-axis. The angle of each IACM is described below:
Φij =
Fig. 2.12. The layout scheme of IACMs, left side for odd number and right side for even number
where i∗ = 2 when i is odd, and i∗ = 1 when i is even. The force and moment of each (i, j) IACM is presented as
Ftbyij = KcΦijsij Ftbzij =−KsΦijsij
(2.124)
and
Mtbyij =−liFtbzij Mtbzij = liFtbyij
(2.125)
where K is the force of each IACM; li is the moment arm of ith circle from the center of gravity of the missile to the location of the ith circle; and sij is defined as if (i, j) IAMC is untapped or used once, sij = 0 and if (i, j) IAMC is opened, sij = 1. The total components of the lateral force and moment are
Ftby =∑10 i=1
∑18 j=1Ftbyij Ftbz =∑10
i=1
∑18 j=1Ftbzij
(2.126)
and
Mtby =∑10 i=1−li
∑18 j=1Ftbzij Mtbz =∑10
i=1li
∑18 j=1Ftbyij
(2.127)
CHAPTER THREE
OUTPUT TRACKING CONTROL FOR A NONLIN-EAR SYSTEM
Similar to system stabilizability analysis and synthesis, the task of output tracking has received considerable attention in both theoretical and practical industry applications [52]-[54]. The objective of output tracking control is to design a feedback law such that the output of a controlled plant can track a desired reference signal. To solve the tracking-control problem effectively, many methods and techniques have been presented. Those include regulator-based approach [55], inversion-based approach [56]-[58], Lyapunov-based approach [59], Takagi-Sugeno (T-S) fuzzy model-based approach [60] and sliding mode control-based (SMC) approach [61]-[63]. In this thesis, we will study the output tracking problem from a blended control viewpoint via the following three techniques: CSMC, TSMC and NTSMC schemes.
3.1 Problem Formulation
Consider a nonlinear control system as described by [62]
˙x = fo(x) + Go(x)u (3.1)
and y = h(x) (3.2)
where x = [x1,· · · , xn]T ∈ IRn, u = [u1,· · · , um]T ∈ IRm , and y = [y1,· · · , yv]T ∈ IRv denote the state variables, control inputs, and system outputs, respectively. The functions fo(x) ∈ IRn, Go(x) = [go1(x),· · · , gom(x)]∈ IRn×m and h(x) = [h1(x),· · · , hv(x)]T ∈ IRv are smooth functions. Our interest is to construct a control input so that the output
approaches the sliding surface and achieves the desired value. For the decoupled input-output system, the new input-output form is obtained from differentiating several times until it is related to the input. That is, differentiating the output yj with respect to time, we obtain
˙
yj =▽hj· ˙x = ▽hj· (fo+ Gou) = Lfohj(x) +
∑m i=1
Lgoihj(x)ui (3.3)
where Lfohj(x) and Lgoihj(x) are the Lie derivatives of hj with respect to fo and goi (for definition, please see e.g., [64]). If Lgihj(x) is equal to zero for i=1,· · · , m, then we have to differentiate the outputs yj repeatedly until input appears. Assume that kj is the smallest integer such that at least one of the inputs appears in y(kj j), then
yj(kj) = Lkfj
ohj(x) +
∑m i=1
LgoiLkfj−1
o hj(x)ui (3.4)
with LgoiLkfj−1
o hj(x)̸= 0 for at least one i in a neighborhood of the point x0. kj is exactly the number of times one has to differentiate yj in order to have the control u explicitly appearing, in which {k1,· · · , kv} is called the relative degree [64] of the system. We impose the following assumption:
Assumption 3.1 The System (3.1)-(3.2) has the following three properties:
I) The distribution △ := span {go1(x),· · · , gom(x)} is involutive.
II) It has relative degree {k1,· · · , kv}, that is, for all x ∈ IRn, LgoiLkfohj(x) = 0 for 1 ≤ i ≤ m, 1 ≤ j ≤ v and 0 ≤ k < kj − 1, while LgoiLkfohj(x) ̸= 0 for 1 ≤ i ≤ m, 1≤ j ≤ v and k = kj− 1.
III) The control inputs u are divided into two parts u1 ∈ IRm1 and u2 ∈ IRm2 where m1 ≥ v and m2 ≥ v.
Performing the above procedure for each output yj yields
y1(k1)
... yv(kv)
= f(x) + G(x)u (3.5)
where Equation (3.5) can also be rewritten as
Note that, we have introduced d in Eq. (3.10) to represent possible model uncertainties, measurement noise and external disturbances. In this study, we call u1 the main inputs which are continuous and u2 (with components u2j for 1≤ j ≤ m2) the auxiliary inputs which are constant during a short time period once it was triggered with the following form:
u2j :=
{ NjK if t ∈ [toj, toj +△tp]
0 elsewhere (3.11)
where K denotes the minimum level of auxiliary control force; |Nj| is an integer which represents the number of actuators in u2j being activated; toj is the time instant that the actuator u2j is triggered; and △tp denotes the time duration of the constant force.
Note that u1 suffer from the output magnitude constraints, while u2 only provide discrete values and the integer Nj, given by Eq. (3.11), satisfies |Nj| ≤ Nu, where Nu is a positive integer, i.e., Nj ∈ {0, ±1, ±2, · · · , ±Nu}. Besides, we assume that the output magnitude of the auxiliary inputs are much larger than those of the main inputs.
Before designing the control law, we have to check if the nonlinear system is minimum phase. The scalar kr = k1 +· · · + kv is called the total relative degree of the nonlinear
system [20]. The necessary and sufficient condition for the existence of a coordinate trans-formation and a feedback that can linearize the system completely from the Input/Output (I/O) point of view is the total relative degree kr being the same as the order of the sys-tem, i.e., kr = n. If kr < n, then, the nonlinear system can only be partially linearized.
In this case, the stability of the nonlinear system given by Eqs. (3.1) and (3.2) depends not only on the I/O linearized system, but also on the stability of the internal dynamics (or zero dynamics).
According to linear algebra theory, G1(x) can be expressed as G1(x) = G1v1(x)G1u1(x) where a diagonal matrix G1v1(x)∈ IRv×v and G1u1(x)∈ IRv×m1 satisfy rank (G1v1(x)) = rank (G1u1(x)) = v. Given a desired v, the minimum norm solution of u1 that satisfies v1 = G1u1(x)u1 is u1 = GT1u1(x)(
G1u1(x)GT1u1(x))−1
v1, that is, u1 is easily constructed if v1 has been designed. Note that, G1(x)u1 = G1v1v1. Therefore, we may assume without lose any generality that G1(x) is a diagonal matrix and G1(x) = diag[g11(x),· · · , g1v(x)].
In the same manner, we also may assume that G2(x) is a diagonal matrix and G2(x) = diag[g21(x),· · · , g2v(x)]. Under the setting, System (3.5) can be rewritten in more simpler form as follows:
y(kj j)= fj(x) + g1j(x)u1j + g2j(x)u2j + dj (3.12)
where fj, g1j, g2j : IRn → IR, and u1j and u2j are the jth component of u1 and u2, respectively, for j = 1,· · · , v. In addition, we define the output errors to be
ej(t) = yj(t)− yjd(t), j = 1,· · · , v (3.13)
where yjd(t) is the desired output trajectory.
The main goal of this thesis is then to design suitable control laws which integrate the
The main goal of this thesis is then to design suitable control laws which integrate the