# Computational Optimal Impact Angle Control Guidance Laws Weighted by Arbitrary Functions

Authors: Qi Chen; Jinguang Shi; Zhongyuan Wang
Journal: International Journal of Aerospace Engineering
 (2023)

---

## Abstract

This paper presents a computational impact angle control guidance law based on the energy cost weighted by arbitrary functions in order to shape the acceleration command as desired. The optimal guidance problem is established in the impact angle frame and is solved by the Gauss orthogonal collocation method. The proposed guidance law is formulated from the generalized optimal control framework, thus a new guidance law that allows achieving a specific guidance goal can be easily obtained in the way of devising a proper weighting function (smooth, piecewise, nonsmooth, or even discontinuous). This property provides more degrees of freedom in the guidance law design to accomplish a specified guidance objective. A hardware experiment is conducted to evaluate the real-time computational capacity of the proposed computational guidance law. Illustrative examples with several types of weighting functions, including smooth, piecewise, nonsmooth, and discontinuous functions, are provided to demonstrate the advantages and capacity of the proposed guidance law.


---

## Body

## 1. Introduction

The main goal of designing a guidance law for an interception missile is to obtain a minimum miss distance with respect to a given target. However, to enhance the effect of kill probability, modern tactical missiles are expected not only to intercept the target but also to achieve a desired terminal impact angle [1–3]. Driven by these practical needs, intercept angle control guidance laws have been growing rapidly and become an important research topic.In literature, a variety of guidance laws based on modern control theory have been devised to impose desired impact angle constraint, such as optimal guidance laws [4–8], sliding mode guidance laws [9–12], three-dimensional guidance laws [13–16], and L2 gain-based guidance law [17]. It is known that the guidance laws based on optimal control theory can not only impose optimality in some particular performance indexes such as control effort and miss distance but also handle the terminal constraints easily. Therefore, a variety of efficient guidance laws have been widely reported using optimal control theory. Kim and Grider [18] proposed a suboptimal guidance law for reentry vehicles with a constraint on the body attitude angle at impact, where the guidance problem was formulated as a linear quadratic problem. By solving the energy minimization optimal control problem, Ryoo et al. [19] presented an optimal guidance law with impact angle constraints for both the lag free and the first-order missile systems. Therein, three methods of calculating the time-to-go using to construct the closed form of guidance command were also proposed. An optimal guidance law for impact angle control was developed in [20] for antitank and antiship missiles against maneuvering targets. To improve the survivability and effectiveness of missiles as much as possible, an optimal control based guidance law was proposed in [21] to control impact time and impact angle. By analyzing the dynamics of time-to-go error and introducing an impact time controller, an impact angle and time control guidance law is proposed along the corresponding time-to-go estimation. Based on flight model construction, Zhu et al. [22] presented a cruising guidance law with impact time and angle control, in which the time-to-go estimation is not required.In the framework of optimal control-based guidance law design, the selection of the cost function is an important issue because the response of state variables and the guidance performance are affected by the cost function. To shape the missile’s trajectory and thus accomplish some specified guidance objectives, various weighted optimal guidance laws with impact angle constraint have been developed. In [23–25], a power of time-to-go function was adopted to weight the energy cost function, by which the robustness of the guidance law to external disturbances at the end of the homing phase was improved. For aerodynamically maneuvering ground-to-air missiles, an exponential function was employed for the weighting function of the energy cost in [26], where the penalty for late maneuvers was higher than for earlier ones, and thus the acceleration demand was distributed to compensate for the loss of aerodynamic maneuvering. Similarly, exponential weighting functions were suggested in [27] to shape guidance commands for surface-to-air missiles against high-speed incoming targets. An optimal guidance law providing some degree of freedom to prevent command saturation was devised by Lee et al. [28]. In this study, a Gaussian function was used as a weighting function to alleviate sensitivity with respect to initial heading error. Lee et al. [29] used sinusoidal functions to weight the energy performance index and proposed new guidance laws that generate a small acceleration at both the initial and final time for a stationary or slowly moving target. Xiong et al. [30] proposed an optimal impact angle control guidance law for a missile with time-varying velocity, where the impact angle costs and miss distance were weighted by the variants of hyperbolic tangent function while the energy cost was weighted by a power of time-to-go. This guidance law can produce a small acceleration command at the initial phase of the engagement.In such previous works, for achieving the specified guidance purpose, particular weighting functions were introduced to shape the missile’s trajectory or to distribute the acceleration demand during the engagement. Recently, generalized optimal impact angle guidance laws have been developed by using various approaches. In [31], the generalized formulations of optimal guidance laws were provided by solving optimal control problems with a generalized performance index, which showed that any positive function can be used as a weighting function. In [32], by taking into account the dynamic lag effect and the velocity variation, a more practical generalized optimal guidance law was proposed. In [33], the generalized guidance command is determined in a closed-loop feedback form by solving a linear quadratic optimal control problem with the energy consumption weighted by an arbitrary positive function. The above methods analytically studied the solution to the generalized optimal guidance problem, but mainly considered smooth weighting functions. If the weighting function can be devised arbitrarily (for instance, piecewise, nonsmooth, and discontinuous), the flexibility of guidance law design can be significantly increased. Consequently, this paper provides an alternative approach to solve optimal control problem with arbitrary weighting functions from a computational perspective. The main contributions of this paper are summarized as follows:
(1)
A computational optimal guidance law that has the command shaping capacity is proposed by using the collocation method. The proposed guidance law is formulated from the generalized optimal control framework and a new guidance law that allows achieving a specific guidance goal can be easily obtained by devising a proper weighting function. In addition, the proposed guidance law can be easily extended to the dynamic lag system. The real-time computational capacity of the proposed approach is evaluated on a microprocessor(2)
The results presented are the first attempts in literature to introduce piecewise weighting functions. Owing to the elegant characteristics of piecewise functions, several guidance operational goals, such as providing small acceleration command at the initial time or at the final time or both at the initial time and final time, can be easily accomplished by using one single weighting function. However, in previous studies, to accomplish the aforementioned guidance goals several distinct weighting functions are required to designed specifically(3)
The proposed guidance has the capacity to easily tackle piecewise, nonsmooth, and discontinuous weighting functions, which can provide more degrees of freedom in the guidance law design to accomplish a specified guidance objective. Therefore, hybrid functions that consist of different forms of functions during different guidance phases can be considered for the weighting function of energy cost, which can take advantages of various weighting functions during the entire guidance process and thus hold the potential to significantly enhance the guidance performanceThis paper is organized as follows. The engagement kinematics and problem formulation of optimal guidance law are described in Section2. The derivation of optimal guidance law is presented in Section 3. Section 4 presents the hardware experiment and numerical simulation results. Conclusions are presented in Section 5.


## 2. Problem Formulation

Consider the engagement geometry for a stationary or a slowly moving target shown in Figure1. In this geometry, M and T denote the missile and target, respectively. Vm, γm, and am denote the velocity, flight path angle, and normal acceleration of the missile, respectively. σ and γf represent the line-of-sight (LOS) angle and desired impact angle, respectively. The coordinate system OXIYI denotes the inertial reference frame, while the coordinate system TXfYf represents the impact angle frame. The frame TXfYf is fixed with the target and is rotated from OXIYI by the desired impact angle γf. γ¯m and σ¯ denote the fight path angle of missile and LOS angle in the impact angle frame TXfYf and are expressed as
(1)γ¯m=γm−γfσ¯=σ−γf.Figure 1
Engagement geometry.It’s seen from Equation (1) that γ¯m can be regarded as the impact angle error. From Figure 1, the equations of motion for the engagement in impact angle frame TXfYf are given by.
(2)ẏt=Vmtsinγ¯mtγ¯̇mt=amt/Vmt,where y represents the lateral distance perpendicular to the desired impact direction. Under the assumption that Vm is constant and γ¯m is small, Equation (2) can be linearized as
(3)ẏt=Vmtγ¯mtγ¯̇mt=amt/Vmt.From Equation (3), the linearized engagement kinematics can be rewritten in a compact form as follows
(4)xt=Fxt+Gutx0=x0,where
(5)x≜yγ¯mT,u≜amF≜0Vm00,G=01.Note that the linearized kinematics (4) has been widely implemented to devise guidance laws with various constraints in many literatures [19, 24, 28, 29, 31, 32].It can be seen from Figure1 that if the missile velocity vector coincides with the axis TXf, namely, y=0 and γ̅_m=0, then the impact angle interception will be guaranteed. Therefore, the terminal constraint for an interception and impact angle control are given by
(6)ytf=0γ¯mtf=0,where tf denotes the final time of engagement. Now, consider the following finite-time optimal control problem: determine the control ut that minimizes the cost function
(7)J=12xTtfSfxtf+12∫t0tfu2τSτdτ,subject to Equation (4).Here,t0 is the initial time of engagement, Sf is a positive definite matrix and is used to control the terminal states of the missile, St is a weighting function. Note that, to guarantee the impact angle interception, the value of the matrix Sf should be selected as large as possible.


