Power system model important parameter group identification method based on trajectory feature clustering

Abstract

The invention provides an electric power system model important parameter group identification method based on trajectory feature clustering, relates to an electric power system simulation technology, and solves the problem that modeling is not performed due to time domain global error index coarse graining and parameter coupling in traditional electric power system important parameter identification. The parameters are used as to-be-identified high-sensitivity parameters; performing time domain difference on each group of simulation tracks and actually measured tracks to obtain difference tracks, and extracting difference track features; the method comprises the following steps: mapping high-dimensional data to a low-dimensional space through a nonlinear dimension reduction method, and executing K-Means clustering on the data in the low-dimensional space to obtain a cluster label of each sample; sensitivity calculation is carried out on all samples in the clusters, and a parameter combination with the highest occurrence frequency and the maximum contribution to dynamic differences in the clusters is screened out to serve as a final dominant parameter set. According to the scheme, the simulation frequency and the calculation burden can be remarkably reduced, and high-precision and interpretable discrimination of important parameters really influencing the system performance is realized.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of power system simulation, in particular to a power system model important parameter group identification method based on trajectory feature clustering. BACKGROUND

[0002] The existing power system important parameter identification mostly depends on time domain global error indicators such as mean square error (MSE), mean absolute error (MAE), correlation coefficient or dynamic time warping (DTW) distance, and the whole time sequence characteristics before and after the fault or in the disturbance process are condensed into one or a few numerical values for comparison. Although such "coarse-grained" evaluation is simple to implement and has good compatibility with existing monitoring platforms, it essentially discards rich transient information such as phase drift, peak overshoot, damping decay rate and other dynamic details, and it is difficult to capture the differential effects of various parameters on system response at different times. When multiple parameter combinations produce similar gains or deviations on the same global indicator, the identification algorithm often cannot accurately distinguish the true driving factors.

[0003] In order to reduce the complexity of the algorithm, the traditional method usually adopts a "one-by-one scanning" or empirical grouping method for discrete processing of the parameters to be identified, and assumes that the parameters do not interfere with each other during the identification process. This "parameter independent" assumption seriously ignores the complex coupling and interaction effects between the generator excitation system, the governor and the network impedance, and the load characteristics. In actual operation, a small adjustment of the governor parameter may be amplified through the excitation system feedback loop, or it may be suppressed by the parallel network impedance, resulting in non-additive single-parameter test results. That is, the coupling relationship between the important parameters is broken during the acquisition of the important parameters, which further leads to the lack of clear and sufficient physical constraints for the identified parameters, and the credibility is seriously insufficient. Even if principal component analysis (PCA) or factor decomposition is introduced to reduce the dimension, it can only partially alleviate the dimension disaster, but it still cannot establish a nonlinear coupling model between the parameters. In addition, with the popularization of PMU / SCADA high-frequency measurement and the continuous expansion of the power grid, the number of parameters to be identified is increasing, and the comparison cost of coarse-grained error calculation itself needs to be performed hundreds of times large-scale simulation, which grows exponentially, seriously restricting the feasibility of online or near real-time identification, and it is difficult to meet the actual needs of power grid for fast fault location and self-healing scheduling. SUMMARY

[0004] In order to solve the technical problems in the prior art that the traditional power system important parameter identification has the problems of coarse-grained time domain global error indicators and unmodeled parameter coupling, which leads to the lack of effectiveness of the identification results, a power system model important parameter group identification method based on trajectory feature clustering, an electronic device and a computer readable storage medium are provided, which significantly reduces the simulation times and computational burden, and realizes high-precision and explainable discrimination of important parameters that truly affect system performance.

[0005] To solve the above technical problems, in a first aspect, according to some embodiments, the present application provides a power system model important parameter group identification method based on trajectory feature clustering, comprising:

[0006] S1, Latin hypercube sampling method is used for the parameters to be identified in the power system to generate n groups of parameter perturbation samples, standard simulation is run for each group of parameter perturbation samples, corresponding multi-channel time sequence simulation trajectories are obtained, a node near a certain fault point is selected as a monitoring point, data of the node are measured, and coarse-grained sensitive parameters are screened out as high-sensitivity parameters to be identified, wherein n is a positive integer greater than or equal to 1;

[0007] S2, time domain difference is made for each group of simulation trajectories and measured trajectories to obtain difference trajectories, and difference trajectory features are extracted;

[0008] S3, for each difference trajectory, the original time domain signal and the STFT spectrum feature are input in parallel through a trained time-frequency encoder, a D-dimensional normalized embedding vector is generated by splicing and projecting, high-dimensional data is mapped to a low-dimensional space through a nonlinear dimension reduction method, K-Means clustering is performed on the data in the low-dimensional space in the manifold space, and a cluster label of each sample is obtained; the average silhouette coefficient under different labels is calculated to determine the optimal cluster number and the confidence judgment standard; the measured difference trajectory is calculated and reduced to a low-dimensional space, the minimum distance between the low-dimensional space and each cluster center is calculated, and if the confidence judgment standard is met, the measured sample is classified into the corresponding cluster, and the determination result of this category is recorded, wherein D is a positive integer greater than or equal to 2;

[0009] S4, read all samples in the cluster, calculate the sensitivity of all samples in the cluster, and calculate the parameter distribution of all perturbation samples in the cluster and the role of each parameter in feature contribution, and select the parameter combination with the highest frequency and the largest contribution to the dynamic difference in the cluster as the final dominant parameter set.

[0010] In some embodiments, the standard simulation is run for each group of parameter perturbation samples to obtain corresponding multi-channel time sequence simulation trajectories, specifically comprising:

[0011] Run power system simulation to obtain multi-channel time sequence trajectories of the target measurement point before and after the fault or during the disturbance (k) Run power system simulation to obtain multi-channel time sequence trajectories of the target measurement point before and after the fault or during the disturbance Wherein i = 1,...,C, C is the number of measurement point channels, and k is a positive integer greater than or equal to 1.

[0012] In some embodiments, the node near a certain fault point is selected as a monitoring point, the data of the node are measured, the coarse-grained sensitive parameters are screened out as high-sensitivity parameters to be identified, and specifically comprising:

