A wind field data modeling method and system for wind resistance testing of small and medium-sized UAVs

Abstract

This invention discloses a wind field data modeling method and system for wind resistance testing of small and medium-sized unmanned aerial vehicles (UAVs). The method includes collecting terrain data, wind speed data, and wind turbine data of the test wind field; determining a preset area and an unmeasured area based on the test wind field and a preset UAV trajectory; obtaining the boundary conditions of the unmeasured area and the wake of the perturbation grid through fluid simulation based on the terrain data and the preset area; performing wind speed vector decomposition on the grid according to the preset wind field boundary conditions; obtaining a preset wind field model based on the perturbation time series, UAV maneuver weights, and the wake of the perturbation grid; obtaining a predicted wind field for the unmeasured area through Gaussian process regression based on the preset wind field model and the boundary conditions of the unmeasured area; and predicting the wind field using different preset UAV trajectories through grid temporal coupling to obtain a target wind field model. This method, through spatiotemporal co-interpolation of Gaussian process regression and fluid simulation, makes the inlet boundary conditions more consistent with fluid dynamics, improving the physical consistency and reliability of the target wind field.

Description

Technical Field

[0001] This invention relates to the field of wind field modeling technology, and in particular to a wind field data modeling method and system for wind resistance testing of small and medium-sized UAVs. Background Technology

[0002] Small and medium-sized UAVs have been widely used in aerial surveying, power line inspection, agricultural plant protection, and emergency rescue due to their advantages of high flexibility, low cost, and convenient deployment. As application scenarios extend to complex meteorological environments, their wind resistance performance has become a core indicator that determines the safety and reliability of missions. Wind resistance testing is a key link in UAV research and development and verification, which requires accurate wind field models to simulate real wind environments. Currently, wind field modeling mainly adopts a combination of computational fluid dynamics (CFD) simulation, wind tunnel testing, and on-site wind speed collection. By constructing a spatiotemporal distribution model of wind speed and direction, input conditions are provided for testing the wind resistance performance of UAVs.

[0003] Existing technologies have limitations in wind field modeling. On the one hand, although CFD simulation can generate three-dimensional wind fields, the deviation from the actual flight environment is large, making it difficult to reproduce dynamic wind characteristics such as turbulence and gusts. Limited by physical space, it cannot reproduce large-scale, multi-scenario wind field changes. On-site data collection is easily affected by sensor accuracy and environmental interference. Moreover, wind field prediction in unmeasured areas often relies on simple interpolation, lacking the ability to quantify the perturbation effect of wind fields on different trajectories under the requirements of wind resistance tests involving multiple UAVs and multiple trajectories. Therefore, how to perform CFD simulation of wind fields simply and efficiently and construct accurate local wind field data has become an urgent problem to be solved. Summary of the Invention

[0004] This invention provides a wind field data modeling method and system for wind resistance testing of small and medium-sized unmanned aerial vehicles (UAVs), aiming to solve at least one of the above-mentioned technical problems.

[0005] Firstly, to achieve the above objectives, the present invention provides a wind field data modeling method for wind resistance testing of small and medium-sized unmanned aerial vehicles (UAVs), comprising:

[0006] Collect terrain data, wind speed data, and wind turbine data of the test wind field, and determine the preset area and unmeasured area based on the test wind field according to the preset trajectory of the UAV;

[0007] Based on the terrain data, the preset area is divided into grids, and the wind speed data, wind turbine data and grid coordinates are aligned to obtain the preset wind field boundary conditions. Based on the preset wind field boundary conditions, the boundary conditions of the unmeasured area and the wake of the disturbed grid are obtained by CFD fluid simulation.

[0008] Based on the preset wind field boundary conditions, the grid is decomposed into wind speed vectors. Based on the wind speed vectors and UAV flight deviations of different frequencies, anti-outlier filtering and phase space reconstruction are performed to obtain disturbance displacement and UAV maneuver weights.

[0009] The disturbance displacement is analyzed over time based on the grid coordinates, and a preset wind field model is obtained based on the disturbance time sequence, UAV maneuver weights, and the wake of the disturbance grid.

[0010] Based on the preset wind field model, the predicted wind field of the undetermined area is obtained through Gaussian process regression according to the boundary conditions of the undetermined area. The wind field is predicted by using different preset UAV trajectories through grid temporal coupling to obtain the target wind field model.

[0011] Further, the method for determining the preset region and the unmeasured region includes:

[0012] Collect terrain data, wind speed data, and wind turbine data of the test wind field. Based on the trajectory coordinates of the UAV's preset trajectory, calculate the minimum horizontal enclosing rectangle and the maximum vertical height difference according to the flight cycle to obtain the trajectory's three-dimensional extreme value envelope. Based on the trajectory's three-dimensional extreme value envelope, fit a preset safety distance according to the UAV's preset wind resistance level and flight speed. If the UAV's preset wind resistance level increases by one level, the safety distance is expanded by 15%. If the UAV's speed doubles per unit time, the safety distance is expanded by half to obtain the preset area.

[0013] The area within one kilometer of the preset zone is designated as the preset buffer zone. The slope change rate and obstacle density within the preset buffer zone are statistically analyzed based on terrain data. The terrain complexity is obtained by weighted summation of the slope change rate and obstacle density. The terrain complexity is then multiplied by the preset buffer zone as an exponential decay coefficient to obtain the actual buffer zone. The test wind field outside the preset zone and the actual buffer zone is designated as the unmeasured zone.

[0014] Furthermore, the method for obtaining the boundary conditions of the unmeasured region and the wake of the perturbed grid includes:

[0015] Based on the terrain data, the preset area is divided into a structured grid. The grid resolution is adjusted according to the spatial density of the trajectory coordinates. For areas with abrupt changes in terrain, the grid is densified by multiples based on the slope difference. The sampling point coordinates of the wind speed data and the installation coordinates of the wind turbine are converted to the same spatial coordinate system as the grid.

