Novel power system multi-scale network construction energy storage planning optimization method based on source network charge clustering

Abstract

The invention provides a novel power system multi-scale network construction energy storage planning optimization method based on source network load clustering. The method comprises the following steps: firstly, processing source network load data by adopting a multi-dimensional adaptive feature extraction and fusion method, and realizing accurate processing of the data through interval adaptive probability density estimation, adaptive kernel tensor decomposition and marginal information entropy maximization; then, a multi-level data clustering method is adopted, and accurate clustering is achieved through hierarchical distance measurement and boundary self-adaptive adjustment; then, a multi-scale energy storage modeling method is adopted, and accurate modeling of the energy storage system is achieved through the dynamic characteristic matrix and self-adaptive topology embedding; and finally, a multi-objective optimization solution method is adopted, and efficient solution of an optimization problem is realized through adaptive hierarchical decomposition and hybrid evolutionary optimization. According to the invention, the accuracy and reliability of energy storage planning can be obviously improved, and the system operation efficiency is improved.

Description

Technical Field

[0001] The present invention relates to an energy storage planning optimization method, in particular to a novel power system multi-scale network energy storage planning optimization method with source-grid-load clustering. Background Art

[0002] Against the backdrop of the rapid development of new power systems, the access of a high proportion of renewable energy has brought huge challenges to the safe and stable operation of power systems. As an important flexible regulation resource, grid-connected energy storage plays a key role in promoting the consumption of renewable energy, improving the system regulation capability, and ensuring the safe operation of the power grid. Reasonable planning of the configuration scheme of grid-connected energy storage is of great significance to improving the economy and reliability of new power systems. Especially in scenarios with deep interaction between sources, grids and loads, the multi-scale optimization configuration of grid-connected energy storage must consider both the regulation needs of different time scales and the layout requirements of different spatial scales, which puts higher demands on energy storage planning.

[0003] At present, grid-type energy storage planning mainly adopts deterministic planning and probabilistic planning methods. Deterministic planning methods are mainly based on typical scenarios for optimization, and the optimal configuration scheme is solved through linear programming or mixed integer programming, but it is difficult to effectively deal with the randomness and volatility of renewable energy. The probabilistic planning method introduces randomness considerations and uses Monte Carlo simulation and other technologies to generate scenario sets, but the computational complexity is high and it is difficult to apply to large-scale systems. In terms of data processing, traditional clustering algorithms (such as K-means and fuzzy C-means) are mainly used to analyze source-grid-load data. In terms of energy storage characteristic modeling, simplified energy storage models are often used, and energy storage is often regarded as an ideal energy storage unit, ignoring the dynamic characteristics and life characteristics of energy storage equipment.

[0004] The existing technical solutions have the following problems: First, in terms of source-grid-load data processing, traditional clustering algorithms are difficult to effectively process spatiotemporal coupling data, and lack in-depth analysis of data temporal correlation and spatial correlation, resulting in the clustering results being unable to accurately reflect system characteristics; second, in terms of energy storage characteristic modeling, the existing models do not adequately describe the multi-time scale characteristics of the energy storage system, fail to accurately express the coupling relationship between energy storage capacity and power, and the dynamic modeling of the aging process of energy storage equipment is not accurate enough; third, in terms of optimization methods, when dealing with multi-objective optimization problems, the existing algorithms lack theoretical support for the determination of objective function weights, do not fully consider the impact of changes in the power grid topology, and the optimization process is computationally inefficient; finally, in terms of model solving, the existing methods are difficult to balance global search capabilities and local refinement capabilities, and lack an effective coordination mechanism for the parallel solution of large-scale optimization problems, resulting in the quality and efficiency of the solution results being difficult to meet actual needs. Summary of the invention

[0005] The purpose of the invention is to provide a new type of power system multi-scale grid construction and energy storage planning optimization method with source-grid-load clustering, in order to solve at least one technical problem existing in the prior art.

[0006] The technical solution is a new type of power system multi-scale network construction energy storage planning optimization method based on source-grid-load clustering, which is characterized by comprising the following steps:

[0007] Step S1, collecting basic data and organizing it into a three-dimensional data structure according to the time dimension, space dimension and variable type dimension, based on the three-dimensional data structure, calculating the local probability density value of each data point and performing data cleaning; performing kernel mapping and feature decomposition on the cleaned data to obtain core tensors and factor matrices, using the core tensors and factor matrices to calculate spatiotemporal correlation features, generating fused feature vectors, calculating information entropy and mutual information for the fused feature vectors, screening out important features, and forming an optimized feature matrix;

[0008] Step S2, input each eigenvector in the optimized feature matrix into a nonlinear basis function to calculate the feature manifold value; calculate the local structure distance for each pair of feature manifold values ​​and construct a structure matrix; based on the structure matrix, respectively calculate the Wasserstein distance, the topological structure distance and the time series similarity distance, and combine these three distances into a multi-level distance through a distance fusion operator; perform clustering operations using the multi-level distances to obtain an initial clustering result; optimize and adjust the boundary area of ​​the clustering result to obtain an optimized clustering result; construct a feature combination matrix based on the optimized clustering result, optimize the combination weight, and output the final feature combination result;

[0009] Step S3, input the final feature combination result into the energy storage dynamic characteristic function, calculate the energy storage characteristic matrix by combining the exponential function and the trigonometric function; calculate the state of charge function value, power characteristic function value, efficiency function value and life characteristic function value of the energy storage device based on the energy storage characteristic matrix; use the energy storage characteristic matrix to construct the time scale mapping and space scale mapping, perform the scale conversion operation to obtain the multi-scale mapping matrix; combine the multi-scale mapping matrix with the network topology characteristics to generate a topology embedding matrix; construct a multi-objective optimization model based on the topology embedding matrix, perform model verification and adjustment, and output the verified optimization model;

[0010] Step S4, decompose the verified optimization model into multiple sub-problems, calculate the coupling relationship matrix between the sub-problems; construct a hybrid evolution operator for each sub-problem, including a crossover operator, an adaptive mutation operator, and a local search operator, and design a selection mechanism; build a parallel solution framework to collaboratively solve the sub-problems; make convergence judgments and multi-dimensional evaluations on the solution results to screen the optimal solution set; verify and optimize the optimal solution set, and output the final optimized solution.

[0011] The beneficial effect can significantly improve the accuracy and reliability of energy storage planning and improve the system operation efficiency. BRIEF DESCRIPTION OF THE DRAWINGS

[0012] Figure 1 It is a flow chart of the present invention.

[0013] Figure 2 It is a flow chart of step S1 of the present invention.

[0014] Figure 3 It is a flow chart of step S2 of the present invention.

[0015] Figure 4 It is a flow chart of step S3 of the present invention.

[0016] Figure 5 It is a flow chart of step S4 of the present invention. DETAILED DESCRIPTION

[0017] like Figure 1 As shown, according to one aspect of the present application, a novel power system multi-scale grid-building energy storage planning optimization method with source-grid-load clustering is provided, comprising the following steps:

[0018] Step S1: Obtain the real-time power data, geographic location data and operating status data of the wind turbine from the power system database, obtain the wind speed data and temperature data of the corresponding geographic location from the meteorological database, obtain the real-time power data, geographic location coordinates and operating status data of the photovoltaic array from the photovoltaic system database, and obtain the power load data from the load monitoring system, and organize these data into a three-dimensional data structure according to the time dimension, space dimension and variable type dimension. Based on the three-dimensional data structure, calculate the local probability density value of each data point and perform data cleaning. Perform kernel mapping and feature decomposition on the cleaned data to obtain the core tensor and factor matrix. Use the core tensor and factor matrix to calculate the spatiotemporal correlation features and generate a fused feature vector. Calculate the information entropy and mutual information of the fused feature vector, screen out important features, and form an optimized feature matrix.

[0019] Step S2: Input each eigenvector in the optimized feature matrix output by step S1 into the nonlinear basis function to calculate the feature manifold value. Calculate the local structural distance for each pair of feature manifold values ​​and construct a structure matrix. Based on the structure matrix, calculate the Wasserstein distance, topological structure distance and time series similarity distance respectively, and combine these three distances into a multi-level distance through a distance fusion operator. Perform clustering operations using multi-level distances to obtain initial clustering results. Optimize and adjust the boundary area of ​​the clustering results to obtain optimized clustering results. Construct a feature combination matrix based on the optimized clustering results, optimize the combination weights, and output the final feature combination results.

[0020] Step S3: Input the final feature combination result output from step S2 into the energy storage dynamic characteristic function, and calculate the energy storage characteristic matrix by combining the exponential function and the trigonometric function. Calculate the state of charge function value, power characteristic function value, efficiency function value and life characteristic function value of the energy storage device based on the energy storage characteristic matrix. Use the energy storage characteristic matrix to construct the time scale mapping and space scale mapping, and perform scale conversion operations to obtain a multi-scale mapping matrix. Combine the multi-scale mapping matrix with the network topology characteristics to generate a topology embedding matrix. Construct a multi-objective optimization model based on the topology embedding matrix, perform model verification and adjustment, and output the verified optimization model.

[0021] Step S4: Decompose the verified optimization model output from step S3 into multiple sub-problems, and calculate the coupling relationship matrix between the sub-problems. Construct a hybrid evolution operator for each sub-problem, including a crossover operator, an adaptive mutation operator, and a local search operator, and design a selection mechanism. Construct a parallel solution framework to collaboratively solve the sub-problems. Perform convergence judgment and multi-dimensional evaluation on the solution results to screen the optimal solution set. Verify and optimize the optimal solution set, and output the final optimized solution.

[0022] In this embodiment, by constructing a new type of power system multi-scale network energy storage planning optimization method with complete source-grid-load clustering, the whole process innovation from data processing to final optimization solution is realized. The scheme realizes the accurate processing of source-grid-load data through multi-dimensional adaptive feature extraction and fusion method, realizes the accurate clustering of data through multi-level data clustering and feature combination method, realizes the accurate modeling of energy storage system through multi-scale network energy storage optimization modeling method, and realizes the efficient solution of optimization problem through multi-objective network energy storage optimization solution method. The organic connection between each step ensures the continuity of data flow and the consistency of processing process. Especially when dealing with energy storage planning problems under the background of high proportion of new energy access, the scheme can effectively process the spatiotemporal characteristics of source-grid-load data, accurately describe the multi-scale characteristics of energy storage system, accurately model the constraints of system, and efficiently solve optimization problems, thereby providing reliable technical support for energy storage planning of power system. The implementation of this scheme can significantly improve the accuracy and reliability of energy storage planning, reduce planning costs, improve system operation efficiency, and provide a strong guarantee for the construction of new power system.

[0023] According to one aspect of the present application, step S1 specifically comprises:

[0024] Step S11: Read the power generation value, geographic location coordinates and operating status parameters of each wind turbine at different time points from the power system database. Read the wind speed data and temperature data corresponding to the geographic location of each wind turbine from the meteorological database. Read the power generation value, geographic location coordinates and operating status parameters of each photovoltaic array at different time points from the photovoltaic system database. Read the power consumption value and power consumption type identification at different time points from the load monitoring system. Arrange the read data into a three-dimensional data structure according to the time point, geographic location and data type, in which each data point contains the corresponding time tag, location tag and type tag.