[0013] Selecting a node near the fault point as a monitoring point, measuring data of the node, the data including active power P and reactive power Q information of the generator;

[0014] For the i-th channel, the energy sensitivity index is defined as shown in Formula I:

[0015]

[0016] Read the sensitivity matrix of the k-th group of experiments And do Min-Max normalization on the parameter dimension of each experiment, set the threshold τ of the sensitive parameter ∈ [0, 1], and screen out the coarse-grained sensitive parameters as high-sensitivity parameters to be identified.

[0017] In some embodiments, the time-domain difference is performed on each group of the simulation trajectory and the measured trajectory to obtain a difference trajectory, and the difference trajectory feature is extracted, specifically including:

[0018] The time-domain difference is performed on each group of simulation trajectories and the measured trajectory x meas (t) to obtain a difference signal Δx (k) ;

[0019] The difference signal Δx (k) is subjected to short-time Fourier transform in a sliding window, and the spectral average power is calculated in each analysis window using Formula II:

[0020]

[0021] Wherein, is the short-time Fourier transform complex spectrum of the n-th window, used to retain the energy distribution in the transient error;

[0022] The transformed difference signal is input into an encoder and decoded by a decoder to obtain the difference trajectory feature;

[0023] Wherein, the encoder is a parallel time-domain subnetwork and a frequency-domain subnetwork, and the structure and data flow of the encoder include a time-domain branch and a frequency-domain branch.

[0024] In some embodiments, after the short-time Fourier transform of the difference signal Δx (k) in a sliding window, before the transformed difference signal is input into an encoder, it further includes:

[0025] The difference signal Δx (k) and the corresponding frequency spectrum sequence p (k) is subjected to random cropping and linear or high-order interpolation to a fixed size, wherein the random cropping includes randomly selecting a starting point and randomly selecting a length.

[0026] In some embodiments, before inputting the transformed difference signal to the encoder, the method further includes:

[0027] Based on the original spectral features and Formula 3, a graph adjacency weight matrix is ​​constructed, as shown in Formula 3 below:

[0028]

[0029] Where p i p j σ represents the spectral power vector of the two different trajectories, and σ represents the control neighborhood range;

[0030] During training, the encoder incorporates a graph regularization loss, as shown in Formula 4:

[0031] L grap =∑ i,j A ij ||z i -z j || 2 (4)

[0032] Where z i and z j The embedding vector obtained by the self-supervised encoder;

[0033] The encoder consists of parallel time-domain and frequency-domain subnetworks. The encoder's structure and data flow include time-domain and frequency-domain branches, specifically including:

[0034] The time-domain branch receives the enhanced difference signal as input. First, the number of channels is expanded to a1 dimension through a one-dimensional convolutional layer, and the activation function used is the ReLU function;

[0035] After several layers of dilated depth residual blocks, each block takes the form shown in Formula 5:

[0036]

[0037] in This represents a depthwise separable convolution with an inflation rate of d, and BN represents batch normalization.

[0038] The temporal length is compressed to 1 using adaptive average pooling to obtain the temporal feature vector z. t ;

[0039] The input is the cropped and interpolated spectral power sequence. Concatenate along the time window to form a vector Then, the frequency domain vector z is obtained through two fully connected layers. f z f ;

[0040] The time-domain vector z t and frequency domain vector z f Concatenate into vector z cat After linear mapping and activation by the ReLU function, the new comprehensive feature vector is shown in Formula 6:

[0041]

[0042] The obtained results are L2 normalized to obtain the processed vector z, such that all embeddings fall on the unit hypersphere.

[0043] Wherein, the decoder header g of the decoder dec The normalized vector z is then reverse-mapped using Formula 7, as follows:

[0044]

[0045] The output of the decoder Compared with the original spectrum Alignment provides a target for the spectral reconstruction loss. W3 and W4 are the weight matrices of the output layer, and b3 and b4 are the bias terms.

[0046] In some embodiments, calculating the average silhouette coefficient under different labels, determining the optimal number of clusters, and the confidence level judgment criteria specifically include:

[0047] The average contour coefficient under different labels is calculated using Formula Fourteen:

[0048]

[0049] Where a (k) a (k) Let b be the average distance between sample k and other samples in the same cluster (cluster compactness). (k) Let k be the average distance between sample k and samples in the nearest neighbor cluster;

[0050] Take the average of all samples choose The sample k corresponding to the maximum value is taken as the optimal number of clusters;

[0051] The center vector of each cluster cc is μ. c Calculate the Euclidean distance from the sample within the cluster to the center. And take the 95th percentile value τ of the cluster distance in the training set. c As a standard for judging confidence level.

[0052] In some embodiments, the reading all samples in the cluster, performing sensitivity calculation on all samples in the cluster, statistically analyzing the parameter distribution of all perturbed samples in the cluster and the role of each parameter in feature contribution, screening out the parameter combination with the highest frequency of occurrence and the largest contribution to the dynamic difference in the cluster as the final dominant parameter set, specifically comprising:

[0053] Reading the difference ratio b of the parameters of all samples in the cluster and the simulation setting parameters, defining the parameter importance index I, which meets the condition I = b x E;

[0054] Obtaining the parameter importance array of the samples in the cluster;

[0055] For the kth sample, the normalized parameter importance vector meets Where m is the total number of parameters, Indicates the normalized contribution degree of the ith parameter;

[0056] According to the descending order of the numerical value, the indexes π1, π1,..., are arranged, and the first r items are accumulated until the formula sixteen is met:

[0057]

[0058] Where CUMth represents the cumulative threshold, and the threshold is defined as 0.9;

[0059] The first r parameter indexes constitute the core parameter set S (k) L2 normalization weight is performed on the parameters in the parameter set, as shown in formula seventeen, to obtain the soft core vector:

[0060]

[0061] Statistical aggregation of the soft core in the cluster, for each sample set K c = {k | c (k) =c}, the average soft core weight is calculated, which is shown in formula eighteen:

[0062]

[0063] The core frequency of occurrence meets formula nineteen:

[0064]

[0065] Where 1(·) is an indicator function;