## 3. Derivation of Optimal Guidance Law

### 3.1. Optimal Solution Using the Gauss Orthogonal Collocation Method

In this section, the Gauss orthogonal collocation method is taken to solve the finite-time optimal control problem as shown in Equation (7). The Hamiltonian function of the optimal control problem is
(8)H=12u2tSt+λtTFxt+Gut,where λt denotes the costate and satisfies the following adjoint equation
(9)λ̇t=−∂H∂xt=−FTλt,with the transversality condition
(10)λtf=Sfxtf.From the minimum principle, it is necessary that(11)∂H∂ut=0,which yields the following optimal control
(12)u∗t=−1St−1GTλt.Substituting Equation (12) into Equation (4), and combining with the adjoint Equation (9) yields the following two-point boundary-value problem (TPBVP)
(13)ẋt=∂H∂λt=Fxt−G1St−1GTλt,xt0=x0λ̇t=−∂H∂xt=−FTλt,λtf=Sfxtf.To implement the Gauss orthogonal collocation method, the time intervalt∈t0,tf in Equation (13) should be transformed to τ∈−1,1. By using the following affine transformation
(14)τ=2tf−t0t−tf+t0tf−t0.Equation (13) can be transformed into
(15)ẋτ=tf−t02Fxτ−G1Sτ−1GTλτ,x−1=x0λ̇τ=−tf−t02FTλτ,λ1=Sfx1.The Gauss orthogonal collocation method to solve the optimal control problem Equation (7) is to discretize and transcribe the two-point boundary-value problem into a set of algebraic equations by approximating the state and costate trajectories via interpolating polynomials. The state xτ is approximated using a basis of N+1 Lagrange interpolating polynomials,
(16)xτ≈Xτ=x−1L0τ+∑k=1NxτkLkτ,where Lkτk=0,⋯,N are defined as
(17)Lkτ=∏i=0,i≠kNτ−τiτk−τi.Note that the boundary point, -1, and theN Gauss points, τk, are used to construct the approximation of xτ. Differentiating the expression in Equation (16) gives
(18)ẋτ≈Ẋτ=x−1L̇0τ+∑k=1NxτkL̇kτ.Because the initial value ofλτ is unknown, whereas the terminal value of λτ is known and given by
(19)λ1=Sfx1.Therefore, theN Gauss points, τk, and the boundary point, 1, are used to construct the approximation of λτ as follows
(20)λτ≈Λτ=∑k=1NλτkLk∗τ+λ1LN+1∗τ,where Lk∗τ (k=1,⋯,N+1) are defined as
(21)Lk∗τ=∏i=1,i≠kN+1τ−τiτk−τi.Differentiating the expression in Equation (20) gives
(22)λ̇τ≈Λ̇τ=∑k=1NλτkL̇k∗τ+λ1L̇N+1∗τ.The terminal value of the state,x1, can be obtained in terms xτk and λτk via the Gauss quadrature [34].
(23)x1=x−1+tf−t02∑k=1NwkFxτk−G1Sτ−1GTλτk,where wk are the Gaussian weights and can be calculated offline. Substituting Equation (23) into Equation (19) yields
(24)λ1=Sfx−1+Sftf−t02∑k=1NwkFxτk−G1Sτ−1GTλτk.DenoteD¯i=L̇0τi, Dik=L̇kτi, D¯i∗=L̇N+1∗τ, and Dik∗=L̇k∗τ; Substituting Equation (24) into Equation (22) and including Equation (18) yields
(25)D¯i∗Sfx0+D¯i∗Sftf−t02∑k=1NwkFxk−G1Sk−1GTλk+∑k=1NλkDik∗=−tf−t02FTλix0D¯i+∑k=1NxkDik=tf−t02Fxi−G1Si−1GTλi,where i=1,⋯,N,xk=xτk,λk=λτk,Sk=Sτk and k=1,⋯,N.As proven in [35], the differential approximation matrices, D∈ℝN×N, D¯∈ℝN, D∗∈ℝN×N, and D¯∗∈ℝN have the following relationships:
(26)D¯i=−∑k=1NDik,Dik∗=−wkwiDkiD¯i∗=−∑k=1NDik∗.,The elements of the differential approximation matrixD can be determined offline by
(27)Dik=L̇iτk=∑l=0N∏j=0,j≠i,lNτk−τj∏j=0,j≠iNτi−τj.Consequently, the elements of the matricesD¯, D∗, and D¯∗ can be obtained offline via Equation (26).DenoteΛ~=λ1T,λ2T,⋯,λNTT∈ℝnN, X~=x1T,x2T,⋯,xNTT∈ℝnN. Equation (25) can be rewritten as
(28)D¯~∗x0+D~+∗Λ~+A~X~=0D¯~x0+D~−X~+B~Λ~=0,where
(29)D¯~∗=D¯1∗Sf⋮D¯N∗Sf∈ℝnN×n,(30)D¯~=D¯1In⋮D¯NIn∈ℝnN×N,(31)D~+∗=D11∗InD12∗In⋯D1N∗InD21∗InD22∗In⋯D2N∗In⋮⋮⋱⋮DN1∗InDN2∗In⋯DNN∗In−tf−t02D¯1∗Sfw1G1S1−1GT−FTD¯1∗Sfw2G1S2−1GT⋯D¯1∗SfwNG1SN−1GTD¯2∗Sfw1G1S1−1GTD¯2∗Sfw2G1S2−1GT−FT⋯D¯2∗SfwNG1SN−1GT⋮⋮⋱⋮D¯N∗Sfw1G1S1−1GTD¯N∗Sfw2G1S2−1GT⋯D¯N∗SfwNG1SN−1GT−FT∈ℝnN×nN,(32)D~−=D11InD12In⋯D1NInD21InD22In⋯D2NIn⋮⋮⋱⋮DN1InDN2In⋯DNNIn−tf−t02F0n⋯0n0nF⋯0n⋮⋮⋱⋮0n0n⋯F∈ℝnN×nN,(33)A~=tf−t02D¯1∗Sfw1FD¯1∗Sfw2F⋯D¯1∗SfwNFD¯2∗Sfw1FD¯2∗Sfw2F⋯D¯2∗SfwNF⋮⋮⋱⋮D¯N∗Sfw1FD¯N∗Sfw2F⋯D¯N∗SfwNF∈ℝnN×nN,(34)B~=tf−t02G1S1−1GT0n⋯0n0nG1S2−1GT⋯0n⋮⋮⋱⋮0n0n⋯G1SN−1GT∈ℝnN×nN,and 0n is an n×n zero matrix, In is an n×n identity matrix. By denoting D¯aug=D¯~∗D¯~∈ℝ2nN×n and Daug=D~+∗A~B~D~−∈ℝ2nN×2nN, Equation (28) can be expressed into the following compact form
(35)D¯augx0+DaugΛ~X~=0.Solving forΛ~ and X~ in Equation (35) gives the solution to the two-point boundary-value problem (15) as
(36)Λ~X~=−Daug‐1D¯augx0.From (36) the costate at the Gauss points can be directly determined from the initial state valuex0. However, the Gauss points do not contain the boundary points, therefore, Equation (36) cannot provide the initial value of the costate, λ−1. Because the terminal value of the costate, λ1, is given in (24), therefore, λ−1 can be calculated via the following Gauss quadrature [34].
(37)λ−1=λ1−tf−t02∑k=1NwkFTλk.According to Equations (12), (36), and (37), the optimal control at the boundary points as well as the Gauss points is obtained as
(38)ui∗=−1Si−1GTλii=0,⋯,N+1.It is very important to note that the preceding derivation does not impose any constraint onSt and the form of the equations are identical for F, G, and St, indicating that the proposed method can easily and directly handle different and complex weighting functions without any modification. This property makes the guidance law design very flexible to accomplish the specified objective. In other words, arbitrary positive function can be used as a weighting function, thus the penalty for maneuvers can be adjusted accordingly by tuning the shape of the weighting function in order to achieve some specific guidance objective. In addition, it should be also noted that in our method, the optimal control is obtained without any explicit integration or construction of transition matrices, but only by solving algebraic equations.


