# Optimization of Aircraft Flight Trajectory Combined with Thinking Navigation Algorithm

Authors: Yifan Yang
Journal: Advances in Multimedia
 (2022)

---

## Abstract

In order to improve the optimization effect of the flight trajectory of the aircraft, this paper combines the thinking navigation algorithm to optimize the flight trajectory of the aircraft and analyzes the flight trajectory of the aircraft through the intelligent model. By processing the original satellite clock error data by the first-order difference method, the modeling data can be more suitable for nonlinear characteristics. Moreover, this paper chooses a simple network structure and uses the MEA to select the optimal initial parameters of the model for the BP neural network, which can avoid the local optimization of the BP neural network results. In addition, this paper conducts experimental analysis on the MEA-BP model through fitting data of different lengths. The simulation test results show that the thinking navigation algorithm proposed in this paper has a very obvious effect on the optimization of the flight trajectory of the aircraft.


---

## Body

## 1. Introduction

The trajectory optimization of the aircraft refers to the solution of parameters such as speed and altitude during the entire flight process under the condition of given flight performance indicators (the shortest flight time, and the least fuel consumption in the whole process). Under the action of the engine thrust, aerodynamic force, and the gravity of the aircraft, we need to study the law of the movement of the center of mass of the aircraft, that is, to optimize the takeoff and landing performance of the aircraft, as well as the climb, cruise, and descent performance. Therefore, when studying the problem of aircraft flight performance, it is necessary to establish the equation of motion of the center of mass of the aircraft. The establishment of the equation of motion for the center of mass of the aircraft is closely related to the flight environment. In addition to this, flight performance is affected by the interaction between the aircraft and the air, including the aerodynamic characteristics of the aircraft, the total weight of the aircraft, and the relationship of atmospheric and environmental conditions to the aircraft. Therefore, accurate calculation of flight performance must be based on reliable atmospheric data, aerodynamic data, motorized data, and aircraft data.The essence of civil aircraft flight management system is trajectory optimization and tracking. Modern aircraft trajectory optimization is based on aircraft kinematics equations, dynamic equations, and multiple constraints under the premise of designing control variables according to flight plan requirements, flight technical indicators, etc., and constructing target functionals. Therefore, the optimal control variables are solved to find an optimal flight trajectory that satisfies a certain performance index. The trajectory optimization design is an important part of the overall design of modern aircraft. In the overall design of the aircraft, trajectory optimization and tracking play an important guiding role in other subsystems such as overall layout analysis, overall parameter determination, aerodynamic design, material and structural strength, control system, and power plant system. It is one of the overall goals pursued by the design of civil aviation aircraft to combine with the dynamic parameter design to make all the subsystems of the aircraft cooperate with each other to achieve the optimal performance indicators. Only by considering multiple aspects can the designed aircraft be more practical and valuable. In the three-stage (climbing, cruising, and descending) trajectory optimization, in the trajectory optimization tasks of each stage, early people focused on local point optimization and segment optimization, and the performance indicators that were concerned were the shortest time, the most fuel-efficient, and the shortest distance climb. In index [1], the methods used mainly include retrieval method, iterative method, dynamic variation method, energy state management method, and singular perturbation method. The advantage of the retrieval method [2] is that it is convenient and concise. The disadvantage is that it only considers the local optimum and requires a large storage space. When the flight task and flight environment change, the accuracy, and accuracy and search efficiency of this method are low. The iterative method in the static optimization method expresses the quantity to be optimized, such as the cost of flight time and fuel cost, as an objective function composed of certain parameters, and then substitutes the objective function into the parameter estimated value to continuously iterate to finally obtain the optimal solution of the objective function. The method is simple and direct, but the disadvantage is that the whole trajectory cannot be optimized. The variational method transforms the trajectory optimization problem into a state equation based on the optimal control theory and uses the minimum value principle or the variational principle to obtain the global optimal solution [3].At present, the research on aircraft trajectory planning mainly focuses on UAVs and military aircraft, such as terrain avoidance, terrain tracking, and UAV formation flight trajectory planning in low-altitude penetration of military aircraft. The main applied planning algorithms are as follows: grid method, artificial potential field method, A∗ search algorithm, sparse A∗ search algorithm, D∗ search algorithm, genetic algorithm, ant colony algorithm, particle swarm algorithm, etc. Scholars like Amin JN discussed several problems of UAV track planning based on genetic algorithm [4]. The literature [5] applied the navigation algorithm to online track optimization, but this method is complicated and difficult to understand. The literature [6] proposed that the distance propagation method is used, which can effectively solve some dynamic trajectory optimization problems, but with the increase of the number of grids, the optimization time increases significantly. The literature [7] uses the threat cost weighting method to construct a Voronoi diagram and analyzes it according to the characteristics of real-time trajectory planning. Two (strategic and tactical) trajectory planning ideas are proposed. The dual-population particle swarm planning algorithm proposed by the literature [8] uses a parallel search strategy to meet the requirements of real-time planning. The literature [9] uses dynamic path algorithm and terrain tracking. The algorithm replans the helicopter’s online track and improves the mission survivability of the aircraft. The literature [10] uses the differential evolution algorithm (DE) to randomly search the space and points out that this algorithm has the advantages of fast search and good robustness. Regarding the study of aircraft avoiding restricted areas such as flight restricted areas and threats, the literature [11] uses genetic algorithm to solve the flight trajectory with constraints, and the literature [12] uses MILP method to transform the trajectory planning problem into a linear programming problem under constraints. In order to obtain the optimal obstacle avoidance trajectory, the literature [13] discusses the trajectory optimization problem of UAV in the presence of flight restricted areas. The goal of trajectory planning is to detect and analyze real-time environmental data in a timely and effective manner when the aircraft reaches the destination, so as to change the preplanned flight path in real time and reduce the probability of the aircraft being threatened by maneuvering. The literature [14] applied the adaptive genetic algorithm and D∗ search algorithm to the offline and online trajectory planning of the aircraft cruise segment on the basis of referring to relevant theories and algorithms, and achieved good planning results.With the complexity of the mathematical model of the aircraft, the multicoupling and aerodynamic nonlinearity phenomena in the flight process make this method almost no longer applicable [15]. With the development of science and technology and the study of nonlinear control theory, many nonlinear control methods have appeared. Judging from the current research work on aircraft track tracking, most scholars focus on neural network control, fuzzy control, LQR control, sliding mode control, and adaptive control. In the literature [16], a method for decoupling the strongly coupled transformation of the aircraft is given. When designing the sliding mode control law, the motion model needs to be linearized. The literature [17] designed a sliding mode control law to stably track the three-dimensional flight path of the aircraft. The literature [18] combined the predictive function control command with the fuzzy inference compensation to achieve the optimal track tracking control of the aircraft problem, and the simulation shows the feasibility of this method. The literature [19] proposed an improved robust adaptive control, which reduces the complexity of the control law by introducing a dynamic surface, and uses the RBF neural network approximation method to control the uncertain parameters inside the model, and the overall realization of the pitch angle, attack angle, and level tail is achieved. Deflected instruction tracking. The literature [20] designed a robust adaptive reverse thrust sliding mode trajectory tracking method, which realized trajectory tracking with state variables and angle of attack as commands, and realized observation compensation for atmospheric disturbance and modeling error by designing an observer. And we introduce the tangent function and attenuation factor to reduce the signal chattering problem.This paper combines the thinking navigation algorithm to optimize the flight trajectory of the aircraft, analyzes the flight trajectory of the aircraft through the intelligent model, and combines the intelligent algorithm to improve the planning effect of the flight trajectory of the aircraft.


