A wind load calculation method based on joint distribution of wind speed and wind direction
By constructing a joint distribution model of wind speed and direction and an improved PI-DeepONet model, the problem of neglecting the correlation between wind speed and direction in wind load calculation is solved, enabling refined calculation and dynamic risk assessment of structural wind load, and improving the accuracy of calculation and engineering applicability.
Patent Information
- Application Number
- CN202611292031.8
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2026-08-25
- Publication Date
- 2026-11-06
- Estimated Expiration
- 2046-08-25
AI Technical Summary
Existing wind load calculation methods ignore the strong correlation between wind speed and wind direction, resulting in calculation results that deviate from reality and fail to truly reflect the joint distribution characteristics of wind loads. Furthermore, they lack dynamic response analysis to changes in the wind environment during the service life of the structure.
A calculation method based on the joint distribution of wind speed and direction is adopted. By constructing a joint probability density function of wind speed and direction, and combining an improved PI-DeepONet model and a time derivative decoupling mechanism, real-time evaluation of structural response and dynamic wind load updates are achieved.
It improves the accuracy and efficiency of wind load calculation, ensures the physical interpretability and engineering applicability of the model, and realizes refined solution of structural wind load and risk management throughout the entire life cycle.
Smart Images

Figure CN122819078B_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of structural wind engineering, probability statistics and artificial intelligence, and in particular to a method for calculating wind loads based on the joint distribution of wind speed and wind direction. Background Technology
[0002] Wind load is one of the main lateral loads acting on building structures, and accurate calculation of wind load is crucial for ensuring the safety of engineering structures such as high-rise buildings and long-span bridges. In existing wind load calculation theories and engineering practices, long-term meteorological observation data are typically analyzed using mathematical statistics methods based on standards such as the "Code for Design of Building Structures." Traditional approaches often treat wind speed and wind direction as independent random variables, fitting separate probability distribution models for wind speed and wind direction, or simply statistically analyzing the frequency of each wind direction and then combining this with the wind pressure coefficient given in the standard for calculation. While this method is simple and easy to implement, it ignores the strong correlation between wind speed and wind direction; that is, a specific wind direction often corresponds to a specific wind speed distribution characteristic. This independent assumption leads to calculated wind load extreme values that often deviate from reality, failing to truly reflect the combined distribution characteristics of wind loads. This can result in overly conservative structural designs leading to waste, or underestimating risks and causing safety hazards. Furthermore, when calculating extreme wind speeds, existing technologies mostly use a single extreme distribution model to fit data across all wind directions, lacking a refined fitting for the tail characteristics of wind speed distribution in different wind direction sectors, resulting in insufficient accuracy of the calculated extreme wind speeds.
[0003] With the development of computer technology, computational fluid dynamics (CFD) has been widely used to determine building shape coefficients. However, existing CFD simulation results are often disconnected from probabilistic statistical models, failing to fully consider the weighted influence of wind direction probability, resulting in insufficiently refined calculated wind pressure standard values. On the other hand, in the field of structural response analysis, traditional numerical simulation methods such as the finite element method, while highly accurate, are computationally time-consuming and cannot meet the needs of real-time risk assessment. While data-driven methods such as neural networks, which have emerged in recent years, offer fast computation speeds, they often lack physical constraints, making them "black box" models. Their predictions may violate fundamental physical laws such as Newton's laws of motion, resulting in limited generalization ability and reliability. More importantly, existing wind load calculation methods are mostly static offline analyses, lacking mechanisms for online correction of model parameters using real-time meteorological data, and failing to reflect the impact of wind environment changes on structural failure risk during structural service. Therefore, how to construct a wind load calculation method that integrates the joint distribution characteristics of wind speed and direction, physical information neural networks, and dynamic update mechanisms to improve calculation accuracy, physical interpretability, and timeliness has become a pressing technical problem in the field of structural wind engineering.
[0004] Therefore, how to provide a wind load calculation method based on the combined distribution of wind speed and wind direction is a problem that urgently needs to be solved by those skilled in the art. Summary of the Invention
[0005] One objective of this invention is to propose a wind load calculation method based on the joint distribution of wind speed and direction. This invention fully integrates key steps such as joint probabilistic modeling of wind speed and direction, extreme wind speed estimation by wind direction, construction of an improved PI-DeepONet model, and dynamic wind load updating. It constructs an intelligent calculation process that accurately captures the probabilistic characteristics of the wind environment, embeds physical constraints of structural response, assesses failure risks in real time, and dynamically corrects design parameters. This enables refined solution of structural wind loads, physical consistency verification, and full life-cycle risk management in complex wind field environments. This invention possesses advantages such as accurate edge distribution fitting, reasonable extreme wind speed estimation, strong physical interpretability of the model, and reliable failure probability prediction. By introducing an improved PI-DeepONet model and embedding a time derivative decoupling mechanism, it effectively solves the problems in existing technologies, such as calculation bias caused by ignoring the correlation between wind speed and direction, low efficiency of traditional numerical simulation, and prediction results that violate dynamic laws due to the lack of physical constraints in purely data-driven models. This significantly improves the accuracy, efficiency, and engineering applicability of wind load calculation.
[0006] A wind load calculation method based on the combined distribution of wind speed and wind direction according to an embodiment of the present invention includes the following steps:
[0007] S1. Collect historical wind speed and direction data and preprocess them. Extract wind speed and wind direction sequences as random variables, calculate statistical feature values, and fit the edge distribution function of wind speed and the edge distribution function of wind direction.
[0008] S2. Using the edge distribution function of wind speed and the edge distribution function of wind direction, construct the joint probability density function of wind speed and wind direction, calculate the correlation coefficient between wind speed and wind direction and quantify the dependency relationship, and output the probability density data.
[0009] S3. Divide the entire wind direction space into N sectors. Based on the probability density data, statistically analyze the wind speed data falling into each sector. Fit the tail data of the wind speed distribution in each sector to estimate the extreme wind speed and output the extreme wind speed dataset by wind direction.
[0010] S4. Obtain the aerodynamic shape parameters of the target building structure, calculate the shape coefficient under different wind angles, and perform weighted calculation with the extreme wind speed dataset of wind direction, and output a refined wind pressure standard value by combining the wind pressure height variation coefficient.
[0011] S5. The refined wind pressure standard value is used as the input load function to input the improved PI-DeepONet model to obtain the failure risk index of the structure during its service life. The improved PI-DeepONet model includes a load feature branch network, a spatiotemporal feature backbone network, a physical information time integration module, and a response statistical output module. The physical information time integration module embeds a time derivative decoupling mechanism.
[0012] S6. Collect real-time wind speed and direction data, correct the parameters of the wind speed edge distribution function and the wind direction edge distribution function, calculate the corresponding wind direction extreme value wind speed dataset and refined wind pressure standard value, update the failure risk index of the structure during its service life, and output the dynamically optimized wind load calculation results.
[0013] Optionally, S1 specifically includes:
[0014] S11. Export historical wind speed and direction records from the database of long-term meteorological observation stations in the target area, remove abnormal data points with values of zero or exceeding the preset physical range threshold, and perform linear interpolation to fill in the gaps to obtain preprocessed historical wind speed and direction data.
[0015] S12. The preprocessed historical wind speed and direction data are split into two data sequences in chronological order. The wind speed value column is extracted as the wind speed sequence, and the wind direction angle column is extracted as the wind direction sequence. The mean and variance of the wind speed sequence and the mean and standard deviation of the wind direction sequence are calculated respectively, and the statistical feature values are output.
[0016] S13. Substitute the mean wind speed and standard deviation of wind speed in the statistical characteristic values into the moment estimation formula, calculate the ratio of the mean wind speed to the standard deviation of wind speed, and refer to the preset comparison table to obtain the initial value of the shape parameter. Substitute the mean wind speed into the algebraic expression containing the initial value of the shape parameter to calculate the initial value of the scale parameter. Substitute the initial value of the shape parameter, the initial value of the scale parameter, and the wind speed sequence into the log-likelihood function, calculate the log probability value of each value in the wind speed sequence, and sum all the log probability values to obtain the log-likelihood function value.
[0017] S14. Iteratively adjust the initial values of the shape parameter and the scale parameter using the gradient ascent algorithm until the log-likelihood function value no longer increases. Substitute the adjusted parameters into the Weibull distribution probability density formula, calculate the shape parameter adjustment value raised to the power of the scale parameter adjustment value, divide the wind speed sequence value by the scale parameter adjustment value and calculate the corresponding shape parameter adjustment value minus one. Multiply the two calculation results by a negative exponent term with the natural constant as the base and the wind speed sequence value divided by the scale parameter adjustment value as the exponent to obtain the wind speed marginal distribution function.
[0018] S15. Substitute the mean wind direction value in the statistical feature value into the mixed von Mises distribution formula with a preset number of components, calculate the initial values of the mixing ratio, mean direction, and concentration parameter of each component, and substitute them with the wind direction sequence into the log-likelihood function formula for cumulative calculation. Execute the expectation step in the expectation maximization algorithm to calculate the posterior probability of the latent variables, and execute the maximization step to update the initial values of the mixing ratio, mean direction, and concentration parameter using the posterior probability until the change in the log-likelihood function value is less than the preset convergence threshold.
[0019] S16. Substitute the updated parameters into the mixed von Mises distribution probability density formula, multiply the updated mixing ratio value by the exponential function value with the product of the natural constant, the concentration parameter update value and the cosine function value as the exponent, to obtain the component probability value. Sum the component probability values over the preset number of components to obtain the wind direction edge distribution function.
[0020] Optionally, S2 specifically includes:
[0021] S21. Substitute the cumulative probability value of wind speed output by the wind speed marginal distribution function and the cumulative probability value of wind direction output by the wind direction marginal distribution function into the expression of the preset Copula function model, calculate the natural logarithm of the product of the cumulative probability value of wind speed and the cumulative probability value of wind direction and square it, and sum the squared values of all samples to obtain the goodness-of-fit statistic.
[0022] S22. Substitute the goodness-of-fit statistic into the correlation coefficient calculation formula, divide the goodness-of-fit statistic by the difference between the total number of samples minus one, take the square root and take the negative value to obtain the correlation coefficient between wind speed and wind direction.
[0023] S23. Substitute the wind speed sequence, wind direction sequence, and correlation coefficient into the likelihood function expression, calculate the probability value of each wind speed sample and wind direction sample occurring simultaneously, and sum the natural logarithms of all probability values to obtain the likelihood function value. Set the step size search range, select multiple candidate dependency parameter values within the search range, substitute each candidate dependency parameter value into the likelihood function expression to calculate the corresponding likelihood function value, compare and select the candidate dependency parameter value with the largest likelihood function value as the optimal dependency parameter, and quantify the dependency relationship between wind speed and wind direction.
[0024] S24. Substitute the optimal dependency parameters into the expression of the preset Copula function model, calculate the bivariate joint probability value of the cumulative probability value of wind speed and the cumulative probability value of wind direction, calculate the partial derivatives of the bivariate joint probability value with respect to the wind speed and wind direction variables respectively, multiply the two partial derivative calculation results to obtain the joint probability density function of wind speed and wind direction, arrange the discrete point values of the joint probability density function on the wind speed axis and the wind direction axis in a grid, and output the probability density data.
[0025] Optionally, S3 specifically includes:
[0026] S31. Starting from the due north direction, divide the entire wind direction space clockwise according to the preset angle interval to obtain N sectors. Read the probability density data, accumulate and normalize the probability density data falling into each sector to obtain the occurrence probability of each sector, and select the sector with the occurrence probability greater than the preset probability threshold as the effective sector.
[0027] S32. Extract wind speed samples falling into each effective sector to construct a sub-sample set. Sort the wind speed data in the sub-sample set in ascending order of numerical value. Select the wind speed value at a preset proportion position in the sorting sequence as the threshold. Extract the wind speed data greater than the threshold as the tail data of the wind speed distribution.
[0028] S33. Substitute the tail data of the wind speed distribution into the log-likelihood function of the generalized Pareto distribution, calculate the partial derivatives of the shape parameter and the scale parameter and set the partial derivatives to zero, solve the system of equations to obtain the shape parameter and the scale parameter, substitute the tail data of the wind speed distribution, the shape parameter and the scale parameter into the cumulative distribution function formula of the generalized Pareto distribution, calculate the power term with the product of the scale parameter plus the difference between the shape parameter and the tail data of the wind speed distribution minus the threshold as the denominator and the product of the shape parameter and the difference between the tail data of the wind speed distribution minus the threshold as the numerator, add one to the negative power term to obtain the tail model of the wind speed distribution for each effective sector;
[0029] S34. Obtain the design reference period of the meteorological observation station, calculate the ratio of the total number of hours in the design reference period to the total number of tail data of the wind speed distribution, substitute the ratio as the exceedance probability into the inverse function of the tail model of the wind speed distribution, and calculate the extreme wind speed of each effective sector in the design reference period.
[0030] S35. Summarize the extreme wind speeds of all valid sectors, set the extreme wind speeds of non-valid sectors to zero, arrange them in order of sector number, and output the extreme wind speed dataset by wind direction.
[0031] Optionally, S4 specifically includes:
[0032] S41. Obtain the building geometric dimension parameters of the target building structure, establish a three-dimensional geometric model, place the three-dimensional geometric model in the wind field of the computational fluid dynamics software, adjust the inlet boundary conditions according to the preset wind direction angle interval, simulate the airflow process, extract the average wind pressure value of the building surface and the reference wind pressure value of the far reference point, calculate the quotient of the average wind pressure value divided by the reference wind pressure value, and obtain the shape coefficient under different wind direction angles.
[0033] S42. Read the extreme wind speed dataset by wind direction, square the extreme wind speed in the extreme wind speed dataset, multiply the squared value of the extreme wind speed by the corresponding body shape coefficient, and multiply the product by half of the air density value to obtain the extreme wind pressure of each sector.
[0034] S43. Obtain the ground roughness index of the target building structure location and the height value of the calculated location. Divide the height value of the calculated location by the preset standard height value to obtain the height ratio. Raise the height ratio to the power of the ground roughness index. Multiply the result of the power operation by the wind pressure coefficient at the preset standard height to obtain the wind pressure height variation coefficient. The preset standard height value is 10 meters and the wind pressure coefficient at the preset standard height is 1.0.
[0035] S44. Read the probability density data, use the values in the probability density data as weights, multiply the extreme wind pressure, wind pressure height variation coefficient and weights to obtain the weighted wind pressure value, sum the weighted wind pressure of all effective sectors and divide by the total weights to output the refined wind pressure standard value.
[0036] Optionally, S5 specifically includes:
[0037] S51. Construct an improved PI-DeepONet model that includes a load feature branch network, a spatiotemporal feature backbone network, a physical information time integration module, and a response statistical output module. Use the refined wind pressure standard value as the input variable of the function space. Construct a fully connected neural network containing an input layer, a hidden layer, and an output layer as the load feature branch network. Initialize the network weight matrix and bias vector. Perform a linear transformation on the refined wind pressure standard value to obtain the input of the first hidden layer and substitute it into the ReLU activation function to obtain the output of the first hidden layer.
[0038] S52. Perform matrix multiplication between the output of the first hidden layer and the weight matrix of the next layer, and add the bias vector to obtain the input of the second hidden layer. Substitute the input of the second hidden layer into the Tanh activation function to obtain the output of the second hidden layer. Iterate in this way until the output layer is reached. The output layer outputs the load encoding vector.
[0039] S53. Discretize the preset time domain according to the preset integration step size to generate a time coordinate sampling point sequence. Construct a multilayer perceptron containing an input layer, a hidden layer and an output layer as the spatiotemporal feature backbone network. Perform a linear transformation on the spatiotemporal coordinate sampling point sequence to obtain the hidden layer input and substitute it into the Sigmoid activation function for nonlinear mapping to obtain the hidden layer output.
[0040] S54. Pass the hidden layer output to the next layer until the output layer, output the spatiotemporal coding vector, calculate the dot product of the load coding vector and the spatiotemporal coding vector, and obtain the scalar form of the structural response prediction value.
[0041] S55. Input the predicted structural response values into the physical information time integration module. Based on the time derivative decoupling mechanism, split the predicted structural response values into a sequence of predicted displacement values and a sequence of predicted velocity values, and calculate the internal force response vector and physical constraint loss function value. Based on the physical constraint loss function value, use the gradient descent algorithm to update the weight matrix and bias vector of the load feature branch network and the spatiotemporal feature backbone network.
[0042] S56. Input the displacement prediction value sequence output by the improved PI-DeepONet model after updating the weight parameters into the response statistics output module, read the preset displacement failure threshold, traverse each value in the displacement prediction value sequence, calculate the absolute value of the value and compare it with the preset displacement failure threshold. When it is greater than the preset displacement failure threshold, increment the value of the failure counter by one. After the traversal is completed, calculate the ratio of the value of the failure counter to the total length of the displacement prediction value sequence, use the ratio as the failure probability value, and output the failure risk index of the structure during its service life.
[0043] Optionally, S55 specifically includes:
[0044] S551. Input the predicted structural response value into the physical information time integration module. Based on the time derivative decoupling mechanism, expand the predicted structural response value into a response sequence according to the time step. Set the element at the odd index position in the response sequence as the displacement prediction value to construct the displacement prediction value sequence, and set the element at the even index position as the velocity prediction value to construct the velocity prediction value sequence.
[0045] S552. Construct a chain for calculating the derivative of the displacement prediction value with respect to the time variable. Calculate the first derivative of the displacement prediction value with respect to the time variable as a physical constraint on the velocity prediction value. Calculate the product of the velocity prediction value and the preset time step to obtain the displacement increment prediction value. Add the displacement increment prediction value to the displacement prediction value at the previous moment to obtain the updated displacement value at the current moment.
[0046] S553. Calculate the second derivative of the displacement prediction value with respect to the time variable as the acceleration prediction value, obtain the mass matrix, damping matrix and stiffness matrix of the target building structure, calculate the matrix product of the mass matrix and the acceleration prediction value to obtain the inertial force vector, calculate the matrix product of the damping matrix and the velocity prediction value to obtain the damping force vector, and calculate the matrix product of the stiffness matrix and the displacement prediction value to obtain the restoring force vector.
[0047] S554. Calculate the sum of corresponding elements of the inertial force vector, damping force vector, and restoring force vector to obtain the structural internal force response vector. Read the refined wind pressure standard value as the external load vector. Calculate the difference between the external load vector and the structural internal force response vector to obtain the dynamic residual vector. Calculate the sum of the squares of each element in the dynamic residual vector and divide it by the total number of samples to obtain the physical constraint loss function value. Use the backpropagation algorithm to calculate the gradient of the physical constraint loss function value with respect to the weight matrix and bias vector in the improved PI-DeepONet model. Use the gradient descent algorithm to update the weight matrix and bias vector of the load feature branch network and the spatiotemporal feature backbone network.
[0048] Optionally, S6 specifically includes:
[0049] S61. Collect the current real-time wind speed and wind direction angle values of the target area, construct real-time sample data pairs, and add them to the database of historical wind speed and wind direction data. Substitute the updated historical wind speed and wind direction data into the wind speed edge distribution function calculation process, update the shape parameter adjustment value and scale parameter adjustment value, and output the updated wind speed edge distribution function.
[0050] S62. Substitute the updated historical wind speed and direction data into the wind direction edge distribution function calculation process, output the updated wind direction edge distribution function, substitute the updated wind speed edge distribution function and the updated wind direction edge distribution function into the joint probability density function calculation process, update the correlation coefficient and optimal dependency parameter between wind speed and wind direction, and output the updated probability density data.
[0051] S63. Substitute the updated probability density data into the extreme wind speed calculation process, output the updated extreme wind speed data by wind direction, substitute the updated extreme wind speed dataset by wind direction into the refined wind pressure standard value calculation process, output the updated refined wind pressure standard value, input the updated refined wind pressure standard value into the improved PI-DeepONet model, and output the dynamically optimized failure risk index.
[0052] S64. Read the updated refined wind pressure standard value as the wind pressure load intensity, read the windward area of the target building structure, calculate the product of the refined wind pressure standard value and the windward area to obtain the concentrated wind load value, calculate the product of the concentrated wind load value and the dynamically optimized failure risk index, add the product result to the preset wind load partial factor, and output the dynamically optimized wind load calculation result. The preset wind load partial factor is 1.4.
[0053] The beneficial effects of this invention are:
[0054] This invention collects historical wind speed and direction data and fits marginal distribution functions to construct a joint probability density function for wind speed and direction. Addressing the issue of neglected correlation between wind speed and direction leading to calculation errors, it quantifies the dependency between the two and outputs probability density data, thereby calculating extreme wind speeds and refined wind pressure standard values for each wind direction. By acquiring the aerodynamic shape parameters and shape coefficients of the target building structure and combining them with the extreme wind speed dataset for each wind direction for weighted calculation, it achieves a refined solution for wind pressure considering wind direction distribution characteristics. Finally, it constructs an improved PI-DeepONet model that includes a load feature branch network, a spatiotemporal feature backbone network, and a physical information time integration module. To address the issues of low computational efficiency in traditional numerical simulations and the lack of physical constraints in purely data-driven models, this paper utilizes a time derivative decoupling mechanism to decompose structural response predictions into sequences of displacement and velocity predictions. By combining the mass, damping, and stiffness matrices, the paper calculates the structural internal force response vector and dynamic residual vector, constructs a physical constraint loss function, and updates the network weights using a backpropagation algorithm, achieving a deep integration of data-driven approaches and physical mechanisms. Furthermore, by collecting real-time wind speed and direction data, the parameters of the wind speed and wind direction edge distribution functions are corrected, updating the failure risk indicators and wind load calculation results, achieving dynamic optimization throughout the entire lifecycle. Ultimately, this approach enables refined calculation of wind loads on building structures, physical consistency verification, and dynamic risk assessment, effectively improving the accuracy of wind load calculations, the reliability of model predictions, and the adaptability of engineering applications. Attached Figure Description
[0055] The accompanying drawings are provided to further illustrate the invention and form part of the specification. They are used in conjunction with embodiments of the invention to explain the invention and do not constitute a limitation thereof. In the drawings:
[0056] Figure 1 This is a flowchart of a wind load calculation method based on the joint distribution of wind speed and wind direction proposed in this invention;
[0057] Figure 2 This is a flowchart of the construction of the joint wind speed and wind direction distribution and the calculation of extreme wind speeds by wind direction proposed in this invention;
[0058] Figure 3 This is a flowchart illustrating the construction of the improved PI-DeepONet model and the calculation of structural failure risk indicators proposed in this invention. Detailed Implementation
[0059] The present invention will now be described in further detail with reference to the accompanying drawings. These drawings are simplified schematic diagrams, illustrating only the basic structure of the invention, and therefore only show the components relevant to the invention.
[0060] refer to Figures 1-3A wind load calculation method based on the combined distribution of wind speed and wind direction includes the following steps:
[0061] S1. Collect historical wind speed and direction data from long-term meteorological observation stations in the target area and preprocess them. Extract wind speed and wind direction sequences as random variables, calculate statistical characteristic values, and fit the edge distribution function of wind speed and the edge distribution function of wind direction.
[0062] S2. Using the edge distribution function of wind speed and the edge distribution function of wind direction, construct the joint probability density function of wind speed and wind direction, calculate the correlation coefficient between wind speed and wind direction, quantify the dependence between the two, and output the probability density data.
[0063] S3. Divide the entire wind direction space into N sectors. Based on the probability density data, statistically analyze the wind speed data falling into each sector. Fit the tail data of the wind speed distribution of each sector to estimate the extreme wind speed of each sector within the design reference period and output the extreme wind speed dataset by wind direction.
[0064] S4. Obtain the aerodynamic shape parameters of the target building structure, calculate the shape coefficient under different wind angles, weight the shape coefficient with the extreme wind speed dataset of each wind direction, and combine it with the wind pressure height variation coefficient to output a refined wind pressure standard value that takes into account the wind direction distribution characteristics.
[0065] S5. The refined wind pressure standard value is used as the input load function to input the improved PI-DeepONet model to obtain the failure risk index of the structure during its service life. The improved PI-DeepONet model includes a load feature branch network, a spatiotemporal feature backbone network, a physical information time integration module, and a response statistical output module. The physical information time integration module embeds a time derivative decoupling mechanism.
[0066] S6. Collect real-time wind speed and direction data, correct the parameters of the wind speed edge distribution function and the wind direction edge distribution function, calculate the corresponding wind direction extreme value wind speed dataset and refined wind pressure standard value, update the failure risk index of the structure during its service life, and output the dynamically optimized wind load calculation results.
[0067] This implementation significantly improves the precision of wind load calculation and the timeliness of structural failure risk assessment. By collecting long-term meteorological observation data and fitting the edge distribution functions of wind speed and direction, a joint probability density function is constructed, effectively solving the problem of traditional methods neglecting the correlation between wind speed and direction, and achieving accurate modeling of wind environment characteristics. In the extreme wind speed estimation, combined with sector tail fitting technology, a dataset of extreme wind speeds by wind direction is output, ensuring the statistical rationality of design parameters. By acquiring the building's aerodynamic shape parameters and weighting them with probability weights, a refined wind pressure standard value considering wind direction distribution characteristics is generated, overcoming the deficiency of the standard method in poor adaptability to complex shapes. Especially in the structural response analysis stage, the wind pressure standard value is input into the improved PI-DeepONet model, and a time derivative decoupling mechanism is embedded in the physical information time integral module. This retains the efficient reasoning ability of the neural network while strictly following the laws of dynamic physics, achieving millisecond-level failure risk index output. Finally, by correcting model parameters with real-time data, a dynamically optimized wind load calculation closed loop was constructed, which significantly enhanced the safety early warning capability and the reliability of the engineering structure throughout its entire life cycle operation and maintenance in complex wind field environments.
[0068] In this embodiment, S1 specifically includes:
[0069] S11. Export historical wind speed and direction records from the database of long-term meteorological observation stations in the target area, remove abnormal data points with values of zero or exceeding the preset physical range threshold, and perform linear interpolation to fill in the gaps to obtain preprocessed historical wind speed and direction data. The preset physical range threshold is a wind speed of 60 m / s and a wind direction of 0 to 360 degrees.
[0070] S12. The preprocessed historical wind speed and direction data are split into two data sequences in chronological order. The wind speed value column is extracted as the wind speed sequence, and the wind direction angle column is extracted as the wind direction sequence. The mean and variance of the wind speed sequence and the mean and standard deviation of the wind direction sequence are calculated respectively, and the statistical feature values are output.
[0071] S13. Substitute the mean wind speed and standard deviation of wind speed in the statistical characteristic values into the moment estimation formula, calculate the ratio of the mean wind speed to the standard deviation of wind speed, and refer to the preset comparison table to obtain the initial value of the shape parameter. Substitute the mean wind speed into the algebraic expression containing the initial value of the shape parameter to calculate the initial value of the scale parameter. Substitute the initial value of the shape parameter, the initial value of the scale parameter, and the wind speed sequence into the log-likelihood function, calculate the log probability value of each value in the wind speed sequence, and sum all the log probability values to obtain the log-likelihood function value.
[0072] S14. Iteratively adjust the initial values of the shape parameter and the scale parameter using the gradient ascent algorithm until the log-likelihood function value no longer increases. Substitute the adjusted parameters into the Weibull distribution probability density formula, calculate the shape parameter adjustment value raised to the power of the scale parameter adjustment value, divide the wind speed sequence value by the scale parameter adjustment value and calculate the corresponding shape parameter adjustment value minus one. Multiply the two calculation results by a negative exponent term with the natural constant as the base and the wind speed sequence value divided by the scale parameter adjustment value as the exponent to obtain the wind speed marginal distribution function.
[0073] S15. Substitute the mean wind direction value in the statistical feature value into the mixed von Mises distribution formula with a preset number of components, calculate the initial values of the mixing ratio, mean direction, and concentration parameter of each component, and substitute them with the wind direction sequence into the log-likelihood function formula for cumulative calculation. Execute the expectation step in the expectation maximization algorithm to calculate the posterior probability of the latent variables, and execute the maximization step to update the initial values of the mixing ratio, mean direction, and concentration parameter using the posterior probability until the change in the log-likelihood function value is less than the preset convergence threshold.
[0074] S16. Substitute the updated parameters into the mixed von Mises distribution probability density formula, multiply the updated mixing ratio value by the exponential function value with the product of the concentration parameter update value and the cosine function value as the exponent, with the natural constant as the base, to obtain the component probability value. Sum the component probability values of a preset number of components to obtain the wind direction edge distribution function. The preset number of components is 3, the preset convergence threshold is 0.001, and the mixed von Mises distribution probability density formula is the weighted sum of the probability density values of each component.
[0075] In this embodiment, S2 specifically includes:
[0076] S21. Substitute the cumulative probability value of wind speed output by the wind speed marginal distribution function and the cumulative probability value of wind direction output by the wind direction marginal distribution function into the expression of the preset Copula function model, calculate the natural logarithm of the product of the cumulative probability value of wind speed and the cumulative probability value of wind direction and square it, and sum the squared values of all samples to obtain the goodness-of-fit statistic.
[0077] S22. Substitute the goodness-of-fit statistic into the correlation coefficient calculation formula, divide the goodness-of-fit statistic by the difference between the total number of samples minus one, take the square root and take the negative value to obtain the correlation coefficient between wind speed and wind direction.
[0078] S23. Substitute the wind speed sequence, wind direction sequence, and correlation coefficient into the likelihood function expression, calculate the probability value of each wind speed sample and wind direction sample occurring simultaneously, and sum the natural logarithms of all probability values to obtain the likelihood function value. Set the step size search range, select multiple candidate dependency parameter values within the search range, substitute each candidate dependency parameter value into the likelihood function expression to calculate the corresponding likelihood function value, compare and select the candidate dependency parameter value with the largest likelihood function value as the optimal dependency parameter, and quantify the dependency relationship between wind speed and wind direction.
[0079] S24. Substitute the optimal dependency parameters into the expression of the preset Copula function model, calculate the bivariate joint probability value of the cumulative probability value of wind speed and the cumulative probability value of wind direction, calculate the partial derivatives of the bivariate joint probability value with respect to the wind speed and wind direction variables respectively, multiply the two partial derivative calculation results to obtain the joint probability density function of wind speed and wind direction, arrange the discrete point values of the joint probability density function on the wind speed axis and the wind direction axis in a grid, and output the probability density data.
[0080] In this embodiment, S3 specifically includes:
[0081] S31. Taking due north as the starting edge, divide the entire wind direction space clockwise according to a preset angle interval to obtain N sectors. Read the probability density data, accumulate and normalize the probability density data falling into each sector to obtain the occurrence probability of each sector, and select sectors with an occurrence probability greater than a preset probability threshold as valid sectors. The preset angle interval is 10 degrees, and the preset probability threshold is 0.01.
[0082] S32. Extract wind speed samples falling within each valid sector to construct a sub-sample set. Sort the wind speed data in the sub-sample set in ascending order of numerical value. Select the wind speed value at a preset proportion position in the sorting sequence as a threshold. Extract wind speed data greater than the threshold as the tail data of the wind speed distribution. The preset proportion position is 90%.
[0083] S33. Substitute the tail data of the wind speed distribution into the log-likelihood function of the generalized Pareto distribution, calculate the partial derivatives of the shape parameter and the scale parameter and set the partial derivatives to zero, solve the system of equations to obtain the shape parameter and the scale parameter, substitute the tail data of the wind speed distribution, the shape parameter and the scale parameter into the cumulative distribution function formula of the generalized Pareto distribution, calculate the power term with the product of the scale parameter plus the difference between the shape parameter and the tail data of the wind speed distribution minus the threshold as the denominator and the product of the shape parameter and the difference between the tail data of the wind speed distribution minus the threshold as the numerator, add one to the negative power term to obtain the tail model of the wind speed distribution for each effective sector;
[0084] S34. Obtain the design reference period of the meteorological observation station, calculate the ratio of the total number of hours in the design reference period to the total number of tail data of the wind speed distribution, substitute the ratio as the exceedance probability into the inverse function of the tail model of the wind speed distribution, and calculate the extreme wind speed of each effective sector in the design reference period.
[0085] S35. Summarize the extreme wind speeds of all valid sectors, set the extreme wind speeds of non-valid sectors to zero, arrange them in order of sector number, and output the extreme wind speed dataset by wind direction.
[0086] In this embodiment, S4 specifically includes:
[0087] S41. Obtain the building geometric dimension parameters of the target building structure, establish a three-dimensional geometric model, place the three-dimensional geometric model in the wind field of the computational fluid dynamics software, adjust the inlet boundary conditions according to the preset wind direction angle interval, simulate the airflow process, extract the average wind pressure value of the building surface and the reference wind pressure value of the far reference point, calculate the quotient of the average wind pressure value divided by the reference wind pressure value, and obtain the shape coefficient under different wind direction angles.
[0088] S42. Read the extreme wind speed dataset by wind direction, square the extreme wind speed in the extreme wind speed dataset, multiply the squared value of the extreme wind speed by the corresponding body shape coefficient, and multiply the product by half of the air density value to obtain the extreme wind pressure of each sector.
[0089] S43. Obtain the ground roughness index of the target building structure location and the height value of the calculated location. Divide the height value of the calculated location by the preset standard height value to obtain the height ratio. Raise the height ratio to the power of the ground roughness index. Multiply the result of the power operation by the wind pressure coefficient at the preset standard height to obtain the wind pressure height variation coefficient. The preset standard height value is 10 meters and the wind pressure coefficient at the preset standard height is 1.0.
[0090] S44. Read the probability density data, use the values in the probability density data as weights, multiply the extreme wind pressure, the wind pressure height variation coefficient and the weights to obtain the weighted wind pressure value, sum the weighted wind pressure of all effective sectors and divide by the total weights to output the refined wind pressure standard value that takes into account the wind direction distribution characteristics.
[0091] In this embodiment, S5 specifically includes:
[0092] S51. Construct an improved PI-DeepONet model that includes a load feature branch network, a spatiotemporal feature backbone network, a physical information time integration module, and a response statistical output module. Use the refined wind pressure standard value as the input variable of the function space. Construct a fully connected neural network containing an input layer, a hidden layer, and an output layer as the load feature branch network. Initialize the network weight matrix and bias vector. Perform a linear transformation on the refined wind pressure standard value to obtain the input of the first hidden layer and substitute it into the ReLU activation function to obtain the output of the first hidden layer.
[0093] S52. Perform matrix multiplication between the output of the first hidden layer and the weight matrix of the next layer, and add the bias vector to obtain the input of the second hidden layer. Substitute the input of the second hidden layer into the Tanh activation function to obtain the output of the second hidden layer. Iterate in this way until the output layer is reached. The output layer outputs the load encoding vector.
[0094] S53. Discretize the preset time domain according to the preset integration step size to generate a time coordinate sampling point sequence. Construct a multilayer perceptron containing an input layer, a hidden layer, and an output layer as the spatiotemporal feature backbone network. Perform a linear transformation on the spatiotemporal coordinate sampling point sequence to obtain the hidden layer input and substitute it into the Sigmoid activation function for nonlinear mapping to obtain the hidden layer output. The preset time domain is 0 to 3600 seconds, and the preset integration step size is 0.1 seconds.
[0095] S54. Pass the hidden layer output to the next layer until the output layer, output the spatiotemporal coding vector, calculate the dot product of the load coding vector and the spatiotemporal coding vector, and obtain the scalar form of the structural response prediction value.
[0096] S55. Input the predicted structural response values into the physical information time integration module. Based on the time derivative decoupling mechanism, split the predicted structural response values into a sequence of predicted displacement values and a sequence of predicted velocity values, and calculate the internal force response vector and physical constraint loss function value. Based on the physical constraint loss function value, use the gradient descent algorithm to update the weight matrix and bias vector of the load feature branch network and the spatiotemporal feature backbone network.
[0097] S56. Input the displacement prediction value sequence output by the improved PI-DeepONet model after updating the weight parameters into the response statistics output module, read the preset displacement failure threshold, traverse each value in the displacement prediction value sequence, calculate the absolute value of the value and compare it with the preset displacement failure threshold. When it is greater than the preset displacement failure threshold, increment the value of the failure counter by one. After the traversal, calculate the ratio of the value of the failure counter to the total length of the displacement prediction value sequence, use the ratio as the failure probability value, and output the failure risk index of the structure during its service life. The preset displacement failure threshold is 0.05 meters.
[0098] In this embodiment, S55 specifically includes:
[0099] S551. Input the predicted structural response value into the physical information time integration module. Based on the time derivative decoupling mechanism, expand the predicted structural response value into a response sequence according to the time step. Set the element at the odd index position in the response sequence as the displacement prediction value to construct the displacement prediction value sequence, and set the element at the even index position as the velocity prediction value to construct the velocity prediction value sequence.
[0100] S552. Construct a chain for calculating the derivative of the displacement prediction value with respect to the time variable. Calculate the first derivative of the displacement prediction value with respect to the time variable as a physical constraint on the velocity prediction value. Calculate the product of the velocity prediction value and the preset time step to obtain the displacement increment prediction value. Add the displacement increment prediction value to the displacement prediction value at the previous moment to obtain the updated displacement value at the current moment.
[0101] S553. Calculate the second derivative of the displacement prediction value with respect to the time variable as the acceleration prediction value, obtain the mass matrix, damping matrix and stiffness matrix of the target building structure, calculate the matrix product of the mass matrix and the acceleration prediction value to obtain the inertial force vector, calculate the matrix product of the damping matrix and the velocity prediction value to obtain the damping force vector, and calculate the matrix product of the stiffness matrix and the displacement prediction value to obtain the restoring force vector.
[0102] S554. Calculate the sum of corresponding elements of the inertial force vector, damping force vector, and restoring force vector to obtain the structural internal force response vector. Read the refined wind pressure standard value as the external load vector. Calculate the difference between the external load vector and the structural internal force response vector to obtain the dynamic residual vector. Calculate the sum of the squares of each element in the dynamic residual vector and divide it by the total number of samples to obtain the physical constraint loss function value. Use the backpropagation algorithm to calculate the gradient of the physical constraint loss function value with respect to the weight matrix and bias vector in the improved PI-DeepONet model. Use the gradient descent algorithm to update the weight matrix and bias vector of the load feature branch network and the spatiotemporal feature backbone network.
[0103] This implementation achieves high-precision prediction of structural wind-induced response and quantification of failure risk by constructing an improved PI-DeepONet model that integrates data-driven and physical mechanisms. The refined standard wind pressure value is encoded by a load feature branch network and mapped to the predicted structural response value through a dot product operation with the time coordinate vector generated by the spatiotemporal feature backbone network. The core lies in the application of the physical information time integration module, which intelligently decomposes the response sequence into displacement and velocity sequences using a time derivative decoupling mechanism. It then calculates the inertial force, damping force, and restoring force vectors by combining mass, damping, and stiffness matrices, constructing a physical constraint loss function based on the dynamic equations. This mechanism forces the model prediction to follow Newton's laws of motion, effectively avoiding the "black box" divergence problem of purely data-driven models. Finally, the response statistics output module traverses the displacement sequence and compares it with a threshold to accurately calculate the failure probability. This invention significantly improves the physical interpretability and computational efficiency of structural response assessment under complex wind field environments, achieving millisecond-level risk indicator output.
[0104] The improved PI-DeepONet model of this invention is similar to the traditional PI-DeepONet model in that both retain the basic architecture of deep operator networks, namely, processing input function parameters through branch networks, processing spatiotemporal coordinate variables through the backbone network, and outputting structural response prediction values through dot product operations. At the same time, both introduce physical information constraints, using the residuals of physical equations to construct loss functions and updating network weights through backpropagation algorithms, thereby ensuring that the model predictions conform to physical laws.
[0105] The difference lies in that this invention breaks away from the conventional processing mode of traditional PI-DeepONet models, which typically use automatic differentiation techniques to calculate the derivative of the network output with respect to the input coordinates to construct the residuals of the physical equations. Instead, it creatively introduces a time derivative decoupling mechanism and explicit numerical integration logic into the physical information time integration module. In traditional models, physical constraints often rely on the direct differentiation of the network output, which can easily lead to gradient instability or accumulated errors when dealing with long-term or complex dynamic systems. In steps S551 and S552 of this invention, firstly, based on the time derivative decoupling mechanism, the structural response prediction values output by the model are explicitly split into displacement prediction sequences and velocity prediction sequences according to their index positions, constructing independent displacement and velocity solution spaces. Subsequently, the velocity physical constraint is constructed using the displacement-time derivative relationship, and the displacement increment is calculated by multiplying the velocity by the time step, thereby obtaining the updated displacement value. In steps S553 and S554, the mass, damping and stiffness matrices of the structure are further combined to explicitly calculate the inertial force, damping force and restoring force vectors, and a dynamic residual vector containing the difference between the external load and internal force response is constructed based on the dynamic equilibrium equation, thereby driving the network optimization.
[0106] Based on the aforementioned improvements, the beneficial effects of this invention lie in effectively solving the numerical instability and lack of physical consistency problems that traditional PI-DeepONet easily encounters when solving second-order dynamic equations through the time derivative decoupling mechanism and the embedding of an explicit integral structure. This design forces the displacement and velocity sequences output by the neural network to satisfy the temporal logical correlation at the calculus level, significantly reducing the optimization difficulty of the physical constraint loss function and enhancing the model's prediction accuracy for long-term dynamic responses. Simultaneously, the explicit matrix operations make the physical constraint process transparent and interpretable, enabling more accurate capture of the nonlinear dynamic characteristics of the structure under extreme wind loads, thus improving the robustness and reliability of failure risk assessment.
[0107] In this embodiment, S6 specifically includes:
[0108] S61. Collect the current real-time wind speed and wind direction angle values of the target area, construct real-time sample data pairs, and add them to the database of historical wind speed and wind direction data. Substitute the updated historical wind speed and wind direction data into the wind speed edge distribution function calculation process, update the shape parameter adjustment value and scale parameter adjustment value, and output the updated wind speed edge distribution function.
[0109] S62. Substitute the updated historical wind speed and direction data into the wind direction edge distribution function calculation process, output the updated wind direction edge distribution function, substitute the updated wind speed edge distribution function and the updated wind direction edge distribution function into the joint probability density function calculation process, update the correlation coefficient and optimal dependency parameter between wind speed and wind direction, and output the updated probability density data.
[0110] S63. Substitute the updated probability density data into the extreme wind speed calculation process, output the updated extreme wind speed data by wind direction, substitute the updated extreme wind speed dataset by wind direction into the refined wind pressure standard value calculation process, output the updated refined wind pressure standard value, input the updated refined wind pressure standard value into the improved PI-DeepONet model, and output the dynamically optimized failure risk index.
[0111] S64. Read the updated refined wind pressure standard value as the wind pressure load intensity, read the windward area of the target building structure, calculate the product of the refined wind pressure standard value and the windward area to obtain the concentrated wind load value, calculate the product of the concentrated wind load value and the dynamically optimized failure risk index, add the product result to the preset wind load partial factor, and output the dynamically optimized wind load calculation result. The preset wind load partial factor is 1.4.
[0112] Example 1: To verify the feasibility of this invention in the field of structural wind load calculation and failure risk assessment, it was applied to the full life-cycle health monitoring system of a super high-rise building in a city center. The area where this building is located is frequently hit by typhoons, with a design reference period of 50 years and a relatively high basic wind pressure setting. The building's exterior features a complex hyperboloid streamlined design, exhibiting significant uneven distribution of local wind pressure on its surface, and its first-order natural vibration period is relatively long, making it extremely sensitive to wind loads. The monitoring system integrates on-site meteorological stations, wind tunnel test databases, and a structural health monitoring sensor network. Traditional wind load calculation methods mainly rely on the constant value coefficient method according to the "Code for Design of Building Structures," and are verified using periodic wind tunnel test reports. However, traditional methods have revealed many limitations in practical applications: First, the shape coefficient and wind pressure height variation coefficient in the specifications are usually based on regular shapes or empirical formulas, which are difficult to accurately describe the unsteady wind pressure distribution characteristics caused by the complex aerodynamic shape of the building; Second, traditional extreme wind speed calculations often ignore the strong correlation between wind speed and wind direction, simply assuming that the two are independent, which leads to the underestimation of wind pressure extremes under the prevailing wind direction and overly conservative calculations under non-prevailing wind directions; Third, in structural dynamic response analysis, the traditional time history analysis method is extremely time-consuming and cannot meet the need for real-time early warning of structural status during typhoons, while the frequency domain analysis method is difficult to handle strongly nonlinear aerodynamic damping problems, resulting in significant biases in failure risk assessment.
[0113] In practical deployment, the method of this invention utilizes nearly 20 years of long-term meteorological observation data collected by an ultrasonic anemometer deployed on the roof of a building as a basic database. During the system initialization phase, historical wind speed and direction data are first cleaned and preprocessed to remove abnormal data points caused by instrument malfunctions during typhoons, and missing records are filled using linear interpolation. The marginal distributions of wind speed and direction are accurately fitted using Weibull and mixed von Mises distributions, respectively, and a joint probability density model of wind speed and direction is constructed using the Copula function. The system divides the entire wind direction space into 36 10-degree sectors, and identifies 18 effective sectors with an occurrence probability greater than 1% based on the joint probability density data. For the tail data of the wind speed distribution within each effective sector, the extreme wind speed of each sector within the design reference period is calculated using a generalized Pareto distribution model. After acquiring the extreme wind speed dataset, the system incorporates computational fluid dynamics simulation technology to simulate the wind pressure distribution on building surfaces under 36 wind direction angles. It calculates the shape coefficient for each sector and, combined with the wind pressure height variation coefficient and probability weights, outputs a refined standard wind pressure value that considers the characteristics of wind direction distribution. This standard wind pressure value is input into an improved PI-DeepONet model, which includes a load feature branch network and a spatiotemporal feature backbone network. Through a time derivative decoupling mechanism, the network's output structural response predictions are decomposed into displacement and velocity sequences. Combining the building structure's mass, stiffness, and damping matrices, a dynamic residual is constructed as a physical constraint loss function. The network weights are updated in real-time using a backpropagation algorithm, enabling rapid and accurate output of structural failure risk indicators without iteratively solving large finite element equation systems. Furthermore, the system integrates real-time meteorological monitoring data streams. Whenever real-time wind speed and direction data are updated, the system automatically corrects the edge distribution parameters and joint probability model, achieving self-evolutionary updates to the wind load calculation model. Table 1 below shows the comparative data of the method of this invention and the traditional standard method in the wind load calculation and failure risk assessment of this super high-rise building during a strong typhoon:
[0114] Table 1. Comparison of wind load calculation and risk assessment performance between the present invention and traditional methods.
[0115] Based on the comparative data shown in Table 1, it can be seen that the wind load calculation method based on the joint distribution of wind speed and wind direction proposed in this invention shows significant performance advantages over traditional standard methods in wind resistance analysis of super high-rise buildings. In particular, it has achieved comprehensive improvement in key indicators such as the accuracy of extreme wind speed estimation, the refinement of wind pressure calculation, the timeliness of risk warning, and the accuracy of displacement prediction.
[0116] In terms of extreme wind speed estimation and wind pressure calculation, this invention fully considers the coupling characteristics of wind speed and wind direction. Taking the downwind condition as an example, traditional methods, based on the independence assumption, ignore the strong correlation between the prevailing typhoon wind direction (southeast) and extreme wind speed in the region, resulting in a calculated extreme wind speed of only 38.5 m / s, significantly lower than the actual situation. This invention, however, quantifies the dependency relationship through a Copula function, estimating an extreme wind speed of 42.3 m / s in the southeast sector, correcting a calculation error of nearly 10%. Correspondingly, the refined wind pressure standard value is also increased from 1.15 to 1.28, effectively avoiding insufficient safety reserves in wind-resistant weak points of structural design, reflecting the core value of "refined" calculation. In crosswind vibration analysis, this invention, through sector-based weighted calculation, more accurately captures the wind direction range with significant vortex-induced shedding effects, significantly improving the agreement between the wind pressure standard value and actual wind tunnel test data.
[0117] In terms of failure risk warning and response timeliness, this invention demonstrates a revolutionary speed advantage. Traditional methods, when conducting dynamic time-history analysis of the entire typhoon process, often require building massive finite element models, with a single analysis taking 4 to 5 hours, unable to keep up with the rapid changes in the typhoon's path, resulting in severely delayed warning information. In contrast, this invention utilizes an improved PI-DeepONet model, embedding physical equations into a neural network, achieving millisecond-level inference speed. During the typhoon's passage, the system's average warning response time is only about 12 to 15 seconds, more than a thousand times faster than traditional methods, truly achieving "real-time" risk assessment. This allows operation and maintenance management teams ample time to adjust damper parameters or evacuate personnel before the typhoon makes landfall, greatly improving emergency management efficiency.
[0118] Regarding prediction accuracy and physical consistency, this invention effectively solves the problem of unstable "black box" prediction in pure data-driven models by embedding a time derivative decoupling mechanism and a physical constraint loss function. Data shows that in controlling peak displacement prediction errors, traditional time history analysis methods are affected by grid partitioning and iterative convergence, resulting in large error fluctuations, averaging around 18%; while the prediction error of this invention is stably controlled within 3%. Particularly noteworthy is the extremely low dynamic residual calculated by this invention (around 0.004), indicating that the model's predicted structural response strictly adheres to Newton's laws of motion, avoiding the physical law violations that often occur in traditional neural network models, significantly enhancing engineers' confidence in the early warning results. This invention also performs excellently in local nodal stress assessment, significantly improving the accuracy of failure risk probability assessment, effectively filling the technical gap in dynamic monitoring of local stress concentration areas using traditional methods, and providing solid technical support for wind-induced disaster prevention and control of large and complex engineering structures.
[0119] The above description is only a preferred embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any equivalent substitutions or modifications made by those skilled in the art within the scope of the technology disclosed in the present invention, based on the technical solution and inventive concept of the present invention, should be covered within the scope of protection of the present invention.
Claims
1. A method for calculating wind loads based on the combined distribution of wind speed and direction, characterized in that, Includes the following steps: S1. Collect historical wind speed and direction data and preprocess them. Extract wind speed and wind direction sequences as random variables, calculate statistical feature values, and fit the edge distribution function of wind speed and the edge distribution function of wind direction. S2. Using the edge distribution function of wind speed and the edge distribution function of wind direction, construct the joint probability density function of wind speed and wind direction, calculate the correlation coefficient between wind speed and wind direction and quantify the dependency relationship, and output the probability density data. S3. Divide the entire wind direction space into N sectors. Based on the probability density data, statistically analyze the wind speed data falling into each sector. Fit the tail data of the wind speed distribution in each sector to estimate the extreme wind speed and output the extreme wind speed dataset by wind direction. S4. Obtain the aerodynamic shape parameters of the target building structure, calculate the shape coefficient under different wind angles, and perform weighted calculation with the extreme wind speed dataset of wind direction, and output a refined wind pressure standard value by combining the wind pressure height variation coefficient. S5. The refined wind pressure standard value is used as the input load function to input the improved PI-DeepONet model to obtain the failure risk index of the structure during its service life. The improved PI-DeepONet model includes a load feature branch network, a spatiotemporal feature backbone network, a physical information time integration module, and a response statistical output module. The physical information time integration module embeds a time derivative decoupling mechanism. S6. Collect real-time wind speed and direction data, correct the parameters of the wind speed edge distribution function and the wind direction edge distribution function, calculate the corresponding wind direction extreme value wind speed dataset and refined wind pressure standard value, update the failure risk index of the structure during its service life, and output the dynamically optimized wind load calculation results. S5 specifically includes: S51. Construct an improved PI-DeepONet model that includes a load feature branch network, a spatiotemporal feature backbone network, a physical information time integration module, and a response statistical output module. Use the refined wind pressure standard value as the input variable of the function space. Construct a fully connected neural network containing an input layer, a hidden layer, and an output layer as the load feature branch network. Initialize the network weight matrix and bias vector. Perform a linear transformation on the refined wind pressure standard value to obtain the input of the first hidden layer and substitute it into the ReLU activation function to obtain the output of the first hidden layer. S52. Perform matrix multiplication between the output of the first hidden layer and the weight matrix of the next layer, and add the bias vector to obtain the input of the second hidden layer. Substitute the input of the second hidden layer into the Tanh activation function to obtain the output of the second hidden layer. Iterate in this way until the output layer is reached. The output layer outputs the load encoding vector. S53. Discretize the preset time domain according to the preset integration step size to generate a time coordinate sampling point sequence. Construct a multilayer perceptron containing an input layer, a hidden layer and an output layer as the spatiotemporal feature backbone network. Perform a linear transformation on the spatiotemporal coordinate sampling point sequence to obtain the hidden layer input and substitute it into the Sigmoid activation function for nonlinear mapping to obtain the hidden layer output. S54. Pass the hidden layer output to the next layer until the output layer, output the spatiotemporal coding vector, calculate the dot product of the load coding vector and the spatiotemporal coding vector, and obtain the scalar form of the structural response prediction value. S55. Input the predicted structural response values into the physical information time integration module. Based on the time derivative decoupling mechanism, split the predicted structural response values into a sequence of predicted displacement values and a sequence of predicted velocity values, and calculate the internal force response vector and physical constraint loss function value. Based on the physical constraint loss function value, use the gradient descent algorithm to update the weight matrix and bias vector of the load feature branch network and the spatiotemporal feature backbone network. S56. Input the displacement prediction value sequence output by the improved PI-DeepONet model after updating the weight parameters into the response statistics output module, read the preset displacement failure threshold, traverse each value in the displacement prediction value sequence, calculate the absolute value of the value and compare it with the preset displacement failure threshold. When it is greater than the preset displacement failure threshold, increment the value of the failure counter by one. After the traversal is completed, calculate the ratio of the value of the failure counter to the total length of the displacement prediction value sequence, use the ratio as the failure probability value, and output the failure risk index of the structure during its service life.
2. The wind load calculation method based on the combined distribution of wind speed and direction according to claim 1, characterized in that, S1 specifically includes: S11. Export historical wind speed and direction records from the database of long-term meteorological observation stations in the target area, remove abnormal data points with values of zero or exceeding the preset physical range threshold, and perform linear interpolation to fill in the gaps to obtain preprocessed historical wind speed and direction data. S12. The preprocessed historical wind speed and direction data are split into two data sequences in chronological order. The wind speed value column is extracted as the wind speed sequence, and the wind direction angle column is extracted as the wind direction sequence. The mean and variance of the wind speed sequence and the mean and standard deviation of the wind direction sequence are calculated respectively, and the statistical feature values are output. S13. Substitute the mean wind speed and standard deviation of wind speed in the statistical characteristic values into the moment estimation formula, calculate the ratio of the mean wind speed to the standard deviation of wind speed, and refer to the preset comparison table to obtain the initial value of the shape parameter. Substitute the mean wind speed into the algebraic expression containing the initial value of the shape parameter to calculate the initial value of the scale parameter. Substitute the initial value of the shape parameter, the initial value of the scale parameter, and the wind speed sequence into the log-likelihood function, calculate the log probability value of each value in the wind speed sequence, and sum all the log probability values to obtain the log-likelihood function value. S14. Iteratively adjust the initial values of the shape parameter and the scale parameter using the gradient ascent algorithm until the log-likelihood function value no longer increases. Substitute the adjusted parameters into the Weibull distribution probability density formula, calculate the shape parameter adjustment value raised to the power of the scale parameter adjustment value, divide the wind speed sequence value by the scale parameter adjustment value and calculate the corresponding shape parameter adjustment value minus one. Multiply the two calculation results by a negative exponent term with the natural constant as the base and the wind speed sequence value divided by the scale parameter adjustment value as the exponent to obtain the wind speed marginal distribution function. S15. Substitute the mean wind direction value in the statistical feature value into the mixed von Mises distribution formula with a preset number of components, calculate the initial values of the mixing ratio, mean direction, and concentration parameter of each component, and substitute them with the wind direction sequence into the log-likelihood function formula for cumulative calculation. Execute the expectation step in the expectation maximization algorithm to calculate the posterior probability of the latent variables, and execute the maximization step to update the initial values of the mixing ratio, mean direction, and concentration parameter using the posterior probability until the change in the log-likelihood function value is less than the preset convergence threshold. S16. Substitute the updated parameters into the mixed von Mises distribution probability density formula, multiply the updated mixing ratio value by the exponential function value with the product of the natural constant, the concentration parameter update value and the cosine function value as the exponent, to obtain the component probability value. Sum the component probability values over the preset number of components to obtain the wind direction edge distribution function.
3. The wind load calculation method based on the combined distribution of wind speed and direction according to claim 1, characterized in that, S2 specifically includes: S21. Substitute the cumulative probability value of wind speed output by the wind speed marginal distribution function and the cumulative probability value of wind direction output by the wind direction marginal distribution function into the expression of the preset Copula function model, calculate the natural logarithm of the product of the cumulative probability value of wind speed and the cumulative probability value of wind direction and square it, and sum the squared values of all samples to obtain the goodness-of-fit statistic. S22. Substitute the goodness-of-fit statistic into the correlation coefficient calculation formula, divide the goodness-of-fit statistic by the difference between the total number of samples minus one, take the square root and take the negative value to obtain the correlation coefficient between wind speed and wind direction. S23. Substitute the wind speed sequence, wind direction sequence, and correlation coefficient into the likelihood function expression, calculate the probability value of each wind speed sample and wind direction sample occurring simultaneously, and sum the natural logarithms of all probability values to obtain the likelihood function value. Set the step size search range, select multiple candidate dependency parameter values within the search range, substitute each candidate dependency parameter value into the likelihood function expression to calculate the corresponding likelihood function value, compare and select the candidate dependency parameter value with the largest likelihood function value as the optimal dependency parameter, and quantify the dependency relationship between wind speed and wind direction. S24. Substitute the optimal dependency parameters into the expression of the preset Copula function model, calculate the bivariate joint probability value of the cumulative probability value of wind speed and the cumulative probability value of wind direction, calculate the partial derivatives of the bivariate joint probability value with respect to the wind speed and wind direction variables respectively, multiply the two partial derivative calculation results to obtain the joint probability density function of wind speed and wind direction, arrange the discrete point values of the joint probability density function on the wind speed axis and the wind direction axis in a grid, and output the probability density data.
4. The wind load calculation method based on the combined distribution of wind speed and direction according to claim 1, characterized in that, S3 specifically includes: S31. Starting from the due north direction, divide the entire wind direction space clockwise according to the preset angle interval to obtain N sectors. Read the probability density data, accumulate and normalize the probability density data falling into each sector to obtain the occurrence probability of each sector, and select the sector with the occurrence probability greater than the preset probability threshold as the effective sector. S32. Extract wind speed samples falling into each effective sector to construct a sub-sample set. Sort the wind speed data in the sub-sample set in ascending order of numerical value. Select the wind speed value at a preset proportion position in the sorting sequence as the threshold. Extract the wind speed data greater than the threshold as the tail data of the wind speed distribution. S33. Substitute the tail data of the wind speed distribution into the log-likelihood function of the generalized Pareto distribution, calculate the partial derivatives of the shape parameter and the scale parameter and set the partial derivatives to zero, solve the system of equations to obtain the shape parameter and the scale parameter, substitute the tail data of the wind speed distribution, the shape parameter and the scale parameter into the cumulative distribution function formula of the generalized Pareto distribution, calculate the power term with the product of the scale parameter plus the difference between the shape parameter and the tail data of the wind speed distribution minus the threshold as the denominator and the product of the shape parameter and the difference between the tail data of the wind speed distribution minus the threshold as the numerator, add one to the negative power term to obtain the tail model of the wind speed distribution for each effective sector; S34. Obtain the design reference period of the meteorological observation station, calculate the ratio of the total number of hours in the design reference period to the total number of tail data of the wind speed distribution, substitute the ratio as the exceedance probability into the inverse function of the tail model of the wind speed distribution, and calculate the extreme wind speed of each effective sector in the design reference period. S35. Summarize the extreme wind speeds of all valid sectors, set the extreme wind speeds of non-valid sectors to zero, arrange them in order of sector number, and output the extreme wind speed dataset by wind direction.
5. The wind load calculation method based on the combined distribution of wind speed and direction according to claim 1, characterized in that, S4 specifically includes: S41. Obtain the building geometric dimension parameters of the target building structure, establish a three-dimensional geometric model, place the three-dimensional geometric model in the wind field of the computational fluid dynamics software, adjust the inlet boundary conditions according to the preset wind direction angle interval, simulate the airflow process, extract the average wind pressure value of the building surface and the reference wind pressure value of the far reference point, calculate the quotient of the average wind pressure value divided by the reference wind pressure value, and obtain the shape coefficient under different wind direction angles. S42. Read the extreme wind speed dataset by wind direction, square the extreme wind speed in the extreme wind speed dataset, multiply the squared value of the extreme wind speed by the corresponding body shape coefficient, and multiply the product by half of the air density value to obtain the extreme wind pressure of each sector. S43. Obtain the ground roughness index of the target building structure location and the height value of the calculated location. Divide the height value of the calculated location by the preset standard height value to obtain the height ratio. Raise the height ratio to the power of the ground roughness index. Multiply the result of the power operation by the wind pressure coefficient at the preset standard height to obtain the wind pressure height variation coefficient. The preset standard height value is 10 meters and the wind pressure coefficient at the preset standard height is 1.
0. S44. Read the probability density data, use the values in the probability density data as weights, multiply the extreme wind pressure, wind pressure height variation coefficient and weights to obtain the weighted wind pressure value, sum the weighted wind pressure of all effective sectors and divide by the total weights to output the refined wind pressure standard value.
6. The wind load calculation method based on the combined distribution of wind speed and direction according to claim 1, characterized in that, Specifically, S55 includes: S551. Input the predicted structural response value into the physical information time integration module. Based on the time derivative decoupling mechanism, expand the predicted structural response value into a response sequence according to the time step. Set the element at the odd index position in the response sequence as the displacement prediction value to construct the displacement prediction value sequence, and set the element at the even index position as the velocity prediction value to construct the velocity prediction value sequence. S552. Construct a chain for calculating the derivative of the displacement prediction value with respect to the time variable. Calculate the first derivative of the displacement prediction value with respect to the time variable as a physical constraint on the velocity prediction value. Calculate the product of the velocity prediction value and the preset time step to obtain the displacement increment prediction value. Add the displacement increment prediction value to the displacement prediction value at the previous moment to obtain the updated displacement value at the current moment. S553. Calculate the second derivative of the displacement prediction value with respect to the time variable as the acceleration prediction value, obtain the mass matrix, damping matrix and stiffness matrix of the target building structure, calculate the matrix product of the mass matrix and the acceleration prediction value to obtain the inertial force vector, calculate the matrix product of the damping matrix and the velocity prediction value to obtain the damping force vector, and calculate the matrix product of the stiffness matrix and the displacement prediction value to obtain the restoring force vector. S554. Calculate the sum of corresponding elements of the inertial force vector, damping force vector, and restoring force vector to obtain the structural internal force response vector. Read the refined wind pressure standard value as the external load vector. Calculate the difference between the external load vector and the structural internal force response vector to obtain the dynamic residual vector. Calculate the sum of the squares of each element in the dynamic residual vector and divide it by the total number of samples to obtain the physical constraint loss function value. Use the backpropagation algorithm to calculate the gradient of the physical constraint loss function value with respect to the weight matrix and bias vector in the improved PI-DeepONet model. Use the gradient descent algorithm to update the weight matrix and bias vector of the load feature branch network and the spatiotemporal feature backbone network.
7. The wind load calculation method based on the combined distribution of wind speed and direction according to claim 1, characterized in that, S6 specifically includes: S61. Collect the current real-time wind speed and wind direction angle values of the target area, construct real-time sample data pairs, and add them to the database of historical wind speed and wind direction data. Substitute the updated historical wind speed and wind direction data into the wind speed edge distribution function calculation process, update the shape parameter adjustment value and scale parameter adjustment value, and output the updated wind speed edge distribution function. S62. Substitute the updated historical wind speed and direction data into the wind direction edge distribution function calculation process, output the updated wind direction edge distribution function, substitute the updated wind speed edge distribution function and the updated wind direction edge distribution function into the joint probability density function calculation process, update the correlation coefficient and optimal dependency parameter between wind speed and wind direction, and output the updated probability density data. S63. Substitute the updated probability density data into the extreme wind speed calculation process, output the updated extreme wind speed data by wind direction, substitute the updated extreme wind speed dataset by wind direction into the refined wind pressure standard value calculation process, output the updated refined wind pressure standard value, input the updated refined wind pressure standard value into the improved PI-DeepONet model, and output the dynamically optimized failure risk index. S64. Read the updated refined wind pressure standard value as the wind pressure load intensity, read the windward area of the target building structure, calculate the product of the refined wind pressure standard value and the windward area to obtain the concentrated wind load value, calculate the product of the concentrated wind load value and the dynamically optimized failure risk index, add the product result to the preset wind load partial factor, and output the dynamically optimized wind load calculation result. The preset wind load partial factor is 1.4.
Citation Information
Patent Citations
Turnout rigidity degradation limit value determination method and system based on time-varying reliability
CN122197279A
Four-footed robot dog sliding mode tracking control method and device for cable inspection, terminal equipment and computer readable storage medium
CN122569402A