### 3.2. Guidance Law Implementation

Note that the optimal solution obtained in Section3.1 is open-loop for the guidance problem. To obtain the closed-loop solution, the receding-horizon technique is implemented. The receding-horizon length is chosen as the time-to-go (i.e., the time remaining to interception, and denoted as tgo), which is approximated by tgo=R/Vm. Consequently, the final time of engagement can be expressed as tf=t0+tgo. The procedure for implementing the proposed guidance law with impact angle constraint can be summarized as follows.
(1)
Initial number of the Gauss points,N, and set t0=0(2)
Solve the sequence of TPBVP(a)
Calculate the time-to-go,tgo(b)
Calculate the matricesD¯aug and D¯aug, and solve Equation (36).(c)
Update the controlui (i=0,⋯,N) using the solution from step 2b and Equation (38).(3)
Apply the control input at the current time at the first discretization point,u0, to the missile’s dynamics and update the missile’s states at the next time(4)
Repeat steps 2-3 until intercepting the targetFor complex weighting functions, the proposed method can easily obtain the optimal control via solving a set of algebraic equations and does not require any explicit propagation techniques (e.g. the backward sweep method). In addition, the derivation of the proposed guidance law does not rely on the concrete form of the weighting function, such that arbitrary positive function can be used as a weighting function. We will show later that complex weighting functions (smooth, piecewise, nonsmooth, and even discontinuous) can be easily employed in the proposed guidance law.


## 3.1. Optimal Solution Using the Gauss Orthogonal Collocation Method

In this section, the Gauss orthogonal collocation method is taken to solve the finite-time optimal control problem as shown in Equation (7). The Hamiltonian function of the optimal control problem is
(8)H=12u2tSt+λtTFxt+Gut,where λt denotes the costate and satisfies the following adjoint equation
(9)λ̇t=−∂H∂xt=−FTλt,with the transversality condition
(10)λtf=Sfxtf.From the minimum principle, it is necessary that(11)∂H∂ut=0,which yields the following optimal control
(12)u∗t=−1St−1GTλt.Substituting Equation (12) into Equation (4), and combining with the adjoint Equation (9) yields the following two-point boundary-value problem (TPBVP)
(13)ẋt=∂H∂λt=Fxt−G1St−1GTλt,xt0=x0λ̇t=−∂H∂xt=−FTλt,λtf=Sfxtf.To implement the Gauss orthogonal collocation method, the time intervalt∈t0,tf in Equation (13) should be transformed to τ∈−1,1. By using the following affine transformation
(14)τ=2tf−t0t−tf+t0tf−t0.Equation (13) can be transformed into
(15)ẋτ=tf−t02Fxτ−G1Sτ−1GTλτ,x−1=x0λ̇τ=−tf−t02FTλτ,λ1=Sfx1.The Gauss orthogonal collocation method to solve the optimal control problem Equation (7) is to discretize and transcribe the two-point boundary-value problem into a set of algebraic equations by approximating the state and costate trajectories via interpolating polynomials. The state xτ is approximated using a basis of N+1 Lagrange interpolating polynomials,
(16)xτ≈Xτ=x−1L0τ+∑k=1NxτkLkτ,where Lkτk=0,⋯,N are defined as
(17)Lkτ=∏i=0,i≠kNτ−τiτk−τi.Note that the boundary point, -1, and theN Gauss points, τk, are used to construct the approximation of xτ. Differentiating the expression in Equation (16) gives
(18)ẋτ≈Ẋτ=x−1L̇0τ+∑k=1NxτkL̇kτ.Because the initial value ofλτ is unknown, whereas the terminal value of λτ is known and given by
(19)λ1=Sfx1.Therefore, theN Gauss points, τk, and the boundary point, 1, are used to construct the approximation of λτ as follows
(20)λτ≈Λτ=∑k=1NλτkLk∗τ+λ1LN+1∗τ,where Lk∗τ (k=1,⋯,N+1) are defined as
(21)Lk∗τ=∏i=1,i≠kN+1τ−τiτk−τi.Differentiating the expression in Equation (20) gives
(22)λ̇τ≈Λ̇τ=∑k=1NλτkL̇k∗τ+λ1L̇N+1∗τ.The terminal value of the state,x1, can be obtained in terms xτk and λτk via the Gauss quadrature [34].
(23)x1=x−1+tf−t02∑k=1NwkFxτk−G1Sτ−1GTλτk,where wk are the Gaussian weights and can be calculated offline. Substituting Equation (23) into Equation (19) yields
(24)λ1=Sfx−1+Sftf−t02∑k=1NwkFxτk−G1Sτ−1GTλτk.DenoteD¯i=L̇0τi, Dik=L̇kτi, D¯i∗=L̇N+1∗τ, and Dik∗=L̇k∗τ; Substituting Equation (24) into Equation (22) and including Equation (18) yields
(25)D¯i∗Sfx0+D¯i∗Sftf−t02∑k=1NwkFxk−G1Sk−1GTλk+∑k=1NλkDik∗=−tf−t02FTλix0D¯i+∑k=1NxkDik=tf−t02Fxi−G1Si−1GTλi,where i=1,⋯,N,xk=xτk,λk=λτk,Sk=Sτk and k=1,⋯,N.As proven in [35], the differential approximation matrices, D∈ℝN×N, D¯∈ℝN, D∗∈ℝN×N, and D¯∗∈ℝN have the following relationships:
(26)D¯i=−∑k=1NDik,Dik∗=−wkwiDkiD¯i∗=−∑k=1NDik∗.,The elements of the differential approximation matrixD can be determined offline by
(27)Dik=L̇iτk=∑l=0N∏j=0,j≠i,lNτk−τj∏j=0,j≠iNτi−τj.Consequently, the elements of the matricesD¯, D∗, and D¯∗ can be obtained offline via Equation (26).DenoteΛ~=λ1T,λ2T,⋯,λNTT∈ℝnN, X~=x1T,x2T,⋯,xNTT∈ℝnN. Equation (25) can be rewritten as
(28)D¯~∗x0+D~+∗Λ~+A~X~=0D¯~x0+D~−X~+B~Λ~=0,where
(29)D¯~∗=D¯1∗Sf⋮D¯N∗Sf∈ℝnN×n,(30)D¯~=D¯1In⋮D¯NIn∈ℝnN×N,(31)D~+∗=D11∗InD12∗In⋯D1N∗InD21∗InD22∗In⋯D2N∗In⋮⋮⋱⋮DN1∗InDN2∗In⋯DNN∗In−tf−t02D¯1∗Sfw1G1S1−1GT−FTD¯1∗Sfw2G1S2−1GT⋯D¯1∗SfwNG1SN−1GTD¯2∗Sfw1G1S1−1GTD¯2∗Sfw2G1S2−1GT−FT⋯D¯2∗SfwNG1SN−1GT⋮⋮⋱⋮D¯N∗Sfw1G1S1−1GTD¯N∗Sfw2G1S2−1GT⋯D¯N∗SfwNG1SN−1GT−FT∈ℝnN×nN,(32)D~−=D11InD12In⋯D1NInD21InD22In⋯D2NIn⋮⋮⋱⋮DN1InDN2In⋯DNNIn−tf−t02F0n⋯0n0nF⋯0n⋮⋮⋱⋮0n0n⋯F∈ℝnN×nN,(33)A~=tf−t02D¯1∗Sfw1FD¯1∗Sfw2F⋯D¯1∗SfwNFD¯2∗Sfw1FD¯2∗Sfw2F⋯D¯2∗SfwNF⋮⋮⋱⋮D¯N∗Sfw1FD¯N∗Sfw2F⋯D¯N∗SfwNF∈ℝnN×nN,(34)B~=tf−t02G1S1−1GT0n⋯0n0nG1S2−1GT⋯0n⋮⋮⋱⋮0n0n⋯G1SN−1GT∈ℝnN×nN,and 0n is an n×n zero matrix, In is an n×n identity matrix. By denoting D¯aug=D¯~∗D¯~∈ℝ2nN×n and Daug=D~+∗A~B~D~−∈ℝ2nN×2nN, Equation (28) can be expressed into the following compact form
(35)D¯augx0+DaugΛ~X~=0.Solving forΛ~ and X~ in Equation (35) gives the solution to the two-point boundary-value problem (15) as
(36)Λ~X~=−Daug‐1D¯augx0.From (36) the costate at the Gauss points can be directly determined from the initial state valuex0. However, the Gauss points do not contain the boundary points, therefore, Equation (36) cannot provide the initial value of the costate, λ−1. Because the terminal value of the costate, λ1, is given in (24), therefore, λ−1 can be calculated via the following Gauss quadrature [34].
(37)λ−1=λ1−tf−t02∑k=1NwkFTλk.According to Equations (12), (36), and (37), the optimal control at the boundary points as well as the Gauss points is obtained as
(38)ui∗=−1Si−1GTλii=0,⋯,N+1.It is very important to note that the preceding derivation does not impose any constraint onSt and the form of the equations are identical for F, G, and St, indicating that the proposed method can easily and directly handle different and complex weighting functions without any modification. This property makes the guidance law design very flexible to accomplish the specified objective. In other words, arbitrary positive function can be used as a weighting function, thus the penalty for maneuvers can be adjusted accordingly by tuning the shape of the weighting function in order to achieve some specific guidance objective. In addition, it should be also noted that in our method, the optimal control is obtained without any explicit integration or construction of transition matrices, but only by solving algebraic equations.