## 2. Thinking Navigation Algorithm

### 2.1. The Basic Theory of Satellite Clock Error Forecasting

Satellite clock error data, also known as phase data, are the time offset data obtained by standard time and atomic clock output time, in nanoseconds or seconds. In satellite navigation and positioning, phase data represent satellite clock error data and receiver clock error data. The raw phase data can be directly used for navigation positioning and satellite clock error prediction. However, the magnitude of the phase data is large, and it is impossible to judge and identify the location of outliers such as gross errors and clock jumps. Therefore, the data need to be preprocessed before the actual clock error forecasting. By making corresponding mathematical changes to the phase data, the corresponding frequency data can be obtained, and by processing the frequency data, the gross errors and clock jumps existing in the clock error can be effectively detected.The relative deviation of the standard output frequencies of different frequencies can obtain frequency data, which is defined as:(1)D=fA−fBf0.In the formula,fA is the output frequency of the frequency standard A, fB is the output frequency of the frequency standard B, f0 is the nominal frequency, and D is the obtained frequency data.Over time, the actual output frequency of the atomic clock will be affected by the internal components and the external environment, and it usually presents a linear trend of increasing or decreasing. When time changes, the frequency of satellite clock also changes, and the rate of change is called frequency drift rate or aging rate. The calculation formula is as follows:(2)zi=yi+1−yiτ0.In the formula,yi is the frequency value of epoch i, yi+1 is the frequency value of epoch i + 1, and τ0 is the sampling interval of adjacent epochs.The frequency of atomic clocks is often not a specific value, but fluctuates within a certain range. It obtains correct information on whether the atomic clock frequency is accurate or not by comparing the actual frequency obtained by measuring the 1BP neural network algorithm or calculation with the nominal frequency. The calculation formula is as follows:(3)A=f−f0f0.In the formula,f is the actual frequency measured or calculated, f0 is the real frequency, that is, the nominal frequency, and A is the calculated frequency accuracy. Since the frequency of atomic clocks fluctuates within a certain range, the average of multiple measurement results is usually taken as the final result.


### 2.2. Preprocessing of Satellite Clock Error Data

