#### Abstract

A mathematical model of the dive phase is an important research content for improving the accuracy of terminal control in the small unmanned aerial vehicle. The acquisition of the diving model poses new challenges, such as the small installation space, ultra-low flying height of small suicide drones, short flight time, strong coupling, less observable measurement, and elastic deformation of the wings during the drone dive phase. Based on the autoregressive moving average method, a multi-input multioutput noise term topology mathematical model is proposed in this paper. Through an improved least squares identification method, the diving model in the flight test is analyzed and verified. The identification results of the diving model obtained by the proposed method are compared with the least squares method dive model. The results indicate that the mathematical model and identification method proposed in this paper can effectively obtain the parameters of the drone dive model.

#### 1. Introduction

With the development of SUAV technology, the function of unmanned aerial vehicles (UAVs) has gradually expanded from a single function of scouting, interference, and damage assessment to the integration of scouting and strikes. During the execution of scouting and searching above a target area, the SUAV can immediately hit the target when a time-sensitive target is found. Due to strong maneuverability and high attack precision (with a controllable range at the meter level), the SUAV can effectively destroy the target while carrying fewer combat damage units, thereby exhibiting high benefits in actual combat. Representative SUAVs are the Switchblade, the Coyote, and the Hero loitering munition.

During attacks, dive modeling is quite important in the end trajectory control. To enhance the controlling precision of SUAV during attacks, the diving model should reflect the SUAV's actual performance. Researchers mainly adopted wind tunnel tests and computational fluid dynamics (CFD) software for the establishment of mathematical models of conventional air vehicles [1]. Both methods require high professional skills and are time-consuming as well as costly. Currently, certain engineering calculation software of aerodynamic parameters exhibit poor computational accuracy, which can thereby affect the precision of the entire mathematical model. By utilizing the system response to specific input single utilization, the system identification method can rapidly and effectively acquire the mathematical model of the entire air vehicle, which has been extensively applied during flight tests of conventional air vehicles [2]. Since the introduction of “identification” by American scholar Zadeh in 1956 [3], system identification, a method of system modeling, has been continuously developed for over 60 years. Following booming, the system identification was rapidly applied to system modeling of air vehicles [4–6]. The aforementioned has become a mature modeling method in the conventional air vehicle field.

Recently, with the development of low-cost microelectromechanical system (MEMS) sensor technology, system identification has been utilized in extensive applications in the SUAV field. With the unscented Kalman filter (UKF), Mararu Naruoka et al. [7] performed system identification on a Mitsubishi MARS07AF. This was an SUAV carrying MEMS devices in 2009. However, the identification results significantly differed from the wind tunnel test data. In 2010, Khaled S. Hatamleh et al. used an inertial measurement unit (IMU) on the UAV for the related model parameters estimation [8]. In 2011, Andrei Dorobantu et al. employed the frequency-domain identification technique for both longitudinal and transverse model identification of a low-sized fixed-wing UAV [9]. In 2012, Aaron Wypyszynski conducted model-identification flight tests on a small UAV and identified its aerodynamic parameters [10]. In 2012, Arnaud Koehl et al. established a 6-DOF model of GLMAV, a coaxial-rotor UAV, and performed parameter identification through Kalman filtering [11]. In 2015, Konstantin identified the four-rotor UAV model by using the recursive algorithm [12]. Kailong Liu et al. employed a recursive algorithm to investigate the battery temperature prediction method with regard to the battery safety problem [13, 14]. In addition, with the improvement of hardware processing capabilities, methods such as machine learning and Kalman filtering are also used in model identification and parameter prediction of complex systems [15–18].