## 3.2. Guidance Law Implementation

Note that the optimal solution obtained in Section3.1 is open-loop for the guidance problem. To obtain the closed-loop solution, the receding-horizon technique is implemented. The receding-horizon length is chosen as the time-to-go (i.e., the time remaining to interception, and denoted as tgo), which is approximated by tgo=R/Vm. Consequently, the final time of engagement can be expressed as tf=t0+tgo. The procedure for implementing the proposed guidance law with impact angle constraint can be summarized as follows.
(1)
Initial number of the Gauss points,N, and set t0=0(2)
Solve the sequence of TPBVP(a)
Calculate the time-to-go,tgo(b)
Calculate the matricesD¯aug and D¯aug, and solve Equation (36).(c)
Update the controlui (i=0,⋯,N) using the solution from step 2b and Equation (38).(3)
Apply the control input at the current time at the first discretization point,u0, to the missile’s dynamics and update the missile’s states at the next time(4)
Repeat steps 2-3 until intercepting the targetFor complex weighting functions, the proposed method can easily obtain the optimal control via solving a set of algebraic equations and does not require any explicit propagation techniques (e.g. the backward sweep method). In addition, the derivation of the proposed guidance law does not rely on the concrete form of the weighting function, such that arbitrary positive function can be used as a weighting function. We will show later that complex weighting functions (smooth, piecewise, nonsmooth, and even discontinuous) can be easily employed in the proposed guidance law.


## 4. Hardware Experiment and Simulation Results

In this section, we will first evaluate the real-time computational capacity of the proposed guidance law, and then demonstrate its guidance performance under various weighting functions. The homing conditions for the simulations are given in Table1.Table 1
Conditions for nonlinear simulations.
ParameterValueMissile initial position, (Xm0,Ym0)(1000 m, 5000 m)Target position, (Xt0,Yt0)(6000 m, 0 m)Missile velocity,Vm300 m/sInitial flight path angle,γm(0)-45°, -20°, 20°, 45°Desired impact angle,γf-70°, -90°


### 4.1. Evaluation of Real-Time Computational Capacity

To assess the time-consuming property and the real-time computational capacity of the proposed guidance law, a hardware experimental platform based on the ARM Cortex-M7 processor is addressed in this section. The guidance law is performed on a STM32 Nucleo-F767ZI development board. The hardware experimental platform is shown in Figure2. The microcontroller (STM32F767ZIT6U) used in the development board is based on the high-performance ARM Cortex-M7 32-bit RISC core operating at up to 216 MHz frequency, and features a floating-point unit (FPU) which supports ARM double-precision and single-precision processing instructions. To implement the guidance law on the hardware platform, the STM32Cube MX software is used to generate the HAL-based embedded project of IAR compiler. All the matrix operations needed in the guidance law are carried out using the corresponding C functions provided by CMSIS DSP library, which is a built-in suite of common signal processing functions.Figure 2
Hardware experiment setup.In the test, the compiler’s code optimization of IAR is turned on and off, respectively, while the mode of single-precision FPU is turned on all the time. Table2 summarizes the computation time of the proposed guidance law in one computation loop. It is seen that the computation time increases as the discretization point increases. The compiler optimization for the embedded code is helpful to reduce the computation time. By turning on the compiler optimization mode of IAR and selecting a proper number of discretization points, the computational time can be reduced to a satisfactory level, which means the proposed guidance law is computationally viable for onboard implementation. In addition, it is expected that the computation time will be further reduced by using other high performance processors (for instance, DSP and FPGA).Table 2
Computation time of the proposed guidance law.
C/C++ compiler optimizationFPU: VFPv5 single precisionComputation time (ms)N=5N=10N=15High level with speed modeTurn ON1615NoneTurn ON31021


### 4.2. Performance of the Proposed Guidance Law

In this subsection, numerical examples are performed to illustrate the performance of the proposed guidance law. The initial conditions for the following simulation scenarios are listed in Table1. In these simulations, five different types of weighting functions, including continuous, piecewise, nonsmooth, and discontinuous functions, are employed to demonstrate the capacity and the performance of the proposed guidance law. The number of the Gauss points, N, is selected to be 15. Sf is chosen to be [1×104, 1×104]T.


#### 4.2.1. Case 1: Gaussian and Sinusoidal Weighting Functions

In this case, we take the Gaussian weighting function and the sinusoidal weighting function, which were originally presented in [28, 29], as an example to demonstrate that the proposed guidance law has the capacity to tackle existing weighting functions easily. The Gaussian weighting function S1t and the sinusoidal weighting function S2t are as follows:
(39)S1t=1et−tf/22/2τσ2,(40)S2t=sinπtftf−t.S1t is adopted to shape the specific guidance command: introducing a small acceleration at the initial time. While S2t is adopted to shape the specific guidance command: introducing a small acceleration at the initial and final time. In this case, we set τσ as tf/3.0 for S1t. The profiles of S1t and S2t are shown in Figure 3. As can be seen, the weighting profile of S1t approaches zero values at the initial time, while the weighting profile of S2t approaches zero values at the initial time and the final time. Detailed descriptions for S1t and S2t can be referred to [28, 29].Figure 3
Profiles of the weighting functionS1t and S2t.For comparison purposes, the impact angle control optimal guidance law called OGL that is presented in [19] is also applied in the simulation. Figures 4 and 5 show the intercept trajectories and guidance commands for various weighting functions. The proposed approach can easily tackle the Gaussian weighting function and the sinusoidal weighting function, and successfully provide zero miss distance with the desired impact angle for all cases as shown in Figures 4(a) and 5(a). The guidance commands are clearly shaped by different weighting functions as shown in Figures 4(b) and 5(b). Specifically, the command of OGL approaches a nonzero value at both the initial and final time. This profile may result in an abrupt change of guidance command at the initial time and a performance degradation at the final time under the presences of external disturbances. However, the Gaussian weighting S1t generates a small acceleration command at the initial time compared with OGL, while the sinusoidal weighting function S2t always generates zero initial command and also converges to a zero value as the missile approaches a target. These properties can improve the guidance command, alleviating sensitivity with respect to initial heading error, or improving the lethality of the warhead and guaranteeing operational margin to cope with external disturbances in the vicinity of target. The above results are consistent with the analyses in [28, 29]. Since the proposed guidance law is formulated from the generalized optimal control framework, the existing weighting functions (e.g., Gaussian and sinusoidal) can be easily incorporated into the proposed approach for achieving various guidance goals.Figure 4
The comparison results of different weighting functions forγm0=−45°, γf=−90°.
(a)
Trajectories(b)
Acceleration commandsFigure 5
The comparison results of different weighting functions forγm0=20°, γf=−70°.
(a)
Trajectories(b)
Acceleration commands


#### 4.2.2. Case 2: Sigmoid-Based Weighting Function