[0016] Based on the grid coordinates, wind speed data and wind turbine data are distributed to each grid cell through a spatial interpolation algorithm to form preset wind field boundary conditions. Based on the preset wind field boundary conditions, a computational domain containing preset areas and unmeasured areas is constructed, and inlet boundaries, outlet boundaries, and wall boundaries are set.

[0017] Based on the computational domain, the fluid motion control equations are solved by vortex simulation to obtain the boundary conditions of the unmeasured area and extract the perturbation grid wake. The perturbation grid wake includes the wake diffusion range, wind speed attenuation law, gust duration, turbulence intensity distribution characteristics and turbulence time scale.

[0018] Furthermore, the method for obtaining the disturbance displacement and UAV maneuver weights is as follows:

[0019] The wind speed vector of each grid is obtained based on the preset wind field boundary conditions. The wind speed vector is decomposed into an axial component along the preset trajectory of the UAV and a radial component perpendicular to the trajectory. The initial state of position, velocity and attitude of the UAV in the preset trajectory is used as the initial value of integration based on the axial component and the radial component. The theoretical deviation of wind disturbance of the preset trajectory under the action of wind field alone is obtained through the preset aerodynamic model based on the initial value of integration.

[0020] Using the initial state as the initial value for simulation, the drone dynamics equations are used to perform windless simulation based on the flight control commands of the preset trajectory to obtain the disturbance displacement components under the individual action of each flight control command. The actual flight deviation is obtained based on the actual flight trajectory of the drone and the preset trajectory of the drone. The abnormal points of sensor noise and electromagnetic interference are removed from the actual flight deviation by anti-outlier filtering. The Pearson correlation coefficient between the filtered deviation and the disturbance displacement components is calculated. If the correlation coefficient is greater than 0.7, the disturbance displacement components are deducted to obtain the actual wind disturbance deviation. If the correlation coefficient is less than 0.3, the filtered deviation is taken as the actual wind disturbance deviation.

[0021] The embedding dimension is determined by the spurious nearest neighbor method based on the actual wind disturbance deviation, and the delay time is determined by the mutual information method. The phase space of the actual wind disturbance deviation is reconstructed based on the embedding dimension and the delay time, and the one-dimensional time series actual wind disturbance deviation is mapped to a high-dimensional disturbance displacement vector. The L2 norm of the high-dimensional disturbance displacement vector and the theoretical wind disturbance deviation is used as the disturbance displacement, and the ratio of the theoretical wind disturbance deviation to the disturbance displacement is used as the UAV maneuver weight.

[0022] Furthermore, the method for obtaining the preset wind field model includes:

[0023] The disturbance displacement is extracted based on grid coordinates and flight data sampling time step to obtain a three-dimensional matrix. The duration of gusts is used as the short-term window length, and the long-term window length is determined according to the turbulence time scale. Based on the three-dimensional matrix, the short-term disturbance standard deviation, maximum modulus, minimum modulus and long-term disturbance mean are extracted by sliding between the short-term window length and the long-term window length. A four-dimensional feature matrix is ​​obtained based on the grid coordinates.

[0024] The disturbance displacement is corrected coordinate-by-coordinate based on the three-dimensional matrix. If the UAV maneuver weight of the grid is greater than 0.7, the product of the disturbance displacement and the UAV maneuver inverse weight is used as the correction value. If the UAV maneuver weight is less than 0.3, it remains unchanged. If the network is located within the wake diffusion range, the wind speed attenuation law value corresponding to the grid coordinates is used as the weight coefficient of the disturbance displacement to obtain the corrected three-dimensional matrix. The corrected three-dimensional matrix is ​​then used to perform inversion through a preset aerodynamic model to obtain the inverted wind speed and inverted wind direction. The dynamic characteristics of the wind field are calculated based on the inverted wind speed, inverted wind direction, and four-dimensional feature matrix. The formula for calculating the dynamic characteristics of the wind field is as follows:

[0025]

[0026] Where γ(t) represents the wind field dynamics at time t, s is the short-term window length, Δt is the flight data sampling time step, h is the time index, and d(h) is the disturbance displacement. Let denoted as the mean displacement of the disturbance within the short-term window, tanh be the hyperbolic tangent function, f be the length of the long-term window, and d be the mean displacement of the disturbance within the short-term window. max d represents the maximum disturbance displacement within the short-term window. min This represents the minimum disturbance displacement within a short-term window.

[0027] Based on the inverted wind speed, inverted wind direction, and dynamic characteristics of the wind field, Kriging interpolation is performed on the preset area and the unmeasured area to obtain the preset wind field model.

[0028] Furthermore, the method for obtaining the predicted wind field includes:

[0029] Based on the boundary conditions of the unmeasured area, the wake diffusion boundary and wind speed attenuation law are extracted. According to the wake diffusion boundary and wind speed attenuation law, the grid data of the unmeasured area is corrected using a preset wind field model to obtain the prior mean of the prediction points. Using the prior mean of the prediction points in the unmeasured area and the grid data of the preset area, the prediction points in the unmeasured area are calculated through Gaussian process regression. The Gaussian process regression formula is as follows:

[0030]

[0031] in For the predicted points in the undetermined area, μ(x) u ,t u ) represents the prior mean of the predicted points, N represents the total number of spatial grids in the preset area, M represents the time step of the preset area, and x represents the value of x. i Let t be the spatial coordinates of the i-th grid in the preset area. k The timestamp of the k-th time step in the preset area. The inversion wind speed is the i-th grid and the k-th time step in the preset area. x is the average wind speed u The coordinates of the predicted points in the undetermined area, t uFor the predicted timestamp, x u The predicted spatial coordinates;