[0025] Step S12: For each data point in the three-dimensional data structure obtained in step S11, extract the data of five position points that are adjacent in space and 12 hours before and after the time position of the point, and calculate the mean and standard deviation of these data. Calculate the probability density value of each data point and the data sample set extracted from it. Calculate the time-varying adjustment coefficient based on the probability density distribution characteristics of the data point 24 hours before and after. Set the product of the time-varying adjustment coefficient and the standard deviation as the judgment threshold of the current data point. Mark the data points that exceed the judgment threshold as abnormal points. For each abnormal point, extract the normal data points in the 24-hour time window before and after it, calculate the weighted average of these normal data points, and replace the data value of the abnormal point with the calculated weighted average.

[0026] Step S13: For each data point in the data set obtained in step S12, calculate the Gaussian kernel value and polynomial kernel value between it and all other data points. Dynamically calculate the variance parameter value of the Gaussian kernel according to the distribution range of the data points, and dynamically determine the order parameter value of the polynomial kernel according to the change trend of the data. Combine the calculated Gaussian kernel value and polynomial kernel value into a composite kernel value through weight coefficients. Construct a kernel matrix with the composite kernel values ​​between all pairs of data points. Perform eigenvalue decomposition on the kernel matrix to extract eigenvectors and eigenvalues. Reconstruct the data based on the eigenvectors and eigenvalues, calculate the reconstruction error value, set the error threshold to filter the eigenvector, and obtain the reduced-dimensional core tensor and factor matrix.

[0027] Step S14: Separate the factor matrix of the time dimension and the factor matrix of the space dimension from the factor matrix obtained in step S13. Calculate the time correlation coefficient for the time dimension factor matrix, and calculate the spatial correlation coefficient for the space dimension factor matrix. Calculate the weight coefficient of the time dimension according to the time-varying characteristics of the data, and calculate the weight coefficient of the space dimension according to the spatial distribution characteristics of the data. Multiply the time correlation coefficient with the time weight coefficient to obtain the time correlation matrix, and multiply the space correlation coefficient with the space weight coefficient to obtain the space correlation matrix. Superimpose the time correlation matrix and the space correlation matrix to form a comprehensive correlation matrix. Multiply the comprehensive correlation matrix with the core tensor to generate a fused feature vector.

[0028] Step S15: For each feature component in the fused feature vector obtained in step S14, count the frequency of occurrence of its value and calculate the probability distribution. Calculate the information entropy value of each feature component according to the probability distribution. Calculate the mutual information value between each feature component and all other feature components. Combine the information entropy value and the mutual information value through the weight parameter to obtain the importance score of each feature component. Sort the feature components in descending order of importance score and select feature components from high to low. Calculate the matrix rank of the selected feature subset, and stop selecting when the matrix rank reaches the preset threshold. The selected feature components are combined into an optimized feature matrix.

[0029] In this embodiment, by establishing a multi-dimensional adaptive feature extraction and fusion method, efficient processing and feature extraction of source-grid-load data are achieved. Specifically, by introducing the interval adaptive probability density estimation method in the data cleaning stage, the identification process of abnormal data can dynamically adjust the judgment criteria according to the data distribution characteristics of different time periods and different regions, avoiding the misjudgment and missed judgment caused by the traditional fixed threshold method; the adaptive kernel tensor decomposition method is adopted in the feature extraction stage, and the ability to capture the nonlinear characteristics of the data is enhanced through the design of the composite kernel function, while the adaptive adjustment of the kernel function parameters ensures the accuracy of feature extraction; the dynamic spatiotemporal correlation matrix method is introduced in the feature fusion stage, and the adaptive fusion of features is realized through the dynamic weight function of the spatiotemporal dimension, effectively retaining the spatiotemporal correlation characteristics of the data; finally, the feature optimization is performed through the marginal information entropy maximization method, which not only ensures the maximization of the information amount of the selected features, but also reduces the redundancy between features through the introduction of mutual information. The combined application of this series of innovative methods enables the preprocessing and feature extraction process of source-grid-load data to better adapt to the variability and complexity of power system data, and provides high-quality feature data support for subsequent clustering analysis and optimization modeling.

[0030] According to one aspect of the present application, step S2 is specifically:

[0031] Step S21: Input each eigenvector in the optimized feature matrix output by step S1 into the hyperbolic tangent function and the sine function to obtain a nonlinear transformation value. Dynamically calculate the scaling parameter of the hyperbolic tangent function and the frequency parameter of the sine function according to the numerical distribution of the eigenvector. Combine the hyperbolic tangent function value with the sine function value through an adaptive weight coefficient to obtain a characteristic manifold value. Calculate the Euclidean distance for each pair of characteristic manifold values, divide the Euclidean distance by the dynamically calculated distance parameter to obtain the local structure distance. The local structure distances between all pairs of characteristic manifold values ​​form a structure matrix.

[0032] Step S22: Extract the distance value between each pair of data points from the structure matrix obtained in step S21. Calculate the Wasserstein distance between the probability distributions of each pair of data points, take the shortest path length of the data point in the network as the topological structure distance, and calculate the dynamic time warping distance of the data point time series as the time series similarity distance. Dynamically calculate the weight coefficients of the three distances according to the data distribution characteristics. Combine the weighted three distances into multi-level distances through the distance fusion operator. Perform hierarchical clustering operations based on the multi-level distances, divide the data points into different clusters, and form an initial clustering result.

[0033] Step S23: For each cluster obtained in step S22, extract the data points in its boundary area. For each boundary data point, calculate its membership with the current cluster center and the adjacent cluster center. Dynamically calculate the variance parameter of the cluster center according to the local distribution characteristics of the data point. Substitute the distance from the data point to the cluster center, the variance parameter and the membership into the boundary function to calculate the boundary function value. Substitute the difference between the actual boundary function value and the ideal boundary function value, the first-order derivative and the second-order derivative of the boundary function into the optimization objective function. Adjust the attribution of the boundary data points by minimizing the optimization objective function to obtain the optimized clustering result.

[0034] Step S24: Input the feature vector of each pair of clusters in the clustering result optimized in step S23 into the nonlinear mapping function, and calculate the sigmoid function value and the hyperbolic tangent function value respectively. Multiply the sigmoid function value with the hyperbolic tangent function value to obtain the nonlinear feature combination value. Construct a feature combination matrix for the feature combination values ​​of all cluster pairs. Calculate the one norm (L1 norm) and 2,1 norm (L2,1 norm) of the feature combination matrix, and substitute the norm value into the sparse regularization term. By minimizing the weighted sum of the square of the difference between the feature combination matrix and the target matrix and the sparse regularization term, a preliminary optimized feature combination weight matrix is ​​obtained. The L1 norm is used to calculate the sum of the absolute values ​​of all elements in the vector, tends to produce sparse solutions, and is often used in feature selection and compressed sensing. L2,1 norm, first calculate the L2 norm for each row of the matrix, and then calculate the L1 norm for the result. Consider group sparsity and feature selection at the same time.

[0035] Ω(W)=||W||1+α||W||2,1; where: -||W||1 is used to promote the sparsity of the weight matrix; -||W||2,1 is used to promote group sparsity, that is, the weights of the entire row or column are zero at the same time. α is a balancing parameter used to adjust the relative importance of the two norms.

[0036] Step S25: For the feature combination weight matrix obtained in step S24, calculate the average distance between the data points within each cluster to obtain the cluster compactness, calculate the minimum distance between different cluster centers to obtain the inter-class separation, calculate the time series continuity of the data points within the cluster to obtain the temporal continuity, and calculate the spatial distribution consistency of the clusters to obtain the spatial consistency. Combine these four indicators into a quality assessment value through weight coefficients. Calculate the entropy value of the quality assessment value, and add the entropy value to the weighted combination value to obtain a feedback signal. Fine-tune the feature combination weight matrix according to the feedback signal to obtain the final feature combination result.

[0037] In this embodiment, by proposing a multi-level data clustering and feature combination method, accurate clustering and feature fusion of source-network-load data are achieved. In the feature space construction stage, an adaptive feature manifold construction method is adopted, and the description ability of the intrinsic structure of the data is enhanced through the combination of nonlinear basis functions and the dynamic adjustment of parameters; in the clustering process, through the application of hierarchical distance metric clustering method, Wasserstein distance, topological structure distance and time series similarity distance are organically integrated, so that the clustering results can simultaneously reflect the distribution characteristics, network structure characteristics and time series change characteristics of the data; in the boundary optimization stage, the boundary adaptive adjustment method is adopted, and the accurate division of cluster boundaries is achieved by constructing boundary functions and optimizing objective functions; in the feature combination stage, the adaptive fusion of source-network-load features is achieved through the dynamic feature combination optimization method; finally, the multi-dimensional quality evaluation method is used to evaluate the clustering quality from multiple dimensions such as clustering compactness, inter-class separation, time series continuity and spatial consistency, and the continuous optimization of clustering results is achieved through the feedback mechanism. The combined use of these methods has significantly improved the accuracy and reliability of source-grid-load data clustering, providing a reliable data basis for optimizing grid energy storage configuration.

[0038] According to one aspect of the present application, step S3 is specifically:

[0039] Step S31: Expand the feature combination results output from step S2 according to the time dimension and the space dimension. For each data point in the time and space dimension, calculate the exponential decay function value and the cosine function value. According to the distribution characteristics of the time and space data, dynamically calculate the attenuation coefficient of the exponential decay function and the frequency and phase parameters of the cosine function. Multiply the exponential decay function value with the cosine function value to generate the energy storage characteristic basis function value. Perform a weighted combination of the energy storage characteristic basis function values ​​of all time and space dimension points to obtain an energy storage characteristic matrix. According to the energy storage characteristic matrix, calculate the state of charge value, charge and discharge power value, conversion efficiency value and remaining life value at each moment to form a set of energy storage constraint conditions.

[0040] Step S32: For the energy storage characteristic matrix obtained in step S31, extract the time dimension data and the space dimension data respectively. Input the time dimension data into the basis function, calculate the mapping value of each time point, and multiply the mapping value with the Gaussian function value to obtain the time scale mapping value. Input the space dimension data into the basis function, calculate the mapping value of each space point, and multiply the mapping value with the Gaussian function value to obtain the space scale mapping value. According to the data distribution characteristics of the time dimension and the space dimension, dynamically calculate the variance parameter of the Gaussian function. Perform tensor product operation on the time scale mapping value and the space scale mapping value to generate a scale mapping function matrix. Perform matrix multiplication operation on the scale mapping function matrix and the energy storage characteristic matrix to obtain a multi-scale mapping matrix.

[0041] Step S33: Extract node data and edge data of the power grid. Input the node data into the nonlinear mapping function, calculate the product of the sigmoid function value and the hyperbolic tangent function value to obtain the node eigenvalue. Input the edge data into the nonlinear mapping function, calculate the product of the sigmoid function value and the hyperbolic tangent function value to obtain the edge eigenvalue. Construct a node feature matrix based on the node eigenvalue, and construct an edge feature matrix based on the edge eigenvalue. Perform a weighted combination of the node feature matrix and the edge feature matrix to obtain a network topology matrix. Perform a tensor product operation on the network topology matrix and the multi-scale mapping matrix obtained in step S32 to generate a topology embedding matrix.