In this case, we constructed a new kind of piecewise weighting function based on the sigmoid function as(41)S3t=A1+B11+e−4.6/wt−t10≤t≤tf2A2+B21+e−4.6/wt−tf+t1tf2<t≤tf,where A1,B1,A2,B2,w,andt1 are designed parameters. According to selections of the designed parameters (i.e., A1,B1,A2,B2,w,andt1), various weighting profiles can be generated. Table 3 provides three kind of combinations of A1,B1,A2,B2,w,andt1, and the corresponding profiles of S3t are shown in Figure 6.Table 3
Different parameters forS3t.
Weighting functionParameter of the sigmoid function based weighting functionA1B1A2B2wt1S31t1.0-1.01.0-1.0tf/2.0tf/4.0S32t0.01.00.01.0tf/2.0tf/4.0S33t0.01.01.0-1.0tf/6.0tf/5.0Figure 6
Profiles of the weighting functionS3t with different designed parameters.As seen in Figure6, S3t can be utilized to shape the specific guidance command: providing small acceleration command at the initial time (S31t) or at the final time (S32t) or both at the initial and final time (S33t) by simply tuning the parameters A1,B1,A2,B2,w,andt1. As can be seen, by adopting a piecewise function, the guidance objectives of the Gaussian weighting function and the sinusoidal weighting function can be easily achieved by using a single weighting function S3t.Figure7 shows the intercept trajectories and guidance commands for S3t with different parameters as listed in Table 3. As shown in Figure 7(a), the proposed guidance law can successfully provide zero miss distance with the desired impact angle for all cases. In Figure 7(b), the guidance command provided by the proposed guidance law is flexibly shaped by simply tuning the designed parameters (i.e., A1,B1,A2,B2,w,andt1) of S3t. S31t requires small acceleration demands at the final time, which is desired for ensuring the acceleration margin to cope with unanticipated circumstances at the terminal phase. S32t requires small acceleration demands at the initial time, which is desired for alleviating the sensitivity against the initial heading error and is helpful for smooth guidance handover from a midcourse phase to a homing phase. While S33t can simultaneously satisfy both recommended guidance properties of S31t and S33t.Figure 7
The comparison results of different weighting functions forγm0=−45°, γf=−90°.
(a)
Trajectories(b)
Acceleration commandsThe intercept trajectories and guidance commands for various initial flight path angles with the desired impact angle asγf=−70° under the weighting functions S31t, S32t, and S33t are depicted in Figures 8–10. The simulation results show that although the proposed method is derived in the linear system equations, it can be applied to nonlinear engagement cases with a wide range of initial flight path angle sets. In addition, different guidance objectives can be clearly observed with different weighting functions.Figure 8
The comparison results for various initial flight path angles with the weighting functionS31t.
(a)
Trajectories(b)
Acceleration commandsFigure 9
The comparison results for various initial flight path angles with the weighting functionS32t.
(a)
Trajectories(b)
Acceleration commandsFigure 10
The comparison results for various initial flight path angles with the weighting functionS33t.
(a)
Trajectories(b)
Acceleration commandsThe profile ofS33t can be further shaped by tuning parameters w and t1 as illustrated in Figure 11. The corresponding intercept trajectories and guidance commands are illustrated in Figure 12. As can be seen, the duration of the achieved small guidance command at the initial and final time decrease as t1 decreases, meanwhile the energy expenditure in the middle portion of the engagement decreases. Therefore, for practical guidance applications there is a tradeoff between the energy expenditure and the requirement of small acceleration command in the initial and final phase of the engagement. In addition, through appropriate choices in the design parameters, w and t1, the missile’s trajectory and command profile can be further shaped as desired.Figure 11
Profiles of the weighting functionS33t with different values of w and t1.Figure 12
The comparison results of weighting functionsS33t with different values of w and t1.
(a)
Trajectories(b)
Acceleration commandsThe above simulations show that the proposed method can provide various shapes of the guidance command according to appropriate choices of the weighting functions. Most importantly, our method does not need to make any modifications for each new weighting function. This advantage provides more degrees of freedom in the guidance law design with specified guidance objectives. We will show that the proposed guidance law can easily consider nonsmooth or even discontinuous weighting functions in the following cases.


#### 4.2.3. Case 3: Nonsmooth Weighting Function

In this case, we construct a nonsmooth weighting function that consists of a sigmoid function and a Gaussian function as(42)S4t=21+e−4.6/wt−tf/20≤t≤tf21et−tf/22/2τσ2tf2<t≤tf.By selecting the designed parameters asw=tf/2.0 and τσ=tf/6.0, the profile of S4t is illustrated in Figure 13. As shown in Figure 13, the weighting function of S4t is expected to provide small acceleration demands at the initial time and the final time, which is similar to S2t and S33t. However, different from S2t and S33t, S4t is nonsmooth at t=tf/2.0. This kind of nonsmooth weighting function has not been discussed in previous studies. The main purpose that we design this weighting function is to assess the capacity of the proposed guidance in solving complex weighting functions.Figure 13
Profile of the weighting functionS4t.Figures14 and 15 show the engagement trajectories and guidance commands for the OGL and the proposed guidance law weighted by S4t. Although the weighting function is nonsmooth, the proposed guidance can also successfully satisfy the terminal constraints for interception as provided in Figures 14(a) and 15(a). From Figures 14(b) and 15(b), the guidance commands of the proposed guidance law start from small values and also converge to small values as the missile approaches the target. These simulation results are similar to those for S2t and S33t, since the profile of S4t is similar to those of S2t and S33t. It can be seen from the above results that the proposed guidance is able to handle nonsmooth weighting functions easily, which provides more degrees of freedom in the design of weighting functions.Figure 14
The comparison results of different weighting functions forγm0=−45°, γf=−90°.
(a)
Trajectories(b)
Acceleration commandsFigure 15
The comparison results of different weighting functions forγm0=20°, γf=−70°.
(a)
Trajectories(b)
Acceleration commands


#### 4.2.4. Case 4: Discontinuous Weighting Function

In this case, we take a discontinuous weighting function as an example to further demonstrate the capacity of the proposed guidance in solving complex weighting functions. The discontinuous weighting function is chosen as(43)S5t=0.10≤t≤tf31.0tf3<t≤2tf30.12tf3<t≤tf.The profile ofS5t is shown in Figure 16. S5t is a switching function and is discontinuous at t=tf/3.0 and t=2tf/3.0. This kind of weighting function has not been considered in previous studies. Simulation results are depicted in Figures 17 and 18. As can be seen, the guidance law can successfully provide zero miss distance with the desired impact angle for all cases. However, due to discontinuities in the weighting function, abrupt command changes are occurred during the engagement. From a practical standpoint, the discontinuous function (43) may not be a good candidate for the design of weighting function. However, this case is mainly used to validate the capacity of the proposed guidance law, and has demonstrated that the proposed guidance can handle the discontinuous weighting function easily. This property makes it possible to design different weighting functions in different engagement phases in order to achieve some specific guidance operational goals, which provides a huge degree of freedom in shaping guidance commands. By using the proposed approach, a new guidance law that allows achieving a specific guidance goal can be easily obtained by devising a proper weighting function (smooth, nonsmooth, or even discontinuous).Figure 16
Profile of the weighting functionS5t.Figure 17
The comparison results of different weighting functions forγm0=−45°, γf=−90°.
(a)
Trajectories(b)
Acceleration commandsFigure 18
The comparison results of different weighting functions forγm0=20°, γf=−70°.
(a)
Trajectories(b)
Acceleration commands


#### 4.2.5. Case 5: Extension to the First-Order Lag System

Although the proposed guidance law is derived for the lag-free missile system in Section3, our approach can be easily extended to handle the guidance problem with dynamic lag effect. Suppose that the missile autopilot model is a first-order delay system:
(44)ȧmt=1tτact−amt,where act and tτ denote the guidance command and autopilot time delay constant, respectively. Consequently, for the first-order lag missile system, the whole derivation in Section 3 does not require any modification but simply set Equation (5) as x≜yγ¯mamT, u≜ac, F≜0Vm0001/Vm00−1/tτ, and G≜001/tτ. Then the corresponding matrix Sf is chosen to be [1×104, 1×104, 0]T. Figure 19 provides the guidance performance of the proposed guidance law to the first-order lag system for the weighting function S31t with the initial flight path angle as γm0=−45° and the desired impact angle as γf=−90°. From this figure, we observe that the proposed guidance law can successfully intercept the target for all cases and there is a time delay between the command and response.Figure 19
Trajectories and commands of the weighting functionS31t for various dynamic lag.
(a)
Trajectories(b)
Acceleration commands


## 4.1. Evaluation of Real-Time Computational Capacity

