Engineering Comprehensive Model of Complex Wind Fields for Flight Simulation

A complex wind field refers to the typical atmospheric disturbance phenomena existing in nature that have a great influence on the flight of aircrafts. Aimed at the issues involving large volume of data, complex computations and a single model in the current wind field simulation approaches for flight environments, based on the essential principles of fluid mechanics, in this paper, wind field models for two kinds of wind shear such as micro-downburst and low-level jet plus three-dimensional atmospheric turbulence are established. The validity of the models is verified by comparing the simulation results from existing wind field models and the measured data. Based on the principle of vector superposition, three wind field models are combined in the ground coordinate system, and a comprehensive model of complex wind fields is established with spatial location as the input and wind velocity as the output. The model is applied to the simulated flight of a rocket projectile, and the change in the rocket projectile’s flight attitude and flight trajectory under different wind fields is analyzed. The results indicate that the comprehensive model established herein can reasonably and efficiently reflect the influence of various complex wind field environments on the flight process of aircrafts, and that the model is simple, extensible, and convenient to use.


Introduction
A complex wind field refers to the typical atmospheric disturbance phenomena in nature, for instance gusts, typhoons, wind shear, turbulence and so on, which are closely related to factors such as weather and terrain. In flying aircrafts, rockets, missiles, and airships, etc., the whole flight process is affected by complex wind fields. For example, the flight stability and safety of aircrafts will reduce, and the firing accuracy of rockets and missiles will become worse. Therefore, establishing a computational model of complex wind fields suitable for flight simulation is of great significance to the design of aircraft and the study of flight control.
At present, there are three types of modeling method for complex wind fields. The first method is the wind field data measurement method. The data obtained from this method are authentic and reliable, but the volume of data is large and the computation cost is high. [1,2], respectively, substituted the exploratory data of 14 Atlantic tropical storms from 1982 to 1989 collected by the Hurricane Research Division (HRD) and the wind vector data from the HY-2 satellite into the wind field data model of the National Center for Environmental Prediction (NCEP) for calculation, analysis and verification. [3] used particle filter and Monte Carlo methods to process the measured data and set up the wind velocity model. [4] reanalyzed the atmospheric circulation data taken from the Black Sea region during the period of 1958 to 2001 and established the distribution of circulation wind field affected by sea temperature. The second method is the numerical method of atmospheric dynamics, which needs to solve the computation-heavy nonlinear differential equations, due to which the methodology is too complicated. [5] used the mesoscale model (MM5) of the National Center for Atmospheric Research of Pennsylvania State University to calculate

