A method and system for evaluating wind-solar resource complementarity based on physical embedded deep learning

Abstract

This invention discloses a method and system for evaluating the complementarity of wind and solar resources based on physically embedded deep learning, belonging to the field of resource assessment. The method includes: acquiring meteorological elements and digital elevation model data for the target area; calculating the elevation change rate between adjacent grids, distinguishing between flat and non-flat areas based on the IGC algorithm, and extracting terrain features separately; generating periodic features of solar altitude angle through radiation data, and concatenating them with terrain features to obtain fused features; constructing a physically embedded NODE-SA-ConvLSTM model, combining a self-attention mechanism to capture the global dependencies of wind and solar elements, and outputting wind speed and radiation prediction results; using a constant error compensation mechanism to perform sliding balance correction on the prediction residuals; and quantifying the complementarity of regional wind and solar resources through a wind-solar complementarity coefficient. This invention, by introducing a physically embedded mechanism and continuous-time modeling, effectively improves the accuracy and stability of wind and solar resource prediction, realizing the complementary evaluation of wind and solar resource endowments.

Description

Technical Field

[0001] This invention relates to the field of new energy resource assessment and intelligent modeling technology, and in particular to a wind and solar resource complementarity assessment method and system based on physical embedded deep learning. Background Technology

[0002] Wind and solar energy, due to their clean, renewable, and sustainable characteristics, have become key supports for energy transition and green, low-carbon development. Wind and solar resources exhibit significant spatiotemporal complementarity: during periods of ample sunshine and weak winds, solar energy can compensate for insufficient wind power; while at night or during cloudy or rainy weather, wind energy can, to some extent, substitute for solar energy, achieving a dynamic balance in energy supply. However, influenced by factors such as terrain complexity, meteorological disturbances, and seasonal climate changes, the spatiotemporal distribution and energy conversion efficiency of wind and solar resources vary significantly, leading to substantial differences in the accuracy of resource endowment assessment and prediction across different regions. This uncertainty is particularly pronounced in areas with dramatic topographic relief, such as mountainous and basin areas.

[0003] Currently, wind and solar resource assessment methods mainly fall into three categories: physical modeling, statistical analysis, and data-driven approaches. Physical modeling methods (such as WRF models and CFD numerical simulations) are based on dynamic equations and energy balance equations, and can reflect large-scale meteorological processes well. However, they are highly sensitive to the accuracy of initial boundary conditions and topographic input, and have extremely high computational complexity, making them unsuitable for large-scale, long-term dynamic assessments. Statistical analysis methods (such as Weibull distribution fitting for wind speed and gamma distribution fitting for radiation intensity) rely on sample distribution assumptions and can achieve rapid estimation, but they struggle to capture the coupling relationship between topographic nonlinear effects and meteorological variables, and are insufficiently adaptable to non-stationary and multi-scale characteristics. Data-driven methods, especially deep learning models such as LSTM and ConvLSTM widely used in recent years, significantly improve short-term prediction accuracy by learning nonlinear mapping relationships through samples; however, their discrete-time step-based state transfer mechanism cannot characterize continuous evolution characteristics, and they are prone to information loss and prediction lag under complex terrain, strong disturbances, and non-stationary time series. Furthermore, most traditional deep models lack physical interpretability and struggle to reveal the physical coupling mechanism between wind speed and radiation.

[0004] On the other hand, the original meteorological sequences of wind and solar resources often simultaneously contain periodic variations (such as diurnal cycles and seasonal variations), local disturbances (such as orographic winds and cumulus albedo), and observational noise, exhibiting typical multi-scale non-stationary characteristics. Single frequency domain analysis methods (such as Fourier transform and wavelet transform) struggle to extract features from both stationary and non-stationary signals, resulting in insufficient time-frequency information representation and affecting the accuracy of subsequent modeling. Wind and solar signals under complex terrain also exhibit significant spatial heterogeneity, with obvious differences in wind speed and radiation characteristics between adjacent areas. Traditional convolutional neural networks often demonstrate insufficient generalization ability when processing such spatially discontinuous features.

[0005] Furthermore, deep learning models suffer from error accumulation and gradient vanishing during long-term iterative training, leading to drift or oscillations in prediction results over time. While some improved methods (such as attention mechanisms and residual structures) can alleviate information loss, they still struggle to accurately model continuous spatiotemporal processes. Particularly in wind and solar resource assessment scenarios, wind speed and solar radiation exhibit significant dynamic coupling, making it difficult for data-driven discrete update models to reveal their underlying physical laws.

[0006] To address the aforementioned issues, a hybrid modeling framework that integrates physical mechanisms with the feature representation capabilities of deep learning is needed, in order to provide a technical system that combines physical interpretability, continuous-time modeling capabilities, and multi-scale time-frequency fusion feature extraction. Summary of the Invention

[0007] To address the problems existing in the prior art, this invention provides a method and system for evaluating the complementarity of wind and solar resources based on physical embedded deep learning. By introducing a physical embedding mechanism and continuous-time modeling, the accuracy and stability of wind and solar resource prediction are effectively improved, and the complementary evaluation of wind and solar resource endowments is realized.

[0008] This invention adopts the following technical solution: a method for evaluating wind and solar resource complementarity based on physically embedded deep learning, comprising the following steps:

[0009] S1. Obtain gridded meteorological data and gridded terrain data for the target area;

[0010] S2. Calculate the elevation change rate between adjacent grids in the gridded terrain data, determine the elevation change rate threshold using the IGC (Iterative Gradient Clustering) algorithm, and distinguish between flat and non-flat areas based on the elevation change rate threshold;

[0011] S3. Perform S-Transform time-frequency analysis on meteorological data in flat areas to extract the instantaneous energy spectrum and main frequency characteristics of core meteorological elements. For meteorological data in non-flat areas, use a combination of VMD and Hilbert transform to extract the instantaneous amplitude and frequency characteristics of core meteorological elements.