To assess the time-consuming property and the real-time computational capacity of the proposed guidance law, a hardware experimental platform based on the ARM Cortex-M7 processor is addressed in this section. The guidance law is performed on a STM32 Nucleo-F767ZI development board. The hardware experimental platform is shown in Figure2. The microcontroller (STM32F767ZIT6U) used in the development board is based on the high-performance ARM Cortex-M7 32-bit RISC core operating at up to 216 MHz frequency, and features a floating-point unit (FPU) which supports ARM double-precision and single-precision processing instructions. To implement the guidance law on the hardware platform, the STM32Cube MX software is used to generate the HAL-based embedded project of IAR compiler. All the matrix operations needed in the guidance law are carried out using the corresponding C functions provided by CMSIS DSP library, which is a built-in suite of common signal processing functions.Figure 2
Hardware experiment setup.In the test, the compiler’s code optimization of IAR is turned on and off, respectively, while the mode of single-precision FPU is turned on all the time. Table2 summarizes the computation time of the proposed guidance law in one computation loop. It is seen that the computation time increases as the discretization point increases. The compiler optimization for the embedded code is helpful to reduce the computation time. By turning on the compiler optimization mode of IAR and selecting a proper number of discretization points, the computational time can be reduced to a satisfactory level, which means the proposed guidance law is computationally viable for onboard implementation. In addition, it is expected that the computation time will be further reduced by using other high performance processors (for instance, DSP and FPGA).Table 2
Computation time of the proposed guidance law.
C/C++ compiler optimizationFPU: VFPv5 single precisionComputation time (ms)N=5N=10N=15High level with speed modeTurn ON1615NoneTurn ON31021


## 4.2. Performance of the Proposed Guidance Law

In this subsection, numerical examples are performed to illustrate the performance of the proposed guidance law. The initial conditions for the following simulation scenarios are listed in Table1. In these simulations, five different types of weighting functions, including continuous, piecewise, nonsmooth, and discontinuous functions, are employed to demonstrate the capacity and the performance of the proposed guidance law. The number of the Gauss points, N, is selected to be 15. Sf is chosen to be [1×104, 1×104]T.


### 4.2.1. Case 1: Gaussian and Sinusoidal Weighting Functions

In this case, we take the Gaussian weighting function and the sinusoidal weighting function, which were originally presented in [28, 29], as an example to demonstrate that the proposed guidance law has the capacity to tackle existing weighting functions easily. The Gaussian weighting function S1t and the sinusoidal weighting function S2t are as follows:
(39)S1t=1et−tf/22/2τσ2,(40)S2t=sinπtftf−t.S1t is adopted to shape the specific guidance command: introducing a small acceleration at the initial time. While S2t is adopted to shape the specific guidance command: introducing a small acceleration at the initial and final time. In this case, we set τσ as tf/3.0 for S1t. The profiles of S1t and S2t are shown in Figure 3. As can be seen, the weighting profile of S1t approaches zero values at the initial time, while the weighting profile of S2t approaches zero values at the initial time and the final time. Detailed descriptions for S1t and S2t can be referred to [28, 29].Figure 3
Profiles of the weighting functionS1t and S2t.For comparison purposes, the impact angle control optimal guidance law called OGL that is presented in [19] is also applied in the simulation. Figures 4 and 5 show the intercept trajectories and guidance commands for various weighting functions. The proposed approach can easily tackle the Gaussian weighting function and the sinusoidal weighting function, and successfully provide zero miss distance with the desired impact angle for all cases as shown in Figures 4(a) and 5(a). The guidance commands are clearly shaped by different weighting functions as shown in Figures 4(b) and 5(b). Specifically, the command of OGL approaches a nonzero value at both the initial and final time. This profile may result in an abrupt change of guidance command at the initial time and a performance degradation at the final time under the presences of external disturbances. However, the Gaussian weighting S1t generates a small acceleration command at the initial time compared with OGL, while the sinusoidal weighting function S2t always generates zero initial command and also converges to a zero value as the missile approaches a target. These properties can improve the guidance command, alleviating sensitivity with respect to initial heading error, or improving the lethality of the warhead and guaranteeing operational margin to cope with external disturbances in the vicinity of target. The above results are consistent with the analyses in [28, 29]. Since the proposed guidance law is formulated from the generalized optimal control framework, the existing weighting functions (e.g., Gaussian and sinusoidal) can be easily incorporated into the proposed approach for achieving various guidance goals.Figure 4
The comparison results of different weighting functions forγm0=−45°, γf=−90°.
(a)
Trajectories(b)
Acceleration commandsFigure 5
The comparison results of different weighting functions forγm0=20°, γf=−70°.
(a)
Trajectories(b)
Acceleration commands


### 4.2.2. Case 2: Sigmoid-Based Weighting Function

In this case, we constructed a new kind of piecewise weighting function based on the sigmoid function as(41)S3t=A1+B11+e−4.6/wt−t10≤t≤tf2A2+B21+e−4.6/wt−tf+t1tf2<t≤tf,where A1,B1,A2,B2,w,andt1 are designed parameters. According to selections of the designed parameters (i.e., A1,B1,A2,B2,w,andt1), various weighting profiles can be generated. Table 3 provides three kind of combinations of A1,B1,A2,B2,w,andt1, and the corresponding profiles of S3t are shown in Figure 6.Table 3
Different parameters forS3t.
Weighting functionParameter of the sigmoid function based weighting functionA1B1A2B2wt1S31t1.0-1.01.0-1.0tf/2.0tf/4.0S32t0.01.00.01.0tf/2.0tf/4.0S33t0.01.01.0-1.0tf/6.0tf/5.0Figure 6
Profiles of the weighting functionS3t with different designed parameters.As seen in Figure6, S3t can be utilized to shape the specific guidance command: providing small acceleration command at the initial time (S31t) or at the final time (S32t) or both at the initial and final time (S33t) by simply tuning the parameters A1,B1,A2,B2,w,andt1. As can be seen, by adopting a piecewise function, the guidance objectives of the Gaussian weighting function and the sinusoidal weighting function can be easily achieved by using a single weighting function S3t.Figure7 shows the intercept trajectories and guidance commands for S3t with different parameters as listed in Table 3. As shown in Figure 7(a), the proposed guidance law can successfully provide zero miss distance with the desired impact angle for all cases. In Figure 7(b), the guidance command provided by the proposed guidance law is flexibly shaped by simply tuning the designed parameters (i.e., A1,B1,A2,B2,w,andt1) of S3t. S31t requires small acceleration demands at the final time, which is desired for ensuring the acceleration margin to cope with unanticipated circumstances at the terminal phase. S32t requires small acceleration demands at the initial time, which is desired for alleviating the sensitivity against the initial heading error and is helpful for smooth guidance handover from a midcourse phase to a homing phase. While S33t can simultaneously satisfy both recommended guidance properties of S31t and S33t.Figure 7
The comparison results of different weighting functions forγm0=−45°, γf=−90°.
(a)
Trajectories(b)
Acceleration commandsThe intercept trajectories and guidance commands for various initial flight path angles with the desired impact angle asγf=−70° under the weighting functions S31t, S32t, and S33t are depicted in Figures 8–10. The simulation results show that although the proposed method is derived in the linear system equations, it can be applied to nonlinear engagement cases with a wide range of initial flight path angle sets. In addition, different guidance objectives can be clearly observed with different weighting functions.Figure 8
The comparison results for various initial flight path angles with the weighting functionS31t.
(a)
Trajectories(b)
Acceleration commandsFigure 9
The comparison results for various initial flight path angles with the weighting functionS32t.
(a)
Trajectories(b)
Acceleration commandsFigure 10
The comparison results for various initial flight path angles with the weighting functionS33t.
(a)
Trajectories(b)
Acceleration commandsThe profile ofS33t can be further shaped by tuning parameters w and t1 as illustrated in Figure 11. The corresponding intercept trajectories and guidance commands are illustrated in Figure 12. As can be seen, the duration of the achieved small guidance command at the initial and final time decrease as t1 decreases, meanwhile the energy expenditure in the middle portion of the engagement decreases. Therefore, for practical guidance applications there is a tradeoff between the energy expenditure and the requirement of small acceleration command in the initial and final phase of the engagement. In addition, through appropriate choices in the design parameters, w and t1, the missile’s trajectory and command profile can be further shaped as desired.Figure 11
Profiles of the weighting functionS33t with different values of w and t1.Figure 12
The comparison results of weighting functionsS33t with different values of w and t1.
(a)
Trajectories(b)
Acceleration commandsThe above simulations show that the proposed method can provide various shapes of the guidance command according to appropriate choices of the weighting functions. Most importantly, our method does not need to make any modifications for each new weighting function. This advantage provides more degrees of freedom in the guidance law design with specified guidance objectives. We will show that the proposed guidance law can easily consider nonsmooth or even discontinuous weighting functions in the following cases.