[0032] Align the set of predicted points in the unmeasured area with the grid and time step of the preset area to obtain the predicted wind field in the unmeasured area.

[0033] Furthermore, the method for obtaining the target wind field model includes:

[0034] Based on different drone preset trajectories, discrete preset trajectory points are obtained according to grid coordinates and flight data sampling time step. The required wind speed and required wind direction are obtained from the preset trajectory points according to the wind resistance test requirements.

[0035] Based on the grid coordinates of preset trajectory points, predicted wind speed and direction are obtained through predicted wind field. Using the predicted wind speed and direction of the preset trajectory points as initial target values, the target wind field loss function is iteratively optimized based on the required wind speed, required wind direction, and initial target values ​​to obtain the target wind field with minimized trajectory deviation. The expression for the target wind field loss function is:

[0036]

[0037] Where L is the loss function, representing the deviation between the predicted wind field and the multi-UAV trajectory, n is the number of discrete points on a single trajectory, and v p,i To determine the required wind speed at point i on the preset trajectory, v t,i To predict the actual wind speed at point i on the trajectory, θ p,i To determine the required wind direction at point i on the preset trajectory, θ t,i To predict the actual wind direction at point i on the trajectory, λ is the smoothing coefficient, and Δt is the time step of the flight data sampling.

[0038] Secondly, in order to achieve the above objectives, the present invention also provides a wind field data modeling system for wind resistance testing of small and medium-sized unmanned aerial vehicles (UAVs), comprising:

[0039] Wind field preset area acquisition module: used to collect terrain data, wind speed data and wind turbine data of the test wind field, and determine the preset area and unmeasured area based on the test wind field and the preset trajectory of the UAV;

[0040] Wind field boundary condition acquisition module: used to divide the preset area into grids based on the terrain data, align the wind speed data and wind turbine data with the grid coordinates, obtain the preset wind field boundary conditions, and obtain the boundary conditions of the unmeasured area and the wake of the disturbed grid through CFD fluid simulation based on the preset wind field boundary conditions;

[0041] Wind field disturbance acquisition module: used to perform wind speed vector decomposition on the grid according to the preset wind field boundary conditions, and perform anti-outlier filtering and phase space reconstruction based on the wind speed vector and UAV flight deviations of different frequencies to obtain disturbance displacement and UAV maneuver weights;

[0042] Preset wind field modeling module: used to perform time-series analysis on the disturbance displacement based on grid coordinates, and obtain a preset wind field model based on the disturbance time sequence, UAV maneuver weights and disturbance grid wake;

[0043] Target wind field modeling module: This module is used to obtain the predicted wind field of the undetermined area based on the boundary conditions of the undetermined area through Gaussian process regression according to the preset wind field model. It uses different preset UAV trajectories to predict the wind field through grid temporal coupling to obtain the target wind field model.

[0044] The beneficial effects of this invention are as follows:

[0045] This invention utilizes a spatiotemporal co-interpolation model based on spatial-temporal Gaussian process regression, combined with terrain data and anemometer data to calculate grid wind speed vectors. It fully leverages prior historical wind speed information and multi-temporal point measured data to improve the prediction accuracy of wind speed vectors in preset wind field areas. Simultaneously, based on UAV attitude data, it inverts wind field interference characteristics, capturing local wind field details around the UAV flight trajectory and enhancing the feature integrity of the three-dimensional wind field. By combining fluid simulation with filtering to eliminate noise and abrupt changes, it makes the inlet boundary conditions more consistent with fluid dynamics laws, thereby improving the physical consistency and reliability of the target wind field. Attached Figure Description

[0046] Figure 1 This is a flowchart illustrating the steps of a wind field data modeling method for wind resistance testing of small and medium-sized unmanned aerial vehicles (UAVs) according to the present invention. Detailed Implementation

[0047] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the invention.

[0048] Reference Figure 1 As shown, this invention provides a wind field data modeling method for wind resistance testing of small and medium-sized unmanned aerial vehicles (UAVs), including:

[0049] Collect terrain data, wind speed data, and wind turbine data of the test wind field, and determine the preset area and unmeasured area based on the test wind field according to the preset trajectory of the UAV;

[0050] Based on the terrain data, the preset area is divided into grids, and the wind speed data, wind turbine data and grid coordinates are aligned to obtain the preset wind field boundary conditions. Based on the preset wind field boundary conditions, the boundary conditions of the unmeasured area and the wake of the disturbed grid are obtained by CFD fluid simulation.

[0051] Based on the preset wind field boundary conditions, the grid is decomposed into wind speed vectors. Based on the wind speed vectors and UAV flight deviations of different frequencies, anti-outlier filtering and phase space reconstruction are performed to obtain disturbance displacement and UAV maneuver weights.

[0052] The disturbance displacement is analyzed over time based on the grid coordinates, and a preset wind field model is obtained based on the disturbance time sequence, UAV maneuver weights, and the wake of the disturbance grid.

[0053] Based on the preset wind field model, the predicted wind field of the undetermined area is obtained through Gaussian process regression according to the boundary conditions of the undetermined area. The wind field is predicted by using different preset UAV trajectories through grid temporal coupling to obtain the target wind field model.