[0012] S4. Encode the radiation data in the meteorological data to generate the periodic feature of the solar altitude angle, and then concatenate it with the main frequency feature and the frequency feature to obtain the fused feature;

[0013] S5. Construct a physical embedded NODE-SA-ConvLSTM (Spatio-Temporal Attention Convolutional Long Short-Term Memory Network) model, introduce the NODE (Neural Ordinary Differential Equation) mechanism into SA-ConvLSTM, transform the discrete hidden state transfer process into the form of solving differential equations in continuous time, and input the fused features into the model;

[0014] S6. Output the prediction results of wind speed and solar radiation from the NODE-SA-ConvLSTM model;

[0015] S7. Use the error constant compensation mechanism to perform residual balance correction on the prediction results;

[0016] S8. Construct the WSCC (Wind-Solar Complementarity Coefficient) according to the corrected prediction results, and output the complementary evaluation results of wind and solar resource endowments.

[0017] Preferably, step S2 includes the following sub-steps:

[0018] S2.1. Process the elevation height elements of adjacent grids in the target area. By performing gradient operations on the elevation height difference Δh of adjacent grid cells, obtain the elevation change rate .

[0019] S2.2. Use the IGC algorithm to obtain the elevation change rate threshold. Through clustering and iterative optimization of the grid terrain gradient data, obtain the converged clustering centers and , and set the threshold D of the elevation change rate as:

[0020] .

[0021] S2.3. When G i < D, determine that the area where the adjacent grids of the i-th group are located is a flat area. When G i > D, determine that the area where the adjacent grids of the i-th group are located is a non-flat area.

[0022] Preferably, in step S3, time-frequency feature extraction is performed on the flat area and the non-flat area respectively, and the method is as follows:

[0023] For the meteorological sequence in the flat area, use the time-frequency analysis algorithm to extract the instantaneous energy spectrum and main frequency features, and obtain the flat area feature vector set , where t is the time variable in the time series;

[0024] The meteorological sequence in the non-flat region was decomposed using the VMD-Hilbert algorithm. VMD was used to decompose the original signal into several eigenmode components with different frequency band characteristics. The decomposed mode results were then input into the Hilbert transform to obtain their instantaneous energy and phase change characteristics, resulting in a feature vector set for the non-flat region. ;

[0025] Align the two sets of feature vectors according to their time indices to form time-frequency feature datasets for two regions.

[0026] Preferably, step S4 includes the following sub-steps:

[0027] S4.1 Calculate the solar altitude angle, solar zenith angle, and solar azimuth angle based on the observation time and geographic coordinates. Construct a radiation periodic feature code based on this, including the diurnal periodic variation of the solar altitude angle, the time-period symmetrical variation of the azimuth angle, and the cosine correction term of the zenith angle. Combine the coding results with the radiation variables in meteorological elements to form a radiation prior feature set containing solar geometric parameters and radiation characteristics. .

[0028] S4.2, respectively for , , Perform Z-score normalization to obtain the standardized feature set. , , .

[0029] S4.3 Under the condition of consistent time index, the time-frequency features of flat areas, time-frequency features of non-flat areas, radiation prior features, and wind and solar element vector E(t) are weighted and fused as four independent scales to obtain multi-scale time-frequency feature fusion data F(t).

[0030] Preferably, in step S5, an SA-ConvLSTM model is constructed, and the NODE (Neural Ordinary Differential Equation) mechanism is introduced into the SA-ConvLSTM to transform its originally discrete hidden state transmission process into a continuous-time differential equation solution.

[0031] The SA-ConvLSTM model consists of 1 layer, each layer comprising an input gate i(t), a forget gate f(t), an output gate o(t), a hidden state h(t), a cell state c(t), and a candidate cell state g(t). The three gating units i(t), f(t), and o(t) jointly regulate the flow and update of information in h(t). After introducing NODE, these gating mechanisms are processed continuously, so that the evolution of the cell state c(t) can be described by ordinary differential equations, thereby achieving smooth modeling of temporal features.

[0032] Preferably, in step S6, the predicted results of wind speed and solar radiation are output from the NODE-SA-ConvLSTM model, including: the meridional component U10 and the zonal component V10 of wind speed, as well as the total solar radiation SSRD and the direct solar radiation FDIR, forming a four-dimensional meteorological element output matrix Y(t).

[0033] Preferably, step S7 includes the following sub-steps:

[0034] S7.1 Design a constant error compensation mechanism to calculate the prediction result y at each time step during model training and prediction. t Compared with measured values The residuals between .

[0035] Set the sliding window length k, and perform a moving average on the residuals over k time steps within the sliding window to obtain the dynamic residual mean. .

[0036] When the change in the mean residual in two consecutive iterations satisfies When, the constant compensation unit is triggered, where, The preset error fluctuation threshold is used;

[0037] S7.2, Constant Compensation Unit After inverse normalization, the result is superimposed on the prediction output of S6 and then normalized again to correct the prediction output.

[0038] Preferably, in step S8, WSCC is used to quantify the complementarity between wind speed and solar radiation, and Spearman rank correlation coefficient ρ is calculated for the time series of wind speed components U and V and solar radiation components SSRD and FDIR, respectively. s Kendall rank correlation coefficient T k The wind-solar complementarity coefficient (WSCC) is obtained by combining two correlation measures using a parameter adaptive weighted fusion method.

[0039] The present invention also provides: a wind-solar resource complementarity assessment system based on physical embedded deep learning, used to implement the above-mentioned wind-solar resource complementarity assessment method, comprising:

[0040] Terrain data processing module: Receives terrain data of the target area, calculates the gradient of the elevation difference between adjacent grid cells, and obtains a gridded terrain gradient field;

[0041] Gradient Iterative Clustering Algorithm Module: Based on terrain gradient data, it performs an iterative clustering process. By setting initial cluster centers and weight coefficients, it calculates the Euclidean distance between sample points and cluster centers to determine the region category. After each iteration, it updates the cluster centers until convergence, obtains the optimal threshold D, and distinguishes between flat and non-flat regions.