The frequency data are obtained from the phase data by mathematical transformation, and the conversion formula is as follows:(4)yi=xi+1−xiτ0.In the formula,xi is the value of the satellite clock error at epoch i, xi+1 is the value of the satellite clock error at epoch i + 1, τ0 is the sampling interval of the clock error data in adjacent epochs, and yi is the obtained clock difference frequency data.After the detection of gross errors and clock jumps is completed for the frequency data, the processed frequency data can be recovered by the following formula to obtain new phase data, that is, the clock error data.(5)xi=∫0tyλdλ.In order to check the characteristics of phase data and frequency data, in this paper, IGS postevent satellite clock error is used for experimental analysis, and one day's GPS satellite clock error data are randomly selected with a sampling rate of 30 s. Taking the PRN10 satellite as an example, its satellite clock error data and corresponding frequency data are given, as shown in Figure1 and Figure 2.Figure 1
PRN10 satellite clock error data.Figure 2
PRN10 satellite frequency data.It can be seen from Figures1 and 2 that the clock error data (phase data) of the No. 10 satellite are in a linear trend as a whole, and its order of magnitude is 10–4. It is difficult to see with the naked eye whether there are outliers such as gross errors in the data. However, the frequency data obtained through mathematical changes are far less than the order of magnitude of the satellite clock error. Since the gross error value is small, it is easier to find outliers such as gross error using frequency data.Most gross error detection methods are processed based on satellite frequency data. Currently, the robust estimation method based on median (MAD) is used more, and the median representation is as follows:(6)MAD=Medianyi−m0.6745.In the formula,yi is the satellite frequency data, and m=Medianyi is the middle number of the frequency data sequence.After the algorithm obtains the median (MAD), the middle numberm of the frequency data is added to the n times of the obtained median, and the sum of the above two items is compared with the magnitude of the absolute value yi of each frequency. According to its relationship to judge gross error, the formula is as follows:(7)yi>m+n⋅MAD.Among them, the integern is selected according to the actual situation, and the other parameters are the same as those in formula (6).When the satellite frequency data satisfy the above formula, the data are judged as gross error and marked. After the gross error detection is carried out on all the frequency data, according to the number of marked epochs, the frequency data of the corresponding epoch are set to 0, and the frequency data are filled by means of mathematical calculation. The method is simple and easy to implement, and has weak sensitivity to gross error size, strong error resistance, and good applicability to frequency data. Using this method, the satellite frequency sequence in Figure2 is processed for gross error detection, and the processed satellite frequency data are shown in Figure 3.Figure 3
Frequency data of PRN10 satellite after removing gross errors.By comparing Figures2 and 3, it can be seen that the frequency data quality after using the median method to remove gross errors has been significantly improved, and the overall fluctuation is relatively stable. The clock error data with better quality are obtained by restoring the processed satellite frequency data. Using these data for modeling prediction will greatly improve the accuracy of the clock error prediction value.In order to judge whether there is a clock jump in the clock difference data, this paper adopts the following methods to detect: the algorithm first processes the satellite frequency data using the median method mentioned above. After the processing is completed, the epochs marked as gross errors are stored, and sometimes adjacent epochs are marked at the same time. There is a high probability of a clock jump in this case. The algorithm sums the adjacent two items of the satellite frequency data to form a new data sequenceyi¯=yi+yi+1 and marks the epoch according to the above stored gross errors. If yi¯ only has a peak value in the marked epoch, the epoch is considered to be a gross error, and if both the marked epoch and its next epoch have peak values, it can be judged that the abnormal value is a clock jump. The specific judgment is shown in Figure 4.Figure 4
Schematic diagram of the difference between gross errors and clock jumps in satellite clock errors.


### 2.3. Commonly Used Clock Error Prediction Models

The quadratic polynomial model of satellite clock error forecast is as follows:(8)xi=a0+a1ti−t0+a2ti−t02+ε.In the formula,i=1,2,…,n;xi is the value of the satellite clock offset at time ti; t0 is the satellite clock reference time; a0,a1,a2 are the satellite phase, frequency offset, and frequency drift at the reference time; and ε is the random error of the model.The linearization formula of this model is as follows:(9)L=Ba+v.In the formula, there are(10)B=1t1−t0t1−t021t2−t0t2−t02⋮⋮⋮1ti−t0ti−t02,L=x1x2⋮xi,v=v1v2⋮vi.When the number of satellite clock errors exceeds three, based on the principle of least squares, the estimation of the parameters in the fitting function can be obtained as follows:(11)a^=a^0a^1a^2T=BTB−1BTL.The fitting of the clock error sequence can be completed by substituting the parameter estimates into formula (8). After the model is determined, the time extrapolation can be used to complete the prediction of the satellite clock error at the next moment(12)x^i+1x^i+2⋮x^i+n=1ti+1−t0ti+1−t021ti+2−t0ti+2−t02⋮⋮⋮1ti+n−t0ti+n−t02a^0a^1a^2.In the formula,n is the number of forecast clock errors.In the actual satellite clock error prediction, the appropriate polynomial order is generally selected according to the actual satellite clock characteristics. When the frequency drift of the satellite clock cannot be ignored, the quadratic polynomial fitting model is often used in the actual clock error prediction; otherwise, the linear model is used.The clock error prediction model with additional periodic term is as follows:(13)xi=a0+a1ti−t0+a2ti−t02+∑k=1pAksin2πfkti−t0+φk+ε.In the formula, the first half is the quadratic polynomial model, and the latter is the periodic term. Among them,Ak,fk,andφk are related to the period term, which are amplitude, frequency, and phase, respectively; p is the number of main periods; and the rest parameters are the same as formula (8).After the parametersp and fk are determined, the unknown parameters in the model are solved by using the clock difference data after substituting into formula (13). For the convenience of calculation, we set as follows:(14)bk=Akcosφk,ck=Aksinφk.Then, the linearized form of formula (13) is as follows:(15)xi=a0+a1ti−t0+a2ti−t02+∑k=1pbksin2πfkti−t0+ckcos2πfkti−t0+ε.Its matrix form is as follows:(16)L=Aa+v.In the formula,L is the n-dimensional satellite clock error data, a=a0a1a2b1c1⋯⋯bpcpT is the parameter to be determined, A is the coefficient matrix, and the dimension is n×2p+3. The specific expression is as follows:(17)A=A1A2,Among them, there are(18)A1=1δt1δt121δt2δt22⋮⋮⋮1δtnδtn2,A2=sin2πf1δt1cos2πf1δt1⋯sin2πfpδt1cos2πfpδt1sin2πf1δt2cos2πf1δt2⋯sin2πfpδt2cos2πfpδt2⋮⋮⋮⋮⋮sin2πf1δtncos2πf1δtn⋯sin2πfpδtncos2πfpδtn.In the formula,δtn=tn−t0.The same polynomial model is solved, and the least squares principle is used as the criterion to obtain the estimated value of the unknown parametera^:(19)a^=a^0a^1a^2b^1c^1⋯⋯b^pc^pT=ATA−1ATL.After obtaining the parameter estimates in the model, we substitute them into formula (13) to complete the model fitting and then perform the clock error prediction at the next moment as follows:(20)L′=A′a^.In the formula,A′=A1′A2′. Among them, there are(21)A1′=1δtn+iδtn+i2,A2′=sin2πf1δtn+icos2πf1δtn+i…sin2πfpδtn+icos2πfpδtn+i.The basic contents of the grey model suitable for satellite clock error forecasting are as follows.We assume that there is a sequence of original clock offsets, as follows:(22)X0=x01,x02,…,x0n.The new clock error sequence obtained by accumulation is as follows:(23)X1=x11,x12,…,x1n,Among them, x1k=∑i=1kx0i,k=1,2,…,n.The differential formula of the grey model is as follows:(24)dxdt+bx=u.Substitute the new data sequence generated by accumulating satellite clock error data into the above formula, we get(25)dX1dt+bX1=u.In the formula,t∈0,+∞.Integrating formula (25) in the interval [k, k+1], we have(26)X1k+1−X1k+b∫kk+1X1tdt=u.In the formula,k=1,2,…,n−1.Taking into accountX0k+1=X1k+1−X1k, and expressing the integral in the form of mean, the above formula can be written as follows:(27)X0k+1=−bZ1k+1+u.In the formula,Z1k+1=X1k+X1k+1/2 is the arithmetic mean of points X1k、X1k+1.The matrix form of formula (27) is as follows:(28)L=Cb+v.In the formula, there are(29)L=X02X03⋮X0n,C=−Z121−Z131⋮⋮−Z1n1,b=au.Based on the principle of least squares, the parameterb^=a^u^=CTC−1CTL of the grey model is obtained.Substituting the obtained parameters into the differential formula, taking into account the initial parameterX11=X01, we get(30)X1k+1=X01−uae−ak+ua.The above formula can be used for prediction, and the obtained prediction value is the new data after one accumulation, and the accumulation and subtraction operation is performed, and then, the prediction value of the original clock error data before accumulation is as follows:(31)X0k+1=X1k+1−X1k=1−eaX01−uae−ak.Formula (31) is the grey model suitable for satellite clock error forecasting. The model can perform modeling prediction with more than 4 data volumes, and its effect is better in the long-term prediction of clock errors.The autoregressive moving average model ARMA (p, q) is expressed in the following form:(32)xt=φ1xt−1+φ2xt−2+⋯+φpxt−p+at−θ1at−1−θ2at−2−⋯−θqat−q.In the formula,xt is the data sequence calculated by the substitution model; φi and θi are the unknown autoregressive (AR) coefficient and moving average (MA) coefficient, respectively; p and q refer to the order of the autoregressive model and the order of the moving average model, respectively; and at is white noise, which conforms to the normal distribution N0,σ2.