[0042] Step S34: Based on the topological embedding matrix obtained in step S33, the sum of the construction cost and operation and maintenance cost of the energy storage system, the system stability index, the reliability index and the life index are calculated respectively. These indexes are constructed into an objective function set. The state of charge constraint, power constraint, efficiency constraint and life constraint are extracted from the energy storage constraint condition set of step S31. The objective function set and the constraint condition set are combined into an optimization model.

[0043] Step S35: For the optimization model constructed in step S34, respectively calculate the accuracy index, computational efficiency index, robustness index and convergence index of the model. Combine these indexes into a verification index vector. Calculate the gradient value of the verification index vector. Adjust the parameters of the optimization model according to the gradient value. Repeat the process of calculating the verification index and adjusting the parameters until the verification index reaches the preset threshold, and output the verified optimization model and the verification index vector.

[0044] In this embodiment, by establishing a multi-scale network energy storage optimization modeling method, accurate modeling and multi-scale optimization of energy storage system characteristics are achieved. In the energy storage characteristic modeling stage, the dynamic characteristic matrix method is adopted to describe the attenuation characteristics and periodic characteristics of the energy storage system through the combination of exponential basis functions and triangular basis functions. The adaptive adjustment of the basis function parameters ensures the accurate description of the energy storage characteristics by the model; in the spatiotemporal scale mapping stage, the multi-dimensional scale adaptive mapping method is used to achieve seamless conversion of different time scales (second response to seasonal regulation) and spatial scales (equipment layer to system layer), solving the limitation that traditional modeling methods are difficult to handle multi-scale problems at the same time; in the network topology embedding stage, the energy storage characteristics are organically combined with the grid topology structure through the adaptive topology embedding method, which improves the adaptability of the model to network constraints; in the optimization model construction stage, the unified description of multiple objectives such as economy, stability, reliability and life is achieved through the multi-objective constraint optimization model method; finally, the dynamic optimization and adjustment of the model is achieved through the adaptive model verification method. The coordinated application of these innovative methods has significantly improved the accuracy and adaptability of energy storage system modeling, laying a solid theoretical foundation for subsequent optimization solutions.

[0045] According to one aspect of the present application, step S4 is specifically:

[0046] Step S41: Group the objective functions and constraints in the optimization model outputted from step S3 according to the coupling relationship. Calculate the weight coefficient for each group of objective functions to generate the objective function of the sub-problem. Calculate the weight coefficient for each group of constraints to generate the constraints of the sub-problem. Combine the objective function and constraints of each sub-problem into a complete sub-problem. Calculate the first-order and second-order derivatives of the objective function and constraints between each pair of sub-problems. Combine the derivative values ​​into a coupling relationship matrix, and output the decomposed sub-problem set and the coupling relationship matrix.

[0047] Step S42: For each sub-problem outputted from step S41, two candidate solutions are selected, a new solution is generated through linear combination, and a time-varying disturbance term is added to the new solution to obtain the output solution of the crossover operator. A random disturbance that obeys a normal distribution is added to the candidate solution, and the disturbance amplitude is dynamically adjusted during the optimization process to obtain the output solution of the mutation operator. The objective function gradient of the candidate solution is calculated, and an adaptive step size is moved in the opposite direction of the gradient to obtain the output solution of the local search operator. The fitness function value is calculated for the output solutions of the three operators, and the selection probability is calculated based on the fitness value and the adaptive temperature parameter, and the three operators and their selection mechanisms are outputted.

[0048] Step S43: Divide the sub-problem set outputted by step S41 into the parallel computing framework. For each parallel process, record the position vector, velocity vector and acceleration vector of the solution. Calculate the synergy coefficient between processes according to the coupling relationship matrix of the sub-problems. Substitute the Euclidean distance between the current solution set and the optimal solution set into the synergy kernel function to obtain the synergy strength between processes. Adjust the evolution direction and step size of each process according to the synergy strength, perform parallel solution operations, and output the parallel solution set.

[0049] Step S44: Extract the solution sets of two adjacent iterations from the parallel solution sets output in step S43. Calculate the relative error and objective function gradient value of the two solution sets. Compare the relative error with the preset convergence threshold, and compare the gradient value with the preset gradient threshold to determine whether the convergence condition is met. For the converged solution set, calculate the objective function value, constraint satisfaction, solution diversity, and solution stability. Filter the solution set according to the multi-dimensional evaluation results and output the optimal solution set.

[0050] Step S45: For the optimal solution set outputted in step S44, calculate the economic index, technical feasibility index, operation stability index and system adaptability index respectively. Combine the calculation results of each index into a verification index vector. Calculate the gradient value of the verification index vector. Determine the optimization direction and step size according to the gradient value, and make local adjustments to each solution in the optimal solution set. Repeat the verification and adjustment process until the verification index meets the preset requirements, and output the final optimized solution and evaluation index.

[0051] In this embodiment, by proposing a multi-objective grid energy storage optimization solution method, efficient solution of complex optimization problems is achieved. In the optimization problem decomposition stage, the large-scale optimization problem is decomposed into multiple sub-problems that can be solved in parallel through the adaptive hierarchical decomposition method, and the correlation between the sub-problems is maintained through the construction of the coupling matrix; in the optimization algorithm construction stage, the crossover operator, adaptive mutation operator and local search operator are organically combined through the hybrid evolutionary optimization method, which enhances the global search ability and local refinement ability of the algorithm; in the parallel solution implementation stage, the dynamic parallel collaborative solution method is used to achieve efficient parallel solution of sub-problems, and the introduction of collaborative kernel functions ensures the consistency of solutions; in the solution convergence and screening stage, the multi-dimensional solution evaluation method is used to evaluate the quality of the solution from multiple dimensions such as objective function value, constraint satisfaction, solution diversity and stability; finally, the optimal solution is finely adjusted through the adaptive solution optimization method to further improve the practicality of the solution. The innovative application of this series of methods has significantly improved the solution efficiency and solution quality of optimization problems, making it possible to solve large-scale grid energy storage optimization problems.

[0052] According to one aspect of the present application, step S11 is specifically:

[0053] Step S111: Extract the 10-minute power data, geographic coordinate data, and operating status data of each wind turbine from the power system database, and extract the wind speed data and temperature data of the corresponding geographic location at 10-minute time intervals from the meteorological database. Align the extracted wind speed data and temperature data with the wind power data according to the timestamp. Calculate the correlation coefficient between wind power and wind speed, and remove data points with correlation coefficients below the threshold. Correct the wind power according to the temperature data to obtain a temperature-corrected wind power data set.

[0054] Step S112: Extract the 10-minute power data, geographic coordinate data, and operating status data of each photovoltaic array from the photovoltaic system database, and extract the irradiance data and component temperature data of the corresponding geographic location at 10-minute time intervals from the meteorological database. Align the irradiance data and component temperature data with the photovoltaic power data according to the timestamp. Calculate the ratio of photovoltaic power to irradiance to obtain the photoelectric conversion efficiency sequence. Based on the photoelectric conversion efficiency sequence and component temperature data, calculate the temperature correction coefficient to obtain the corrected photovoltaic power data set.

[0055] Step S113: Extract 10-minute load data, power type identification and geographic location data of each power node from the load monitoring system. Divide the load data into industrial load, commercial load and residential load according to the power type identification. Calculate the power factor for each type of load data, and mark the data points with a power factor lower than the preset threshold. For the marked data points, extract the valid data 12 hours before and after to calculate the correction coefficient, and multiply the correction coefficient by the original load data to obtain the corrected load data. Merge the three types of corrected load data according to timestamps and geographic locations to obtain a complete load data set.

[0056] Step S114: Align the wind power data set output by step S111, the photovoltaic power data set output by step S112, and the load data set output by step S113 according to the time dimension. Calculate the mean and standard deviation for each data set. Subtract the mean of the corresponding data set from each data point and divide it by the standard deviation to obtain the standardized data value. Perform interval mapping on the standardized data values ​​and adjust the numerical range to the [0,1] interval. The mapped data are recorded as the standardized wind power data set, the standardized photovoltaic data set, and the standardized load data set.

[0057] Step S115: Sort the three standardized data sets outputted from step S114 according to timestamp and geographic location. Construct a three-dimensional array structure, in which the first dimension represents the time series, the second dimension represents the geographic location, and the third dimension represents the data type. For each time slice in the three-dimensional array, calculate the correlation coefficient between different geographic locations. For each geographic location slice, calculate the autocorrelation coefficient at different time points. For each data type slice, calculate the spatiotemporal cross-correlation coefficient. Use the calculated correlation coefficient as a weight coefficient to perform weighted processing on the three-dimensional array to obtain the final three-dimensional data structure.

[0058] In this embodiment, high-quality preliminary processing of source-grid-load data is achieved through a refined data reading and preprocessing process. In the wind power data processing stage, by introducing data screening and temperature correction mechanisms based on correlation coefficients, the abnormal points of the wind speed-power curve are effectively eliminated, and the reliability of wind power data is improved; in the photovoltaic data processing stage, the actual operating characteristics of the photovoltaic system are accurately reflected through the dynamic calculation of the photoelectric conversion efficiency sequence and the real-time update of the temperature correction coefficient; in the load data processing stage, the accuracy of the load data is improved through classification processing and power factor correction; in the data standardization stage, the relative distribution characteristics of the data are maintained through adaptive interval mapping; in the data structure construction stage, the temporal and spatial correlation characteristics of the data are strengthened by introducing correlation coefficient weights. This refined data preprocessing scheme significantly improves the basic data quality of subsequent analysis through multi-level data quality control and feature retention mechanisms, and provides reliable data support for the entire optimization process.

[0059] According to one aspect of the present application, step S13 is specifically:

[0060] Step S131: Calculate the kurtosis and skewness of the data distribution for the cleaned data set output from step S12. Dynamically determine the variance parameter of the Gaussian kernel function based on the kurtosis value, and dynamically determine the order parameter of the polynomial kernel function based on the skewness value. Perform block processing on the data set and calculate the local statistical features of each data block. Construct an adaptive weight coefficient based on the local statistical features, multiply the Gaussian kernel parameters and the polynomial kernel parameters with the weight coefficient, and obtain a modified kernel function parameter set.

[0061] Step S132: Using the kernel function parameter set outputted in step S131, calculate the Gaussian kernel value and polynomial kernel value between each pair of data points in the data set. Construct a distance decay function according to the temporal and spatial position relationship of the data points, and multiply the distance decay function value by the kernel function value to obtain a corrected kernel value. Construct a symmetric matrix for the corrected kernel value, and check the semi-positive definiteness of the matrix. Perform eigenvalue correction on the kernel matrix that does not meet the semi-positive definiteness to obtain a corrected kernel matrix.

[0062] Step S133: Perform eigenvalue decomposition on the kernel matrix output from step S132. Calculate the cumulative contribution rate of the decomposed eigenvalue sequence. Extract the eigenvector subset corresponding to the main eigenvalue from the eigenvector matrix. Calculate the projection value of each eigenvector in the eigenvector subset with the original data. The projection values ​​are combined into a feature projection matrix to obtain a low-dimensional representation of the data.