[0042] S-Transform Feature Extraction Module: Performs S-Transform time-frequency transformation on meteorological sequences in flat areas to extract the instantaneous energy spectrum and dominant frequency features of wind speed and radiation sequences, forming a stationary time-frequency feature matrix, and obtains the time-varying characteristics of meteorological elements in flat areas in a time-frequency adaptive manner;

[0043] VMD-Hilbert feature extraction module: Performs VMD decomposition on meteorological sequences in non-flat areas to obtain multiple narrowband intrinsic mode functions, and then obtains the instantaneous amplitude and instantaneous frequency of each mode through Hilbert transformation, constructing a non-stationary time-frequency feature matrix to reflect the non-stationary evolution law of meteorological signals under complex terrain;

[0044] NODE-SA-ConvLSTM module: Constructs NODE-SA-ConvLSTM model, inputs fused features, models spatial dependencies through self-attention mechanism, and describes the continuous time evolution of hidden state and cell state in the form of ordinary differential equations, outputting the prediction results of wind speed component and radiation.

[0045] Error constant compensation module: During the model training and prediction phase, the prediction residuals at each time step are subjected to sliding statistics, the dynamic residual mean is calculated, and compensation is triggered when the threshold condition is met. The compensation unit performs inverse normalization superposition and re-normalization correction on the residuals to achieve error constant control and prediction drift suppression.

[0046] Wind-Solar Complementarity Assessment Module: Calculates the wind-solar complementarity coefficient (WSCC), obtains the Spearman rank correlation coefficient and Kendall rank correlation coefficient for wind speed and radiation components respectively, and obtains the comprehensive complementarity index through adaptive weighted fusion to quantify the complementary relationship between wind and solar elements and output the assessment results.

[0047] Compared with the prior art, the present invention, employing the above technical solution, has the following technical effects:

[0048] 1. This invention proposes a wind-solar resource complementarity assessment framework based on physical embedded deep learning. It introduces the Neuron Regular Differential Equation (NODE) ​​into a spatiotemporal attention convolutional long short-term memory network (SA-ConvLSTM), transforming the traditional discrete time-transfer process into a continuous evolution modeling process. This enables the model to have continuous solution capability in the time dimension, thereby more accurately depicting the dynamic changes in wind speed and solar radiation, and improving the temporal consistency and physical interpretability of the prediction.

[0049] 2. This invention extracts stationary and non-stationary time-frequency features of meteorological sequences by combining S-Transform and VMD-Hilbert algorithms, and introduces periodic feature encoding of solar altitude angle and azimuth angle to make the input features more consistent with the natural radiation variation law. In the process of multi-scale feature fusion, the model's ability to perceive diurnal and seasonal changes is enhanced, and high-precision spatiotemporal fusion modeling of wind and light elements is realized.

[0050] 3. In the model output stage, the present invention designs a constant error compensation mechanism to perform moving average and dynamic balance correction on the prediction residual. When the error fluctuation is lower than the threshold, the compensation unit is automatically triggered to achieve constant control of the prediction error, effectively eliminate prediction drift caused by complex terrain or data non-stationarity, and significantly improve the convergence stability and robustness of the model.

[0051] 4. This invention constructs a wind-solar complementarity coefficient (WSCC), which quantitatively characterizes the complementary relationship between wind speed and solar radiation by adaptively fusing Spearman and Kendall correlation coefficients, thereby realizing a spatial and quantitative assessment of the complementary potential of wind and solar resources and providing a scientific basis for the coordinated allocation of new energy sources and the optimized utilization of wind-solar complementarity. Attached Figure Description

[0052] Figure 1 Flowchart of the wind-solar resource complementarity assessment method of this invention;

[0053] Figure 2 This is a diagram of the overall network model structure of the wind-solar resource complementarity assessment method of the present invention;

[0054] Figure 3 The flat region features, non-flat region features, and fused feature map are shown in the method of this invention.

[0055] Figure 4 This is a network model structure diagram of the physical embedding module of the method of the present invention;

[0056] Figure 5 This is a line graph comparing the prediction results, correction results, and actual values ​​of each variable at the 100th time step of the method of the present invention.

[0057] Figure 6 This is a visual comparison chart of the prediction results of the method of the present invention, the prediction results of the comparison method, and the actual values;

[0058] Figure 7 This is a line graph comparing the prediction results of the method of the present invention with the prediction results of the comparison method and the actual values. Detailed Implementation

[0059] To make the objectives, technical solutions, and advantages of this invention clearer, the technical solutions of the application will be further described in detail below with reference to the accompanying drawings. The described embodiments are only a part of the embodiments involved in this invention. All non-innovative embodiments based on these embodiments by other researchers in the art are within the protection scope of this invention. Furthermore, the step numbers in the embodiments of this invention are only set for ease of explanation and do not limit the order of the steps. The execution order of each step in the embodiments can be adaptively adjusted according to the understanding of those skilled in the art.

[0060] In one embodiment of the present invention, a wind-solar resource complementarity assessment method based on physically embedded deep learning is provided, such as... Figure 1 As shown, it includes the following steps:

[0061] Step 1: Select ERA5 wind and solar data and DEM data within the target area.

[0062] Step 1.1: ERA5 is the fifth-generation global atmospheric reanalysis product released by the European Centre for Medium-Range Weather Forecasts (ECMWF), which has high spatiotemporal resolution and physical consistency. This dataset is based on the four-dimensional variational assimilation (4D-Var) method, which integrates satellite remote sensing observations, ground station data, and numerical weather prediction results, and can simultaneously reflect large-scale circulation structures and local meteorological changes.

[0063] In this embodiment, the selected ERA5 meteorological elements include the longitude and zonal components (U10, V10) of the 10 m surface wind speed, as well as the surface incident shortwave radiation (SSRD) and direct radiation (FDIR), which can comprehensively describe the spatiotemporal distribution characteristics of wind and solar energy. The temporal resolution of the data is 1 hour, and the spatial resolution is 0.25°×0.25°.