### 2.4. BP Neural Network Algorithm

Neural networks are effective in fitting and predicting nonlinear time series, in which BP neural network is the most widely used. Through iteration, the output value is gradually approached to the expected value of the output.The BP neural network structure can set up multiple hidden layers. However, too many hidden layers will increase the computational workload. In order to shorten the calculation amount and speed up the calculation processing, this paper sets the hidden layer as one layer. During the training process, the input layer of each layer node is calculated as follows:(33)Inn,j=∑j=1nωijxi+θj.In the formula,n is the number of network layers, ωij is the weight of the n-th layer of the network, θj is the threshold of the n-th layer, xi is the input data, which corresponds to the input layer in Figure 1, and Inn,j=∑j=1nωijxi+θj in the formula (33) is the intermediate parameter of the model.The activation function often uses the sigmoid function(34)fx=11+e−x.The function is continuously differentiable, and it can project data from−∞,+∞ into the interval [0, 1].After the input layer is calculated, the corresponding output layer is calculated as follows:(35)Outn,j=fInn,j=11+e−∑j=1nωijxi+θj.In the formula,Outn,j is the output value obtained by the BP neural network.According to the gradient descent algorithm, the estimated valueOutn,j of the network output is compared with the real value of the product released by IGS, and the error value is obtained. Through the error back propagation, the BP neural network completes the update of the weight and threshold parameters, and the update formula is as follows:(36)ωt+1=ωt+ηεtyt,θt+1=θt+ηεtyt.In the formula,t is the number of iterations; εt is the error term, that is, the difference between the expected output value and the actual output result; η is the learning efficiency of the neural network; and y(t) is the output value of the neuron.


### 2.5. Thinking Evolution Algorithm