[0054] In the actual evaluation, a test area was a 5km x 3km mountainous area, including valleys and steep slopes. The maximum takeoff weight of the small UAV was 3.5kg, and the designed wind resistance level was 12m / s. DEM data with a resolution of 1m was collected for the test area, along with wind speed data from 5 points around the test environment. Two small wind turbines were used in the area, and their start / stop status and wind speed intervention range were recorded. Three preset UAV trajectories were planned, each 5km long, with a flight altitude of 60m, a distance of 8m from the route, and a flight cycle of 10 minutes, including 3 gust simulations. Based on the three-dimensional extreme value envelope of the preset trajectories, a preset buffer zone of 3.45km x 0.58km x 30m was obtained. The terrain complexity within the preset buffer zone was calculated, with areas having a slope change rate greater than 15° accounting for 30%, and an obstacle density of 2 / km. 2 The terrain complexity was calculated to be 0.6 by weighted summation. 0.6 was used as the exponential decay coefficient and multiplied by the preset buffer. The preset buffer range was expanded by 1.2km to obtain the actual buffer. The test wind field outside the preset buffer and the actual buffer was regarded as the unmeasured area.

[0055] In the actual assessment, the preset area was divided into structured grids with a grid resolution of 5m squares. Grids with a slope difference greater than 20° were densified to 2m, resulting in a total of 1.2 million grids. The anemometer sampling points and wind turbine installation coordinates were transformed to the grid coordinate system. Wind speed data was distributed to the grids using Kriging interpolation to obtain preset wind field boundary conditions: an inlet wind speed of 8m / s, free flow at the outlet, and wind speed in the wind turbine area with superimposed intervention values. CFD fluid simulation was performed using these preset wind field boundary conditions to construct a computational domain containing the preset area and the unmeasured area. Inlet, outlet, and wall boundary conditions were set. The fluid motion control equations were solved using DES to obtain boundary conditions for the unmeasured area: a wake diffusion range of 200m, a wind speed attenuation of 30% after 100m from the wind turbine, a gust duration of 5 to 7s, and a turbulence intensity of 15% in areas with complex terrain. The wake diffusion range and wind speed attenuation law in the wake of the disturbed grid were also extracted.

[0056] In the actual evaluation, the grid wind speed vector was decomposed according to the preset wind field boundary conditions of the test area. Using the initial position (0,0,0), velocity (8 m / s), and horizontal attitude of the UAV's preset trajectory as the initial integration values, the theoretical deviation of wind disturbance under the sole action of the wind field was calculated using an aerodynamic model with preset lift and drag coefficients. The trajectory offset was ±0.5 m, and the wind direction deviation was ±5°. Using the initial state as the simulation initial value, the flight control commands were simulated using the UAV's dynamic equations to obtain the disturbance displacement components under the sole action of the flight control, with a trajectory offset of ±0.2 m. The actual flight deviation was filtered to remove anomalies caused by electromagnetic interference. For constant points, the Pearson correlation coefficient between the filtered deviation and the flight control disturbance displacement component is calculated. For the portion with a correlation coefficient greater than 0.7 in the gust section, the flight control disturbance displacement component is subtracted to obtain the actual wind disturbance deviation as a trajectory offset of ±0.8m. For the portion with a correlation coefficient less than 0.3 in the stable wind section, the trajectory offset of ±0.3m is directly used. For the actual wind disturbance deviation, the embedding dimension is determined to be 3 using the spurious nearest neighbor method, and the delay time is determined to be 0.2s using the mutual information method. Phase space reconstruction is performed on the actual wind disturbance deviation to obtain a disturbance displacement of 1.2 in the gust section and 0.5 in the stable section. The UAV maneuver weight is 0.7 in the gust section and 0.3 in the stable section.

[0057] In the actual assessment, based on the three-dimensional matrix of disturbance displacement, a short-term window length of 5s and a long-term window length of 20s were set. The standard deviation of the short-term disturbance and the mean of the long-term disturbance were extracted using a sliding motion to form a four-dimensional feature matrix. For the gust core area grid with a maneuver weight greater than 0.7, the displacement and inverse weight product was used for correction. For the wake diffusion range grid 200m downstream of the wind turbine, the wind speed attenuation law was used for correction. Through aerodynamic model inversion, the inverted wind speed of 12m / s in the gust segment and 8m / s in the stable segment were obtained. The inverted wind direction with an angle of 30° to the vertical component of the trajectory was also obtained. Combined with the four-dimensional feature matrix, the dynamic characteristics of the wind field were calculated, where γ(t) was 0.8 in the gust segment and 0.3 in the stable segment. The grid data of the unmeasured area were corrected using a preset wind field model, and a Gaussian process regression was employed. The predicted wind speed in the unmeasured area is calculated to obtain the predicted wind field for the unmeasured area. In areas with complex terrain, the wind speed is reduced by 15%, and the influence range of gusts is expanded by 10%. Three preset trajectories are discretized to obtain 6000 trajectory points with an accuracy of 10 min / 0.1 s. The required wind speed for the stable segment is extracted as 8 m / s ± 1 m / s, and for the gust segment as 12 m / s ± 2 m / s. The required wind direction is the angle between the wind speed and the trajectory is ±10°. The predicted wind speed for the stable segment is obtained as 7.8 m / s, and for the gust segment as 11.5 m / s, using the predicted wind field as the target initial value. The trajectory deviation is minimized by iterative optimization using a loss function with a smoothing coefficient of 0.5. The deviation between the target wind field and the predicted wind field is 0.6 for the gust segment and 0.3 for the stable segment. The target wind field is obtained based on the deviation.

[0058] In this embodiment, the method for determining the preset area and the unmeasured area includes:

[0059] Collect terrain data, wind speed data, and wind turbine data of the test wind field. Based on the trajectory coordinates of the UAV's preset trajectory, calculate the minimum horizontal enclosing rectangle and the maximum vertical height difference according to the flight cycle to obtain the trajectory's three-dimensional extreme value envelope. Based on the trajectory's three-dimensional extreme value envelope, fit a preset safety distance according to the UAV's preset wind resistance level and flight speed. If the UAV's preset wind resistance level increases by one level, the safety distance is expanded by 15%. If the UAV's speed doubles per unit time, the safety distance is expanded by half to obtain the preset area.