Micro-Downburst
A micro-downburst (abbreviated as MD in this section) is a type of low-level wind shear associated with convective weather [17,18]. Figure 1 shows a schematic diagram of airflow distribution in the MD wind field. It can be noticed that the MD manifests itself as a local vertically downward airflow in the strong convective cloud cluster. After the airflow sinks and touches the ground, it diverges in all directions and curls up to form an area of vortex ring above the ground. In view of such flow characteristics, a bunch of ground-symmetric vortex rings were constructed in the vertical direction of the horizontal plane to simulate the vertical airflow generation [19]. The coordinate system of the model is illustrated in Figure 2.
Considering the ground coordinate system oxyz as the datum and a point O P above the plane xoz as the center, a closed vortex ring of radius R is set, named the main vortex ring, and its curve equation is: where  denotes the intensity of the vortex ring, determined by vertical velocity  Considering the ground coordinate system oxyz as the datum and a point P O above the plane xoz as the center, a closed vortex ring of radius R is set, named the main vortex ring, and its curve equation is: The circulation line equation of the main vortex is: where  denotes the intensity of the vortex ring, determined by vertical velocity where Γ denotes the intensity of the vortex ring, determined by vertical velocity V z (0) of the vortex ring center and the vortex ring radius R: r max and r min are the maximum and the minimum distances of the main vortex ring from an arbitrary point O M (x M , y M , z M ) in the ground coordinate system, and F(k) is the elliptic integral function, where: Based on the theories of higher functions [7], when 0 ≤ k ≤ 1, F(k) is approximated as: Based on Equation (2), the induced radial and axial velocities of the main vortex ring can be calculated by partial derivatives: where r P is distance of the point O M from the axis of vortex ring (r P = (x M − x P ) 2 + (y M − y P ) 2 ). By decomposing the induced velocities calculated by Equation (6) along the ox and oz axes of the ground system, the velocity components can be written as: In the actual circumstances, the vertical velocity component on the ground should be zero after the vertical airflow of the vortex ring center reaches the ground. For that reason, by setting the method of mirror vortex ring for the symmetry of main vortex ring about the plane xoy, the vertical induced velocities on the ground are reversed with equal value and thus cancel each other out. The signs of the streamline equations for the two vortex rings are opposite ψ I = −ψ P . Given the mirrored vortex ring center O I (x p , y p , −z p ), its connection line O P O I with the main vortex ring center O P is vertical to the ground. According to the streamline equation of the mirrored vortex ring, analogous to the derivation of the induced velocity of the main vortex ring, the induced velocities at the spatial points of the mirrored vortex ring can be calculated. Velocity superposition is performed by combining Equations (6) and (7), and the resultant velocity at point O M can be obtained as: Then, the streamline equation of point O M is: where r max and r min denote the maximum and the minimum distances of the mirrored vortex ring from the spatial point O M , and similarly k = r max −r min r max +r min .
The main parameters [20] of the vortex ring model are shown in Table 1. Table 1. Parameters of micro-downburst model.

Parameter Value
The simulation results of the model are shown in Figure 3. It can be noted that the simulated wind field distribution of the vortex ring model is consistent with the measured data. The simulation results of the model are shown in Figure 3. It can be noted that the simulated wind field distribution of the vortex ring model is consistent with the measured data.

Low-Level Jet
Low-level jet (abbreviated as LLJ in this section) refers to the wind velocity zone in the lower troposphere, which is significantly affected by the mesoscale weather system. It is a surface inversion phenomenon occurring in the stable surface boundary layer [23]. According to the principle of the plane wall jet, when a plane free jet with a large width flows through a narrow slit, the velocity distribution of the jet near the wall is similar to that of the LLJ [24]. For that reason, the plane wall jet is employed to simulate the wind velocity distribution of the LLJ wind field, as shown in Figure 4.
According to the principles of fluid dynamics, the relationship between the horizontal velocity component u(x, H) of a free jet and its maximum jet velocity u m (x) is: Low-level jet (abbreviated as LLJ in this section) refers to the wind velocity zone in the lower troposphere, which is significantly affected by the mesoscale weather system. It is a surface inversion phenomenon occurring in the stable surface boundary layer [23]. According to the principle of the plane wall jet, when a plane free jet with a large width flows through a narrow slit, the velocity distribution of the jet near the wall is similar to that of the LLJ [24]. For that reason, the plane wall jet is employed to simulate the wind velocity distribution of the LLJ wind field, as shown in Figure 4.  According to the principles of fluid dynamics, the relationship between the horizontal velocity component   where k 1 denotes a scale factor and th refers to hyperbolic tangent function. To simplify the model, the basic flow characteristics of LLJ are retained, while the wind velocity is evenly distributed in the horizontal direction x, that is, where H s denotes the height of maximum velocity of the free jet with symmetrical distribution, C s is employed to describe the relationship between H s and the width B of the jet in the vertical direction; and the width of jet denoting the velocity range determined by 7% of the maximum velocity. The expression is By superimposing the velocity distribution of the free jet onto the exponential model of the mean wind in the boundary layer, the vertical velocity profile of LLJ on the surface boundary layer can be obtained: In the formula above, u R is the wind velocity corresponding to the height of H R , where the distribution index m p can be calculated in keeping with the relevant meteorological data [25]: where Z 0 denotes terrain roughness; analogous to Equation (13), the deviation of wind direction between height H and height H 0 is derived as: where H 0 , H L and H T denote reference height, height of maximum deviation of wind direction in the jet stream layer and height of top of the jet stream layer, respectively; α H 0 , α H L and α H T are the angles included between the wind direction and the geostrophic wind at the three corresponding heights, respectively; w 0 and w L are wind velocities at heights H 0 and H L , respectively. Thus, the wind velocity components in the ground coordinate system are calculated as: The model parameters [26] are configured as shown in Table 2 and the simulation results are shown in Figure 5.  It can be seen from Figure 5a,b that the wall jet model can accurately simulate the wind velocity distribution of the actual LLJ wind field.