Model identification of small UAVs research is currently mainly focused on longitudinal or horizontal models in the cruise phase [19–23]. The research on the diving phase of SUAVs is mainly focused on control rate, guidance rate, and recognition of the diving model. In 2017, Zhou [24] conducted an integrated design for the guidance and controller of the suicide UAV. Moreover, the author verified the effectiveness of the system. He [25] proposed a hypersonic gliding vehicle spiral dive maneuver penetration concept and designed an adaptive proportional guidance rate for virtual portrait targets. Karpenko [26] investigated the evaluation method of attack aircraft's dive characteristics and analyzed the time and height loss of attack aircraft in different dive and pull modes. Wang [27] studied the guidance and control laws of antitank missiles with curved firing and self-seeking. Moreover, the author established a guidance and control law suitable for ballistic control. With the rapid expansion of military applications of suicide attacks by small UAVs and to improve the accuracy of UAV strikes, the problem of accurate modeling and identification of the diving model needs to be solved urgently. As an SUAV gradually approaches the target during diving, the position and attitude should be constantly corrected. This phase has strong nonlinear characteristics and is a multi-input multioutput system [8, 9]. In general, due to installation space and landing weight limitations, testing equipment such as attack angle and sideslip angle sensors cannot be installed. Moreover, only a certain amount of data can be acquired. In addition, since the SUAV always flies at a low velocity, the flight is easily disturbed by many factors, such as gusts.

In this paper, a topological mathematical model of multi-input and multioutput noise terms for SUAV is first constructed. Then, an improved multivariable augmented recursive least squares identification method is proposed, and the feasibility of the method is verified by the simulation model. Finally, the identification methods are compared and verified through flight tests of SUAV.

#### 2. Mathematical Model

The diving phase of SUAV is a very complex and highly coupled nonlinear dynamic system. Its six-degree-of-freedom motion equations can be expressed as [28, 29]where , , and are body velocities along with the *x*-axis, *y*-axis, and *z*-axis. Parameters *p, q*, and *r* are angular rates along the *x*-axis, *y*-axis, and *z*-axis. Parameters , , and are Euler angles along the *x*-axis, *y*-axis, and *z*-axis. Coefficients , , , , , and are aerodynamic force and moment coefficients. is the coordinate transformation matrix, and *m*, *S*, and *b* are the quality, reference area, and wingspan of the small UAV, respectively. Parameter *q* is the air density, is the inertia matrix, with being the rolling moment of inertia, being the pitching moment of inertia, being the yawing moment of inertia, and and are inertia products. The parameter is the airspeed, while , , and are airspeed vector components along the *x*-axis, *y*-axis, and *z*-axis in the body frame. Parameter *n* is the rotational speed of the propeller, *D* is the diameter of the propeller, and is the thrust coefficient.

The aerodynamic force of SUAV has a very complicated nonlinear relationship with many factors such as the geometry, flight speed, aerodynamic angle, and even atmospheric temperature of the aircraft. During the flight test, factors such as aerodynamic angle and wing deformation do not have the measurement conditions. If all these factors are simultaneously taken into account, model identification and solution will be relatively difficult. Limited by the research funding, the abovementioned complex parameters can be hardly obtained. Therefore, the ARMA model is used in this paper to solve this problem. This method treats the model as a gray box, which can fit the model according to the input and output data of the system. Moreover, it can predict the output features of the model with sufficient accuracy, thereby reducing the amount of calculation and the difficulty of the identification test.

The ARMA model is a common and high-precision time-series short-term prediction method. It exhibits extensive applications in many domains such as fault detection, controller model identification, sensor correction, and structural testing and identification [30–33]. The dynamic discrete system input and the output sampling value series can be expressed as the following linear difference equation:where and are the coefficients of autoregression items.

The multi-input multioutput system of the SUAV, consisting of *r*-dimensional inputs and *m*-dimensional outputs, is displayed in Figure 1. The control input matrix and system output matrix are represented as and , respectively, while the sustained colored noise matrix can be expressed as .

According to the analysis of the characteristics during diving, the structure of the diving model can be expressed aswhere denotes a 22 matrix to be identified, is also a 22 matrix to be identified, and denotes a 22 colored noise matrix in the model indicating interference factors such as wind gusts and temperature during the test. Parameters and are the control input matrix and the state output matrix, respectively. and can be written aswhere and are the *k*-th PWM control input data in pitching and rolling channels, while and denote the *k*-th pitching angular rate and the rolling angular rate corresponding to control inputs, respectively.

To determine the structural parameters and , the one-step-ahead prediction value of the model output can be written as