[0066] Combining And f i (c) For each internal parameter of the cluster, the parameters are sorted according to the comprehensive score, and the top several parameters with the highest score are defined as the key parameter group of the cluster.

[0067] In a second aspect, the embodiments of the present application further provide an electronic device, comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor implements the steps of the method of any one of the first aspect when executing the program.

[0068] In a third aspect, the embodiments of the present application further provide a computer readable storage medium, having stored thereon a computer program, wherein the computer program is executable on a processor to implement the steps of the method of any one of the first aspect.

[0069] The above technical solution of the present application has at least the following beneficial technical effects: the present application adopts the full-parameter space uniform perturbation and sensitivity fast compression technology based on Latin hypercube sampling (LHS), calculates the Min-Max normalized parameter sensitivity by one-time LHS sampling of all to-be-identified parameters in combination with the trajectory energy index, can quickly filter out the most significant parameter subset affecting the system dynamics on the basis of one simulation result, directly compresses the parameter scale from hundreds or even thousands to tens of levels, and greatly reduces the calculation amount of subsequent analysis; and the present application realizes adaptive optimal cluster number determination in combination with UMAP dimension reduction and contour coefficient method, preferentially projects the high-dimensional time-frequency embedding into a low-dimensional manifold space, then performs K-Means clustering on the possible cluster number range and calculates the average contour coefficient, and autonomously selects the cluster number at the highest contour coefficient; subsequently, the confidence threshold is constructed based on the quantile value of the distance from each cluster sample to the cluster center, and the automation and quantification of trajectory category division and unknown mode / abnormality judgment are realized. BRIEF DESCRIPTION OF DRAWINGS

[0070] In order to more clearly illustrate the technical solutions in the embodiments of the present application or in the prior art, the drawings needed in the embodiments will be briefly introduced as follows. Obviously, the drawings in the following description are only some embodiments of the present application, and other drawings can be obtained by those skilled in the art without creative effort on the basis of these drawings.

[0071] Figure 1 is a power system model important parameter group identification method based on trajectory feature clustering provided by the embodiments of the present application.

[0072] Figure 2 is a 10-machine 39-node system for example verification provided by the embodiments of the present application.

[0073] Figure 3 is a parameter sensitivity statistical chart provided by the embodiments of the present application.

[0074] Figure 4 is a cluster number schematic diagram of clustering after dimension reduction and feature extraction of the encoded provided by the embodiments of the present application.

[0075] Figure 5 is a simulation result schematic diagram before and after parameter correction provided by an embodiment of the present application.

[0076] Figure 6 is a schematic block diagram of an electronic device provided by an embodiment of the present application. DETAILED DESCRIPTION

[0077] The technical solutions in the embodiments of the present application will be clearly and completely described below with reference to the drawings in the embodiments of the present application. Obviously, the described embodiments are only part of the embodiments of the present application, rather than all the embodiments of the present application. Based on the embodiments in the present application, all other embodiments obtained by a person of ordinary skill in the art without creative work fall within the protection scope of the present application.

[0078] If the description of "first", "second" and the like is involved in the embodiments of the present application, the description of "first", "second" and the like is only for the purpose of description, and cannot be understood as indicating or implying the relative importance of the indicated technical features or implicitly indicating the number of the indicated technical features. Therefore, the features limited by "first", "second" can explicitly or implicitly include at least one of the features; the technical solutions of various embodiments can be combined with each other, based on the realization of a person of ordinary skill in the art.

[0079] In addition, in the following description, the description of well-known structures and technologies is omitted to avoid unnecessary confusion of the concept of the present application.

[0080] It should be noted that the sequence numbers of the sequence mentioned in the present application do not necessarily represent the strict sequence numbers in the actual specific implementation process. The sequence numbers are used to distinguish each step, so as to prevent confusion.

[0081] In the drawings, a layer structure schematic diagram according to an embodiment of the present application is shown. These drawings are not drawn to scale, in which some details are enlarged for the purpose of clarity, and some details can be omitted. The shapes of various regions, layers and their relative sizes and positional relationships shown in the drawings are only exemplary, and in actuality, there can be deviations due to manufacturing tolerances or technical limitations, and a person skilled in the art can additionally design regions / layers with different shapes, sizes and relative positions according to actual needs.

[0082] In addition, the technical features involved in different embodiments of the present application described below can be combined with each other as long as there is no conflict.

[0083] At present, the identification of important parameters of the existing power system mostly depends on time-domain global error indicators such as mean square error (MSE), mean absolute error (MAE), correlation coefficient or dynamic time warping (DTW) distance, and the whole time sequence characteristics before and after the fault or in the disturbance process are condensed into one or a few numerical values for comparison. Although such "coarse-grained" evaluation is simple to implement and has good compatibility with existing monitoring platforms, it essentially discards rich transient information such as phase drift, peak overshoot, damping decay rate and other dynamic details, and it is difficult to capture the differential influence of each parameter on the system response at different times.

[0084] In order to reduce the complexity of the algorithm, the traditional method usually adopts a "one-by-one scanning" or empirical grouping method for discrete processing of the parameters to be identified, and it is assumed by default that each parameter does not interfere with each other in the identification process. This "parameter independent" assumption seriously ignores the complex coupling and interaction between the generator excitation system, the speed regulator and the network impedance, and the load characteristics.

[0085] The present application aims to solve the technical problems of the traditional power system important parameter identification, such as the coarse-grained time-domain global error indicators and the unmodeled parameter coupling, which lead to the lack of effectiveness of the identification results. To this end, the present application proposes an important parameter intelligent identification method based on time-frequency feature fusion and soft core clustering, which realizes accurate extraction of fine-grained dynamic features in multi-channel system variable trajectories through parallel coding of time-domain convolution network and frequency-domain STFT-MLP; transforms the traditional single parameter non-mechanism identification into an interpretable parameter group identification with clear physical meaning, and fully considers the coupling relationship between parameters. In addition, a "soft core" contribution matrix is further constructed to quantify the cumulative influence between parameters, and UMAP dimensionality reduction and clustering analysis are combined to mine the parameter coupling mode; finally, based on the clustering distance threshold of the training set, confidence evaluation and important parameter group screening are carried out, so as to significantly reduce the simulation times and computational burden, and realize high-precision and interpretable discrimination of important parameters that truly affect system performance.