Atmospheric Turbulence
Atmospheric turbulence is the irregular and uneven random eddy motion of the Earth's atmosphere. The Dryden model is a classical form of atmospheric turbulence, wherein the longitudinal and transverse correlation functions [28] are expressed as: The spatial correlation function [29] of one-dimensional turbulence is: For two (or three)-dimensional turbulence, the spatial correlation function is: It can be seen from Figure 5a,b that the wall jet model can accurately simulate the wind velocity distribution of the actual LLJ wind field.

Atmospheric Turbulence
Atmospheric turbulence is the irregular and uneven random eddy motion of the Earth's atmosphere. The Dryden model is a classical form of atmospheric turbulence, wherein the longitudinal and transverse correlation functions [28] are expressed as: The spatial correlation function [29] of one-dimensional turbulence is: For two (or three)-dimensional turbulence, the spatial correlation function is: where ξ denotes the location difference between the two spatial points, ξ 1 , ξ 2 , ξ 3 are the components of ξ in the ground coordinate system, σ is the turbulence intensity parameter, L is the turbulence scale parameter, and the subscripts u, v, w denote the components of turbulent velocity in the ground coordinate system. For a one-dimensional turbulent sequence, the recursion equation is: where W(x) denotes the turbulent velocity along the one-dimensional x direction, r(x) is the one-dimensional white noise sequence, h is the simulation step size, and a, b are the recursion parameters to be solved. According to the correlation function definitions of the Dryden model, Combining Equations (18) and (21), Thus, the one-dimensional atmospheric turbulent velocity can be calculated. Similarly to Equation (20), the recursion equation of two-dimensional turbulence can be written as: The correlation function values are The equation obtained by expanding Equation (24) is Substitute R 00 , R 01 , R 10 , R 11 calculated by Equation (19) into Equation (25), where R ij = R u (i * h, j * h), so that the parameters a 1 , a 2 , a 3 , b can be obtained. Thus, the twodimensional atmospheric turbulent velocity can be calculated according to the recursive Equation (23).
For the three-dimensional space turbulence, the construction idea is as follows: taking one-dimensional turbulence as the boundary value to figure out the two-dimensional plane turbulence sequence, and subsequently taking the two-dimensional plane turbulence sequence as the boundary value to figure out the three-dimensional space turbulence, as shown in Figure 6.  Analogous to Equations (20) and (23), the three-dimensional turbulence recursion equation is formulated as w

x y z a w x h y h z h a w x h y h z a w x y h z h a w x h y z h a w x y z h a w x h y h z a w x h y z br x y z
The correlation function of the three-dimensional space turbulence can be expressed as Substituting Equation (19) into Equation (27) By solving Equation (28), the recursive parameters  Analogous to Equations (20) and (23), the three-dimensional turbulence recursion equation is formulated as The correlation function of the three-dimensional space turbulence can be expressed as Substituting Equation (19) into Equation (27), where R ijk = R u (i * h, j * h, k * h), the following equations are obtained R 000 = a 1 R 111 + a 2 R 110 + a 3 R 011 + a 4 R 101 + a 5 R 001 + a 6 R 010 + a 7 R 100 + σ 2 R 001 = a 1 R 110 + a 2 R 111 + a 3 R 010 + a 4 R 100 + a 5 R 000 + a 6 R 011 + a 7 R 101 By solving Equation (28), the recursive parameters a 1 ∼ a 7 and b can be obtained, and then the velocity recursion equation of three-dimensional space turbulence w u (x, y, z) is obtained. By using the method above, w v (x, y, z) and w w (x, y, z) can also be obtained. AT refers to the atmospheric turbulence herein.
According to the parameters in Table 3, the three-dimensional atmospheric turbulence was simulated in space. The simulation results are shown in Figures 7 and 8. From Figure 8, the random variation of the atmospheric turbulent wind velocity in space is observed. As presented in Figure 8, the theoretical value is calculated according to correlation function (Equation (19)) of the Dryden model. It can be seen that the variation trend in the correlation value of turbulent velocity simulated by simulation is actually the same as that obtained from the Dryden model, which proves the rationality and effectiveness of the established model. Table 3. Parameters of atmospheric turbulence model [30].