### 4.2.3. Case 3: Nonsmooth Weighting Function

In this case, we construct a nonsmooth weighting function that consists of a sigmoid function and a Gaussian function as(42)S4t=21+e−4.6/wt−tf/20≤t≤tf21et−tf/22/2τσ2tf2<t≤tf.By selecting the designed parameters asw=tf/2.0 and τσ=tf/6.0, the profile of S4t is illustrated in Figure 13. As shown in Figure 13, the weighting function of S4t is expected to provide small acceleration demands at the initial time and the final time, which is similar to S2t and S33t. However, different from S2t and S33t, S4t is nonsmooth at t=tf/2.0. This kind of nonsmooth weighting function has not been discussed in previous studies. The main purpose that we design this weighting function is to assess the capacity of the proposed guidance in solving complex weighting functions.Figure 13
Profile of the weighting functionS4t.Figures14 and 15 show the engagement trajectories and guidance commands for the OGL and the proposed guidance law weighted by S4t. Although the weighting function is nonsmooth, the proposed guidance can also successfully satisfy the terminal constraints for interception as provided in Figures 14(a) and 15(a). From Figures 14(b) and 15(b), the guidance commands of the proposed guidance law start from small values and also converge to small values as the missile approaches the target. These simulation results are similar to those for S2t and S33t, since the profile of S4t is similar to those of S2t and S33t. It can be seen from the above results that the proposed guidance is able to handle nonsmooth weighting functions easily, which provides more degrees of freedom in the design of weighting functions.Figure 14
The comparison results of different weighting functions forγm0=−45°, γf=−90°.
(a)
Trajectories(b)
Acceleration commandsFigure 15
The comparison results of different weighting functions forγm0=20°, γf=−70°.
(a)
Trajectories(b)
Acceleration commands


### 4.2.4. Case 4: Discontinuous Weighting Function

In this case, we take a discontinuous weighting function as an example to further demonstrate the capacity of the proposed guidance in solving complex weighting functions. The discontinuous weighting function is chosen as(43)S5t=0.10≤t≤tf31.0tf3<t≤2tf30.12tf3<t≤tf.The profile ofS5t is shown in Figure 16. S5t is a switching function and is discontinuous at t=tf/3.0 and t=2tf/3.0. This kind of weighting function has not been considered in previous studies. Simulation results are depicted in Figures 17 and 18. As can be seen, the guidance law can successfully provide zero miss distance with the desired impact angle for all cases. However, due to discontinuities in the weighting function, abrupt command changes are occurred during the engagement. From a practical standpoint, the discontinuous function (43) may not be a good candidate for the design of weighting function. However, this case is mainly used to validate the capacity of the proposed guidance law, and has demonstrated that the proposed guidance can handle the discontinuous weighting function easily. This property makes it possible to design different weighting functions in different engagement phases in order to achieve some specific guidance operational goals, which provides a huge degree of freedom in shaping guidance commands. By using the proposed approach, a new guidance law that allows achieving a specific guidance goal can be easily obtained by devising a proper weighting function (smooth, nonsmooth, or even discontinuous).Figure 16
Profile of the weighting functionS5t.Figure 17
The comparison results of different weighting functions forγm0=−45°, γf=−90°.
(a)
Trajectories(b)
Acceleration commandsFigure 18
The comparison results of different weighting functions forγm0=20°, γf=−70°.
(a)
Trajectories(b)
Acceleration commands


### 4.2.5. Case 5: Extension to the First-Order Lag System

Although the proposed guidance law is derived for the lag-free missile system in Section3, our approach can be easily extended to handle the guidance problem with dynamic lag effect. Suppose that the missile autopilot model is a first-order delay system:
(44)ȧmt=1tτact−amt,where act and tτ denote the guidance command and autopilot time delay constant, respectively. Consequently, for the first-order lag missile system, the whole derivation in Section 3 does not require any modification but simply set Equation (5) as x≜yγ¯mamT, u≜ac, F≜0Vm0001/Vm00−1/tτ, and G≜001/tτ. Then the corresponding matrix Sf is chosen to be [1×104, 1×104, 0]T. Figure 19 provides the guidance performance of the proposed guidance law to the first-order lag system for the weighting function S31t with the initial flight path angle as γm0=−45° and the desired impact angle as γf=−90°. From this figure, we observe that the proposed guidance law can successfully intercept the target for all cases and there is a time delay between the command and response.Figure 19
Trajectories and commands of the weighting functionS31t for various dynamic lag.
(a)
Trajectories(b)
Acceleration commands


## 4.2.1. Case 1: Gaussian and Sinusoidal Weighting Functions

In this case, we take the Gaussian weighting function and the sinusoidal weighting function, which were originally presented in [28, 29], as an example to demonstrate that the proposed guidance law has the capacity to tackle existing weighting functions easily. The Gaussian weighting function S1t and the sinusoidal weighting function S2t are as follows:
(39)S1t=1et−tf/22/2τσ2,(40)S2t=sinπtftf−t.S1t is adopted to shape the specific guidance command: introducing a small acceleration at the initial time. While S2t is adopted to shape the specific guidance command: introducing a small acceleration at the initial and final time. In this case, we set τσ as tf/3.0 for S1t. The profiles of S1t and S2t are shown in Figure 3. As can be seen, the weighting profile of S1t approaches zero values at the initial time, while the weighting profile of S2t approaches zero values at the initial time and the final time. Detailed descriptions for S1t and S2t can be referred to [28, 29].Figure 3
Profiles of the weighting functionS1t and S2t.For comparison purposes, the impact angle control optimal guidance law called OGL that is presented in [19] is also applied in the simulation. Figures 4 and 5 show the intercept trajectories and guidance commands for various weighting functions. The proposed approach can easily tackle the Gaussian weighting function and the sinusoidal weighting function, and successfully provide zero miss distance with the desired impact angle for all cases as shown in Figures 4(a) and 5(a). The guidance commands are clearly shaped by different weighting functions as shown in Figures 4(b) and 5(b). Specifically, the command of OGL approaches a nonzero value at both the initial and final time. This profile may result in an abrupt change of guidance command at the initial time and a performance degradation at the final time under the presences of external disturbances. However, the Gaussian weighting S1t generates a small acceleration command at the initial time compared with OGL, while the sinusoidal weighting function S2t always generates zero initial command and also converges to a zero value as the missile approaches a target. These properties can improve the guidance command, alleviating sensitivity with respect to initial heading error, or improving the lethality of the warhead and guaranteeing operational margin to cope with external disturbances in the vicinity of target. The above results are consistent with the analyses in [28, 29]. Since the proposed guidance law is formulated from the generalized optimal control framework, the existing weighting functions (e.g., Gaussian and sinusoidal) can be easily incorporated into the proposed approach for achieving various guidance goals.Figure 4
The comparison results of different weighting functions forγm0=−45°, γf=−90°.
(a)
Trajectories(b)
Acceleration commandsFigure 5
The comparison results of different weighting functions forγm0=20°, γf=−70°.
(a)
Trajectories(b)
Acceleration commands


## 4.2.2. Case 2: Sigmoid-Based Weighting Function