[0086] The present application realizes time-domain-frequency domain double coding on the multi-channel measurement trajectory before and after the fault, generates normalized embedded vectors combined with adaptive projection, and then uses UMAP dimensionality reduction and K-means clustering algorithm to mine the sample cluster structure, and introduces a "soft core" cumulative contribution threshold screening strategy to normalize the parameter weight matrix, so as to quantify the sensitivity and importance of each physical parameter and its parameter combination to the simulation error. This method can quickly locate the high-impact parameter set in the power system model, and is suitable for simulation correction of system dynamic response after system operation parameter drift and device aging.

[0087] The following is illustrated by specific embodiments, as shown in the following Figure 1 The present application provides a power system model important parameter group identification method based on trajectory feature clustering, comprising:

[0088] S1, Latin hypercube sampling method is adopted for parameters to be identified in the power system, n groups of parameter perturbation samples are generated, standard simulation is run for each group of parameter perturbation samples, corresponding multi-channel time sequence simulation trajectories are obtained, a node near a fault point is selected as a monitoring point, data of the node are measured, coarse-grained sensitive parameters are screened out as high-sensitivity parameters to be identified, wherein n is a positive integer greater than or equal to 1.

[0089] S2, time domain difference is made for each group of simulation trajectories and measured trajectories to obtain difference trajectories, and difference trajectory features are extracted.

[0090] S3, for each difference trajectory, the original time domain signal and the STFT spectrum feature are input in parallel through a trained time-frequency encoder, a D-dimensional normalized embedding vector is generated by splicing and projecting, high-dimensional data is mapped to a low-dimensional space through a nonlinear dimension reduction method, K-Means clustering is performed on data in the low-dimensional space in a manifold space, and a cluster label of each sample is obtained; the average silhouette coefficient under different labels is calculated to determine an optimal cluster number and a confidence judgment standard; the measured difference trajectory is calculated and reduced to a low-dimensional space, and the minimum distance between the low-dimensional space and each cluster center is calculated, if the confidence judgment standard is met, the measured sample is classified into the corresponding cluster, and the determination result of the category is recorded, wherein D is a positive integer greater than or equal to 2.

[0091] S4, all samples in the cluster are read, sensitivity calculation is performed on all samples in the cluster, parameter distribution of all perturbation samples in the cluster and the role of each parameter in feature contribution are counted, and the parameter combination with the highest frequency and the largest contribution to dynamic difference in the cluster is selected as the final dominant parameter set.

[0092] The identification idea of the model important parameters of the application is as follows:

[0093] (1) sensitive parameter screening based on LHS method of coarse granularity

[0094] Latin hypercube sampling (LHS) is adopted for all parameters to be identified in the power system to generate a plurality of parameter combinations, and standard simulation is run under each parameter to obtain corresponding multi-channel time sequence trajectories. The global error index of each trajectory and the reference simulation result is calculated, and the low-sensitivity parameter combination with a very small contribution to the event response is directly removed, thereby quickly reducing the dimension from a large number of parameters to obtain a small-scale candidate parameter subset.

[0095] (2) secondary perturbation feature extraction and clustering of parameters

[0096] The LHS perturbation is performed on the screened parameter subset again, and n groups of perturbed parameters and their corresponding simulation trajectories are obtained. The time-domain difference between each perturbed trajectory and the reference trajectory is calculated to obtain the difference sequence; then the fine-grained dynamic features are extracted, including phase shift amplitude, peak overshoot, damping decay rate, and multi-dimensional descriptions such as oscillation amplitude-frequency distribution, and these features are input into the clustering algorithm for classification. The clustering results make the perturbed samples in the same cluster highly consistent in dynamic response mechanism, revealing the coupling effect of different parameter combinations.

[0097] (3) Mapping of measured trajectory and identification of important parameters

[0098] The difference sequence between the measured trajectory and the reference simulation is calculated, and the same dynamic features are extracted, which are mapped to the trained clustering model to determine the cluster to which they belong. Then the parameter distribution of all perturbed samples in the cluster and the role of each parameter in feature contribution are counted, and the parameter combination with the highest frequency and the largest contribution to the dynamic difference in the cluster is selected as the final dominant parameter set. The elements in this set are the important parameters that truly drive the dynamic response of the system in the current event.

[0099] The specific application process is as follows:

[0100] 1. First screening of key parameters based on trajectory energy:

[0101] The Latin Hypercube Sampling (LHS) method is used on the target parameter set {p1,...,p m n groups of parameter perturbation samples are generated to ensure uniform coverage of the value range in each parameter dimension, avoiding the aggregation or void problem that may be caused by random sampling.

[0102] For each group of perturbed parameters p (k) Run the power system simulation to obtain the multi-channel time sequence trajectory of the target measurement point before and after the fault or during the perturbation where i = 1,...,C, and C is the number of measurement channels.

[0103] Select the nodes near the fault point as the monitoring points and measure the data of the nodes. For example, the active power P and reactive power Q of the generator. For the i-th channel, the energy sensitivity index is defined as

[0104]

[0105] Read the sensitivity matrix of the k-th experiment and perform Min-Max normalization on each experiment in the parameter dimension. Set the threshold τ of the sensitive parameter to be in the range [0, 1], and then screen the coarse-grained sensitive parameters as high-sensitivity parameters to be identified.

[0106] 2. Feature extraction of difference trajectory:

[0107] First, for each group of LHS perturbation, the simulated trajectory is compared with the measured trajectory x meas to get the difference signal Δx (k) = x (k) - x (k) (t) in time domain.

[0108]

[0109] where p is the STFT complex spectrum of the nth window, which is used to preserve the energy distribution in transient errors and provide rich local frequency domain information for subsequent self-supervised learning.

[0110] To enable the model to reliably extract features under different truncation lengths or local missing conditions, the difference signal Δx (k) and the corresponding frequency spectrum sequence p (k) are randomly cropped (randomly select the starting point and length) and linearly or high-order interpolated to a fixed size. This strategy produces multiple "views" that not only maintain the original dynamic information but also improve the robustness to time shift and length changes during training.

[0111] To ensure that similar samples in the frequency domain power sequence space can also maintain proximity in the embedding space, a graph adjacency weight matrix is constructed based on the original spectral features:

[0112]

[0113] where p i and p j represent the spectral power vectors of two difference trajectories, and σ controls the neighborhood range.

[0114] Graph regularization loss is added during training:

[0115] L grapg =∑ i,k A ij ||z i -z j || 2 (4)

[0116] where z i and z j are the embedding vectors obtained by the self-supervised encoder, and this term encourages samples that are close in the frequency spectrum to also be close in the low-dimensional embedding space z, which helps to maintain the manifold structure of the frequency domain itself, so that downstream clustering or graph methods based on embedding can directly utilize this local continuity.

[0117] In the present application, the encoder is designed as a set of parallel time-domain and frequency-domain sub-networks, whose structure and data flow include time-domain branch and frequency-domain branch.

[0118] In the time-domain branch, the input is the enhanced difference signal First, the channel number is expanded to a1 dimension through a one-dimensional convolution layer, and the activation function used here is the ReLU function.

[0119] Then, several layers of dilated depth residual blocks (DW) are passed through, each block is in the form of:

[0120]

[0121] Among them , represents the depth separable convolution with dilated rate d, and BN represents batch normalization. In this way, early phase drift, peak overshoot and decay process can be captured under different receptive fields. Finally, the adaptive average pooling is used to compress the time series length to 1, and the time-domain feature vector z t .

[0122] The input is the spectrum power sequence after cutting and interpolation First, it is spliced into a vector along the time window Then, two fully connected layers are used to obtain the frequency domain vector z f z f This process can ensure that the spectral information is compressed without losing its overall energy distribution characteristics.

[0123] The time-domain vector z t and the frequency-domain vector z f are spliced into a vector z cat , and after linear mapping and ReLU activation, the new comprehensive feature vector is:

[0124]

[0125] The result is L2 normalized to obtain the processed vector z, so that all embeddings fall on the unit hypersphere, which is convenient for distance or similarity calculation.

[0126] In order to force to retain the frequency domain details in training, the decoder head g dec is designed to perform reverse mapping on the normalized vector z:

[0127]

[0128] The output of this reconstruction head is aligned with the original spectrum , providing a target for the spectrum reconstruction loss. Among them, W3 and W4 are the weight matrices of the output layer, and b3 and b4 are the bias terms.

[0129] To fully utilize the above structure and data augmentation strategies, we construct a joint training objective consisting of five losses to guide the model to learn high-quality embeddings.

[0130] The loss function includes the following five points, among which points (1), (2), and (3) are commonly used default losses, and points (4) and (5) are self-defined losses of the present application to significantly improve the frequency domain interpretability of the training results:

[0131] (1) Contrastive loss:

[0132] In the form of NT-Xent, two augmented views are generated for each sample in a batch The other views constitute negative pairs. The contrastive loss is defined as:

[0133]

[0134] where τ is the temperature coefficient, <·,·> represents the inner product, and N represents the number of sample pairs in a training batch. This loss encourages different views of the same track to be closely clustered in the embedding space, while different track views are pushed away, thereby improving the discrimination ability.

[0135] (2) Variance loss:

[0136] Prevent all samples from "collapsing" to a constant in some dimensions, and ensure that the standard deviation Std(z :,d ) of each dimension is not less than the variance threshold γ:

[0137]

[0138] z :,d represents the value vector of all samples in the embedding dimension d in the batch, and maz(0, γ-Std(z :,d )) only generates a loss when the standard deviation is lower than γ, otherwise the dimension loss is 0.

[0139] (3) Covariance loss:

[0140] Ensure that there is no redundant information between different embedding dimensions, and minimize the covariance of each pair of dimensions:

[0141]

[0142] where z :,i ,z :,j represent the value vectors of all samples in the i-th and j-th embedding dimensions in the batch.

[0143] Cov(z :,i ,z :,j) represents the covariance of the two vectors, measuring their linear correlation, further improving feature diversity by reducing the linear correlation between the dimensions of the embedding.

[0144] (4) Spectral reconstruction loss:

[0145] Reconstruction head with the original spectrum to minimize reconstruction error:

[0146]

[0147] With this loss, the frequency domain information is explicitly preserved in the embedding, avoiding information loss in the fusion or normalization process.

[0148] (5) Graph regularity loss:

[0149] Based on the spectral adjacency matrix A ij , the local structure of the spectral space is kept smooth in the embedding space:

[0150] L graph =∑ i,j A ij ||z i -z j || 2 (12)

[0151] This term ensures that samples similar in the spectrum also maintain proximity in the embedding space, which helps subsequent clustering and local structure analysis.

[0152] The five loss terms are combined with weights:

[0153] L=λ contrast L contrast +λ var L var +λ cov L cov +λ spec L spec +λ graph L graph (13)

[0154] where λ represents the weight.

[0155] 3. Trajectory clustering and cluster identification:

[0156] For each difference trajectory (LHS perturbed sample and measured-simulation difference sample), we first input the original time-domain signal and STFT spectral features in parallel through the trained time-frequency encoder, concatenate and project to generate a D-dimensional normalized embedding vector z (k) .

[0157] The embedding set {z(k)}Use UMAP (Uniform Manifold Approximation and Projection) nonlinear dimensionality reduction algorithm to map high-dimensional data to low-dimensional space y (k) , where d is the dimension after dimensionality reduction, which ensures that the local neighborhood structure is preserved.

[0158] In the manifold space, apply K-Means to {y (k)}Need to perform K-Means clustering on the set of possible cluster numbers K respectively, and get the cluster label c of each sample (k) Calculate the average silhouette coefficient under different labels:

[0159]

[0160] Where a (k) is the average distance between sample k and other samples in the same cluster (intra-cluster closeness), and b (k) is the average distance between sample k and the nearest neighbor cluster (inter-cluster separation). Take the average of all samples Select The corresponding K when the maximum is the optimal cluster number.

[0161] The center vector of each cluster c is μ c , and the Euclidean distance between the samples in the cluster and the center is calculated And take the 95% quantile value τ c of the distance of the cluster in the training set as the confidence judgment standard.

[0162] The measured difference trajectory is also embedded and reduced to y (meas) , and the minimum distance between it and each cluster center μ c is calculated; if:

[0163]

[0164] The measured sample is classified into the corresponding cluster c, and the determination result of this category is recorded.

[0165] 4. Clustering and acquisition of important parameters:

[0166] Read all samples in the cluster, and calculate the sensitivity (trajectory energy E) of all samples in the cluster. In addition, read the difference ratio b of all samples in the cluster and the simulation setting parameters. Define the parameter importance index I, I = b × E. Get the parameter importance array of the samples in the cluster.

[0167] For the kth sample, the normalized parameter importance vector is Where m is the total number of parameters, denotes the normalized contribution of the i-th parameter. The indices are arranged in descending order of their values, and the first r terms are accumulated until the following condition is met:

[0168]

[0169] where CUMth denotes the cumulative threshold, and the threshold recommendation is defined as 0.9. The first r parameter indices constitute the core parameter set S (k) . The parameters are then L2 normalized to obtain the soft kernel vector:

[0170]

[0171] For intra-cluster soft kernel aggregation statistics, for each cluster c, the sample set K c = {k | c (k) = c}, the average soft kernel weight is calculated as:

[0172]

[0173] The core occurrence frequency is:

[0174]

[0175] where 1(·) is the indicator function.

[0176] Combining and f i (c) , the internal parameters of each cluster are ranked by the comprehensive score, and the top several parameters with the highest scores are defined as the key parameter group of the cluster.

[0177] That is, the parameter group composed of these parameters has a direct impact on the accuracy of the simulation. When there is an error between the simulation result and the measured result, the parameter group identified can be adjusted to correct the trajectory.

[0178] In order to confirm the effect of the scheme of the present application, the following example is used to verify the scheme of the present application.

[0179] According to the IEEE standard 39-node (as shown in Figure 2 ), a wind power system is connected to node 15. For parameter identification, the synchronous machine uses a 4th order model, the excitation system uses a 2nd order model, and the load model uses a ZIP model. A three-phase short circuit fault is set at node 2, and the fault starts at 1s and lasts for 0.1s. The data sampling time is 5s.

[0180] The measured data is the measured data when node 2 fails under the standard data. The sampling points are set to the active and reactive power of the generator 10 (G10). After the random perturbation parameter, the simulation trajectory mismatch phenomenon caused by the parameter is simulated to verify the rationality and effectiveness of the technical scheme of the application.

[0181] Figure 3 The parameter sensitivity distribution after the LHS perturbation of the parameter is shown, including the sensitivity of the active and reactive power trajectories. The vertical coordinate axis represents the number of parameter perturbations, a total of 3000 times, and the horizontal axis band represents the sensitivity of the parameter, a total of 130 parameters. The energy threshold value of the sensitive parameter is 0.1, that is, when the sensitivity calculation result is greater than 0.1, the band corresponding to the perturbation position is black, indicating that the parameter plays an important role in this parameter combination.

[0182] Finally, after the first parameter screening, the parameters that do not reflect sensitivity at all in this event can be excluded, and finally 21 parameters are screened, which are the inertia time constants of generators 2, 3, 4, 5, 7, 8, 9 and 10; The direct-axis sub-transient reactance of generators 1, 2, 5, 8 and 9; The direct-axis reactance of generators 2, 3, 7 and 8; The cross-axis sub-transient reactance of generators 2, 7, 8 and 10.

[0183] By performing secondary perturbation on the selected parameters, the perturbation number is 10000 times, and the obtained trajectory data is subtracted from the simulation trajectory to obtain the difference trajectory. The short-time Fourier transform is performed on the difference trajectory to obtain the spectral average power in each window.

[0184] The data is input into the constructed model, the model is trained, and the encoded results are dimensionally reduced and feature extracted. According to the data processing results, the number of clusters is Figure 4 Show.

[0185] According to the calculation results of the silhouette coefficient, when the cluster number is 7, the silhouette coefficient is 0.632. That is, when the cluster number is 7, the data features have the optimal distinguishability. Further calculate the cluster to which the difference trajectory of the simulation and measured data belongs. The calculation results show that the center distance of the difference trajectory of the simulation and measured data belonging to the 2nd cluster is 1.88, and the threshold value is 13.34. The membership degree of the trajectory belonging to the 2nd cluster is significant.

[0186] Next, the distribution of the parameters in cluster 2 is analyzed. The cluster has 2028 groups of data. The cumulative threshold is defined as 0.9, and the distribution of the parameters is analyzed by the clustering algorithm. Two types of parameter arrangements are obtained. Arrangement one is: the inertia time constant of generators 2, 8, 9, and 10, the direct-axis sub-transient reactance of generators 5, 8, and 9, and the cross-axis sub-transient reactance of generator 8, with a total sample size of 1786. Arrangement two is: the inertia time constant of generators 5, 8, 7, and 10, and the direct-axis sub-transient reactance of generators 8 and 9, with a total sample size of 242. Two important parameter groups are obtained. In order to demonstrate the rationality of the two identified parameter groups, the trajectories are corrected using these parameters. In addition, the traditional sensitive parameter identification method is used, that is, the sensitive parameter coarse-grained selection without considering the parameter coupling relationship, and the correction result is shown in Figure 5

[0187] The following effect analysis is performed:

[0188] The present application combines Latin hypercube sampling with coarse-grained energy screening. First, the parameter space is quickly compressed to 21 high-sensitivity parameters, so that subsequent fine identification only needs to analyze about 16% of the original parameter amount, and the simulation times are significantly reduced.

[0189] In the clustering identification stage, the self-supervised time-frequency encoder and multiple loss training are introduced. The embedding space not only has time-domain transient and frequency-domain modal information, but also accurately divides 7 typical response modes through UMAP and contour coefficient adaptive clustering. The measured-simulation difference trajectory is reliably classified into "cluster 2", and the experimental center distance is much smaller than the 95% percentile threshold, verifying the high confidence of the method in unknown sample judgment and pattern recognition.

[0190] Further "soft core" clustering analysis reveals two types of parameter groups in cluster 2. After correcting the trajectories using the two groups of parameters, the correction error is reduced by more than 90% on average compared to before correction. The correction of parameter group 1 and parameter group 2 can be directly observed. The error between the corrected trajectory and the measured trajectory is significantly reduced. The calculation results show that after correcting the trajectories using the two groups of parameters, the correction error is reduced by more than 90% on average compared to before correction. The error reduction rate of the sensitive parameters selected by the traditional method is 68%, which is lower than the effect of the present technical solution, showing the effectiveness of the present method.

[0191] From the test results above, it can be seen that the error between the corrected trajectory and the measured trajectory is significantly reduced, effectively solving the problems of coarse-grained time-domain global error indicators and unmodeled parameter coupling in traditional power system important parameter identification.

[0192] ​According to the test results above, in combination with the scheme of the present application, the beneficial effects are as follows:

[0193] (1) The present application first proposes a full-parameter space uniform perturbation and sensitivity fast compression technology based on Latin hypercube sampling (LHS). By one-time LHS sampling of all to-be-identified parameters and combining with trajectory energy index calculation Min-Max normalized parameter sensitivity, the most significant parameter subset affecting system dynamics can be quickly screened out on the basis of one simulation result, and the parameter scale is directly compressed from hundreds to thousands to tens of levels, thereby greatly reducing the calculation amount of subsequent analysis.

[0194] (2) The present application innovatively designs a self-supervised time-frequency feature fusion encoder, and introduces multiple consistency constraints of contrast loss, variance loss, covariance loss, spectral reconstruction loss and graph regularization loss. The encoder deploys time domain and frequency domain branches in parallel, which can capture time domain transient features such as peak overshoot and damping attenuation, and also retain spectral energy distribution and modal information, so as to train a high-quality, distinguishable and interpretable comprehensive embedding representation in a label-free scene.

[0195] (3) The present application realizes adaptive optimal cluster number determination in combination with UMAP dimension reduction and contour coefficient method. First, the high-dimensional time-frequency embedding is projected into a low-dimensional manifold space, and then K-Means clustering is performed on the possible cluster number range and the average contour coefficient is calculated, and the cluster number at the highest contour coefficient is selected automatically; subsequently, a confidence threshold is constructed based on the quantile value of the distance from each cluster sample to the cluster center, realizing the automation and quantification of trajectory category division and unknown mode / abnormality judgment.

[0196] (4) The present application constructs a “soft core” contribution matrix to deeply mine parameter coupling effect. For each perturbation sample, the normalized response degree is truncated according to the cumulative threshold, only the core parameters are retained, and the contribution degree is L2 normalized to generate a high-dimensional soft core weight vector; then the average weight and frequency of the soft core vector of all samples in the cluster are counted, which comprehensively reflects the typical influence and stability of the parameters, so as to obtain the key parameter combination that really drives the system in each response mode.

[0197] The embodiment of the present application also provides an electronic device 600, as shown in the figure, comprising a memory 601, a processor 602, and a computer program stored in the memory 601 and executable on the processor, wherein the processor 602 executes the program to realize the steps of the method of any one of the above embodiments. Figure 6

[0198] The embodiment of the present application also provides a computer readable storage medium, which stores a computer program, wherein the computer program is executed by a processor to realize the steps of the method of any one of the above embodiments. ​

[0199] Obviously, the described embodiments are only part of the embodiments of the present application, rather than all the embodiments. All other embodiments obtained by those of ordinary skill in the art based on the embodiments of the present application without creative labor fall within the protection scope of the present application.

[0200] It should be understood that the above specific embodiments of the present application are only used for illustrative or explanatory purposes of the principles of the present application, and do not constitute a limitation on the present application. Therefore, any modification, equivalent replacement, improvement, etc. made without departing from the spirit and scope of the present application shall be included in the protection scope of the present application. In addition, the appended claims of the present application are intended to cover all variations and modifications falling within the scope and boundary of the appended claims, or the equivalent forms of such scope and boundary.

Claims

1. A method for identifying important parameter sets in power system models based on trajectory feature clustering, characterized in that, include: S1. The Latin hypercube sampling method is used to generate n sets of parameter disturbance samples for the parameters to be identified in the power system. Standard simulation is run for each set of parameter disturbance samples to obtain the corresponding multi-channel time-series simulation trajectory. A node near a certain fault point is selected as a monitoring point. The data of the node is measured to screen out coarse-grained sensitive parameters as high-sensitivity parameters to be identified, where n is a positive integer greater than or equal to 1. S2. Perform time-domain difference analysis on each group of simulated trajectories and measured trajectories to obtain the difference trajectories, and extract the features of the difference trajectories; S3. For each of the aforementioned difference trajectories, the original time-domain signal and STFT spectral features are input in parallel through a trained time-frequency encoder, spliced ​​and projected to generate a D-dimensional normalized embedding vector. The high-dimensional data is mapped to a low-dimensional space using a nonlinear dimensionality reduction method. K-Means clustering is performed on the data in the low-dimensional space in the manifold space to obtain the cluster label for each sample. The average silhouette coefficient under different labels is calculated to determine the optimal number of clusters and the confidence judgment criterion. The measured difference trajectory is calculated and reduced to a low-dimensional space. The minimum distance between the trajectory and the center of each cluster is calculated. If the confidence judgment criterion is met, the measured sample is assigned to the corresponding cluster and recorded as the judgment result of that category, where D is a positive integer greater than or equal to 2. S4. Read all samples within the cluster, perform sensitivity calculations on all samples within the cluster, statistically analyze the parameter distribution of all perturbed samples within the cluster and the role of each parameter in the feature contribution, and select the parameter combination with the highest frequency of occurrence and the largest contribution to the dynamic differences within the cluster as the final dominant parameter set.

2. The method according to claim 1, characterized in that, The step of running a standard simulation on each set of parameter perturbation samples to obtain the corresponding multi-channel time-series simulation trajectory specifically includes: For each set of disturbance parameters p (k) Run power system simulation to obtain multi-channel time-series trajectories of target measurement points before and after faults or during disturbances. Where i = 1, ..., C, C is the number of measurement point channels, and k is a positive integer greater than or equal to 1.

3. The method according to claim 1, characterized in that, The process of selecting nodes near a specific fault point as monitoring points, measuring the data of these nodes, and filtering out coarse-grained sensitive parameters as high-sensitivity parameters to be identified specifically includes: Nodes near the fault point are selected as monitoring points, and the data of the nodes are measured. The data includes the generator's active power P and reactive power Q information. For the i-th channel, its energy sensitivity index is defined as shown in Formula 1: Read the sensitivity matrix of the k-th experiment For each experiment, Min–Max normalization was performed on the parameter dimension, and the threshold τ∈[0,1] of the sensitive parameters was set to filter out coarse-grained sensitive parameters as high-sensitivity parameters to be identified.

4. The method according to claim 1, characterized in that, The step of performing time-domain difference analysis on each group of simulated trajectories and measured trajectories to obtain the difference trajectories, and extracting the features of the difference trajectories, specifically includes: For each set of simulation trajectories Compared with the measured trajectory x meas (t) Perform time-domain difference to obtain the difference signal Δx (k) ; The difference signal Δx (k) Perform a short-time Fourier transform using a sliding window, and calculate the spectral average power within each analysis window using Formula 2: in, It is the short-time Fourier transform complex spectrum of the nth window, used to preserve the energy distribution in the transient error; The transformed difference signal is input into the encoder and decoded by the decoder to obtain the difference trajectory features; The encoder consists of parallel time-domain and frequency-domain sub-networks, and its structure and data flow include time-domain and frequency-domain branches.

5. The method according to claim 4, characterized in that, For the difference signal Δx (k) After performing a short-time Fourier transform using a sliding window, and before inputting the transformed difference signal into the encoder, the process also includes: The difference signal Δx (k) and the corresponding spectral sequence p (k) The method employs random cropping and linear or higher-order interpolation to a fixed size, wherein the random cropping includes randomly selecting a starting point and randomly selecting a length.

6. The method according to claim 4, characterized in that, Before inputting the transformed difference signal into the encoder, the process also includes: Based on the original spectral features and Formula 3, a graph adjacency weight matrix is ​​constructed, as shown in Formula 3 below: Where p i p j σ represents the spectral power vector of the two different trajectories, and σ represents the control neighborhood range; During training, the encoder incorporates a graph regularization loss, as shown in Formula 4: L grap =∑ i,j A ij ||z i -z j || 2 (4) Where z i and z j The embedding vector obtained by the self-supervised encoder; The encoder consists of parallel time-domain and frequency-domain subnetworks. The encoder's structure and data flow include time-domain and frequency-domain branches, specifically including: The time-domain branch receives the enhanced difference signal as input. First, the number of channels is expanded to a1 dimension through a one-dimensional convolutional layer, and the activation function used is the ReLU function; After several layers of dilated depth residual blocks, each block takes the form shown in Formula 5: in This represents a depthwise separable convolution with an inflation rate of d, and BN represents batch normalization. The temporal length is compressed to 1 using adaptive average pooling to obtain the temporal feature vector z. t ; The input is the cropped and interpolated spectral power sequence. Concatenate along the time window to form a vector Then, the frequency domain vector z is obtained through two fully connected layers. f z f ; The time-domain vector z t and frequency domain vector z f Concatenate into vector z cat After linear mapping and activation by the ReLU function, the new comprehensive feature vector is shown in Formula 6: The obtained results are L2 normalized to obtain the processed vector z, such that all embeddings fall on the unit hypersphere. Wherein, the decoder header g of the decoder dec The normalized vector z is then reverse-mapped using Formula 7, as follows: The output of the decoder Compared with the original spectrum Alignment provides a target for the spectral reconstruction loss. W3 and W4 are the weight matrices of the output layer, and b3 and b4 are the bias terms.

7. The method according to claim 1, characterized in that, The calculation of the average silhouette coefficient under different labels, and the determination of the optimal number of clusters and confidence level criteria, specifically include: The average contour coefficient under different labels is calculated using Formula Fourteen: Where a (k) Let b be the average distance between sample k and other samples in the same cluster (cluster compactness). (k) The distance between sample k and samples in the nearest neighbor cluster is the average distance between them. Take the average of all samples choose The sample k corresponding to the maximum value is taken as the optimal number of clusters; The center vector of each cluster c is μ. c Calculate the Euclidean distance from the sample within the cluster to the center. And take the 95th percentile value τ of the cluster distance in the training set. c As a standard for judging confidence level.

8. The method according to claim 1, characterized in that, The process involves reading all samples within a cluster, calculating sensitivity for all samples within the cluster, statistically analyzing the parameter distribution of all perturbed samples within the cluster and the role of each parameter in feature contribution, and selecting the parameter combination with the highest frequency and the largest contribution to dynamic differences within the cluster as the final dominant parameter set. Specifically, this includes: Read the difference ratio b between the parameters of all samples in the cluster and the simulation settings parameters, and define the parameter importance index I. The parameter importance index I meets the condition I = b × E. Obtain the parameter importance array for samples within the cluster; For the k-th sample, its normalized parameter importance vector conforms to Where m is the total number of parameters, This represents the normalized contribution of the i-th parameter; Arrange the indices π1, π1, ... in descending order of their values, and sum the first r terms until the conditions shown in Formula 16 are met: Where CUMth represents the cumulative threshold, which is defined as 0.9; The first r parameter indices constitute the core parameter set S. (k) The parameters in the parameter set are subjected to L2 normalized weights, as shown in Formula 17, to obtain the soft kernel vector: For cluster-wide soft-core aggregation statistics, for each cluster c, the sample set K c ={k|c (k) =c}, calculate the average soft core weight using Formula 18: The frequency of core occurrences conforms to Formula Nineteen: Where 1(·) is an indicator function, Combination and f i (c) The parameters within each cluster are sorted according to their comprehensive scores, and the parameters with the highest scores are identified as the key parameter group for that cluster.

9. An electronic device comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that, When the processor executes the program, it implements the steps of the method according to any one of claims 1-8.

10. A computer-readable storage medium having a computer program stored thereon, characterized in that: When the computer program is executed by a processor, it implements the steps of the method according to any one of claims 1-8.

Images