[0060] The area within one kilometer of the preset zone is designated as the preset buffer zone. The slope change rate and obstacle density within the preset buffer zone are statistically analyzed based on terrain data. The terrain complexity is obtained by weighted summation of the slope change rate and obstacle density. The terrain complexity is then multiplied by the preset buffer zone as an exponential decay coefficient to obtain the actual buffer zone. The test wind field outside the preset zone and the actual buffer zone is designated as the unmeasured zone.

[0061] In this embodiment, the method for obtaining the boundary conditions of the unmeasured region and the wake of the perturbed grid includes:

[0062] Based on the terrain data, the preset area is divided into a structured grid. The grid resolution is adjusted according to the spatial density of the trajectory coordinates. For areas with abrupt changes in terrain, the grid is densified by multiples based on the slope difference. The sampling point coordinates of the wind speed data and the installation coordinates of the wind turbine are converted to the same spatial coordinate system as the grid.

[0063] Based on the grid coordinates, wind speed data and wind turbine data are distributed to each grid cell through a spatial interpolation algorithm to form preset wind field boundary conditions. Based on the preset wind field boundary conditions, a computational domain containing preset areas and unmeasured areas is constructed, and inlet boundaries, outlet boundaries, and wall boundaries are set.

[0064] Based on the computational domain, the fluid motion control equations are solved by vortex simulation to obtain the boundary conditions of the unmeasured area and extract the perturbation grid wake. The perturbation grid wake includes the wake diffusion range, wind speed attenuation law, gust duration, turbulence intensity distribution characteristics and turbulence time scale.

[0065] In this embodiment, the method for obtaining the disturbance displacement and UAV maneuver weights is as follows:

[0066] The wind speed vector of each grid is obtained based on the preset wind field boundary conditions. The wind speed vector is decomposed into an axial component along the preset trajectory of the UAV and a radial component perpendicular to the trajectory. The initial state of position, velocity and attitude of the UAV in the preset trajectory is used as the initial value of integration based on the axial component and the radial component. The theoretical deviation of wind disturbance of the preset trajectory under the action of wind field alone is obtained through the preset aerodynamic model based on the initial value of integration.

[0067] Using the initial state as the initial value for simulation, the drone dynamics equations are used to perform windless simulation based on the flight control commands of the preset trajectory to obtain the disturbance displacement components under the individual action of each flight control command. The actual flight deviation is obtained based on the actual flight trajectory of the drone and the preset trajectory of the drone. The abnormal points of sensor noise and electromagnetic interference are removed from the actual flight deviation by anti-outlier filtering. The Pearson correlation coefficient between the filtered deviation and the disturbance displacement components is calculated. If the correlation coefficient is greater than 0.7, the disturbance displacement components are deducted to obtain the actual wind disturbance deviation. If the correlation coefficient is less than 0.3, the filtered deviation is taken as the actual wind disturbance deviation.

[0068] The embedding dimension is determined by the spurious nearest neighbor method based on the actual wind disturbance deviation, and the delay time is determined by the mutual information method. The phase space of the actual wind disturbance deviation is reconstructed based on the embedding dimension and the delay time, and the one-dimensional time series actual wind disturbance deviation is mapped to a high-dimensional disturbance displacement vector. The L2 norm of the high-dimensional disturbance displacement vector and the theoretical wind disturbance deviation is used as the disturbance displacement, and the ratio of the theoretical wind disturbance deviation to the disturbance displacement is used as the UAV maneuver weight.

[0069] In this embodiment, the method for obtaining the preset wind field model includes:

[0070] The disturbance displacement is extracted based on grid coordinates and flight data sampling time step to obtain a three-dimensional matrix. The duration of gusts is used as the short-term window length, and the long-term window length is determined according to the turbulence time scale. Based on the three-dimensional matrix, the short-term disturbance standard deviation, maximum modulus, minimum modulus and long-term disturbance mean are extracted by sliding between the short-term window length and the long-term window length. A four-dimensional feature matrix is ​​obtained based on the grid coordinates.

[0071] The disturbance displacement is corrected coordinate-by-coordinate based on the three-dimensional matrix. If the UAV maneuver weight of the grid is greater than 0.7, the product of the disturbance displacement and the UAV maneuver inverse weight is used as the correction value. If the UAV maneuver weight is less than 0.3, it remains unchanged. If the network is located within the wake diffusion range, the wind speed attenuation law value corresponding to the grid coordinates is used as the weight coefficient of the disturbance displacement to obtain the corrected three-dimensional matrix. The corrected three-dimensional matrix is ​​then used to perform inversion through a preset aerodynamic model to obtain the inverted wind speed and inverted wind direction. The dynamic characteristics of the wind field are calculated based on the inverted wind speed, inverted wind direction, and four-dimensional feature matrix. The formula for calculating the dynamic characteristics of the wind field is as follows:

[0072]

[0073] Where γ(t) represents the wind field dynamics at time t, s is the short-term window length, Δt is the flight data sampling time step, h is the time index, and d(h) is the disturbance displacement. Let denoted as the mean displacement of the disturbance within the short-term window, tanh be the hyperbolic tangent function, f be the length of the long-term window, and d be the mean displacement of the disturbance within the short-term window. max d represents the maximum disturbance displacement within the short-term window. min This represents the minimum disturbance displacement within a short-term window.

[0074] Based on the inverted wind speed, inverted wind direction, and dynamic characteristics of the wind field, Kriging interpolation is performed on the preset area and the unmeasured area to obtain the preset wind field model.

[0075] In this embodiment, the method for obtaining the predicted wind field includes:

[0076] Based on the boundary conditions of the unmeasured area, the wake diffusion boundary and wind speed attenuation law are extracted. According to the wake diffusion boundary and wind speed attenuation law, the grid data of the unmeasured area is corrected using a preset wind field model to obtain the prior mean of the prediction points. Using the prior mean of the prediction points in the unmeasured area and the grid data of the preset area, the prediction points in the unmeasured area are calculated through Gaussian process regression. The Gaussian process regression formula is as follows:

[0077]

[0078] in For the predicted points in the undetermined area, μ(x) u ,t u) represents the prior mean of the predicted points, N represents the total number of spatial grids in the preset area, M represents the time step of the preset area, and x represents the value of x. i Let t be the spatial coordinates of the i-th grid in the preset area. k The timestamp of the k-th time step in the preset area. The inversion wind speed is the i-th grid and the k-th time step in the preset area. x is the average wind speed u The coordinates of the predicted points in the undetermined area, t u For the predicted timestamp, x u The predicted spatial coordinates;

[0079] Align the set of predicted points in the unmeasured area with the grid and time step of the preset area to obtain the predicted wind field in the unmeasured area.

[0080] In this embodiment, the method for obtaining the target wind field model includes:

[0081] Based on different drone preset trajectories, discrete preset trajectory points are obtained according to grid coordinates and flight data sampling time step. The required wind speed and required wind direction are obtained from the preset trajectory points according to the wind resistance test requirements.

[0082] Based on the grid coordinates of preset trajectory points, predicted wind speed and direction are obtained through predicted wind field. Using the predicted wind speed and direction of the preset trajectory points as initial target values, the target wind field loss function is iteratively optimized based on the required wind speed, required wind direction, and initial target values ​​to obtain the target wind field with minimized trajectory deviation. The expression for the target wind field loss function is:

[0083]

[0084] Where L is the loss function, representing the deviation between the predicted wind field and the multi-UAV trajectory, n is the number of discrete points on a single trajectory, and v p,i To determine the required wind speed at point i on the preset trajectory, v t,i To predict the actual wind speed at point i on the trajectory, θ p,i To determine the required wind direction at point i on the preset trajectory, θ t,i To predict the actual wind direction at point i on the trajectory, λ is the smoothing coefficient, and Δt is the time step of the flight data sampling.

[0085] Secondly, in order to achieve the above objectives, the present invention also provides a wind field data modeling system for wind resistance testing of small and medium-sized unmanned aerial vehicles (UAVs), comprising:

[0086] Wind field preset area acquisition module: used to collect terrain data, wind speed data and wind turbine data of the test wind field, and determine the preset area and unmeasured area based on the test wind field and the preset trajectory of the UAV;

[0087] Wind field boundary condition acquisition module: used to divide the preset area into grids based on the terrain data, align the wind speed data and wind turbine data with the grid coordinates, obtain the preset wind field boundary conditions, and obtain the boundary conditions of the unmeasured area and the wake of the disturbed grid through CFD fluid simulation based on the preset wind field boundary conditions;

[0088] Wind field disturbance acquisition module: used to perform wind speed vector decomposition on the grid according to the preset wind field boundary conditions, and perform anti-outlier filtering and phase space reconstruction based on the wind speed vector and UAV flight deviations of different frequencies to obtain disturbance displacement and UAV maneuver weights;

[0089] Preset wind field modeling module: used to perform time-series analysis on the disturbance displacement based on grid coordinates, and obtain a preset wind field model based on the disturbance time sequence, UAV maneuver weights and disturbance grid wake;

[0090] Target wind field modeling module: This module is used to obtain the predicted wind field of the undetermined area based on the boundary conditions of the undetermined area through Gaussian process regression according to the preset wind field model. It uses different preset UAV trajectories to predict the wind field through grid temporal coupling to obtain the target wind field model.

[0091] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

Claims

1. A wind field data modeling method for wind resistance test of small and medium-sized unmanned aerial vehicles, characterized in that, The method comprises the following steps: Collecting terrain data, wind speed data and wind turbine data of a test wind field, determining a preset area and an unmeasured area according to a preset track of a UAV based on the test wind field; Based on the terrain data, the preset area is divided into grids, the wind speed data and the wind turbine data are aligned with the grid coordinates, the preset wind field boundary conditions are obtained, and the unmeasured area boundary conditions and the disturbed grid wake are obtained through CFD fluid simulation according to the preset wind field boundary conditions; According to the preset wind field boundary conditions, the grid is subjected to wind speed vector decomposition, and the anti-outlier filtering and phase space reconstruction are carried out based on the wind speed vector and the flight deviation of the UAV with different frequencies, so as to obtain the disturbed displacement and the UAV maneuvering weight; The disturbed displacement is subjected to time sequence analysis according to the grid coordinates, and the preset wind field model is obtained according to the disturbance time sequence, the UAV maneuvering weight and the disturbed grid wake; Based on the preset wind field model, the predicted wind field of the unmeasured area is obtained through Gaussian process regression according to the unmeasured area boundary conditions, and the target wind field model is obtained through grid time sequence coupling prediction wind field using different UAV preset tracks; The method for obtaining the preset wind field model comprises: Based on the grid coordinates and the flight data sampling time step, the disturbed displacement is extracted, a three-dimensional matrix is obtained, the gust duration is taken as a short-term window length, the long-term window length is determined according to the turbulence time scale, the short-term disturbance standard deviation, the maximum modulus, the minimum modulus and the long-term disturbance mean are extracted through the short-term window length and the long-term window length sliding based on the three-dimensional matrix, and a four-dimensional feature matrix is obtained according to the network coordinates; The disturbed displacement is corrected according to the three-dimensional matrix, if the UAV maneuvering weight of the grid is greater than 0.7, the product of the disturbed displacement and the inverse weight of the UAV maneuvering is taken as the correction value, if the UAV maneuvering weight is less than 0.3, it remains unchanged, if the network is located in the wake diffusion range, the wind speed attenuation law value corresponding to the grid coordinates is taken as the weight coefficient of the disturbed displacement, a corrected three-dimensional matrix is obtained, the corrected three-dimensional matrix is used to carry out inversion through a preset aerodynamic model, the inversion wind speed and the inversion wind direction are obtained, and the wind field dynamic characteristics are calculated according to the inversion wind speed, the inversion wind direction and the four-dimensional feature matrix, wherein the wind field dynamic characteristics calculation formula is: ； wherein is the wind field dynamic characteristic at time is the short-term window length, is the flight data sampling time step, is the time index, is the perturbation displacement, is the mean of the perturbation displacement within the short-term window, is the hyperbolic tangent function, is the long-term window length, is the maximum perturbation displacement within the short-term window, is the minimum perturbation displacement within the short-term window;​ The preset area and the unmeasured area are subjected to Kriging interpolation based on the inversion wind speed, the inversion wind direction and the wind field dynamic characteristics, and the preset wind field model is obtained.