[0064] Step 1.2: Digital Elevation Model (DEM) data, as geomorphic feature data, provides precise height information of the Earth's surface. By analyzing the impact of different terrain features on wind speed distribution and photovoltaic radiation, the accuracy of wind speed prediction and photovoltaic power generation prediction can be effectively enhanced.

[0065] In this embodiment, the selected DEM data has a spatial resolution of 30m.

[0066] Step 2: Calculate the elevation change rate of adjacent grids in the gridded terrain data, determine the change rate threshold using the IGCA algorithm, and distinguish between flat and non-flat areas based on the threshold.

[0067] Step 2.1: Process the elevation features in the DEM data. Obtain the elevation change rate G by performing gradient calculations on the elevation difference Δh between adjacent grid cells.

[0068] ;

[0069] where d is the adjacent grid spacing, and x and y are the number of rows and columns of the grid cell respectively, represents the altitude value of the grid coordinates (x, y); represents the horizontal distance between adjacent grid cells, represents the distance between adjacent grid coordinates.

[0070] For the adjacent grid spacing d of the i-th group i , calculate the corresponding altitude change rate G i :

[0071] .

[0072] Step 2.2: Adopt the IGC algorithm, and through clustering iteration optimization of the grid terrain gradient data, determine the optimal threshold. Set the initial clustering centers , :

[0073] ;

[0074] where min(G) and max(G) are the maximum and minimum values of G respectively, , are the weight coefficients of the two clustering centers respectively.

[0075] Based on the initial clustering centers for clustering, determine the range to which the sample points belong by calculating the Euclidean distance between each sample point and the clustering centers. After the division, recalculate the means of the two regions, and set the means as the new clustering centers, and so on. Iterate n rounds until the clustering centers no longer change, and obtain and .

[0076] Set the threshold D as:

[0077] ;

[0078] According to the threshold D, when G i < D, it is determined as a flat area, and when G i > D, it is determined as a non-flat area.

[0079] Step 3: According to the different terrain features, extract the stationary and non-stationary features of the meteorological sequences in the flat area and non-flat area respectively.

[0080] Step 3.1: The wind speed and radiation in the flat area change relatively smoothly. Adopt the S-Transform algorithm to extract its instantaneous energy spectrum and dominant frequency features.

[0081] S-Transform is a time-frequency analysis method developed based on short-time Fourier transform and wavelet transform. It has time-frequency adaptability and is defined as follows:

[0082] ;

[0083] in, Let x(t) be the time variable in the time series, x(t) be the time signal, and b be the time shift. f For frequency, This represents the Fourier transform.

[0084] w(tb,f) is a Gaussian window function, expressed as:

[0085] ;

[0086] The magnitude S(b,f) of the S-Transform characterizes the energy distribution of the signal in the time-frequency domain, representing the instantaneous energy spectrum. The calculation formula is as follows:

[0087] ;

[0088] The dominant frequency characteristic is obtained by energy-weighted averaging:

[0089] ;

[0090] The final output time-frequency feature matrix is:

[0091] .

[0092] Step 3.2: Due to the complex terrain and strong non-stationarity of signals in non-flat areas, a method combining VMD and Hilbert transform is used to extract instantaneous amplitude and instantaneous frequency features.

[0093] VMD decomposes the signal into several narrowband intrinsic mode functions (IMFs) by solving a variational problem, with each IMF having a center frequency. The definition is as follows:

[0094] ;

[0095] The constraints are satisfied:

[0096] ;

[0097] in, For the k-th modal component, * denotes convolution operation. Represents the unit impulse function. It is the time derivative; The transform kernel represents the singular kernel function required for signal analyticalization; For center frequency The demodulation factor, based on the standard, converts the IMF component into a baseband signal;

[0098] A set of stationary sub-signals is obtained by solving the problem using the alternating direction multiplier method (ADMM). {u k (t)} Each modal component is transformed by Hilbert to obtain an analytical signal:

[0099] ;

[0100] in, [] represents the Hilbert operator.

[0101] The instantaneous amplitude is obtained from this. and instantaneous frequency :

[0102] ;

[0103] Finally, the time-frequency characteristic matrix of the non-flat region is obtained:

[0104] ;

[0105] Step 3.3: Calculate the mean and variance of the two sets of feature vectors:

[0106] ;

[0107] ;

[0108] in, μ flat and μ rough These are the means of the two feature vector sets, respectively. σ flat and σ rought These are the standard deviations of the two feature vector sets, respectively.

[0109] Z-Score standardization is performed on the two sets of feature vectors:

[0110] ;

[0111] Obtain the standardized feature set F' flat (t) and F' rough (t).

[0112] Step 4: Standardize and dimension-match the time-frequency characteristics of flat and non-flat areas, and integrate the wind speed components U10 and V10 with solar radiation SSRD and FDIR from the ERA5 reanalysis data.

[0113] Before fusion, solar elevation angle was introduced into the solar radiation-related variables. With azimuth The periodic feature encoding yields a radiation prior feature set containing solar geometric parameters and radiation characteristics. :

[0114] ;

[0115] Normalization yields the standardized radiation feature set. ,and , , Multi-scale time-frequency feature fusion data is obtained after fusion. :

[0116] ;

[0117] in, This represents the landscape element vector in ERA5 data. , , , For adjustable weight parameters, satisfying .

[0118] In this embodiment, the features of flat regions, the features of non-flat regions, and the fused feature map are as follows: Figure 3 As shown.

[0119] Step 5: Construct a physically embedded NODE-SA-ConvLSTM model, such as... Figure 4 As shown, the fused features are input into the model.

[0120] Step 5.1: Construct the SA-ConvLSTM model, where the SA module acts as an attention memory unit to model global spatial dependencies. This process can be represented as follows:

[0121] ;

[0122] in, Functions represent memory modules. This represents the input feature tensor at the current time step. It is a feature aggregated through memory modules. The hidden state one step above the current time step. Features obtained from the input gate in the gating mechanism. For the Gate of Oblivion Candidate cell state, In cellular state, This represents the cell state at the previous moment. This represents element-wise multiplication. For output gate, These are the bias parameters for the input, forget, candidate state, and output gates. The parameters of the convolution kernel for input, forget, candidate states, and output gate are applied to the adjusted... ; The corresponding convolution kernel parameters are applied to the adjusted hidden state. , , , For different activation functions.

[0123] At each time step, the location attention module selectively aggregates the input features at each location by weighted summing of features across all locations, thereby capturing the global spatial dependencies of the wind and light time-series data. An SA self-attention module is proposed to construct a memory unit. Representing spatiotemporal information with a global perspective, the memory module receives two inputs: the input features at the current time step. And the memory of the final step The entire process is divided into three parts: feature aggregation to obtain global context information, memory update, and output, using weight W. mk and W mv 1×1 convolution, memory mapping as key Sum .

[0124] in, and Corresponding to and respectively and The number of channels, It is the height multiplied by the width of the feature map, representing the total number of positions after the two-dimensional feature map is flattened.

[0125] By querying Q h and key K m Matrix multiplication is used to calculate the similarity score between the input and the memory:

[0126] ;

[0127] The superscript T indicates transpose.

[0128] By applying the combined ReLU and sigmoid function along each column of the matrix, the weights of the clustered features are obtained. :

[0129] ;

[0130] Feature Y m The grid point at position i passes through value V mWeighted sum calculation of N positions in the middle:

[0131] ;

[0132] in, is the j-th column of the memory storage, and Y is the aggregation feature.

[0133] Furthermore, a gating mechanism is employed to update the memory storage. This allows SA to capture long-range dependencies in both the spatial and temporal domains, aggregating features Y and the original input. Used to generate input gate i' t and fusion feature g' t The value of the forget gate is replaced with (1-i'). t ) / i' t With fewer parameters, the update progress can be represented as:

[0134] ;

[0135] in, express The weights of the input gate and fusion gate in Y, express middle The weights of the input gate and the fusion gate, This represents the bias parameters for each gate.

[0136] storage unit Through convolution operations and aggregated feature Y t Update to obtain the global spatial correlation of the landscape, and Contains information about the past global spacetime, a memory storage unit Through convolution operations and aggregated feature Y t Update to obtain the global spatial correlation of the landscape, and Includes global spatiotemporal information from time 1 to t-1, and the output characteristics of the memory module. It is an output gate and updated memory The dot product between them is given by the following formula:

[0137] ;

[0138] in, These represent the weight parameters and bias parameters for the aggregated features and the original input, respectively.

[0139] Step 5.2: Introduce the NODE mechanism into SA-ConvLSTM, transforming the hidden state h(t) and cell state c(t) of SA-ConvLSTM into a set of coupled ordinary differential equations. The input gate i(t), forget gate f(t), output gate o(t), and candidate cell state g(t) are all transformed into the continuous-time evolution of cell state c(t).

[0140] ;

[0141] in, Time-scale parameters representing cell state This indicates element-wise multiplication, where the gated signal and the cell state are multiplied element-wise along the same dimension.

[0142] The hidden state h(t) is jointly determined by the output gate o(t) and the attention mechanism, and its evolution equation is:

[0143] ;

[0144] in, It is the timescale of the hidden state. This represents the global context information extracted by the self-attention module, and β is a parameter that controls the intensity of the attention drive.

[0145] In implementation, the evolution of the joint state vector z(t) = (h(t), c(t)) is determined by the vector field. describe:

[0146] ;

[0147] Here, the parameter θ includes the convolution kernel, self-attention weights, and time scale hyperparameters. Given the initial conditions z(t) k ), in the time interval [t k , t k+1 The next round of conditions z(t) is obtained by integrating using the ODE (Ordinary Differential Equation) solver in NODE. k+1 ), and so on, for a total of iterations. wheel.

[0148] Step 6: Output the predicted wind speed and solar radiation from the model.

[0149] Specifically, after training the physically embedded NODE-SA-ConvLSTM model, the fused feature data is input into the model's prediction stage. The model's output layer transforms the hidden state vector into meteorological elements in the physical quantity space through linear mapping. The output results include the meridional component U10 and zonal component V10 of wind speed, as well as the total solar radiation SSRD and direct solar radiation FDIR, forming a four-dimensional meteorological element output matrix.

[0150] The output Y(t) is presented in the form of a time series matrix:

[0151] ;

[0152] in, Y(t) The dimensions are [T, N, 4], where T represents the prediction step size, N represents the number of spatial grid points, and 4 represents the number of output variables.

[0153] Each time step corresponds to a complete spatial prediction field, which can simultaneously describe the changes in wind speed and solar radiation in both time and space dimensions.

[0154] After inverse normalization, the model output results are restored to m / s for wind speed and W / m² for radiation. The results are saved in NetCDF format, containing latitude and longitude coordinates, time index, and corresponding meteorological variable values. The spatial resolution of the data is consistent with the input ERA5 grid (0.25°×0.25°).

[0155] Step 7: Use the constant error compensation mechanism to perform residual balance correction on the prediction results.

[0156] Step 7.1: Design a constant error compensation mechanism. During model training and prediction, for each time step of the prediction result y... t Compared with measured values Calculate residuals :

[0157] ;

[0158] Set the sliding window length k, and perform a moving average on the residuals over k consecutive time steps to obtain the dynamic residual mean. :

[0159] ;

[0160] When the change in the mean residual in two consecutive iterations satisfies When, the constant compensation unit is triggered, where, This is a preset error fluctuation threshold.

[0161] Step 7.2: The compensation unit performs inverse normalization on the residuals and then adds them to the prediction output. The correction formula is as follows:

[0162] ;

[0163] Where α is the compensation coefficient, satisfying 0 < α < 1, and Norm() and DeNorm() represent the normalization and denormalization operations, respectively. The corrected prediction output. The model is re-inputted for the next iteration of training until the early termination condition is met. .

[0164] The final corrected prediction results are obtained. :

[0165] ;

[0166] Step 8: Construct the wind-solar complementarity coefficient based on the corrected prediction results and output the wind-solar endowment complementarity assessment results.

[0167] Specifically, WSCC is used to quantify the complementarity between wind speed and solar radiation, and the time series of wind speed components U10 and V10 and solar radiation SSRD and FDIR are analyzed. Calculate the Spearman rank correlation coefficient ρ respectively. s Kendall rank correlation coefficient T k :

[0168] ;

[0169] in, d i The sample rank difference, where I and J are the number of consistent and inconsistent logs, respectively, and n is the sample size.

[0170] The wind-solar complementarity coefficient (WSCC) is obtained by combining two correlation measures using a parameter adaptive weighted fusion method.

[0171] ;

[0172] in, The weights of the Spearman coefficients in WSCC, The weights of the Kendall coefficients in WSCC, The optimal wind-solar complementarity weighting factor satisfies Its optimal value is determined by the Bayesian optimization method.

[0173] In this embodiment, the specific WSCC of each meteorological element is shown in Table 1:

[0174] Table 1. Comparison of WSCC Correlation Among Variables

[0175]

[0176] Furthermore, this embodiment also provides a wind-solar resource complementarity assessment system based on physically embedded deep learning, with the structure as follows: Figure 2 As shown, it consists of seven modules: Module 1 Terrain Data Processing (Input), Module 2 Gradient Iterative Clustering Algorithm Module, Module 3 S-Transform Feature Extraction Module, Module 4 VMD-Hilbert Feature Extraction Module, Module 5 NODE-SA-ConvLSTM Module, Module 6 Error Constant Compensation Module, and Module 7 Wind-Solar Complementarity Evaluation Module (Output).

[0177] The specific build process is as follows:

[0178] First, the terrain data of the target area is acquired. Then, the gradient of the elevation difference between adjacent grid cells is calculated in the terrain data processing module to obtain the gridded terrain gradient field.

[0179] Secondly, the initial cluster centers and weight coefficients are set through the gradient iterative clustering algorithm module to determine the sample region. The Euclidean distance between the sample point and the cluster center is calculated to determine the cluster to which the sample point belongs. After each iteration, the cluster centers are updated until convergence, the optimal threshold D is obtained, and the region is divided accordingly to distinguish between flat and non-flat regions.

[0180] Then, in the S-Transform feature extraction module, the S-Transform time-frequency transformation is performed on the meteorological sequence of the flat area to extract the instantaneous energy spectrum and dominant frequency features of the wind speed and radiation sequences, forming a stationary time-frequency feature matrix;

[0181] Meanwhile, in the VMD-Hilbert feature extraction module, the meteorological sequence in the non-flat area is decomposed by VMD to obtain multiple narrowband intrinsic mode functions. Then, the instantaneous amplitude and instantaneous frequency of each mode are obtained by Hilbert transformation to construct a non-stationary time-frequency feature matrix.

[0182] Next, the features are fused and input into the NODE-SA-ConvLSTM module. Spatial dependencies are modeled through a self-attention mechanism, physical embedding is performed, and the continuous time evolution of the hidden state and cell state is characterized in the form of ordinary differential equations. The predicted results of wind speed components and radiation are output.

[0183] Furthermore, the error constant compensation module performs sliding statistics on the prediction residuals at each time step during the model training and prediction phases, calculates the mean of the dynamic residuals, and triggers compensation when the threshold condition is met. The compensation unit performs inverse normalization superposition and re-normalization correction on the residuals to achieve error constant control and prediction drift suppression, and outputs the correction results.

[0184] In this embodiment, a line graph comparing the predicted and corrected results of each variable at the 100th time step with the actual values ​​is shown, as follows: Figure 5 As shown.

[0185] Finally, the wind-solar complementarity coefficient is calculated, and the Spearman rank correlation coefficient and Kendall rank correlation coefficient of the wind speed component and radiation component in the correction result are obtained respectively. The comprehensive complementarity index is obtained through adaptive weighted fusion. The complementarity relationship between wind and solar elements is quantified in the wind-solar complementarity assessment module and the assessment results are output.

[0186] Furthermore, to verify the effectiveness of the method of this invention in wind and solar prediction tasks, representative deep learning models were selected as reference objects, including COA-CNN-LSTM, CNN-GRU-LSTM, Bi-A-ConvLSTM, and Mamba. These models have certain advantages in terms of combining convolutional and recurrent structures, gating and temporal dependency modeling, introducing attention mechanisms, and long sequence modeling, respectively.

[0187] By comparison, the performance differences and improvements of the NODE-SA-ConvLSTM proposed in this invention in wind speed and radiation prediction tasks can be clearly seen. Specific experimental results are shown in Table 2. Figure 6 and Figure 7 As shown.

[0188] Table 2 Comparative Experiments

[0189]

[0190] To investigate the roles of the components of the NODE-SA-ConvLSTM model in the overall prediction process, ablation experiments were conducted by sequentially removing the ConvLSTM, NODE, Attention, and physical constraints to construct different variant models for testing. This verified the contribution of each module to the final prediction accuracy. By comparing these experimental results, the importance and irreplaceability of each mechanism in wind speed and radiation prediction tasks can be revealed. The experimental results are shown in Table 3.

[0191] Table 3 Ablation Experiment

[0192]

[0193] The above description is only a preferred embodiment of the present invention. It should be noted that for those skilled in the art, several improvements and modifications can be made without departing from the principle of the present invention, and these improvements and modifications should also be considered within the scope of protection of the present invention.

Claims