In this case, we constructed a new kind of piecewise weighting function based on the sigmoid function as(41)S3t=A1+B11+e−4.6/wt−t10≤t≤tf2A2+B21+e−4.6/wt−tf+t1tf2<t≤tf,where A1,B1,A2,B2,w,andt1 are designed parameters. According to selections of the designed parameters (i.e., A1,B1,A2,B2,w,andt1), various weighting profiles can be generated. Table 3 provides three kind of combinations of A1,B1,A2,B2,w,andt1, and the corresponding profiles of S3t are shown in Figure 6.Table 3
Different parameters forS3t.
Weighting functionParameter of the sigmoid function based weighting functionA1B1A2B2wt1S31t1.0-1.01.0-1.0tf/2.0tf/4.0S32t0.01.00.01.0tf/2.0tf/4.0S33t0.01.01.0-1.0tf/6.0tf/5.0Figure 6
Profiles of the weighting functionS3t with different designed parameters.As seen in Figure6, S3t can be utilized to shape the specific guidance command: providing small acceleration command at the initial time (S31t) or at the final time (S32t) or both at the initial and final time (S33t) by simply tuning the parameters A1,B1,A2,B2,w,andt1. As can be seen, by adopting a piecewise function, the guidance objectives of the Gaussian weighting function and the sinusoidal weighting function can be easily achieved by using a single weighting function S3t.Figure7 shows the intercept trajectories and guidance commands for S3t with different parameters as listed in Table 3. As shown in Figure 7(a), the proposed guidance law can successfully provide zero miss distance with the desired impact angle for all cases. In Figure 7(b), the guidance command provided by the proposed guidance law is flexibly shaped by simply tuning the designed parameters (i.e., A1,B1,A2,B2,w,andt1) of S3t. S31t requires small acceleration demands at the final time, which is desired for ensuring the acceleration margin to cope with unanticipated circumstances at the terminal phase. S32t requires small acceleration demands at the initial time, which is desired for alleviating the sensitivity against the initial heading error and is helpful for smooth guidance handover from a midcourse phase to a homing phase. While S33t can simultaneously satisfy both recommended guidance properties of S31t and S33t.Figure 7
The comparison results of different weighting functions forγm0=−45°, γf=−90°.
(a)
Trajectories(b)
Acceleration commandsThe intercept trajectories and guidance commands for various initial flight path angles with the desired impact angle asγf=−70° under the weighting functions S31t, S32t, and S33t are depicted in Figures 8–10. The simulation results show that although the proposed method is derived in the linear system equations, it can be applied to nonlinear engagement cases with a wide range of initial flight path angle sets. In addition, different guidance objectives can be clearly observed with different weighting functions.Figure 8
The comparison results for various initial flight path angles with the weighting functionS31t.
(a)
Trajectories(b)
Acceleration commandsFigure 9
The comparison results for various initial flight path angles with the weighting functionS32t.
(a)
Trajectories(b)
Acceleration commandsFigure 10
The comparison results for various initial flight path angles with the weighting functionS33t.
(a)
Trajectories(b)
Acceleration commandsThe profile ofS33t can be further shaped by tuning parameters w and t1 as illustrated in Figure 11. The corresponding intercept trajectories and guidance commands are illustrated in Figure 12. As can be seen, the duration of the achieved small guidance command at the initial and final time decrease as t1 decreases, meanwhile the energy expenditure in the middle portion of the engagement decreases. Therefore, for practical guidance applications there is a tradeoff between the energy expenditure and the requirement of small acceleration command in the initial and final phase of the engagement. In addition, through appropriate choices in the design parameters, w and t1, the missile’s trajectory and command profile can be further shaped as desired.Figure 11
Profiles of the weighting functionS33t with different values of w and t1.Figure 12
The comparison results of weighting functionsS33t with different values of w and t1.
(a)
Trajectories(b)
Acceleration commandsThe above simulations show that the proposed method can provide various shapes of the guidance command according to appropriate choices of the weighting functions. Most importantly, our method does not need to make any modifications for each new weighting function. This advantage provides more degrees of freedom in the guidance law design with specified guidance objectives. We will show that the proposed guidance law can easily consider nonsmooth or even discontinuous weighting functions in the following cases.


## 4.2.3. Case 3: Nonsmooth Weighting Function

In this case, we construct a nonsmooth weighting function that consists of a sigmoid function and a Gaussian function as(42)S4t=21+e−4.6/wt−tf/20≤t≤tf21et−tf/22/2τσ2tf2<t≤tf.By selecting the designed parameters asw=tf/2.0 and τσ=tf/6.0, the profile of S4t is illustrated in Figure 13. As shown in Figure 13, the weighting function of S4t is expected to provide small acceleration demands at the initial time and the final time, which is similar to S2t and S33t. However, different from S2t and S33t, S4t is nonsmooth at t=tf/2.0. This kind of nonsmooth weighting function has not been discussed in previous studies. The main purpose that we design this weighting function is to assess the capacity of the proposed guidance in solving complex weighting functions.Figure 13
Profile of the weighting functionS4t.Figures14 and 15 show the engagement trajectories and guidance commands for the OGL and the proposed guidance law weighted by S4t. Although the weighting function is nonsmooth, the proposed guidance can also successfully satisfy the terminal constraints for interception as provided in Figures 14(a) and 15(a). From Figures 14(b) and 15(b), the guidance commands of the proposed guidance law start from small values and also converge to small values as the missile approaches the target. These simulation results are similar to those for S2t and S33t, since the profile of S4t is similar to those of S2t and S33t. It can be seen from the above results that the proposed guidance is able to handle nonsmooth weighting functions easily, which provides more degrees of freedom in the design of weighting functions.Figure 14
The comparison results of different weighting functions forγm0=−45°, γf=−90°.
(a)
Trajectories(b)
Acceleration commandsFigure 15
The comparison results of different weighting functions forγm0=20°, γf=−70°.
(a)
Trajectories(b)
Acceleration commands


## 4.2.4. Case 4: Discontinuous Weighting Function

In this case, we take a discontinuous weighting function as an example to further demonstrate the capacity of the proposed guidance in solving complex weighting functions. The discontinuous weighting function is chosen as(43)S5t=0.10≤t≤tf31.0tf3<t≤2tf30.12tf3<t≤tf.The profile ofS5t is shown in Figure 16. S5t is a switching function and is discontinuous at t=tf/3.0 and t=2tf/3.0. This kind of weighting function has not been considered in previous studies. Simulation results are depicted in Figures 17 and 18. As can be seen, the guidance law can successfully provide zero miss distance with the desired impact angle for all cases. However, due to discontinuities in the weighting function, abrupt command changes are occurred during the engagement. From a practical standpoint, the discontinuous function (43) may not be a good candidate for the design of weighting function. However, this case is mainly used to validate the capacity of the proposed guidance law, and has demonstrated that the proposed guidance can handle the discontinuous weighting function easily. This property makes it possible to design different weighting functions in different engagement phases in order to achieve some specific guidance operational goals, which provides a huge degree of freedom in shaping guidance commands. By using the proposed approach, a new guidance law that allows achieving a specific guidance goal can be easily obtained by devising a proper weighting function (smooth, nonsmooth, or even discontinuous).Figure 16
Profile of the weighting functionS5t.Figure 17
The comparison results of different weighting functions forγm0=−45°, γf=−90°.
(a)
Trajectories(b)
Acceleration commandsFigure 18
The comparison results of different weighting functions forγm0=20°, γf=−70°.
(a)
Trajectories(b)
Acceleration commands


## 4.2.5. Case 5: Extension to the First-Order Lag System

Although the proposed guidance law is derived for the lag-free missile system in Section3, our approach can be easily extended to handle the guidance problem with dynamic lag effect. Suppose that the missile autopilot model is a first-order delay system:
(44)ȧmt=1tτact−amt,where act and tτ denote the guidance command and autopilot time delay constant, respectively. Consequently, for the first-order lag missile system, the whole derivation in Section 3 does not require any modification but simply set Equation (5) as x≜yγ¯mamT, u≜ac, F≜0Vm0001/Vm00−1/tτ, and G≜001/tτ. Then the corresponding matrix Sf is chosen to be [1×104, 1×104, 0]T. Figure 19 provides the guidance performance of the proposed guidance law to the first-order lag system for the weighting function S31t with the initial flight path angle as γm0=−45° and the desired impact angle as γf=−90°. From this figure, we observe that the proposed guidance law can successfully intercept the target for all cases and there is a time delay between the command and response.Figure 19
Trajectories and commands of the weighting functionS31t for various dynamic lag.
(a)
Trajectories(b)
Acceleration commands


## 5. Conclusions

In this paper, a new computational impact angle control homing guidance law is investigated based on the energy cost weighted by arbitrary functions (smooth, nonsmooth, or even discontinuous), which has the command shaping capacity as desired. By virtue of the Gauss orthogonal collocation method, the optimal control is obtained by solving a set of algebraic equations. The real-time computational capacity of the proposed guidance law is evaluated on a microprocessor. Several types of weighting functions, including smooth, piecewise, nonsmooth, and discontinuous functions, have been tested to evaluate the capacity of the proposed guidance law. Simulation results show that the proposed guidance law can provide various shapes of the guidance command according to appropriate choices of the weighting functions. In addition, profiting from the adoption of a piecewise function, the guidance command shape that provides small acceleration demands at the initial phase, or at the terminal phase, or both at the initial and terminal phases can be easily achieved by using the proposed approach with just one single weighting function. The proposed guidance law is formulated from the generalized optimal control framework and a new guidance law that allows achieving a specific guidance goal can be easily obtained by devising a proper weighting function (smooth, nonsmooth, or even discontinuous), providing more degrees of freedom in the guidance law design with specified guidance goals.



---

*Source: 8747790-2023-01-05.xml*