2. The wind field data modeling method for wind resistance test of a small or medium-sized unmanned aerial vehicle according to claim 1, characterized in that, The method for determining the preset area and the unmeasured area comprises: Collecting terrain data, wind speed data and wind turbine data of a test wind field, determining a preset area and an unmeasured area according to a preset track of a UAV based on the test wind field; Based on the terrain data, the preset area is divided into grids, the wind speed data and the wind turbine data are aligned with the grid coordinates, the preset wind field boundary conditions are obtained, and the unmeasured area boundary conditions and the disturbed grid wake are obtained through CFD fluid simulation according to the preset wind field boundary conditions; According to the preset wind field boundary conditions, the grid is subjected to wind speed vector decomposition, and the anti-outlier filtering and phase space reconstruction are carried out based on the wind speed vector and the flight deviation of the UAV with different frequencies, so as to obtain the disturbed displacement and the UAV maneuvering weight; The disturbed displacement is subjected to time sequence analysis according to the grid coordinates, and the preset wind field model is obtained according to the disturbance time sequence, the UAV maneuvering weight and the disturbed grid wake; Based on the preset wind field model, the predicted wind field of the unmeasured area is obtained through Gaussian process regression according to the unmeasured area boundary conditions, and the target wind field model is obtained through grid time sequence coupling prediction wind field using different UAV preset tracks; The method for obtaining the preset wind field model comprises: Based on the grid coordinates and the flight data sampling time step, the disturbed displacement is extracted, a three-dimensional matrix is obtained, the gust duration is taken as a short-term window length, the long-term window length is determined according to the turbulence time scale, the short-term disturbance standard deviation, the maximum modulus, the minimum modulus and the long-term disturbance mean are extracted through the short-term window length and the long-term window length sliding based on the three-dimensional matrix, and a four-dimensional feature matrix is obtained according to the network coordinates; The disturbed displacement is corrected according to the three-dimensional matrix, if the UAV maneuvering weight of the grid is greater than 0.7, the product of the disturbed displacement and the inverse weight of the UAV maneuvering is taken as the correction value, if the UAV maneuvering weight is less than 0.3, it remains unchanged, if the network is located in the wake diffusion range, the wind speed attenuation law value corresponding to the grid coordinates is taken as the weight coefficient of the disturbed displacement, a corrected three-dimensional matrix is obtained, the corrected three-dimensional matrix is used to carry out inversion through a preset aerodynamic model, the inversion wind speed and the inversion wind direction are obtained, and the wind field dynamic characteristics are calculated according to the inversion wind speed, the inversion wind direction and the four-dimensional feature matrix, wherein the wind field dynamic characteristics calculation formula is: The preset area and the unmeasured area are subjected to Kriging interpolation based on the inversion wind speed, the inversion wind direction and the wind field dynamic characteristics, and the preset wind field model is obtained. A one-kilometer range around the preset area is taken as a preset buffer zone, slope change rates and obstacle densities in the preset buffer zone are counted according to terrain data, terrain complexity is obtained by weighting and adding the slope change rates and the obstacle densities, the terrain complexity is taken as an exponential decay coefficient to multiply the preset buffer zone to obtain an actual buffer zone, and a test wind field outside the preset area and the actual buffer zone is taken as an unmeasured area.

3. The wind field data modeling method for wind resistance test of small and medium-sized unmanned aerial vehicles according to claim 1, characterized in that, A method for obtaining boundary conditions and a disturbance grid wake of the unmeasured area comprises: The preset area is divided by a structured grid based on the terrain data, the grid resolution is adjusted according to the spatial density of the trajectory coordinates, and the grid resolution is multiplied in the terrain mutation area according to the slope difference, the sampling point coordinates of the wind speed data and the installation coordinates of the wind turbine are converted to the same spatial coordinate system as the grid; The wind speed data and the wind turbine data are distributed to each grid cell by a spatial interpolation algorithm according to the grid coordinates to form preset wind field boundary conditions, a calculation domain containing the preset area and the unmeasured area is constructed according to the preset wind field boundary conditions, and an inlet boundary, an outlet boundary and a wall boundary are set; The fluid motion control equation is solved by detached eddy simulation based on the calculation domain to obtain the boundary conditions of the unmeasured area, and the disturbance grid wake is extracted, the disturbance grid wake includes a wake diffusion range, a wind speed decay law, a gust duration, a turbulent intensity distribution characteristic and a turbulent time scale.