The parameter is calculated through least squares, while the output prediction error vector can be defined as

The loss function of the model can be defined as

#### 3. Algorithm

The least squares method was proposed by Gauss to determine the orbit of Ceres. Within the method, the optimal functional matching was searched for by employing the sum minimization of squares of the errors. Currently, the least squares method is a common method for parameter estimation on the test data that exhibits extensive applications and popularization in many domains. Despite simplicity and practicability, the least squares method is no longer an unbiased estimation for the system including colored noise. Therefore, given that an SUAV is a multi-input multioutput coupling system during the diving phase, while the colored noise imposed considerable disturbance, a multivariable recursive augmented least squares identification method is derived in this paper.

##### 3.1. Standard Augmented Recursive Least Squares Method

The standard recursive least squares identification algorithm can be expressed as

During the actual small aircraft test process, due to the influence of gust and other environmental effects, the actual noise is not white noise. Therefore, it is necessary to identify the noise model and take it into account. Noise model parameters are jointly identified. The noise model can be expressed aswhere represents random noise with the nonzero mean, represents random noise with the mean value of 0, and represent the parameters of the noise model that require identification. The model can be expressed as

Noise sequence can be replaced by the estimation of :where

The noise sequence can be supplemented as the identified object. Parameters to be reformulated are as follows:

The model can be expressed as

The calculation expression of the augmented least squares recursive method can be formulated as

##### 3.2. Improved Multivariable Augmented Recursive Least Squares Method

The multivariable identification model can be defined aswhere

To improve the model's ability to deal with the noise interference in the experiment, a noise topology model is added, and parameter is augmented. Therefore, (16) can be written aswhere

The *j*-th column in the model can be rewritten as

In a matrix form, (20) can be further rewritten aswhere

If the noise is a zero-mean unideal random series, the consistency and unbiased estimated value of can be described as

If , the estimated values of the parameters in various rows, denoted as , can be acquired. Therefore, the estimated parameter value of the system can be acquired.

The inverse calculation of a matrix demands heavy calculation, whereas a singular matrix easily appears during the test data processing. Therefore, the following expressions can be set:

The multivariable recursive augmented least squares identification algorithm process is shown in Figure 2 and can be consequently described as

For the identification performance analysis, the following multi-input multioutput simulation model was established:where denotes the white noise that follows a normal distribution in the diving phase. A four-displacement register with an input amplitude of 1 is used for the signal input.

#### 4. Simulation Verification

To verify the effectiveness of the algorithm and analyze the identification performance of the modified multivariate augmented recursive least squares method, a simulation model is taken as an example. Figure 3 represents the comparison of the model noise identification results. The noise prediction of the system by the multivariate augmented recursive least squares method is closer to the real noise after approximation. The first channel noise prediction model characteristic parameters are within the range of , and the second channel noise prediction model characteristic parameters are within the range of , . The identification results are shown in Figure 4 and Table 1.

**(a)**

**(b)**

**(a)**

**(b)**

#### 5. Flight Test

An SUAV generally follows a single-wing and V-shaped stern layout. Pitching and yawing are controlled by the V-shaped stern control surface while rolling is manipulated by the aileron. Electrically powered two-blade propellers are used for power supply, and the electric motor is installed behind the fuselage, i.e., a push-type SUAV is employed. As presented in Figure 5, a single ejection device is used for the launch. Table 2 lists the basic parameters of this SUAV.

After launching the ejection device, the SUAV is controlled by the autopilot and a user in the loop to ensure that the SUAV flies in accordance with the preset trajectory at a given velocity and height. During the flight, the target region can be monitored by an airborne camera. When a target is found, the SUAV dives, attacks, and destroys the target.

Many measurement devices including GPS receiver, triaxial MEMS accelerometer, triaxial angular velocity gyroscope, air speedometer, and barometer are integrated into the navigation, guidance, and control system of SUAV. More specifically, longitude, latitude, height, and the heading information of the SUAV are provided by the GPS receiver. The accelerated velocity and angular velocity information are measured by the triaxial accelerometer and triaxial angular velocity gyroscope, respectively. The airspeed is measured by the airborne airspeed sensor.