1. A wind-solar resource complementation evaluation method based on physically embedded deep learning, characterized in that, It comprises the following steps: S1, obtaining the gridded meteorological data and gridded terrain data of the target area; S2, calculating the altitude change rate of adjacent grids of the gridded terrain data, determining the altitude change rate threshold value through the IGC algorithm, and distinguishing the flat area from the non-flat area according to the altitude change rate threshold value; S3, performing S-Transform time-frequency analysis on the meteorological data of the flat area to extract the instantaneous energy spectrum and frequency main feature of the core meteorological element, and using the combination of VMD and Hilbert transform on the meteorological data of the non-flat area to extract the instantaneous amplitude and frequency feature of the core meteorological element; S4, encoding the radiation data in the meteorological data to generate the periodic feature of the solar elevation angle, and splicing the frequency main feature, the frequency feature and the radiation prior feature to obtain the fusion feature; S5, constructing a physical embedded NODE-SA-ConvLSTM model, introducing the NODE mechanism into the SA-ConvLSTM to convert the discrete hidden state transmission process into a continuous time differential equation solution form, and inputting the fusion feature into the model; The NODE-SA-ConvLSTM model comprises layers, each layer consisting of an input gate i(t), a forget gate f(t), an output gate o(t), a hidden state h(t), a cell state c(t), and a candidate cell state g(t); Through the three gating units i(t), f(t) and o(t) in SA-ConvLSTM, the flow and update of information in h(t) are adjusted, the gating mechanism is continuously processed by NODE, the evolution process of the cell state c(t) from discrete to continuous is described by ordinary differential equation, the smooth modeling of time series feature is realized, and the evolution formula is as follows: ； wherein, a time scale parameter indicative of a cell state, denotes an element-wise multiplication; The hidden state h(t) is determined by the output gate o(t) and the attention mechanism, and the evolution equation is: ； wherein, is a time scale of the hidden state, represents global context information extracted by a self-attention module, is a parameter for controlling the attention driving strength; Evolved state vector By vector field Description: ； where the parameter θ contains all the convolutional kernels, self-attention weights, and time-scale hyperparameters in SA-ConvLSTM, given the initial conditions is integrated by the ODE solver in NODE at the time interval to obtain the next round of conditions and so on, for a total of iterations rounds. S6, outputting the prediction results of wind speed and solar radiation from the NODE-SA-ConvLSTM model; S7, using the error constant compensation mechanism to balance and correct the residual error of the prediction results; S8, constructing the wind-solar complementary coefficient according to the corrected prediction results, and outputting the wind-solar resource endowment complementary evaluation results.

2. The physical embedding-based deep learning complementary evaluation method for wind and light resources according to claim 1, characterized in that, The S2 comprises the following sub-steps: S2.1, processing the elevation height elements of the adjacent grids in the target area, obtaining the rate of elevation change by gradient operation on the elevation height difference Ah of the adjacent grid cells : ； where d i is the adjacent grid spacing of the ith group; S2.2, the IGC algorithm is used to obtain the altitude change rate threshold, and the converged clustering center is obtained by clustering and iterative optimization of the grid terrain gradient data With , the threshold D of the altitude change rate is set as ； S2.3, when G i <D, the region where the adjacent grid of the i-th group is determined as a flat region, when G i >D, the region where the adjacent grid of the i-th group is determined as a non-flat region.

3. The physical embedding-based deep learning complementary evaluation method for wind and light resources according to claim 2, characterized in that, The S2.2 performs clustering iterative optimization on the grid terrain gradient data, and the method is as follows: S2.2.1, setting two initial cluster centers , : ； wherein, represents the rate of change of altitude, min(G), max(G) are the maximum and minimum values of G, respectively, , are the weight coefficients of the two cluster centers, respectively; S2.2.2, the first clustering based on the initial clustering center, determining the cluster to which each sample point belongs by calculating the Euclidean distance between the sample point and the clustering center, and recalculating the elevation change rate G of the two regions after the division i Mean, set as the new clustering center S2.2.3, and so on, iterate n rounds until the clustering centers converge, obtaining converged clustering centers with .

4. The physical embedding-based deep learning complementary evaluation method for wind and light resources according to claim 1, characterized in that, The S3 comprises the following sub-steps: S3.1, performing S-Transform time-frequency analysis on the meteorological data of the flat area to extract the instantaneous energy spectrum and frequency main feature of the core meteorological element, and the method is as follows: The S-Transform algorithm is used to extract the instantaneous energy spectrum and main frequency feature of the meteorological data of the flat area, and the formula is as follows: ； wherein is a time variable in a time series, is a time signal, is a time shift, is a frequency, is a Gaussian window function, denotes a Fourier transform; The energy distribution of the signal in the time-frequency domain is represented by the modulus S(b,f) of the S-Transform, and an instantaneous energy spectrum is obtained The dominant frequency feature is obtained by energy weighted average The time-frequency feature matrix is output The flat zone feature vector set is obtained : ； S3.2, using the combination of VMD and Hilbert transform on the meteorological data of the non-flat area to extract the instantaneous amplitude and frequency feature of the core meteorological element, and the method is as follows: The VMD decomposes a signal into several narrow-band intrinsic mode functions (IMFs) by solving a variational problem, each IMF having a center frequency defined as: ； wherein, is the kth modal component, * denotes a convolution operation, denotes a unit impulse function, is a time derivative; is a transform kernel, representing a singular kernel function required for signal analyticization; is a center frequency is a demodulation factor with respect to a reference, making the IMF component into a baseband signal; The constraint condition is satisfied: ； A set of stationary sub-signals is obtained by solving the problem with the alternating direction method of multipliers Each modal component is transformed into an analytical signal by Hilbert transform : ； wherein is the Hilbert operator, j is the imaginary unit, yielding the instantaneous amplitude and the instantaneous frequency : ； wherein, denotes the instantaneous phase term of the k-th IMF modal component, denotes the derivative of the instantaneous phase with respect to time, used to obtain the instantaneous frequency of the modal component, resulting in a time-frequency feature matrix of the non-flat region: ; S3.3, the with The two sets of feature vectors are aligned by time index to form two sets of time-frequency feature data for the regions.