Mind evolutionary algorithm (MEA) deals with optimization problems in the way of evolution of human thinking. The algorithm first uses the convergence operation to optimize the individuals in the subgroup, and then, the mature subgroups compete globally through the alienation operation, which greatly improves the optimization efficiency. Figure5 is a schematic diagram of the structure of the MEA.Figure 5
:Schematic diagram of the structure of the MEA.We assume that there is a set of satellite clock offset datax1,x2,…,xn, and use the data at these n times to model and predict the clock offset value at future times. Usually, the number of output layers is the same as the number of output data types, and the output value in this paper is the clock difference. Therefore, the output node is set to 1, and the neural network is used to establish the mapping relationship between x1,x2,…,xN and xn+1 (N is the number of input nodes), and the idea of sliding window is adopted. On the premise of keeping the number of samples unchanged, the new prediction data are constantly used to replace the previous known data, and the multi-epoch satellite clock error prediction is realized.According to formula (37) and formula (38), the algorithm first calculates the mean square error E and then takes the reciprocal of it to obtain the score function score of various groups and individuals.(37)E=1m∑i=1my^i−yi2,(38)score=1E.In the formula,y^i is the output value after iteration, yi is the sample value, and m is the total number of samples.The algorithm flow is shown in Figure6.Figure 6
MEA to optimize BP neural network process.After completing the neural network training using n-dimensional modeling data, according to the neural network structure, the number of input nodes is N, and the number of output nodes is 1. The algorithm adopts the idea of sliding window to ensure that the number of samples remains unchanged and so on to achieve multistep forecasting.We assume thatX=x1,x2,…,xn is a set of n-dimensional satellite clock error sequences, where xi,i=1,2,…,n is the clock error data of different epochs. The algorithm divides the clock difference data into two groups and makes the difference in turn to obtain a new first-order difference data sequence(39)ΔX=Δx1,Δx2,…,Δxn−1.Among them, Δxi=xi+1−xi.After the original clock error is processed by one time difference, different from the conventional operation, the data used to solve the parameters in the modeling process are the data sequence after one time difference. When determining the model, a neural network is used to model the formed data sequenceΔX, and the satellite clock error of the next m epochs is predicted: the algorithm models the n-1-dimensional data, and after the modeling process is completed, the m-dimensional data sequence is predicted, ΔX′=Δxn,Δxn+1,…,Δxm+n−1. The algorithm recovers the clock error of the epoch to be determined by using the predicted first-order difference and the modeling last clock error data x(n), that is,(40)xk=xn+∑nk−1Δxk.


## 2.1. The Basic Theory of Satellite Clock Error Forecasting

Satellite clock error data, also known as phase data, are the time offset data obtained by standard time and atomic clock output time, in nanoseconds or seconds. In satellite navigation and positioning, phase data represent satellite clock error data and receiver clock error data. The raw phase data can be directly used for navigation positioning and satellite clock error prediction. However, the magnitude of the phase data is large, and it is impossible to judge and identify the location of outliers such as gross errors and clock jumps. Therefore, the data need to be preprocessed before the actual clock error forecasting. By making corresponding mathematical changes to the phase data, the corresponding frequency data can be obtained, and by processing the frequency data, the gross errors and clock jumps existing in the clock error can be effectively detected.The relative deviation of the standard output frequencies of different frequencies can obtain frequency data, which is defined as:(1)D=fA−fBf0.In the formula,fA is the output frequency of the frequency standard A, fB is the output frequency of the frequency standard B, f0 is the nominal frequency, and D is the obtained frequency data.Over time, the actual output frequency of the atomic clock will be affected by the internal components and the external environment, and it usually presents a linear trend of increasing or decreasing. When time changes, the frequency of satellite clock also changes, and the rate of change is called frequency drift rate or aging rate. The calculation formula is as follows:(2)zi=yi+1−yiτ0.In the formula,yi is the frequency value of epoch i, yi+1 is the frequency value of epoch i + 1, and τ0 is the sampling interval of adjacent epochs.The frequency of atomic clocks is often not a specific value, but fluctuates within a certain range. It obtains correct information on whether the atomic clock frequency is accurate or not by comparing the actual frequency obtained by measuring the 1BP neural network algorithm or calculation with the nominal frequency. The calculation formula is as follows:(3)A=f−f0f0.In the formula,f is the actual frequency measured or calculated, f0 is the real frequency, that is, the nominal frequency, and A is the calculated frequency accuracy. Since the frequency of atomic clocks fluctuates within a certain range, the average of multiple measurement results is usually taken as the final result.


## 2.2. Preprocessing of Satellite Clock Error Data

The frequency data are obtained from the phase data by mathematical transformation, and the conversion formula is as follows:(4)yi=xi+1−xiτ0.In the formula,xi is the value of the satellite clock error at epoch i, xi+1 is the value of the satellite clock error at epoch i + 1, τ0 is the sampling interval of the clock error data in adjacent epochs, and yi is the obtained clock difference frequency data.After the detection of gross errors and clock jumps is completed for the frequency data, the processed frequency data can be recovered by the following formula to obtain new phase data, that is, the clock error data.(5)xi=∫0tyλdλ.In order to check the characteristics of phase data and frequency data, in this paper, IGS postevent satellite clock error is used for experimental analysis, and one day's GPS satellite clock error data are randomly selected with a sampling rate of 30 s. Taking the PRN10 satellite as an example, its satellite clock error data and corresponding frequency data are given, as shown in Figure1 and Figure 2.Figure 1
PRN10 satellite clock error data.Figure 2
PRN10 satellite frequency data.It can be seen from Figures1 and 2 that the clock error data (phase data) of the No. 10 satellite are in a linear trend as a whole, and its order of magnitude is 10–4. It is difficult to see with the naked eye whether there are outliers such as gross errors in the data. However, the frequency data obtained through mathematical changes are far less than the order of magnitude of the satellite clock error. Since the gross error value is small, it is easier to find outliers such as gross error using frequency data.Most gross error detection methods are processed based on satellite frequency data. Currently, the robust estimation method based on median (MAD) is used more, and the median representation is as follows:(6)MAD=Medianyi−m0.6745.In the formula,yi is the satellite frequency data, and m=Medianyi is the middle number of the frequency data sequence.After the algorithm obtains the median (MAD), the middle numberm of the frequency data is added to the n times of the obtained median, and the sum of the above two items is compared with the magnitude of the absolute value yi of each frequency. According to its relationship to judge gross error, the formula is as follows:(7)yi>m+n⋅MAD.Among them, the integern is selected according to the actual situation, and the other parameters are the same as those in formula (6).When the satellite frequency data satisfy the above formula, the data are judged as gross error and marked. After the gross error detection is carried out on all the frequency data, according to the number of marked epochs, the frequency data of the corresponding epoch are set to 0, and the frequency data are filled by means of mathematical calculation. The method is simple and easy to implement, and has weak sensitivity to gross error size, strong error resistance, and good applicability to frequency data. Using this method, the satellite frequency sequence in Figure2 is processed for gross error detection, and the processed satellite frequency data are shown in Figure 3.Figure 3
Frequency data of PRN10 satellite after removing gross errors.By comparing Figures2 and 3, it can be seen that the frequency data quality after using the median method to remove gross errors has been significantly improved, and the overall fluctuation is relatively stable. The clock error data with better quality are obtained by restoring the processed satellite frequency data. Using these data for modeling prediction will greatly improve the accuracy of the clock error prediction value.In order to judge whether there is a clock jump in the clock difference data, this paper adopts the following methods to detect: the algorithm first processes the satellite frequency data using the median method mentioned above. After the processing is completed, the epochs marked as gross errors are stored, and sometimes adjacent epochs are marked at the same time. There is a high probability of a clock jump in this case. The algorithm sums the adjacent two items of the satellite frequency data to form a new data sequenceyi¯=yi+yi+1 and marks the epoch according to the above stored gross errors. If yi¯ only has a peak value in the marked epoch, the epoch is considered to be a gross error, and if both the marked epoch and its next epoch have peak values, it can be judged that the abnormal value is a clock jump. The specific judgment is shown in Figure 4.Figure 4
Schematic diagram of the difference between gross errors and clock jumps in satellite clock errors.