Parameter
Value According to the parameters in Table 3, the three-dimensional atmospheric turbulence was simulated in space. The simulation results are shown in Figures 7 and 8. From Figure 8, the random variation of the atmospheric turbulent wind velocity in space is observed. As presented in Figure 8, the theoretical value is calculated according to correlation function (Equation (19)) of the Dryden model. It can be seen that the variation trend in the correlation value of turbulent velocity simulated by simulation is actually the same as that obtained from the Dryden model, which proves the rationality and effectiveness of the established model. Table 3. Parameters of atmospheric turbulence model [30].

Comprehensive Wind Field Model
Since wind velocity is a vector including magnitude and direction, it meets the principle of vector superposition. For any spatial point P , assuming that wind field A induces the wind velocity vector

Comprehensive Wind Field Model
Since wind velocity is a vector including magnitude and direction, it meets the principle of vector superposition. For any spatial point P, assuming that wind field A induces the wind velocity vector → W P A at point P and wind field B induces the wind velocity vector → W P B at point P, the total wind velocity vector → W P at point P induced by the wind fields A and B can be obtained by the following equation: Based on Equation (30), in this section, a superposition method is adopted to synthesize the established typical wind field models. The input and output parameters of the three typical wind field models are listed in Table 4. Table 4. Input and output parameters of the three typical wind field models established.

Wind Field
Input Output It can be seen that the inputs of the three models are the spatial positions, and the outputs are the wind velocity components. Therefore, the inputs and outputs are combined, and the ground coordinate system is used as the frame of reference to establish a comprehensive model. The application structure of this model is shown in Figure 9. As can be seen from Figure 9, the calculation process of the comprehensive wind field model consists of four key steps: • Input spatial location parameters; • Select the wind field models according to simulation requirements; • Calculate the total wind velocity value; • Substitute the wind velocity into the flight simulation.
To ensure the scalability of the model, an expansion module (yellow box in the chart) is set up to add other wind field models such as gust and mountain flow.

Flight Simulations under Different Wind Field Conditions
In order to verify the calculation effect on the whole model, a rocket projectile is As can be seen from Figure 9, the calculation process of the comprehensive wind field model consists of four key steps: • Input spatial location parameters; • Select the wind field models according to simulation requirements; • Calculate the total wind velocity value; • Substitute the wind velocity into the flight simulation.
To ensure the scalability of the model, an expansion module (yellow box in the chart) is set up to add other wind field models such as gust and mountain flow.