[0063] Step S134: Reconstruct the low-dimensional feature representation output by step S133 into the original data dimension through reverse mapping. Calculate the mean square error and relative error between the reconstructed data and the original data. Calculate the statistical characteristic value of the reconstructed error sequence. Construct an error evaluation function based on statistical characteristics, and determine the optimal number of features according to the evaluation function value.

[0064] Step S135: Based on the optimal number of features determined in step S134, select a corresponding number of eigenvectors from the eigenvector matrix. Calculate the orthogonality index of the selected eigenvectors. Construct a feature importance evaluation function that comprehensively considers the orthogonality of the eigenvectors and the contribution rate of the corresponding eigenvalues. Sort and select the eigenvectors according to the evaluation results, and output the final core tensor and factor matrix.

[0065] In this embodiment, efficient dimensionality reduction and feature extraction of data are achieved through kernel mapping and feature decomposition methods. In the kernel function parameter calculation stage, the adaptability of the kernel function to the data distribution characteristics is improved by introducing a dynamic parameter adjustment mechanism based on kurtosis and skewness; in the kernel matrix construction stage, the effectiveness of the kernel matrix is ​​ensured by introducing a distance decay function and a semi-positive definiteness test; in the feature decomposition process, the dynamic calculation of the cumulative contribution rate and the analysis of the projection value are used to achieve effective feature extraction; in the reconstruction error calculation stage, the accuracy of the dimensionality reduction process is ensured through a comprehensive evaluation of statistical features; in the feature screening stage, the optimal selection of features is achieved through a combination of orthogonality indicators and importance evaluation. This feature extraction method, through the innovative design of multiple links, not only ensures the accuracy of feature extraction, but also improves the calculation efficiency, providing high-quality feature representation for subsequent clustering analysis.

[0066] According to one aspect of the present application, step S23 is specifically:

[0067] Step S231: Extract the data point set of each cluster from the initial clustering result output in step S22. Calculate the coordinates and distribution range of the center point of each cluster. Determine the search radius according to the distribution range, and extract data points within the range with the cluster center as the center and the search radius as the radius. Calculate the distance between the extracted data point and the cluster center, and arrange the distance values ​​in ascending order. Set a distance threshold, mark the data points with a distance greater than the threshold as boundary candidate points, and obtain a boundary candidate point set.

[0068] Step S232: For each data point in the set of boundary candidate points outputted in step S231, extract the center point of the cluster to which it currently belongs and the center points of the adjacent clusters. Calculate the Euclidean distance between the data point and each cluster center. Construct a radial basis function based on the distance value and calculate the membership of the data point to each cluster center. Adjust the membership value according to the local density distribution characteristics of the data point. The adjusted membership values ​​are combined into a membership matrix.

[0069] Step S233: Use the membership matrix output from step S232 to construct a boundary function. Extract data points within a predetermined range around each boundary candidate point and calculate the local data density. Construct an adaptive variance function based on the local data density. Substitute the membership value and the adaptive variance into the exponential function to calculate the boundary strength value. Normalize the boundary strength value to obtain a boundary function value sequence.

[0070] Step S234: Calculate the spatial gradient and time gradient of the boundary function value sequence output in step S233. Construct a gradient-based regularization term, and use the absolute value of the first-order derivative and the sum of the squares of the second-order derivative as a penalty term. Use the square of the difference between the actual boundary function value and the ideal boundary function value as a fitting term. Dynamically adjust the weight coefficients of the fitting term and the penalty term according to the data distribution characteristics to construct a complete optimization objective function.

[0071] Step S235: Perform minimization operation on the optimization objective function constructed in step S234. Adjust the clustering of the boundary candidate points according to the optimization results. Recalculate the cluster center and distribution range of the adjusted clustering results. Check the compactness and separation of the adjusted clusters. When the compactness and separation meet the preset conditions, output the optimized clustering results.

[0072] In this embodiment, accurate division and optimization of clustering results are achieved through an innovative boundary optimization process. In the boundary point extraction stage, the boundary candidate points are accurately identified by setting the adaptive search radius and distance threshold; in the membership calculation stage, the accuracy of the membership calculation is improved by introducing a dynamic adjustment mechanism of the local density distribution characteristics; in the boundary function construction stage, the adaptability of the boundary description is enhanced by introducing an adaptive variance function; in the optimization target construction stage, an effective optimization target is formed through comprehensive consideration of gradient information and regularization terms; in the boundary adjustment stage, accurate boundary division is achieved through dynamic evaluation of compactness and separation. This boundary optimization method significantly improves the accuracy and reliability of clustering results through multi-level boundary feature analysis and optimization processing, and provides an accurate clustering basis for energy storage configuration optimization.

[0073] According to one aspect of the present application, step S31 is specifically:

[0074] Step S311: Expand the final feature combination result output from step S2 according to time series and spatial position. Construct exponential basis functions and triangular basis functions. The exponential basis functions are used to describe the energy storage attenuation characteristics, and the triangular basis functions are used to describe the energy storage cycle characteristics. Construct dynamic basis functions according to the time distribution characteristics of the data, and modify the original basis functions. Output the modified basis function set as the basic description function group of energy storage characteristics.

[0075] Step S312: For the basis function group outputted in step S311, extract the variation characteristics of the data at different time scales and spatial scales. Calculate the data volatility at the time scale and the distribution discreteness at the spatial scale. Dynamically adjust the attenuation coefficient of the exponential basis function according to the volatility, and dynamically adjust the frequency parameter and phase parameter of the triangular basis function according to the discreteness. Construct a parameter adaptive mapping function to realize the dynamic update of the basis function parameters, and output the updated basis function parameter set.

[0076] Step S313: Calculate the characteristic function value of each spatiotemporal data point using the basis function parameter set output from step S312. Construct a weighted combination of the characteristic function values, and dynamically assign weight coefficients based on the importance of the data. Arrange the weighted characteristic function values ​​into a matrix form. Calculate the singular values ​​of the matrix, and perform noise reduction on the characteristic matrix based on the singular value distribution characteristics. Perform a standardization operation on the processed characteristic matrix, and output the energy storage characteristic matrix.

[0077] Step S314: extract time series data from the energy storage characteristic matrix outputted in step S313. Calculate the state of charge change rate of the energy storage device at different times to obtain the state of charge function. Calculate the charge and discharge power change characteristics based on the state of charge function to obtain the power characteristic function. Calculate the energy conversion efficiency based on the power characteristic function to obtain the efficiency function. Construct a life prediction model based on historical operation data and calculate the remaining life function. These function values ​​are combined into a constraint parameter sequence.

[0078] Step S315: Classify and sort the constraint parameter sequence output from step S314. Construct the state of charge constraint interval, set the power change constraint range, determine the efficiency constraint threshold, and establish the life constraint boundary. Normalize these constraints. Construct the membership function of the constraint conditions and calculate the degree of constraint satisfaction. The constraints and their membership functions form a complete set of constraints.

[0079] Through the innovative energy storage characteristic modeling method, the accurate description of the multi-dimensional characteristics of the energy storage system is achieved. In the basis function construction stage, the dynamic characteristics of the energy storage system are accurately characterized through the organic combination of exponential decay and periodic characteristics; in the parameter dynamic adjustment stage, the adaptive optimization of the basis function parameters is achieved through the real-time calculation of volatility and distribution discreteness; in the characteristic matrix generation stage, the stability of the characteristic description is improved through singular value analysis and noise reduction processing; in the constraint condition calculation stage, the real-time validity of the constraint conditions is ensured through the dynamic update of multi-dimensional constraint parameters; in the constraint set construction stage, the flexible processing of the constraint conditions is achieved through the introduction of membership functions. This energy storage characteristic modeling method, through the organic combination of multiple innovative links, significantly improves the accuracy and adaptability of energy storage system modeling, and provides a reliable model basis for optimal configuration.

[0080] According to one aspect of the present application, step S34 is specifically:

[0081] Step S341: Extract cost features, stability features, reliability features, and life features from the topology embedding matrix output from step S33. Convert construction cost and operation and maintenance cost data into annualized cost series. Construct stability evaluation series based on system status data. Calculate reliability index series based on failure rate and repair time. Generate life prediction series using equipment aging characteristics. Combine these series to construct a multi-objective function matrix.

[0082] Step S342: Extract various constraints from the energy storage constraint set outputted in step S31. Convert the state of charge constraint into a state variable inequality. Convert the power constraint into a power change rate inequality. Convert the efficiency constraint into an energy conversion efficiency inequality. Convert the life constraint into a cycle count inequality. Normalize all inequality constraints and construct a constraint matrix.

[0083] Step S343: Normalize the multi-objective function matrix output from step S341. Calculate the value range and distribution characteristics of each objective function. Construct a normalization function based on the distribution characteristics to convert the objective functions of different dimensions to a unified interval. Calculate the correlation coefficient matrix for the converted objective functions. Construct the objective function weight coefficient based on the correlation coefficient to obtain a normalized multi-objective function set.

[0084] Step S344: Perform normalization transformation on the constraint matrix outputted in step S342. Extract the boundary value and critical value of the constraint. Construct slack variables of the constraint, and transform the hard constraint into a soft constraint with penalty term. Calculate the compatibility index between the constraint. Adjust the constraint according to the compatibility index, and output a normalized constraint set.

[0085] Step S345: Integrate the normalized multi-objective function set output from step S343 and the normalized constraint condition set output from step S344. Construct an association matrix of the objective function and the constraint condition. Calculate the coupling index of the model based on the association matrix. Decouple the highly coupled constraint conditions. Assemble the processed objective function and constraint conditions into a complete optimization model, and construct a priority sequence for model solution.

[0086] Through the systematic optimization model construction process, the accurate description of multi-objective optimization problems is achieved. In the objective function construction stage, a comprehensive evaluation system is formed through the dynamic combination of multi-dimensional indicators; in the constraint condition conversion stage, the flexibility of constraint condition processing is improved through normalization processing and the introduction of slack variables; in the objective function normalization stage, the reasonable balance of multiple objectives is achieved through distribution feature analysis and dynamic adjustment of weights; in the constraint condition normalization stage, the feasibility of model solution is enhanced through compatibility analysis and soft constraint conversion; in the model integration stage, the overall performance of the model is improved through coupling degree analysis and decoupling processing. This optimization model construction method significantly improves the accuracy and practicality of the optimization model through multi-level model optimization and adjustment, laying a solid foundation for the final solution.

[0087] According to one aspect of the present application, step S42 is specifically:

[0088] Step S421: Select two candidate solutions from the subproblem set output from step S41. Calculate the adaptive weight coefficients based on the objective function values ​​of the candidate solutions. Construct a weighted linear combination function to generate a preliminary crossover solution. Extract the local distribution features of the candidate solutions and construct a time-varying perturbation function. Add the time-varying perturbation value to the preliminary crossover solution to obtain the final crossover solution. Perform a feasibility check on the crossover solution to ensure that it meets the basic constraints.

[0089] Step S422: Construct a mutation operator for the crossover solution output from step S421. Extract the neighborhood distribution characteristics of the solution and calculate the local probability density function. Construct an adaptive mutation amplitude based on the probability density function to generate a random perturbation sequence that obeys a specific distribution. Synthesize the perturbation sequence with the crossover solution to obtain a mutation solution. Dynamically adjust the mutation probability according to the stage characteristics of the optimization process and output the adjusted mutation solution sequence.