5. The physical embedding-based deep learning complementary evaluation method of wind and light resources according to claim 4, characterized in that, The S4 comprises the following sub-steps: S4.1、based on the observation time and geographic coordinates, calculating the solar elevation angle, solar zenith angle and solar azimuth angle, constructing the radiation period characteristic code, including: the diurnal variation of the solar elevation angle, the time period symmetry variation of the azimuth angle and the cosine correction term of the zenith angle, combining the code result with the radiation variables in the meteorological elements to form the radiation prior feature set containing the solar geometric parameters and radiation characteristics ; S4.2, respectively, to , , Z-score normalization, to obtain a standardized feature set , , ; S4.3, under the condition of consistent time index, the flat area time-frequency feature, the non-flat area time-frequency feature, the radiation prior feature and the wind-solar element vector E(t) are weighted and fused as four independent scales to obtain the multi-scale time-frequency feature fusion data F(t): ； wherein , , , is an adjustable weight parameter satisfying .

6. The physical embedding-based deep learning complementary evaluation method of wind and light resources according to claim 1, characterized in that, The predicted results of the wind speed and the solar radiation output from the NODE-SA-ConvLSTM model in the S6 include: the zonal component of the wind speed U 10 and the meridional component V 10 and the total solar radiation SSRD and the direct radiation FDIR to form a four-dimensional meteorological element output matrix : ； wherein, The dimension is [T, N, 4], T represents a prediction step, N represents a number of spatial grid points, and 4 is a number of output variables; each time step corresponds to a complete spatial prediction field, and simultaneously describes the change process of wind speed and solar radiation in the time and spatial dimensions.

7. The physical embedding-based deep learning complementary evaluation method of wind and light resources according to claim 1, characterized in that, The S7 comprises the following sub-steps: S7.1, design error constant compensation mechanism, in the model training and prediction process, calculate the prediction result of each time step y t the residual error between the measured value :​ ； Set the sliding window length k, and slide average the residual errors of the k time steps in the sliding window to obtain a dynamic residual mean : ； When the variation of the residual mean in two consecutive iterations satisfies , a constant compensation unit is triggered, wherein, is a preset error fluctuation threshold. S7.2, constant compensation unit pair After inverse normalization, the result is superimposed to the prediction output of S6 and re-normalized, which realizes the correction of the prediction output, and the correction formula is: ； wherein, is a compensation coefficient, satisfying , and respectively represent normalization and inverse normalization operations, and the modified prediction output is re-input into the NODE-SA-ConvLSTM model for the next round of iterative training until the early termination condition is met wherein, is the current minimum residual value. 8.The physical embedding-based deep learning method for wind-solar resource complementation evaluation according to claim 1, wherein, In the S8, the wind-solar complementary coefficient is used to quantify the complementary degree of wind speed and solar radiation, and the method is as follows: wind speed component U 10、 V 10 with solar radiation component SSRD 、 FDIR Spearman rank correlation coefficient and Kendall rank correlation coefficient : ； wherein, denotes the sample rank difference, and are the number of consistent and inconsistent logs, respectively, is the number of samples; The wind-solar complementary coefficient WSCC is obtained by combining the two correlation measures using a parameter adaptive weighted fusion method: ； wherein, , are the weights of Spearman and Kendall coefficients in WSCC, respectively, is the wind-solar complementary weight factor, satisfying The optimal value is determined by the Bayesian optimization method.

9. A physical-embedded deep learning-based wind-solar resource complementation assessment system, configured to implement the wind-solar resource complementation assessment method according to any one of claims 1 to 8. Comprise: a terrain data processing module configured to receive terrain data of a target region, calculate a gradient of an elevation difference between adjacent grid cells, and obtain a grid terrain gradient field; a gradient iterative clustering algorithm module configured to perform an iterative clustering process based on the terrain gradient data, calculate a Euclidean distance between a sample point and a cluster center to determine a region category by setting an initial cluster center and a weight coefficient, and update the cluster center after each iteration until convergence to obtain an optimal threshold D to distinguish between flat and non-flat regions; an S-Transform feature extraction module configured to perform an S-Transform time-frequency transformation on a meteorological sequence of the flat region, extract instantaneous energy spectrum and dominant frequency features of wind speed and radiation sequences, and form a stationary time-frequency feature matrix to obtain time-varying characteristics of meteorological elements in the flat region in a time-frequency adaptive manner; a VMD-Hilbert feature extraction module configured to perform VMD decomposition on a meteorological sequence of the non-flat region, obtain a plurality of narrowband intrinsic mode functions, and then perform Hilbert transformation to obtain instantaneous amplitude and instantaneous frequency of each mode to construct a non-stationary time-frequency feature matrix to reflect the non-stationary evolution law of the meteorological signal under complex terrain; a NODE-SA-ConvLSTM module configured to build a NODE-SA-ConvLSTM model, input fused features, model spatial dependency relationships through a self-attention mechanism, and describe continuous time evolution of hidden states and cell states in the form of ordinary differential equations to output prediction results of wind speed components and radiation; an error constant compensation module configured to perform sliding statistics on prediction residuals at each time step during model training and prediction, calculate a dynamic residual mean, and trigger compensation when a threshold condition is met, wherein a compensation unit performs inverse normalization stacking and re-normalization correction on the residuals to achieve error constant control and prediction drift suppression; a wind-solar complementary evaluation module configured to calculate the wind-solar complementary coefficient WSCC, calculate Spearman rank correlation coefficients and Kendall rank correlation coefficients of wind speed components and radiation components respectively, obtain a comprehensive complementary degree index through adaptive weighted fusion, quantify the complementary relationship between wind and solar elements, and output evaluation results.

Images