Flight Simulations under Different Wind Field Conditions
In order to verify the calculation effect on the whole model, a rocket projectile is taken as an example and the six degrees-of-freedom rigid body trajectory equation of the rocket projectile is considered as the flight simulation model. The basic parameters [31] of the rocket projectile are listed in Table 5. For the flight state of a rocket projectile, the main concern is its trajectory and flight attitude. The flight trajectory can be directly analyzed by calculating the three-dimensional trajectory curve of the rocket projectile, while the flight attitude needs to be reflected by the flight attack angle. Figure 10 shows the conceptual schematic diagram of the attack angle of rocket, where the coordinate system Oξηζ is taken as the reference coordinate system, O is the center of rocket mass, the Oξ axis coincides with the axis of the rocket, the Oη axis points vertically upward, and the Oζ axis is determined according to the right hand rule. The red vector denotes the velocity v of the rocket's centroid, the blue vector represents the three components v ξ , v η and v ζ of v in the coordinate system Oξηζ, δ represents the angle included between the rocket axis Oξ and the velocity v, and is called the total attack angle, δ 1 is called the pitch attack angle, and δ 2 is called the direction attack angle. The attack angles δ, δ 1 and δ 2 of the rocket describe the positional relationship between the projectile axis and the velocity direction during the flight of the rocket. Through the curve of the attack angle, the attitude changes and the stability of the rocket during the flight can be seen. The specific flight simulation calculation steps are shown in Figure 11.  Figure 11.  While using the integral method to solve the trajectory, the coordinates of the rocket's position within each calculation step are substituted into the wind field model to obtain the wind velocity, which is then substituted into the trajectory equations to calculate the next trajectory parameters. These trajectory parameters include the position coordinates (x, y, z) in addition to attack angles δ 1 , δ 2 of the rocket projectile. The detailed solution of the trajectory equations can be found in the literature [31] (pp. 141-143). Table 6 lists two wind field conditions employed for flight simulation of a rocket projectile. The flight processes under the two wind field conditions are simulated, and are compared with the flight state of the rocket projectile under windless condition.  Figure 12 shows the change curves of the attack angles of the rocket projectile with flight time under the influence of the two wind fields listed in Table 6. The one on the lefthand side of the figure is the attack angle of the whole trajectory, while the small one on the right-hand side is the local attack angle of the first 5 s of the trajectory. It can be noticed that in a windless environment, the amplitude of attack angle of rocket projectile persistently decreases with the increase in the flight time and converges to nearly 0 degrees, indicating that the flight attitude of the rocket projectile gradually tends to become stable. Under the influence of the low-level jet in condition 1, the rocket projectile has a low velocity and weak anti-interference capability in the initial stages of launch. So, the amplitude of attack angle increases suddenly under the action of transient but strong airflow. As the rocket projectile continues to fly, the attack angle converges in a similar way as it would in a windless condition. Under the influence of the MD in condition 2, the attack angle of the rocket projectile also increases suddenly, similar to that observed in condition 1. However, due to the influence of the atmospheric turbulence at the same time, the amplitude of attack angle cannot converge further, fluctuates continuously, and the flight attitude is not stable anymore. Figure 13 shows the three-dimensional curve of flight trajectory of rocket projectile under wind field condition 1, wind field condition 2 and the windless environment. As can be seen from Figure 13, under the influence of conditions 1 and 2, the maximum trajectory height of the rocket projectile decreases, the range decreases, and the whole flight trajectory deviates significantly. In particular, the lateral deviation of the landing point under the influence of condition 1 increases much more than that under condition 2, which reflects the difference in the influence of the two wind conditions.    Figure 14 shows the dispersion of the rocket projectile under the influence of condition 2 using the Monte-Carlo method.
According to the flight simulation results shown in Figures 12-14, it can be concluded that the comprehensive model of complex wind field established herein, can reflect the law of influence of different wind fields on the aircraft when applied in flight simulation, which reflects the practicability and rationality of the model. According to the flight simulation results shown in Figures 12-14, it can be concluded that the comprehensive model of complex wind field established herein, can reflect the law of influence of different wind fields on the aircraft when applied in flight simulation, which reflects the practicability and rationality of the model.