[0090] Step S423: Calculate the objective function gradient using the sequence of variant solutions output from step S422. Extract the directional and amplitude features of the gradient vector. Construct an adaptive step function, and dynamically adjust the step value as the optimization progresses. Move in the opposite direction of the gradient with an adaptive step size to obtain a local search solution. Perform a boundary check on the local search solution, map the out-of-bounds solution back to the feasible domain, and output the corrected local search solution set.

[0091] Step S424: Calculate the fitness values ​​of each solution set output by steps S421, S422 and S423. Extract the fitness distribution characteristics of the population and construct a temperature function based on the distribution characteristics. Substitute the fitness value and temperature parameters into the selection probability calculation function to obtain the selection probability of each solution. Construct a roulette selection operator, filter the solution set according to the selection probability, and output the filtered solution set.

[0092] Step S425: Group the filtered solution sets output from step S424 by operator type. Calculate the statistical characteristics of each solution set and construct the distance matrix between groups. Calculate the diversity index of the solution set based on the distance matrix. Dynamically adjust the use probability of different operators according to the diversity index. The adjusted operators and their use probabilities form a hybrid evolution strategy and output a complete evolution operator system.

[0093] Through the design of hybrid evolutionary operators, efficient solutions of optimization algorithms are achieved. In the crossover operator design stage, the exploration ability of the population is enhanced by the introduction of adaptive weights and time-varying disturbances; in the mutation operator design stage, the ability of the algorithm to jump out of the local optimum is improved by the application of probability density functions and dynamic mutation probabilities; in the local search operator design stage, the refined optimization of the solution is ensured by the adaptive step size and boundary checking mechanism; in the selection mechanism design stage, the reasonable screening of the solution is achieved by combining the temperature function and roulette wheel selection; in the operator integration stage, the overall effect of the evolutionary strategy is optimized through diversity analysis and dynamic adjustment. This hybrid evolutionary operator design significantly improves the convergence efficiency of the optimization algorithm and the quality of the solution through the synergy of multiple innovative mechanisms, providing an efficient tool for solving the optimization problem of grid-connected energy storage.

[0094] In a certain scenario, the specific processing process is as follows:

[0095] S1: Data preprocessing and feature extraction

[0096] S11: Raw data collection and preliminary processing

[0097] Read the wind power data set Dw(t,l), photovoltaic data set Ds(t,l), and load data set Dl(t,l), where t represents the time dimension and l represents the spatial position coordinate;

[0098] Construct a three-dimensional tensor X(t,l,v) where v represents the variable type (v=1 wind power, v=2 photovoltaic, v=3 load);

[0099] Get the original three-dimensional tensor X∈R^(T×L×V);

[0100] S12: Adaptive data cleaning reads the original 3D tensor X

[0101] (1) Calculate the dynamic probability density function: P(x) = ∑(wi × Ki(x-xi) / hi) where wi is the adaptive weight, Ki is the kernel function, and hi is the variable bandwidth; (2) Construct an adaptive confidence interval: CI(t) = μ(t) ± k(t) × σ(t) where k(t) is the time-varying adjustment coefficient, which is dynamically adjusted based on the local data distribution characteristics;

[0102] Get the cleaned tensor X'.

[0103] S13: Tensor decomposition

[0104] Read the cleaned tensor X';