## 2.3. Commonly Used Clock Error Prediction Models

The quadratic polynomial model of satellite clock error forecast is as follows:(8)xi=a0+a1ti−t0+a2ti−t02+ε.In the formula,i=1,2,…,n;xi is the value of the satellite clock offset at time ti; t0 is the satellite clock reference time; a0,a1,a2 are the satellite phase, frequency offset, and frequency drift at the reference time; and ε is the random error of the model.The linearization formula of this model is as follows:(9)L=Ba+v.In the formula, there are(10)B=1t1−t0t1−t021t2−t0t2−t02⋮⋮⋮1ti−t0ti−t02,L=x1x2⋮xi,v=v1v2⋮vi.When the number of satellite clock errors exceeds three, based on the principle of least squares, the estimation of the parameters in the fitting function can be obtained as follows:(11)a^=a^0a^1a^2T=BTB−1BTL.The fitting of the clock error sequence can be completed by substituting the parameter estimates into formula (8). After the model is determined, the time extrapolation can be used to complete the prediction of the satellite clock error at the next moment(12)x^i+1x^i+2⋮x^i+n=1ti+1−t0ti+1−t021ti+2−t0ti+2−t02⋮⋮⋮1ti+n−t0ti+n−t02a^0a^1a^2.In the formula,n is the number of forecast clock errors.In the actual satellite clock error prediction, the appropriate polynomial order is generally selected according to the actual satellite clock characteristics. When the frequency drift of the satellite clock cannot be ignored, the quadratic polynomial fitting model is often used in the actual clock error prediction; otherwise, the linear model is used.The clock error prediction model with additional periodic term is as follows:(13)xi=a0+a1ti−t0+a2ti−t02+∑k=1pAksin2πfkti−t0+φk+ε.In the formula, the first half is the quadratic polynomial model, and the latter is the periodic term. Among them,Ak,fk,andφk are related to the period term, which are amplitude, frequency, and phase, respectively; p is the number of main periods; and the rest parameters are the same as formula (8).After the parametersp and fk are determined, the unknown parameters in the model are solved by using the clock difference data after substituting into formula (13). For the convenience of calculation, we set as follows:(14)bk=Akcosφk,ck=Aksinφk.Then, the linearized form of formula (13) is as follows:(15)xi=a0+a1ti−t0+a2ti−t02+∑k=1pbksin2πfkti−t0+ckcos2πfkti−t0+ε.Its matrix form is as follows:(16)L=Aa+v.In the formula,L is the n-dimensional satellite clock error data, a=a0a1a2b1c1⋯⋯bpcpT is the parameter to be determined, A is the coefficient matrix, and the dimension is n×2p+3. The specific expression is as follows:(17)A=A1A2,Among them, there are(18)A1=1δt1δt121δt2δt22⋮⋮⋮1δtnδtn2,A2=sin2πf1δt1cos2πf1δt1⋯sin2πfpδt1cos2πfpδt1sin2πf1δt2cos2πf1δt2⋯sin2πfpδt2cos2πfpδt2⋮⋮⋮⋮⋮sin2πf1δtncos2πf1δtn⋯sin2πfpδtncos2πfpδtn.In the formula,δtn=tn−t0.The same polynomial model is solved, and the least squares principle is used as the criterion to obtain the estimated value of the unknown parametera^:(19)a^=a^0a^1a^2b^1c^1⋯⋯b^pc^pT=ATA−1ATL.After obtaining the parameter estimates in the model, we substitute them into formula (13) to complete the model fitting and then perform the clock error prediction at the next moment as follows:(20)L′=A′a^.In the formula,A′=A1′A2′. Among them, there are(21)A1′=1δtn+iδtn+i2,A2′=sin2πf1δtn+icos2πf1δtn+i…sin2πfpδtn+icos2πfpδtn+i.The basic contents of the grey model suitable for satellite clock error forecasting are as follows.We assume that there is a sequence of original clock offsets, as follows:(22)X0=x01,x02,…,x0n.The new clock error sequence obtained by accumulation is as follows:(23)X1=x11,x12,…,x1n,Among them, x1k=∑i=1kx0i,k=1,2,…,n.The differential formula of the grey model is as follows:(24)dxdt+bx=u.Substitute the new data sequence generated by accumulating satellite clock error data into the above formula, we get(25)dX1dt+bX1=u.In the formula,t∈0,+∞.Integrating formula (25) in the interval [k, k+1], we have(26)X1k+1−X1k+b∫kk+1X1tdt=u.In the formula,k=1,2,…,n−1.Taking into accountX0k+1=X1k+1−X1k, and expressing the integral in the form of mean, the above formula can be written as follows:(27)X0k+1=−bZ1k+1+u.In the formula,Z1k+1=X1k+X1k+1/2 is the arithmetic mean of points X1k、X1k+1.The matrix form of formula (27) is as follows:(28)L=Cb+v.In the formula, there are(29)L=X02X03⋮X0n,C=−Z121−Z131⋮⋮−Z1n1,b=au.Based on the principle of least squares, the parameterb^=a^u^=CTC−1CTL of the grey model is obtained.Substituting the obtained parameters into the differential formula, taking into account the initial parameterX11=X01, we get(30)X1k+1=X01−uae−ak+ua.The above formula can be used for prediction, and the obtained prediction value is the new data after one accumulation, and the accumulation and subtraction operation is performed, and then, the prediction value of the original clock error data before accumulation is as follows:(31)X0k+1=X1k+1−X1k=1−eaX01−uae−ak.Formula (31) is the grey model suitable for satellite clock error forecasting. The model can perform modeling prediction with more than 4 data volumes, and its effect is better in the long-term prediction of clock errors.The autoregressive moving average model ARMA (p, q) is expressed in the following form:(32)xt=φ1xt−1+φ2xt−2+⋯+φpxt−p+at−θ1at−1−θ2at−2−⋯−θqat−q.In the formula,xt is the data sequence calculated by the substitution model; φi and θi are the unknown autoregressive (AR) coefficient and moving average (MA) coefficient, respectively; p and q refer to the order of the autoregressive model and the order of the moving average model, respectively; and at is white noise, which conforms to the normal distribution N0,σ2.