4. The wind field data modeling method for wind resistance test of small and medium-sized unmanned aerial vehicles according to claim 1, characterized in that, A method for obtaining the disturbance displacement and the unmanned aerial vehicle maneuvering weight comprises: The wind speed vector of each grid is obtained based on the preset wind field boundary conditions, the wind speed vector is decomposed into an axial component along the preset trajectory of the unmanned aerial vehicle and a radial component perpendicular to the trajectory, the initial state of the position, the speed and the attitude in the preset trajectory of the unmanned aerial vehicle is taken as an integral initial value according to the axial component and the radial component, and the wind disturbance theoretical deviation of the preset trajectory under the action of the wind field alone is obtained by a preset aerodynamic model according to the integral initial value; The initial state is taken as a simulation initial value, the disturbance displacement component under the action of each flight control instruction alone is obtained by performing a windless simulation on the flight control instruction of the preset trajectory through unmanned aerial vehicle maneuvering dynamics equation, the actual flight deviation is obtained according to the actual flight trajectory of the unmanned aerial vehicle and the preset trajectory of the unmanned aerial vehicle, the abnormal points of sensor noise and electromagnetic interference are removed through anti-outlier filtering on the actual flight deviation, the Pearson correlation coefficient of the filtered deviation and the disturbance displacement component is calculated, if the correlation coefficient is greater than 0.7, the disturbance displacement component is deducted, the actual wind disturbance deviation is obtained, and if the correlation coefficient is less than 0.3, the filtered deviation is taken as the actual wind disturbance deviation; The embedding dimension is determined by the false nearest neighbor method according to the actual wind disturbance deviation, the delay time is determined by the mutual information method, the actual wind disturbance deviation is reconstructed in a phase space according to the embedding dimension and the delay time, the one-dimensional time series actual wind disturbance deviation is mapped into a high-dimensional disturbance displacement vector, the L2 norm of the high-dimensional disturbance displacement vector and the wind disturbance theoretical deviation is taken as the disturbance displacement, and the ratio of the wind disturbance theoretical deviation to the disturbance displacement is taken as the unmanned aerial vehicle maneuvering weight.

5. The wind field data modeling method for wind resistance test of small and medium-sized unmanned aerial vehicles according to claim 1, characterized in that, A method for obtaining the predicted wind field comprises: The tail flow diffusion boundary and the wind speed decay law are extracted based on the boundary condition of the unmeasured area, the grid data of the unmeasured area is corrected by the preset wind field model according to the tail flow diffusion boundary and the wind speed decay law, the prior mean value of the prediction point is obtained, and the prediction point of the unmeasured area is calculated by using the prior mean value of the prediction point of the unmeasured area and the grid data of the preset area through Gaussian process regression, and the Gaussian process regression formula is: ； in Predicted points for unmeasured areas. Let be the prior mean of the predicted points. The total number of spatial grids in the preset area. The number of time steps for the preset area. For the preset area Spatial coordinates of each grid, For the preset area The timestamp of each time step For the preset area The grid, the first The inversion wind speed at each time step This is the average wind speed. The coordinates of the predicted points in the undetermined area are as follows. For the predicted timestamp, The predicted spatial coordinates; The point set of the prediction point of the unmeasured area is aligned with the grid and time step of the preset area, and the prediction wind field of the unmeasured area is obtained.

6. The wind field data modeling method for wind resistance test of small and medium-sized unmanned aerial vehicles according to claim 1, characterized in that, The method for obtaining the target wind field model comprises: According to different unmanned aerial vehicle preset trajectories, discrete preset trajectory points are obtained according to grid coordinates and flight data sampling time steps, and demand wind speed and demand wind direction are obtained according to anti-wind experiment requirements; According to the grid coordinates of the preset trajectory points, the prediction wind speed and the prediction wind direction are obtained through the prediction wind field, the prediction wind speed and the prediction wind direction of the preset trajectory points are taken as target initial values, and the target wind field loss function is iteratively optimized according to the demand wind speed, the demand wind direction and the target initial values, so that the target wind field with minimized trajectory deviation is obtained, and the target wind field loss function expression is: ； wherein is a loss function, representing the deviation of the predicted wind field from the multi-UAV trajectory, is the number of discrete points of a single trajectory, is the required wind speed at the point of the preset trajectory, is the actual wind speed of the predicted wind field at the point of the trajectory, is the required wind direction at the point of the preset trajectory, is the actual wind direction of the predicted wind field at the point of the trajectory, is a smoothing coefficient, is the flight data sampling time step.

7. A wind field data modeling system for wind resistance test of small and medium-sized unmanned aerial vehicles, used to perform the wind field data modeling method for wind resistance test of small and medium-sized unmanned aerial vehicles according to any one of claims 1-6, characterized in that, It comprises: The wind field preset area acquisition module is used for collecting the terrain data, the wind speed data and the fan data of the test wind field, determining the preset area and the unmeasured area based on the unmanned aerial vehicle preset trajectory of the test wind field; The wind field boundary condition acquisition module is used for carrying out grid division on the preset area based on the terrain data, aligning the wind speed data, the fan data and the grid coordinates, obtaining the preset wind field boundary condition, and obtaining the unmeasured area boundary condition and the disturbed grid tail flow through the CFD fluid simulation according to the preset wind field boundary condition; The wind field disturbance acquisition module is used for carrying out wind speed vector decomposition on the grid according to the preset wind field boundary condition, carrying out anti-wild value filtering and phase space reconstruction based on the wind speed vector and the unmanned aerial vehicle flight deviation of different frequencies, and obtaining the disturbance displacement and the unmanned aerial vehicle maneuvering weight; The preset wind field modeling module is used for carrying out time sequence analysis on the disturbance displacement according to the grid coordinates, and obtaining the preset wind field model according to the disturbance time sequence, the unmanned aerial vehicle maneuvering weight and the disturbed grid tail flow; The target wind field modeling module is used for obtaining the prediction wind field of the unmeasured area based on the preset wind field model according to the unmeasured area boundary condition, and obtaining the target wind field model by coupling the prediction wind field through grid time sequence according to different unmanned aerial vehicle preset trajectories.

Images