Analysis on Influence of Different Model Parameters on Flight Process
For a certain kind of typical wind field, it is necessary to study the effect of different model parameters on the flight process of aircrafts. Taking the low-level jet wind field(condition 1 in   Tables 7 and 8, it can be seen that the ballistic characteristics changed under the influence of low-level jet. The changes are specifically manifested in the reduction in flight time and terminal velocity, the decrease in down range and the increase in cross range (the order of magnitude from 10 to 10 3 ). The simulation results also show that the influence degree of low-level jet on the flight process is positively correlated to the strength parameter w L and the scale parameter H T .
From the attack angle curves in Figure 15, the attitude changes of the rocket projectile under different wind field model parameters can be directly observed. It is shown that the low-level jet causes a sudden change and a convergence process of the attack angle of the attack angle δ 2 , which result in the large change in the cross range (in Tables 7 and 8) of the rocket.  Tables 7 and 8, it can be seen that the ballistic characteristics changed under the influence of low-level jet. The changes are specifically manifested in the reduction in flight time and terminal velocity, the decrease in down range and the increase in cross range (the order of magnitude from 10 to 10 3 ). The simulation results also show that the influence degree of low-level jet on the flight process is positively correlated to the strength parameter

Discussions of Use Conditions of the Comprehensive Model
In fact, the real wind field in nature changes not only in space but also in time, which can significantly influence the accuracy of flight simulation. It is necessary to give some conditions of use for the established comprehensive model. According to the atmospheric dynamics [25], the wind field disturbance can be expressed as: Equation (31) indicates that the total disturbance W is composed of the mean wind w and the stochastic term ∆w. In the conditions of small time and spatial scales, w can be described as a certain wind field, such as micro-downburst, and the wind profile, ∆w, can be described as the atmospheric turbulence. The diagram is shown in Figure 16.
In fact, the real wind field in nature changes not only in space but also in time, which can significantly influence the accuracy of flight simulation. It is necessary to give some conditions of use for the established comprehensive model. According to the atmospheric dynamics [25], the wind field disturbance can be expressed as: Equation (31) indicates that the total disturbance W is composed of the mean wind w and the stochastic term w  . In the conditions of small time and spatial scales, w can be described as a certain wind field, such as micro-downburst, and the wind profile, w  , can be described as the atmospheric turbulence. The diagram is shown in Figure 16.

Certain wind field
Wind profiles Atmospheric turbulence

Stochastic term
Total wind disturbance + Figure 16. Composition of the wind disturbance in small time and spatial scales.
Based on Figure 16, the conditions of use for the established comprehensive model are given as: • When the flight simulation is in a small time scale (in condition 1 of Section 4.1, the rocket went through the low-level jet area within 2 s) or in a small spatial scale (the low-level jet area has a height of 800 m and the rocket has a flight altitude of 10 km), the established comprehensive model can be used to obtain some reasonable results.

•
When the flight simulation is in a large time scale (or in a large spatial scale), such as the total flight process of a long-range missile, an airplane or an airship, the comprehensive model might cause significant errors.

Conclusions
Built on the basic principles of fluid mechanics, the engineering wind velocity calculation models of three typical wind fields, namely micro-downburst, low-level jet and atmospheric turbulence, are herein established. By combining the three models in the ground coordinate system and unifying the input and output parameters, a comprehensive model of complex wind field is established. The wind field model Based on Figure 16, the conditions of use for the established comprehensive model are given as:

•
When the flight simulation is in a small time scale (in condition 1 of Section 4.1, the rocket went through the low-level jet area within 2 s) or in a small spatial scale (the low-level jet area has a height of 800 m and the rocket has a flight altitude of 10 km), the established comprehensive model can be used to obtain some reasonable results.

•
When the flight simulation is in a large time scale (or in a large spatial scale), such as the total flight process of a long-range missile, an airplane or an airship, the comprehensive model might cause significant errors.

Conclusions
Built on the basic principles of fluid mechanics, the engineering wind velocity calculation models of three typical wind fields, namely micro-downburst, low-level jet and atmospheric turbulence, are herein established. By combining the three models in the ground coordinate system and unifying the input and output parameters, a comprehensive model of complex wind field is established. The wind field model simulation and application simulation indicate that the comprehensive model can reasonably and effectively describe the flow characteristics of the relevant wind field, and has the characteristics of simple calculation and model scalability. Additionally, as shown in Section 4.2, the comprehensive model can be used not only to analyze the influence of different wind fields on the flight process, but also to study the influence of different model parameters on the flight simulations.
As an important factor of the simulation accuracy of the wind field model, the wind field data collection methods in different time and space scales of the simulated flights are the main points of our further research. Meanwhile, developing a reasonable wind field model for large time and space scales with the least amount of collection data is also an extension of the work in this paper.