## 2.4. BP Neural Network Algorithm

Neural networks are effective in fitting and predicting nonlinear time series, in which BP neural network is the most widely used. Through iteration, the output value is gradually approached to the expected value of the output.The BP neural network structure can set up multiple hidden layers. However, too many hidden layers will increase the computational workload. In order to shorten the calculation amount and speed up the calculation processing, this paper sets the hidden layer as one layer. During the training process, the input layer of each layer node is calculated as follows:(33)Inn,j=∑j=1nωijxi+θj.In the formula,n is the number of network layers, ωij is the weight of the n-th layer of the network, θj is the threshold of the n-th layer, xi is the input data, which corresponds to the input layer in Figure 1, and Inn,j=∑j=1nωijxi+θj in the formula (33) is the intermediate parameter of the model.The activation function often uses the sigmoid function(34)fx=11+e−x.The function is continuously differentiable, and it can project data from−∞,+∞ into the interval [0, 1].After the input layer is calculated, the corresponding output layer is calculated as follows:(35)Outn,j=fInn,j=11+e−∑j=1nωijxi+θj.In the formula,Outn,j is the output value obtained by the BP neural network.According to the gradient descent algorithm, the estimated valueOutn,j of the network output is compared with the real value of the product released by IGS, and the error value is obtained. Through the error back propagation, the BP neural network completes the update of the weight and threshold parameters, and the update formula is as follows:(36)ωt+1=ωt+ηεtyt,θt+1=θt+ηεtyt.In the formula,t is the number of iterations; εt is the error term, that is, the difference between the expected output value and the actual output result; η is the learning efficiency of the neural network; and y(t) is the output value of the neuron.


## 2.5. Thinking Evolution Algorithm

Mind evolutionary algorithm (MEA) deals with optimization problems in the way of evolution of human thinking. The algorithm first uses the convergence operation to optimize the individuals in the subgroup, and then, the mature subgroups compete globally through the alienation operation, which greatly improves the optimization efficiency. Figure5 is a schematic diagram of the structure of the MEA.Figure 5
:Schematic diagram of the structure of the MEA.We assume that there is a set of satellite clock offset datax1,x2,…,xn, and use the data at these n times to model and predict the clock offset value at future times. Usually, the number of output layers is the same as the number of output data types, and the output value in this paper is the clock difference. Therefore, the output node is set to 1, and the neural network is used to establish the mapping relationship between x1,x2,…,xN and xn+1 (N is the number of input nodes), and the idea of sliding window is adopted. On the premise of keeping the number of samples unchanged, the new prediction data are constantly used to replace the previous known data, and the multi-epoch satellite clock error prediction is realized.According to formula (37) and formula (38), the algorithm first calculates the mean square error E and then takes the reciprocal of it to obtain the score function score of various groups and individuals.(37)E=1m∑i=1my^i−yi2,(38)score=1E.In the formula,y^i is the output value after iteration, yi is the sample value, and m is the total number of samples.The algorithm flow is shown in Figure6.Figure 6
MEA to optimize BP neural network process.After completing the neural network training using n-dimensional modeling data, according to the neural network structure, the number of input nodes is N, and the number of output nodes is 1. The algorithm adopts the idea of sliding window to ensure that the number of samples remains unchanged and so on to achieve multistep forecasting.We assume thatX=x1,x2,…,xn is a set of n-dimensional satellite clock error sequences, where xi,i=1,2,…,n is the clock error data of different epochs. The algorithm divides the clock difference data into two groups and makes the difference in turn to obtain a new first-order difference data sequence(39)ΔX=Δx1,Δx2,…,Δxn−1.Among them, Δxi=xi+1−xi.After the original clock error is processed by one time difference, different from the conventional operation, the data used to solve the parameters in the modeling process are the data sequence after one time difference. When determining the model, a neural network is used to model the formed data sequenceΔX, and the satellite clock error of the next m epochs is predicted: the algorithm models the n-1-dimensional data, and after the modeling process is completed, the m-dimensional data sequence is predicted, ΔX′=Δxn,Δxn+1,…,Δxm+n−1. The algorithm recovers the clock error of the epoch to be determined by using the predicted first-order difference and the modeling last clock error data x(n), that is,(40)xk=xn+∑nk−1Δxk.