Before the flight test, the UAV sensor and navigation equipment should be calibrated to reduce the noise impact. The flight test should be carried out in sunny weather and with low wind presence. Wind speed measurement curve of the flight test is displayed in Figure 6. The initial state of the drone was adjusted to the same settings before each dive test. To improve the validity of the data and reduce the effect of noise in the test data, data processing such as outlier processing, sampling time correction, data interpolation encryption, data coordination analysis, and data reconstruction were carried out on the preidentification data of the dive test [34, 35].

The flight path of the SUAV is displayed in Figure 7. During testing, when the SUAV finds the target, it adjusts the attacking posture and consequently dives to attack. The SUAV performed five diving attacks in total. During the first four diving attacks, when the SUAV dived and approached the target, it was rapidly pulled up by the remote control handle to achieve attitude adjustment. During the fifth diving attack, the SUAV directly stroke the target.

The adjusted pulse width outputs of each control manipulate channel during the simulated diving attack are displayed in Figure 8. The manipulated variables in the pitching and rolling channels exhibited the most drastic changes. Furthermore, these were consequently regarded as the main controlled variables that also constitute the coupled control input. The throttle channel during diving was controlled by the remote device under a constant input control pattern.

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

#### 6. Results and Discussion

Through the model loss functions corresponding to different parameter combinations comparison, the model structure of the SUAV in the diving phase can be written aswhere, , , , , , , and denote the parameter matrices to be identified.

The first four sets of test data in the diving flight test were used as the identification data, while the fifth set of test data was used as the validation data. Table 3 lists the identification results.

Figure 9 presents the convergence curves of the identification parameters of the matrices , , , , , , and . At the beginning of iteration, the identification parameters have significantly changed, which also increased the predicted model noise, as presented in Figure 10. Following approximately 100 iterations, the parameters exhibited gentle changes and slight adjustment.

**(a)**

**(b)**

**(a)**

**(b)**

Figure 11 is a comparison curve of the identification results of a small drone dive model. The identification method proposed in this paper still predicts the changes in the pitch and roll angular rates well. The average errors of the pitch and roll angular rate are 0.0143 rad/s and 0.00868 rad/s, respectively, which is better than conventional least squares identification results (average error of 0.0251 rad/s and 0.0125 rad/s). During the dive, the UAV's longitudinal motion range is larger, and the degree of nonlinearity is greater than that of the lateral direction, which makes the longitudinal identification accuracy slightly worse than that of the lateral direction.

Figure 12 presents the validation data of the diving model. The established model well predicted the change of , with a mean validation error of 0.0356 rad/s; by contrast, the *p* was poorly validated, with a mean validation error of 0.0481 rad/s. The prediction results are better than the conventional identification method (the average error is 0.1054 rad/s and 0.0902 rad/s).

#### 7. Conclusions

Based on the flight test of the SUAV, the mathematical model of the dive phase was identified in this paper. A topological mathematical model of multi-input and multioutput noise terms for SUAV was constructed. Moreover, an improved multivariable augmented recursive least squares method was proposed. The above method was applied and analyzed via the identification test of a small UAV dive model. By comparing the results with the least squares method, it was shown that the proposed mathematical model and identification method can effectively obtain the parameters of the drone dive model. Furthermore, they can achieve high-precision identification of the drone dive model. This method is suitable for building nonmechanical modeling of the system when experimental data are available. In addition, the proposed method does not support the change mechanism of the internal physical parameters of the model. Future research work will focus on the application of tilt-rotor, trans-medium aircraft, and other new configuration aircraft tests. In addition, it will generate relevant test data, study the test modeling technology of special flight phases, and expand the scope of application of the proposed method. In terms of algorithm research, the results of more identification methods will be compared to continuously improve the identification performance of the model algorithm.

#### Data Availability

Some or all data, models, or codes generated or used during the study are proprietary or confidential in nature and may only be provided with restrictions (e.g., anonymized data).

#### Conflicts of Interest

The authors declare that they have no conflicts of interest.