[0105] (1) Construct the kernel mapping function: Φ(X') = ∑(αi × K(X', Xi)) where K is the new composite kernel function: K(a, b) = θ1exp(-||ab|| 2 / σ1 2 )+θ2(a·b+c) 2 θ1, θ2 are adaptive weight coefficients; (2) Perform core tensor decomposition: Φ(X') = G×1U1×2U2×3U3, where G is the core tensor and Ui is the factor matrix;

[0106] Get the decomposed core tensor G and factor matrix {U1,U2,U3};

[0107] S14: Spatiotemporal feature fusion

[0108] Read the core tensor G and factor matrix {U1,U2,U3};

[0109] (1) Constructing the spatiotemporal correlation matrix: R(t,l)=β(t)×U1U1 T +γ(l)×U2U2 T Where β(t) and γ(l) are dynamic weight functions in the spatiotemporal dimension. (2) Calculate the comprehensive feature vector: F = R(t,l) × vec(G) where vec(·) represents the tensor vectorization operation. The fused feature matrix F is obtained.

[0110] S15: Feature Optimization and Dimensionality Reduction

[0111] Read the fused feature matrix F;

[0112] (1) Calculate the marginal information entropy: H(fi) = -∑p(fi)log(p(fi))+λ×MI(fi;F / fi), where MI is the mutual information and λ is the balance parameter; (2) Feature selection: max∑H(fi)strank(Fs)≥r||Fs||0≤k, where Fs is the selected feature subset;

[0113] The optimized feature matrix Fs is obtained.

[0114] S2: Data clustering and feature combination

[0115] S21: Dynamic feature space construction

[0116] Read the optimized feature matrix Fs output by S1;

[0117] (1) Constructing feature manifolds: in, is a nonlinear basis function: is the adaptive parameter;

[0118] (2) Calculate the local structure matrix: L(i,j) = exp(-||M(i)-M(j)||2 / η(i,j)) where η(i,j) is the adaptive distance parameter;

[0119] The characteristic manifold M and structure matrix L are obtained.

[0120] S22: Multi-dimensional Clustering

[0121] Read the feature manifold M and structure matrix L;

[0122] (1) Construct a multi-level distance metric: D(x,y) = θ1d1(x,y)fuθ2d2(x,y)fuθ3d3(x,y) where: d1 is the Wasserstein distance, d2 is the topological structure distance, d3 is the temporal similarity distance, and fu is a custom distance fusion operator;

[0123] (2) Perform hierarchical clustering: C = argmin∑∑(D(xi,xj)×δ(ci,cj)) where δ(ci,cj) is the cluster indicator function;

[0124] Get the initial clustering result C;

[0125] S23: Clustering result optimization

[0126] Read the initial clustering result C;

[0127] (1) Calculate the cluster boundary function: B(x) = ∑(pi × exp(-||x-μi||2 / σi 2 )) Where μi is the cluster center and σi is the variable variance;

[0128] (2) Boundary optimization: min x ∑||B(x)-B*(x)|| 2+λR(x) where B*(x) is the ideal boundary function and R(x) is the regularization term: R(x)=||gradB(x)||1+||grad2B(x)||2;

[0129] Get the optimized clustering result C';

[0130] S24: Source-grid-load characteristic combination

[0131] Read the optimized clustering result C';

[0132] (1) Construct feature combination matrix: Where fu is the tensor product operation ψ, is a nonlinear mapping function: ψ(x) = σ(Wx+b) × tanh(Ux+v);

[0133] (2) Optimize the combination weight: W* = argmin∑||G(W)-G*||2+μΩ(W) where Ω(W) is a new type of sparse regularization term: Ω(W) = ||W||1+α||W||2,1;

[0134] Get the feature combination matrix G and weight matrix W;

[0135] S25: Quality Assessment and Feedback

[0136] Read the feature combination matrix G and weight matrix W

[0137] (1) Calculate the multi-dimensional evaluation index: Q = [Q1, Q2, Q3, Q4] where: Q1 = cluster compactness Q2 = inter-class separation Q3 = temporal continuity Q4 = spatial consistency;

[0138] (2) Construct a feedback mechanism: F(Q) = ∑(νi × Qi) + ξ × H(Q) where H(Q) is the entropy function and νi is the weight;

[0139] Obtain a quality assessment vector Q and a feedback signal F;

[0140] S3: Multiscale Energy Storage Modeling

[0141] S31: Energy Storage Characteristics Modeling

[0142] Read the feature combination matrix G* and quality assessment vector Q output by S2

[0143] (1) Construct the energy storage dynamic characteristic matrix: E(t,s) = ∑(ρi × Ψi(t,s)) where Ψi(t,s) is the energy storage characteristic basis function: is the time dimension, s is the space dimension ρi,Ai,λi,ωi, is the adaptive parameter;

[0144] (2) Calculation characteristic constraints: C(E) = {SOC(t), P(t), η(t), L(t)} where: SOC(t) is the state of charge function, P(t) is the power characteristic function, η(t) is the efficiency function, and L(t) is the life characteristic function;

[0145] Obtain the energy storage characteristic matrix E and constraint set C(E);

[0146] S32: Spatiotemporal Scale Mapping

[0147] Read the energy storage characteristic matrix E and constraint set C(E);

[0148] (1) Construct the space-time mapping function: M(t,s) = Θ(t)fuΦ(s) where: Θ(t) is the time scale mapping: Θ(t) = ∑(αi×Bi(t))·exp(-||t-ti|| 2 / σi 2 )Φ(s) is the spatial scale mapping: Φ(s)=∑(βj×Dj(s))·exp(-||s-sj|| 2 / τj 2 );

[0149] (2) Perform scale conversion: T(E) = M(t,s)·E·M(t,s) T ;

[0150] Get the multi-scale mapping matrix T(E);

[0151] S33: Network topology embedding

[0152] Read the multi-scale mapping matrix T(E);

[0153] (1) Construct the network topology matrix: N(v,e) = Γ(v)fuΩ(e) where: Γ(v) is the node feature function: Γ(v) = σ(Wv·v+bv)·tanh(Uv·v+cv)Ω(e) is the edge feature function: Ω(e) = σ(We·e+be)·tanh(Ue·e+ce); (2) Calculate the topological embedding: H = N(v,e)fu T(E);

[0154] Get the topological embedding matrix H;

[0155] S34: Energy storage configuration model construction

[0156] Read the topological embedding matrix H;

[0157] (1) Construct the objective function set: F = {f1, f2, f3, f4} where: f1 = economic goal: min(CAPEX+OPEX) f2 = stability goal: max(System Stability) f3 = reliability goal: max(Reliability) f4 = life goal: max(Lifecycle); (2) Construct the constraint conditions: g1(x): SOCmin≤SOC(t)≤SOCmax g2(x): -Pmax≤P(t)≤Pmax g3(x): η(t)≥ηmin g4(x): L(t)≥Lmin;

[0158] Get the optimized model {F, g(x)};

[0159] S35: Model verification and adjustment

[0160] Read the optimization model {F, g(x)};

[0161] (1) Construct verification indicators: V = [V1, V2, V3, V4] where: V1 = model accuracy V2 = computational efficiency V3 = robustness V4 = convergence;

[0162] (2) Perform model adjustment: δM = χ(V)·gradM, where χ(V) is the adjustment function and gradM is the model gradient;

[0163] Get the verified optimization model {F', g'(x)} and verification index V;

[0164] S4: Multi-objective optimization solution

[0165] S41: Optimization Problem Decomposition

[0166] Read the optimization model {F', g'(x)} and verification index V output by S3;

[0167] (1) Construct decomposition structure: Π(F',g') = {π1,π2,...,π k} Among them, π i For the sub-problem: π i ={F i (x),g i (x)} Decomposition function: F i (x) = ∑(μij × fj(x))g i (x) = ∑(νij × g'j(x)); (2) Calculate the coupling matrix: Get the subproblem set Π and coupling matrix R;

[0168] S42: Optimization algorithm construction

[0169] Read the subproblem set Π and the coupling matrix R;

[0170] (1) Construct a hybrid evolution operator: E(x) = ω1O1(x) fu ω2O2(x) fu ω3O3(x) where: O1 is the crossover operator: O1(x) = λx1+(1-λ)x2+ε(t) O2 is the adaptive mutation operator: O2(x) = x+N(0,σ(t))·D(x) O3 is the local search operator: O3(x) = x-α(t)·gradf(x); (2) Design a selection mechanism: S(x) = exp(Q(x) / T(t)) / ∑exp(Q(x) / T(t)) where Q(x) is the fitness function and T(t) is the adaptive temperature parameter;

[0171] The optimization operator E and the selection function S are obtained; fu is a direct sum operator / fusion operator.

[0172] S43: Parallel solution implementation

[0173] Read the optimization operator E and the selection function S;

[0174] (1) Construct a parallel framework: P(t) = {p1(t), p2(t), ..., pm(t)}, where pi(t) is a parallel process: pi(t) = {xi(t), vi(t), ai(t)} xi is the position vector, vi is the velocity vector, and ai is the acceleration vector;

[0175] (2) Design a collaborative mechanism: C(P) = ∑(θi × pi) K(P) where K(P) is the collaborative kernel function: K(P) = exp(-||PP*||2 / η(t));

[0176] Get the parallel solution set P*;

[0177] S44: Convergence and screening of solutions

[0178] Read the parallel solution set P*;

[0179] (1) Construct the convergence criterion: ε(t) = ||P*(t)-P*(t-1)|| / ||P*(t)|| Convergence condition: ε(t)≤ε0 and |gradF(P*)|≤δ0;

[0180] (2) Design solution screening criteria: Z(x) = ∑(γi × zi(x)) where: z1 = objective function value z2 = constraint satisfaction z3 = solution diversity z4 = solution stability;

[0181] Get the optimal solution set X*;

[0182] S45: Verification and optimization of solutions;

[0183] Read the optimal solution set X*;

[0184] (1) Construct verification indicators: W = [w1, w2, w3, w4] where: w1 = economic indicator w2 = technical feasibility w3 = operational stability w4 = system adaptability; (2) Perform solution optimization: X** = X* + β(t)·gradW(X*) where β(t) is the adaptive step size and grad is the gradient operator.

[0185] The final optimized solution X** and evaluation index W are obtained.

[0186] In some embodiments, the relevant formulas and parameters are as follows:

[0187] Correct the wind power according to the temperature data. P_corrected = P_original × [1 + α (T-T_ref)];

[0188] Where: P_corrected is the corrected power; P_original is the original power; α is the temperature correction coefficient (usually -0.004 / ℃); T is the actual temperature; T_ref is the reference temperature (usually 25℃);

[0189] The calculation of the judgment threshold is specifically: threshold(t) = μ(t) ± k(t) × σ(t);

[0190] The calculation formula of k(t) is: k(t) = k_base × [1 + β × log (ρ(t))]; k_base is the base coefficient (usually 2.5); β is the adjustment parameter (value range [0.1, 0.5]); ρ(t) is the data density at time t;

[0191] Composite kernel function, specifically: K(x,y)=w1×KGaussian(x,y)+w2×KPolynomial(x,y);

[0192] Where: KGaussian(x,y)=exp(-||xy|| 2 / σ 2 ); KPolynomial(x,y)=(x·y+c)^d;

[0193] Parameter setting: w1, w2 are weight coefficients, w1+w2=1;

[0194] σ is the Gaussian kernel bandwidth; c is the polynomial kernel bias term; d is the polynomial kernel degree;

[0195] The fusion operator can be: D(x,y)=θ1d1(x,y)fuθ2d2(x,y)fuθ3d3(x,y);

[0196] The specific implementation is: D(x,y)=[θ1d1^p+θ2d2^p+θ3d3^p]^(1 / p);

[0197] Where: p is the distance fusion parameter (p≥1); θi is the adaptive weight, satisfying ∑θi=1;

[0198] Weight update rule: θi = exp(-λdi) / ∑exp(-λdi); λ is the temperature parameter, which is dynamically adjusted during the iteration process;

[0199] Boundary function B(x), specifically: B(x) = ∑(pi × exp(-||x-μi||2 / σi 2 ));

[0200] Optimization objective function: J(x) = ||B(x)B*(x)||2+λ1||gradB(x)||1+λ2||grad2B(x)||2;

[0201] Where: pi is the weight of the cluster center; μi is the cluster center; σi is the variable variance parameter; λ1,λ2 are regularization parameters;

[0202] SOC, power and efficiency functions can be: SOC(t) = SOC(t-1) + [η_c × P_c(t) P_d(t) / η_d] × Δt / E;

[0203] P(t) constraint: |P(t)|≤Pmax×[1α×(1-SOC(t))];

[0204] Efficiency function: η(t) = η_0 × [1β × |P(t) / Pmax|];

[0205] Among them: η_c, η_d are charge and discharge efficiency; P_c, P_d are charge and discharge power; E is rated capacity; α, β are attenuation coefficients;

[0206] The space-time scale mapping function is as follows: Θ(t) = ∑αi × Bi(t) × exp(-||t-ti|| 2 / σi 2 ); Bi(t)=sin(ωit)×exp(-λit);

[0207] Space-time mapping matrix: M(t,s)=Θ(t)fuΦ(s);=[Θ1(t)Φ1(s)...Θn(t)Φm(s)];

[0208] Hybrid evolution operator, specifically:

[0209] 1. Improved crossover operator: O1(x) = λx1 + (1-λ)x2 + ε(t); ε(t) = ε_0 × exp(-t / τ) × N(0,1);

[0210] 2. Adaptive mutation operator: O2(x)=x+N(0,σ(t))×D(x);σ(t)=σ_0×(1-t / T)^η;

[0211] 3. Local search operator: O3(x) = xα(t) × gradf(x); α(t) = α_0 / (1+γt);

[0212] The collaborative kernel function is: K(P) = exp(-||PP*|| 2 / n(t));n(t)=n_0×(1-t / T)^β;

[0213] Synergistic effect strength: S(i,j)=K(Pi,Pj)×exp(-dij / r); where:

[0214] dij is the distance between processes i and j; r is the coordination range parameter;

[0215] Multi-dimensional evaluation indicators, including: Z(x) = ∑(γi × zi(x));

[0216] Where: z1 = (f(x)f_min) / (f_max f_min); z2 = ∑max(0,gi(x)) / n; z3 = ∑min(||x-xi||) / |P|;

[0217] z4=1var(f(x)) / f_mean;

[0218] Kernel function parameters: Gaussian kernel bandwidth σ is estimated using the median method; polynomial kernel degree d is 2-3; the initial value of the weight coefficient is w1=w2=0.5;

[0219] Distance fusion parameters: p is 2 (Euclidean fusion) or 1 (Manhattan fusion); initial temperature parameter λ is 1.0; temperature attenuation rate is 0.95; energy storage characteristic parameters: charge and discharge efficiency η_c, η_d is 0.9-0.95; power attenuation coefficient α is 0.1-0.2; efficiency attenuation coefficient β is 0.05-0.1;

[0220] Optimization algorithm parameters: crossover probability is 0.8-0.9; mutation probability is 0.1-0.2; population size is 50-100; maximum number of iterations is 1000.

[0221] The preferred embodiments of the present invention are described in detail above; however, the present invention is not limited to the specific details in the above embodiments. Within the technical concept of the present invention, various equivalent transformations can be made to the technical solutions of the present invention, and these equivalent transformations all belong to the protection scope of the present invention.

Claims

1. A new type of power system multi-scale grid-connected energy storage planning optimization method based on source-grid-load clustering, characterized in that: The steps include: Step S1, collecting basic data and organizing it into a three-dimensional data structure according to the time dimension, space dimension and variable type dimension, based on the three-dimensional data structure, calculating the local probability density value of each data point and performing data cleaning; performing kernel mapping and feature decomposition on the cleaned data to obtain core tensors and factor matrices, using the core tensors and factor matrices to calculate spatiotemporal correlation features, generating fused feature vectors, calculating information entropy and mutual information for the fused feature vectors, screening out important features, and forming an optimized feature matrix; Step S2, input each eigenvector in the optimized feature matrix into a nonlinear basis function to calculate the eigenmanifold value; Calculate the local structural distance for each pair of feature manifold values ​​and construct a structural matrix; based on the structural matrix, calculate the Wasserstein distance, topological structure distance and time series similarity distance respectively, and combine these three distances into a multi-level distance through a distance fusion operator; perform clustering operations using multi-level distances to obtain initial clustering results; optimize and adjust the boundary area of ​​the clustering results to obtain optimized clustering results; construct a feature combination matrix based on the optimized clustering results, optimize the combination weights, and output the final feature combination results; Step S3, input the final feature combination result into the energy storage dynamic characteristic function, calculate the energy storage characteristic matrix by combining the exponential function and the trigonometric function; calculate the state of charge function value, power characteristic function value, efficiency function value and life characteristic function value of the energy storage device based on the energy storage characteristic matrix; construct the time scale mapping and space scale mapping using the energy storage characteristic matrix, and perform the scale conversion operation to obtain the multi-scale mapping matrix; Combine the multi-scale mapping matrix with the network topology features to generate a topology embedding matrix; Build a multi-objective optimization model based on the topological embedding matrix, perform model verification and adjustment, and output the verified optimization model; Step S4, decompose the verified optimization model into multiple sub-problems, calculate the coupling relationship matrix between the sub-problems; construct a hybrid evolution operator for each sub-problem, including a crossover operator, an adaptive mutation operator, and a local search operator, and design a selection mechanism; build a parallel solution framework to collaboratively solve the sub-problems; make convergence judgments and multi-dimensional evaluations on the solution results to screen the optimal solution set; verify and optimize the optimal solution set, and output the final optimized solution.

2. The novel multi-scale grid-building energy storage planning optimization method for power system with source-grid-load clustering according to claim 1 is characterized in that: The step S1 is specifically as follows: Step S11, reading the power generation value, geographical location coordinates and operating status parameters of each wind turbine at different time points from the power system database; reading the wind speed data and temperature data corresponding to the geographical location of each wind turbine from the meteorological database; reading the power generation value, geographical location coordinates and operating status parameters of each photovoltaic array at different time points from the photovoltaic system database; reading the power consumption value and power consumption type identification at different time points from the load monitoring system; organizing the read data into a three-dimensional data structure according to the time point, geographical location and data type, wherein each data point contains a corresponding time tag, location tag and type tag; Step S12: for each data point in the three-dimensional data structure obtained in step S11, extract the data of five position points that are adjacent in space and 12 hours before and after the time position of the point, and calculate the mean and standard deviation of these data; calculate the probability density value of each data point and the data sample set extracted from it; calculate the time-varying adjustment coefficient based on the probability density distribution characteristics of the data point 24 hours before and after; set the product of the time-varying adjustment coefficient and the standard deviation as the determination threshold of the current data point; mark the data points that exceed the determination threshold as abnormal points; for each abnormal point, extract the normal data points in the 24-hour time window before and after it, calculate the weighted average of these normal data points, and replace the data value of the abnormal point with the calculated weighted average; Step S13, for each data point in the data set obtained in step S12, respectively calculate the Gaussian kernel value and the polynomial kernel value between it and all other data points; dynamically calculate the variance parameter value of the Gaussian kernel according to the distribution range of the data points, and dynamically determine the order parameter value of the polynomial kernel according to the change trend of the data; combine the calculated Gaussian kernel value and polynomial kernel value into a composite kernel value through a weight coefficient; Construct a kernel matrix with the composite kernel values ​​between all pairs of data points; perform eigenvalue decomposition on the kernel matrix to extract eigenvectors and eigenvalues; Reconstruct data based on eigenvectors and eigenvalues, calculate reconstruction error values, set error thresholds to filter eigenvectors, and obtain the reduced-dimensional core tensor and factor matrix; Step S14, separating the factor matrix of the time dimension and the factor matrix of the space dimension from the factor matrix obtained in step S13; calculating the time correlation coefficient for the time dimension factor matrix, and calculating the space correlation coefficient for the space dimension factor matrix; calculating the weight coefficient of the time dimension according to the time-varying characteristics of the data, and calculating the weight coefficient of the space dimension according to the spatial distribution characteristics of the data; multiplying the time correlation coefficient with the time weight coefficient to obtain the time correlation matrix, and multiplying the space correlation coefficient with the space weight coefficient to obtain the space correlation matrix; superimposing the time correlation matrix and the space correlation matrix to form a comprehensive correlation matrix; multiplying the comprehensive correlation matrix with the core tensor to generate a fusion feature vector; Step S15: for each feature component in the fused feature vector obtained in step S14, count the occurrence frequency of its value and calculate the probability distribution; Calculate the information entropy value of each feature component according to the probability distribution; calculate the mutual information value between each feature component and all other feature components; combine the information entropy value and the mutual information value through the weight parameter to obtain the importance score of each feature component; The feature components are sorted in descending order according to their importance scores, and the feature components are selected from high to low; the matrix rank of the selected feature subset is calculated, and the selection is stopped when the matrix rank reaches a preset threshold; the selected feature components are combined into an optimized feature matrix.

3. The novel multi-scale grid-building energy storage planning optimization method for power system with source-grid-load clustering according to claim 1 is characterized in that: The step S2 is specifically as follows: Step S21, input each eigenvector in the optimized feature matrix output by step S1 into the hyperbolic tangent function and the sine function to obtain a nonlinear transformation value; dynamically calculate the scaling parameter of the hyperbolic tangent function and the frequency parameter of the sine function according to the numerical distribution of the eigenvector; combine the hyperbolic tangent function value with the sine function value through an adaptive weight coefficient to obtain a characteristic manifold value; calculate the Euclidean distance for each pair of characteristic manifold values, divide the Euclidean distance by the dynamically calculated distance parameter to obtain a local structural distance; form a structural matrix with the local structural distances between all pairs of characteristic manifold values; Step S22, extracting the distance value between each pair of data points from the structure matrix obtained in step S21; For each pair of data points, the Wasserstein distance between their probability distributions is calculated, the shortest path length of the data point in the network is used as the topological distance, and the dynamic time warping distance of the data point time series is calculated as the time series similarity distance; the weight coefficients of the three distances are dynamically calculated according to the data distribution characteristics; the weighted three distances are combined into multi-level distances through the distance fusion operator; hierarchical clustering operations are performed based on the multi-level distances to divide the data points into different clusters to form the initial clustering results; Step S23, for each cluster obtained in step S22, extracting data points in its boundary area; For each boundary data point, calculate its membership with the current cluster center and the adjacent cluster center; dynamically calculate the variance parameter of the cluster center according to the local distribution characteristics of the data point; substitute the distance from the data point to the cluster center, the variance parameter and the membership into the boundary function to calculate the boundary function value; substitute the difference between the actual boundary function value and the ideal boundary function value, the first-order derivative and the second-order derivative of the boundary function into the optimization objective function; adjust the attribution of the boundary data points by minimizing the optimization objective function to obtain the optimized clustering result; Step S24, inputting the feature vector of each pair of clusters in the clustering result optimized in step S23 into a nonlinear mapping function, respectively calculating the sigmoid function value and the hyperbolic tangent function value; multiplying the sigmoid function value and the hyperbolic tangent function value to obtain a nonlinear feature combination value; Construct a feature combination matrix for the feature combination values ​​of all cluster pairs; Calculate the one-norm and 2,1-norm of the feature combination matrix, and substitute the norm value into the sparse regularization term; obtain the initially optimized feature combination weight matrix by minimizing the sum of squares of the difference between the feature combination matrix and the target matrix and the weighted sum of the sparse regularization term; Step S25, for the feature combination weight matrix obtained in step S24, the average distance between the data points in each cluster is calculated to obtain the cluster compactness, the minimum distance between different cluster centers is calculated to obtain the inter-class separation, the time series continuity of the data points in the cluster is calculated to obtain the temporal continuity, and the spatial distribution consistency of the clusters is calculated to obtain the spatial consistency; These four indicators are combined into a quality assessment value through weight coefficients; the entropy value of the quality assessment value is calculated, and the entropy value is added to the weighted combination value to obtain a feedback signal; the feature combination weight matrix is ​​fine-tuned according to the feedback signal to obtain the final feature combination result.

4. The novel multi-scale grid-building energy storage planning optimization method for power system with source-grid-load clustering according to claim 1 is characterized in that: The step S3 is specifically as follows: Step S31, expand the feature combination result output from step S2 according to the time dimension and the space dimension; for each data point in the time and space dimension, calculate the exponential decay function value and the cosine function value; dynamically calculate the attenuation coefficient of the exponential decay function and the frequency and phase parameters of the cosine function according to the distribution characteristics of the time and space data; multiply the exponential decay function value with the cosine function value to generate the energy storage characteristic base function value; perform weighted combination of the energy storage characteristic base function values ​​of all time and space dimension points to obtain the energy storage characteristic matrix; calculate the state of charge value, charge and discharge power value, conversion efficiency value and remaining life value at each moment according to the energy storage characteristic matrix to form a set of energy storage constraint conditions; Step S32, extracting the time dimension data and the space dimension data from the energy storage characteristic matrix obtained in step S31 respectively; inputting the time dimension data into the basis function, calculating the mapping value of each time point, and multiplying the mapping value with the Gaussian function value to obtain the time scale mapping value; inputting the space dimension data into the basis function, calculating the mapping value of each space point, and multiplying the mapping value with the Gaussian function value to obtain the space scale mapping value; dynamically calculating the variance parameter of the Gaussian function according to the data distribution characteristics of the time dimension and the space dimension; performing tensor product operation on the time scale mapping value and the space scale mapping value to generate a scale mapping function matrix; performing matrix multiplication operation on the scale mapping function matrix and the energy storage characteristic matrix to obtain a multi-scale mapping matrix; Step S33, extracting node data and edge data of the power grid; inputting the node data into a nonlinear mapping function, calculating the product of the sigmoid function value and the hyperbolic tangent function value to obtain the node eigenvalue; inputting the edge data into a nonlinear mapping function, calculating the product of the sigmoid function value and the hyperbolic tangent function value to obtain the edge eigenvalue; Constructing a node feature matrix according to the node eigenvalues, and constructing an edge feature matrix according to the edge eigenvalues; performing a weighted combination of the node feature matrix and the edge feature matrix to obtain a network topology matrix; performing a tensor product operation on the network topology matrix and the multi-scale mapping matrix obtained in step S32 to generate a topology embedding matrix; Step S34, based on the topological embedding matrix obtained in step S33, respectively calculate the sum of the construction cost and operation and maintenance cost of the energy storage system, the system stability index, the reliability index and the life index; construct these indexes into an objective function set; extract the state of charge constraint, power constraint, efficiency constraint and life constraint from the energy storage constraint condition set in step S31; combine the objective function set and the constraint condition set into an optimization model; Step S35, for the optimization model constructed in step S34, respectively calculate the accuracy index, computational efficiency index, robustness index and convergence index of the model; combine these indexes into a verification index vector; calculate the gradient value of the verification index vector; adjust the parameters of the optimization model according to the gradient value; repeat the process of calculating the verification index and adjusting the parameters until the verification index reaches a preset threshold, and output the verified optimization model and the verification index vector.

5. The novel multi-scale grid-building energy storage planning optimization method for power system with source-grid-load clustering according to claim 1 is characterized in that: The step S4 is specifically as follows: Step S41, grouping the objective functions and constraints in the optimization model outputted from step S3 according to coupling relationships; calculating weight coefficients for each group of objective functions to generate objective functions of subproblems; calculating weight coefficients for each group of constraints to generate constraints of subproblems; combining the objective functions and constraints of each subproblem into a complete subproblem; calculating first-order and second-order derivatives of the objective functions and constraints between each pair of subproblems; combining the derivative values ​​into a coupling relationship matrix, and outputting the decomposed subproblem set and the coupling relationship matrix; Step S42: for each sub-problem outputted in step S41, two candidate solutions are selected, a new solution is generated by linear combination, and a time-varying disturbance term is added to the new solution to obtain an output solution of the crossover operator; Add random perturbations that obey normal distribution to the candidate solutions, and adjust the perturbation amplitude dynamically during the optimization process to obtain the output solution of the mutation operator. Calculate the objective function gradient of the candidate solution, move in the opposite direction of the gradient with an adaptive step size, and obtain the output solution of the local search operator; calculate the fitness function value of the output solution of the three operators, calculate the selection probability according to the fitness value and the adaptive temperature parameter, and output the three operators and their selection mechanism; Step S43, dividing the sub-problem set outputted in step S41 into a parallel computing framework; For each parallel process, the position vector, velocity vector and acceleration vector of the solution are recorded; the synergy coefficient between processes is calculated based on the coupling relationship matrix of the sub-problems; the Euclidean distance between the current solution set and the optimal solution set is substituted into the synergy kernel function to obtain the synergy intensity between processes; the evolution direction and step size of each process are adjusted according to the synergy intensity, parallel solution operations are performed, and the parallel solution set is output; Step S44, extracting the solution sets of two adjacent iterations from the parallel solution sets outputted in step S43; Calculate the relative error of two solution sets and the gradient value of the objective function; Compare the relative error with the preset convergence threshold, and compare the gradient value with the preset gradient threshold to determine whether the convergence condition is met; For the converged solution set, the objective function value, constraint satisfaction degree, solution diversity and solution stability are calculated respectively; The solution set is screened according to the multi-dimensional evaluation results, and the optimal solution set is output; Step S45, for the optimal solution set output by step S44, respectively calculate the economic index, technical feasibility index, operation stability index and system adaptability index; combine the calculation results of various indicators into a verification index vector; calculate the gradient value of the verification index vector; determine the optimization direction and step size according to the gradient value, and make local adjustments to each solution in the optimal solution set; repeat the verification and adjustment process until the verification index meets the preset requirements, and output the final optimized solution and evaluation index.

6. The novel multi-scale grid-connected energy storage planning optimization method for power system with source-grid-load clustering according to claim 2 is characterized in that: The step S11 is specifically as follows: Step S111, extracting 10-minute power data, geographic coordinate data and operating status data of each wind turbine from the power system database, and extracting wind speed data and temperature data at 10-minute time intervals corresponding to the geographic location from the meteorological database; aligning the extracted wind speed data and temperature data with the wind power data according to the timestamp; calculating the correlation coefficient between wind power and wind speed, and eliminating data points with correlation coefficients below a threshold; correcting the wind power according to the temperature data to obtain a temperature-corrected wind power data set; Step S112, extracting 10-minute power data, geographic coordinate data and operating status data of each photovoltaic array from the photovoltaic system database, and extracting irradiance data and component temperature data at 10-minute time intervals corresponding to the geographic location from the meteorological database; aligning the irradiance data and component temperature data with the photovoltaic power data according to the timestamp; Calculate the ratio of photovoltaic power and irradiance to obtain the photoelectric conversion efficiency sequence; calculate the temperature correction coefficient based on the photoelectric conversion efficiency sequence and component temperature data to obtain the corrected photovoltaic power data set; Step S113, extracting 10-minute load data, power consumption type identification and geographic location data of each power consumption node from the load monitoring system; The load data is divided into industrial load, commercial load and residential load according to the power consumption type identification; the power factor is calculated for each type of load data, and the data points with power factor lower than the preset threshold are marked; For the marked data points, extract the valid data of 12 hours before and after to calculate the correction coefficient, and multiply the correction coefficient with the original load data to obtain the corrected load data; merge the three types of corrected load data according to timestamp and geographic location to obtain a complete load data set; Step S114, aligning the wind power data set outputted from step S111, the photovoltaic power data set outputted from step S112, and the load data set outputted from step S113 according to the time dimension; calculating the mean and standard deviation for each data set respectively; subtracting the mean of the corresponding data set from each data point, and then dividing by the standard deviation to obtain a standardized data value; Perform interval mapping on the standardized data values ​​and adjust the value range to the [0,1] interval; record the mapped data as a standardized wind power data set, a standardized photovoltaic data set, and a standardized load data set; Step S115, sorting the three standardized data sets outputted from step S114 according to timestamps and geographic locations; constructing a three-dimensional array structure, wherein the first dimension represents the time series, the second dimension represents the geographic location, and the third dimension represents the data type; For each time slice in the three-dimensional array, calculate the correlation coefficient between different geographical locations; for each geographical location slice, calculate the autocorrelation coefficient at different time points; For each data type slice, calculate the spatiotemporal cross-correlation coefficient; The calculated correlation coefficient is used as the weight coefficient to perform weighted processing on the three-dimensional array to obtain the final three-dimensional data structure.

7. The novel multi-scale grid-building energy storage planning optimization method for power system with source-grid-load clustering according to claim 2 is characterized in that: The step S13 is specifically as follows: Step S131, for the cleaned data set outputted in step S12, calculate the kurtosis value and skewness value of the data distribution; dynamically determine the variance parameter of the Gaussian kernel function based on the kurtosis value, and dynamically determine the order parameter of the polynomial kernel function based on the skewness value; Perform block processing on the data set and calculate the local statistical features of each data block; construct an adaptive weight coefficient based on the local statistical features, multiply the Gaussian kernel parameters and the polynomial kernel parameters with the weight coefficient to obtain a modified kernel function parameter set; Step S132: using the kernel function parameter set outputted in step S131, calculate the Gaussian kernel value and the polynomial kernel value between each pair of data points in the data set; construct a distance decay function according to the spatiotemporal position relationship of the data points, and multiply the distance decay function value by the kernel function value to obtain a modified kernel value; Construct a symmetric matrix for the modified kernel value and test the semi-positive definiteness of the matrix; Perform eigenvalue correction on the kernel matrix that does not satisfy semi-positive definiteness to obtain the corrected kernel matrix; Step S133, performing eigenvalue decomposition on the kernel matrix outputted in step S132; calculating the cumulative contribution rate of the eigenvalue sequence obtained by the decomposition; extracting a subset of eigenvectors corresponding to the dominant eigenvalues ​​from the eigenvector matrix; For each eigenvector in the eigenvector subset, the projection value between the eigenvector and the original data is calculated; the projection values ​​are combined into a feature projection matrix to obtain a low-dimensional representation of the data; Step S134, reconstructing the low-dimensional feature representation outputted in step S133 into the original data dimension through reverse mapping; calculating the mean square error and relative error between the reconstructed data and the original data; calculating the statistical characteristic value of the reconstructed error sequence; constructing an error evaluation function based on statistical characteristics, and determining the optimal number of features according to the evaluation function value; Step S135: Based on the optimal number of features determined in step S134, select a corresponding number of eigenvectors from the eigenvector matrix; calculate the orthogonality index of the selected eigenvectors; construct a feature importance evaluation function, comprehensively considering the orthogonality of the eigenvectors and the contribution rate of the corresponding eigenvalues; sort and screen the eigenvectors according to the evaluation results, and output the final core tensor and factor matrix.

8. The novel power system multi-scale network construction and energy storage planning optimization method based on source-grid-load clustering according to claim 3 is characterized in that: The step S23 is specifically as follows: Step S231, extracting a set of data points of each cluster from the initial clustering result outputted in step S22; calculating the coordinates of the center point and the distribution range of each cluster; determining the search radius according to the distribution range, and extracting data points within a range with the cluster center as the circle center and the search radius as the radius; Calculate the distance between the extracted data point and the cluster center, and arrange the distance values ​​in ascending order; set a distance threshold, mark the data points with a distance greater than the threshold as boundary candidate points, and obtain a set of boundary candidate points; Step S232: for each data point in the set of boundary candidate points outputted in step S231, extract the center point of the cluster to which it currently belongs and the center points of the adjacent clusters; calculate the Euclidean distance between the data point and the centers of each cluster; construct a radial basis function based on the distance value, and calculate the membership of the data point to each cluster center; adjust the membership value according to the local density distribution characteristics of the data point; and form a membership matrix with the adjusted membership values; Step S233, constructing a boundary function using the membership matrix outputted in step S232; extracting data points within a predetermined range around each boundary candidate point, and calculating local data density; An adaptive variance function is constructed based on local data density; the membership value and the adaptive variance are substituted into the exponential function to calculate the boundary strength value; the boundary strength value is normalized to obtain a boundary function value sequence; Step S234, calculating the spatial gradient and the temporal gradient of the boundary function value sequence outputted in step S233; constructing a gradient-based regularization term, taking the absolute value of the first-order derivative and the sum of the squares of the second-order derivative as a penalty term; The square of the difference between the actual boundary function value and the ideal boundary function value is used as the fitting term; the weight coefficients of the fitting term and the penalty term are dynamically adjusted according to the data distribution characteristics to construct a complete optimization objective function; Step S235, perform minimization operation on the optimization objective function constructed in step S234; adjust the clustering affiliation of the boundary candidate points according to the optimization results; recalculate the cluster center and distribution range of the adjusted clustering results; check the compactness and separation of the adjusted clustering, and when the compactness and separation meet the preset conditions, output the optimized clustering results.

9. The novel multi-scale grid-building energy storage planning optimization method for power system with source-grid-load clustering according to claim 4 is characterized in that: The step S31 is specifically as follows: Step S311, expand the final feature combination result outputted from step S2 according to time series and spatial position; construct exponential basis function and triangular basis function, the exponential basis function is used to describe the energy storage attenuation characteristics, and the triangular basis function is used to describe the energy storage cycle characteristics; construct a dynamic basis function according to the time distribution characteristics of the data, and modify the original basis function; output the modified basis function set as the basic description function group of the energy storage characteristics; Step S312: extracting the variation characteristics of the data at different time scales and spatial scales from the basis function group outputted in step S311; Calculate the data volatility on the time scale and the distribution discreteness on the spatial scale; dynamically adjust the attenuation coefficient of the exponential basis function according to the volatility, and dynamically adjust the frequency parameters and phase parameters of the triangular basis function according to the discreteness; construct a parameter adaptive mapping function to realize the dynamic update of the basis function parameters, and output the updated basis function parameter set; Step S313, using the basis function parameter set outputted in step S312, calculate the characteristic function value of each spatiotemporal data point; construct a weighted combination of the characteristic function values, and dynamically allocate the weight coefficient according to the importance of the data; Arrange the weighted characteristic function values ​​into a matrix form; calculate the singular values ​​of the matrix, and perform noise reduction on the characteristic matrix according to the distribution characteristics of the singular values; Performing a standardization operation on the processed characteristic matrix and outputting an energy storage characteristic matrix; Step S314, extracting time series data from the energy storage characteristic matrix outputted in step S313; Calculate the charge state change rate of the energy storage device at different times to obtain the charge state function; calculate the charge and discharge power change characteristics based on the charge state function to obtain the power characteristic function; calculate the energy conversion efficiency based on the power characteristic function to obtain the efficiency function; build a life prediction model based on historical operation data and calculate the remaining life function; combine these function values ​​into a constraint parameter sequence; Step S315, classify and organize the constraint parameter sequence output from step S314; construct a state of charge constraint interval, set a power change constraint range, determine an efficiency constraint threshold, and establish a life constraint boundary; normalize these constraints; construct a membership function of the constraint conditions, and calculate the degree of constraint satisfaction; and combine the constraint conditions and their membership functions into a complete constraint set.

10. The novel multi-scale grid-building energy storage planning optimization method for power system with source-grid-load clustering according to claim 4 is characterized in that: The step S34 is specifically as follows: Step S341, extracting cost features, stability features, reliability features and life features from the topological embedding matrix output in step S33; converting the construction cost and operation and maintenance cost data into an annualized cost series; Construct a stability evaluation sequence based on system status data; calculate a reliability index sequence based on failure rate and repair time; generate a life prediction sequence using equipment aging characteristics; combine these sequences to construct a multi-objective function matrix; Step S342, extracting various constraints from the energy storage constraint set outputted in step S31; converting the state of charge constraint into a state variable inequality; converting the power constraint into a power change rate inequality; converting the efficiency constraint into an energy conversion efficiency inequality; converting the life constraint into a cycle number inequality; normalizing all inequality constraints and constructing a constraint matrix; Step S343, normalizing the multi-objective function matrix outputted in step S341; calculating the value range and distribution characteristics of each objective function; A normalization function is constructed based on the distribution characteristics to convert the objective functions of different dimensions into a unified interval; the correlation coefficient matrix is ​​calculated for the converted objective functions; The objective function weight coefficient is constructed according to the correlation coefficient to obtain a normalized multi-objective function set; Step S344, performing normalization transformation on the constraint matrix outputted in step S342; extracting the boundary value and critical value of the constraint; constructing slack variables of the constraint, converting the hard constraint into a soft constraint with a penalty term; calculating the compatibility index between the constraint; adjusting the constraint according to the compatibility index, and outputting a normalized constraint set; Step S345, integrate the normalized multi-objective function set output by step S343 and the normalized constraint condition set output by step S344; construct an association matrix of the objective function and the constraint conditions; calculate the coupling index of the model based on the association matrix; decouple the constraints with high coupling degree; assemble the processed objective functions and constraints into a complete optimization model, and construct a priority sequence for solving the model.

Images