## 3. Comparative Test Analysis

In order to verify the optimization effect of the thinking navigation algorithm on the flight trajectory of the aircraft, the effect of the method in this paper is studied by means of comparative experiments.The modeling uses 12 hours of data, the satellite clock error data sampling interval is 30 s, and there are 1440 data in total. In this paper, the input layer is set to 1, and the number of hidden layer nodes is in the range of [1 ,10]. This paper takes the G02 satellite as an example and uses its clock error data to conduct experiments. Different mean square errors (RMSE) can be obtained by selecting different hidden layer nodes. The relationship between the two is shown in Figure 7.Figure 7
The relationship between the number of hidden layer nodes and the model accuracy.It can be seen from the analysis in Figure7 that the overall fitting accuracy of the model is relatively stable, fluctuating up and down at 0.3 ns. When the number of hidden layer nodes is 5, the prediction accuracy suddenly deteriorates, and the curve fluctuates greatly. The reason is that when there are too many hidden layer nodes, the BP neural network algorithm is prone to “overfitting.” That is, although the model has a high fitting accuracy at this time, the prediction accuracy suddenly deteriorates. If the number of hidden layer nodes is not selected properly, it will seriously affect the prediction accuracy of the model. According to the experimental results and fully considering the fitting and prediction accuracy, this paper selects 3 as the number of nodes in the hidden layer of the BP neural network.In summary, in this paper, the input layer node is set to 1, the output layer node is set to 1, and the hidden layer node is calculated to be 3.The random selection of initial parameters in BP neural network will increase the convergence time and workload in the training process, which may cause abnormal results in the model solving process and affect the stability of the model. Moreover, the results of the algorithm are easy to fall into the local minimum, which reduces the performance of the model to obtain the global optimum. In this paper, the MEA is used to select the initial parameters required by the BP neural network. In the MEA, the population size is set to 1440, and 5 winning subpopulations and 5 temporary subpopulations are selected. In order to make the algorithm fully find the optimal individual, the number of iterations is set to 200.The prediction performance of the BP model and the MEA-BP model for satellite clock errors is compared. In order to fully compare the forecasting effects, the clock offset data of the first 12 hours are used to independently predict the clock offset values of the next 3 hours, 6 hours, and 12 hours for 10 times. Figure8 shows the changes of the RMS values predicted by the two models.Figure 8
Variation of 10 times RMS of MEA-BP and BP forecast (a) (b) (c).
(a)(b)(c)In this experiment, the MEA-BP model and the traditional BP model are used to forecast the No. 2 satellite in different time periods. It can be seen from Figure8 that under the same network structure, the multiple prediction accuracy of the MEA-BP model and the BP model does not fluctuate significantly, which proves the feasibility of the network structure. When the forecast time increases, the MEA-BP model obtains higher accuracy and has better forecast performance than the traditional BP model. In the clock error prediction of the three time periods, the prediction accuracy of the MEA-BP model is better than that of the BP model. It shows that the initial parameters (weights and thresholds) required by the BP neural network are better selected by MEA, so as to avoid the local optimal results obtained by the BP model, and effectively improve the accuracy of the prediction results. Moreover, it shows the feasibility of the MEA-BP model for clock error prediction, so it is a relatively stable aircraft flight trajectory optimization algorithm.On the basis of the above research, the role of the thinking navigation proposed in this paper in the optimization of the flight trajectory of the aircraft is analyzed through multiple sets of simulation experiments, and the simulation evaluation results shown in Table1 are obtained.Table 1
Evaluation of aircraft flight trajectory optimization.
NumTrajectory optimizationNumTrajectory optimization188.481091.48288.491190.77388.051291.57489.391389.74590.291490.46691.141591.02789.401691.06889.071791.11988.041889.94From the above experimental research, it can be seen that the thinking navigation algorithm proposed in this paper has a very obvious effect on the optimization of the flight trajectory of the aircraft.


## 4. Conclusion

Track tracking is based on the precompleted flight track, and it performs flight track planning in real time according to the external interference, environment, and other internal and external factors in the actual flight. When completing a certain task or operation, it calculates and tracks the required optimal flight path, and generates navigation and control instructions. Correspondingly, it controls the force and attitude angle of the aircraft (including engine thrust, attack angle, roll angle, and heading angle) to make the aircraft fly along the planned flight path. This paper combines the thinking navigation algorithm to optimize the flight trajectory of the aircraft and analyzes the flight trajectory of the aircraft through the intelligent model. The simulation test study shows that the thinking navigation algorithm proposed in this paper has a very obvious effect on the optimization of the flight trajectory of the aircraft.



---

*Source: 6169331-2022-08-03.xml*