A power grid reliability evaluation method based on successive smoothing estimation of optimal importance sampling density function

Abstract

This method relates to the field of large power grid reliability assessment, specifically a power grid reliability assessment method based on successive smoothing estimation of the optimal importance sampling density function. By directly estimating the optimal zero-variance IS-PDF stepwise, a more efficient IS-PDF is obtained, improving the efficiency of the importance sampling method. This includes: constructing a unified importance sampling framework for all random variables in a multivariate standard normal space using a normal transformation method; building upon this importance sampling framework, constructing an intermediate smoothed joint IS-PDF sequence, using a truncated Gaussian distribution as its initial distribution, and sequentially estimating the subsequent distributions in the intermediate smoothed joint IS-PDF sequence stepwise, while simultaneously using the estimated smoothed joint IS-PDF for sampling and power grid reliability index assessment; until the variance coefficient of the reliability index obtained at the Kth estimation meets the preset accuracy, the current power grid reliability index is used as the final assessment index.

Description

Technical Field

[0001] This method relates to the field of large power grid reliability assessment, specifically a power grid reliability assessment method based on successive smoothing estimation of the optimal importance sampling density function. Background Technology

[0002] How to effectively accelerate Monte Carlo simulation and improve the efficiency of reliability assessment has long been a focus of research in the field of large power grid reliability assessment. Among numerous variance reduction techniques, the Cross-entropy-based Importance Sampling (CE-IS) method has attracted much attention due to its excellent acceleration effect. The traditional CE-IS method, as one of the most effective methods for improving reliability assessment efficiency, is widely used in large power grid reliability assessment. The CE-IS method achieves iterative optimization estimation of the parameters of the Importance Sampling Robability Density Function (IS-PDF) through cross-entropy optimization, thus providing rigorous theoretical rules and a generalized implementation approach for effectively obtaining the IS-PDF.

[0003] Traditional cross-entropy importance sampling methods also have some intractable drawbacks.

[0004] First, the traditional cross-entropy importance sampling method aims at the theoretically optimal zero-variance IS-PDF by performing cross-entropy iterative optimization of the parameters of the actual IS-PDF. However, only samples in the fault domain can participate in the iterative optimization calculation of IS-PDF parameters. Therefore, the parameter update of the traditional cross-entropy importance sampling method does not fully utilize the samples. When the system reliability is low, obtaining fault samples is relatively easy, so the traditional cross-entropy importance sampling method can obtain enough fault samples for the iterative optimization calculation of IS-PDF parameters within a limited sample capacity. However, when the system reliability is high, obtaining fault samples is very difficult. Therefore, the fault samples that the traditional cross-entropy importance sampling method can obtain within a limited sample capacity will become very scarce, or even non-existent, leading to reduced efficiency or even failure in the iterative optimization of IS-PDF parameters.

[0005] Secondly, traditional cross-entropy importance sampling methods are essentially parameter estimation methods, requiring prior selection of the IS-PDF type for various random variables before IS-PDF parameter estimation can be performed. However, there are no strict theoretical rules to guide the selection of IS-PDF types. Typically, the choice of IS-PDF type is based on considerations of parameter update formula derivation and subjective experience. Since the desired IS-PDF g(x|Q) appears in the factor ln(g(x|Q)), choosing an exponential family distribution as the form of g(x|Q) can reduce the difficulty of deriving the analytical formula for parameter updates. In power system reliability assessment, for discrete variables in the system, i.e., system component state variables, their original PDF belongs to the Bernoulli distribution of the exponential family. Therefore, the IS-PDF for these variables usually adopts the same functional form, and the component unavailability rate is used as the parameter to be updated iteratively. For continuous variables, due to their complexity and variability, the selection of IS-PDF for these variables mainly aims at the feasibility of deriving the IS-PDF parameter update formula. For example, when considering system load, exponential family distributions such as normal, log-normal, and GMM are typically chosen as the IS-PDF. Normal and log-normal distributions are unimodal distributions with fixed shapes, so they may not effectively highlight the critical regions of the system load. While GMM can construct IS-PDFs of arbitrary shapes, its parameter updates are cumbersome, potentially reducing the computational efficiency of traditional cross-entropy importance sampling methods. Therefore, choosing a suitable IS-PDF type is difficult, and traditional cross-entropy importance sampling methods may experience significant performance degradation due to inappropriate IS-PDF type selection, even leading to serious errors in reliability assessment results.

[0006] In summary, the shortcomings of traditional cross-entropy important sampling methods lie in the fact that the iterative optimization of its IS-PDF parameters can only use faulty samples and the lack of rigorous theoretical rules as a basis for IS-PDF type selection. These shortcomings will affect the efficiency and performance of traditional cross-entropy methods to varying degrees, leading to serious problems in some situations. Summary of the Invention

[0007] The present invention proposes an important sampling method based on stepwise estimation of intermediate smooth distribution. By directly estimating the optimal zero variance IS-PDF stepwise, a more efficient IS-PDF is obtained, which further improves the efficiency of the important sampling method.

[0008] The power grid reliability assessment method based on successive smoothing estimation of the optimal importance sampling density function in this invention includes:

[0009] Step 1 transforms the mixed variable space of the original discrete and continuous variables into a multivariate standard normal space through the normal transformation method, and on this basis, constructs a unified importance sampling framework for all random variables in the multivariate standard normal space.

[0010] Step 2: Based on the important sampling framework, construct an intermediate smooth joint IS-PDF sequence, using a truncated Gaussian distribution as its initial distribution, and sequentially estimate the subsequent distributions in the intermediate smooth joint IS-PDF sequence. Simultaneously, use the estimated smooth joint IS-PDF for sampling and power grid reliability index evaluation. Until the variance coefficient (Cov) of the reliability index obtained at the Kth estimation meets the preset accuracy, the current power grid reliability index is used as the final evaluation index.

[0011] To avoid the shortcomings of traditional CE-IS and iCEM methods in IS-PDF selection, this invention proposes an importance sampling method based on stepwise estimation of intermediate smooth distributions. This method achieves an IS-PDF similar to the theoretically optimal one and with higher efficiency without relying on cross-entropy optimization. First, a unified importance sampling framework is established for all random variables in a multivariate standard normal space, realizing the continuity of the original distribution f(x|P). Then, within this importance sampling framework, a set of smooth intermediate joint IS-PDF sequences is constructed using a smoothing function. These intermediate joint IS-PDF sequences are then stepwise estimated using sample resampling and conditional MH sampling. As the estimation progresses, the intermediate joint IS-PDF sequences gradually and smoothly approximate the theoretically optimal IS-PDF, achieving a direct approximation of the theoretically optimal IS-PDF and obtaining the joint IS-PDF for all random variables. Since this method does not require selecting the IS-PDF type and does not rely on cross-entropy optimization for parameter updates, it avoids the shortcomings of cross-entropy-based importance sampling methods and improves the efficiency of reliability assessment. Numerical analysis shows that the importance sampling method based on stepwise estimation of intermediate smooth distribution achieves an approximate estimate of the theoretically optimal IS-PDF, obtaining a more reasonable IS-PDF, and has superior performance and efficiency compared to the cross-entropy importance sampling method.

[0012] In a further embodiment, a sampling method based on sample resampling operation and conditional MH sampling technique free from the curse of dimensionality is proposed to efficiently sample high-dimensional joint IS-PDF.

[0013] In a further embodiment, a comprehensive reliability metric is designed to fully utilize samples drawn from all intermediate joint IS-PDFs to further improve the efficiency of SSAM. Attached Figure Description

[0014] Figure 1 A schematic diagram illustrating the variation of the edge IS-PDF of system load in the SSAM method of this invention embodiment.

[0015] Figure 2 This is a schematic diagram comparing the original PDF of system load and IS-PDFs obtained by different methods in an embodiment of the present invention.

[0016] Figure 3 In the SSAM method of this embodiment of the invention, different Below, the smoothing parameter σ k and normalized parameter p k A diagram illustrating the changes.

[0017] Figure 4 This is a schematic diagram comparing the acceptance rates of conditional MH sampling and traditional MH sampling in the SSAM method of this invention. Detailed Implementation

[0018] The power grid reliability assessment method (SSAM) based on successive smoothing estimation of the optimal importance sampling density function in this embodiment includes the following steps:

[0019] 1. A unified importance sampling framework for all random variables in a multivariate standard normal space.

[0020] Directly to the optimal zero variance IS-PDFg op (x|Q op )=H(x)f(x|P) / p F Sampling estimation is not feasible because estimating unknown normalized parameters is quite difficult, especially when the system being evaluated has high reliability. Furthermore, even if p is obtained... F Even with approximate estimates, it is difficult to use optimal zero variance IS-PDF for Monte Carlo sampling because not only does the indicator function H(x) lack an analytical form, but the factor f(x|P) is essentially a complex distribution containing both discrete and continuous variables. Therefore, even using the MH sampling method, which can obtain samples from complex or potentially unanalytical distributions, it is difficult to sample optimal zero variance IS-PDF, since the MH sampling method is only suitable for obtaining samples from continuous distributions.

[0021] To address the aforementioned issues, it is necessary to make H(x) and f(x|P) continuous. This example proposes a smoothing function Φ(-S(x) / σ), which replaces the original indicator function H(x) with a smoothing parameter σ that has a smaller smoothing parameter. K The smoothing function achieves smoothing of the indicator function and obtains the optimal zero-variance IS-PDF for smoothing:

[0022]

[0023] Therefore, this example will focus on solving the problem of making f(x|P) continuous, and will use a normal transformation to convert the original probability space following f(x|P) into a space following f(x|P). We have established an independent standard normal space for all random variables and a unified importance sampling framework for all random variables in the multivariate standard normal space, as detailed below:

[0024] 1) System component state vector x c Each discrete variable in All can be converted into their corresponding standard normal variables by equation (2).

[0025]

[0026] Where, δ j =Φ -1 (ξ j ξ is the threshold between the fault state and the normal state. j Φ is the initial unavailability rate of the j-th component. -1 (·) is the inverse function of the standard normal cumulative distribution function. It is an indicator function: if but otherwise,

[0027] 2) System load variable x L Transform it into the corresponding standard normal variable u through equal probability transformation. L :

[0028]

[0029] Among them, F L -1 () is the original cumulative distribution function The inverse function of .

[0030] 3) Correlated wind speed vectors It can be converted into the corresponding independent standard normal random vector through Nataf transformation.

[0031] First, transform x using equal probability. w Transform into a vector with zero mean and a correlation matrix C z The normally distributed random vector z, i.e. in It is the cumulative distribution function The inverse function of .

[0032] Let C x yes The correlation matrix, ρ x (i,j) and ρ z (i,j) are C x and C z The element in the i-th row and j-th column. ρ z (i,j) can be obtained by solving the following formula:

[0033]

[0034] Where, μ i and μ j They are and The mean, σ i and σ j They are and The standard deviation. It has a correlation coefficient C z The bivariate standard normal probability density function of (i,j). After obtaining C... z Then, C was calculated using Cholesky decomposition. z The lower triangular matrix L z C z =L z L z T Then, a correlated normal random vector z can be derived from an independent standard normal random vector u. w This means that z = L z u w Therefore, the wind speed vectors with correlation The corresponding independent standard normal random vector The conversion relationship between them can be expressed by equation (5):

[0035]

[0036] In summary, the above normal transformation converts the original x-space, which contains discrete random variables, continuous variables, and correlated random variables, into a u-space containing only independent standard normal random variables. This transforms the discrete and complex original distribution f(x|P) into a simple, continuous multivariate standard normal distribution. x = [x c ,x w ,x L ] and u=[u c ,u w ,u L The mapping relationship between x and u can be simply represented as x = T(u), which means that each sample drawn from u space can be mapped to a sample in the original x space.

[0037] Step 2: Based on the aforementioned important sampling framework, an intermediate smooth joint IS-PDF sequence is constructed, using a truncated Gaussian distribution as its initial distribution, and it is progressively estimated. As the estimation proceeds, this intermediate distribution sequence will gradually approach the optimal zero-variance IS-PDF, thereby achieving an approximate estimate of the optimal zero-variance IS-PDF, and finally obtaining a joint IS-PDF that is highly similar to the optimal zero-variance IS-PDF in functional form, as detailed below:

[0038] Based on the aforementioned sampling framework, the optimal IS-PDF g in the discrete x-space is... op (x|Q op ) is transformed into the optimal g in the continuous u space. op (u|Q op ):

[0039]

[0040] in Let n represent the n-dimensional standard normal probability density function, where n = D + W + 1. and Let represent the D-dimensional, W-dimensional, and one-dimensional standard normal probability density functions, respectively. S(u) is the system performance function defined in the u-space, which can be determined by power flow or optimal power flow analysis of the system state x = T(u) in the x-space.

[0041] Combining formula (6) with the smoothing function, we can obtain the optimal zero variance IS-PDF for smoothing:

[0042]

[0043] Where, σ k It is a sufficiently small smoothing parameter, p k It is the normalized parameter, Q k ={σ k ,p k The following section focuses on the principle and steps of an important sampling method based on stepwise estimation of intermediate smooth distribution sequences. Although the optimal zero variance IS-PDF g(u|Q) has been constructed in formula (7), k However, accurately estimating the smoothing parameter σ is difficult. k and normalized parameter p k It is very difficult, therefore, how to achieve g(u|Q) kEstimating the value of u is a challenging problem. Therefore, this invention proposes an important sampling method (SSAM) based on stepwise smoothing estimation using intermediate smoothing distributions to address this issue. The basic idea of ​​SSAM is to establish a set of intermediate smoothing joint ISPDF sequences {g0(u|Q0), g1(u|Q1), ..., g...}. k (u|Q k )}, where the initial distribution g0(u|Q0) is a distribution with a known specific form and is easy to sample.

[0044] Then, using the distribution g obtained from the previous one... k-1 (u|Q k-1 The sample drawn from ) is used to apply to the next unknown distribution g. k (u|Q k The parameter Q k ={σ k ,p k Estimate the smoothing parameter. As the estimation progresses, the smoothing parameter will continuously decrease, i.e., {∞=σ0>σ1>…σ K >0}, while the obtained intermediate smoothed joint IS-PDF sequence {g k (u|Q k The expression g(k) = 0, 1, ..., K will gradually approach the optimal zero-variance IS-PDF, i.e., g(k) = 0, 1, ..., K. k (u|Q k )≈g op (u|Q op ).

[0045] Note that the implementation of SSAM requires solving two key problems: 1) How to base the distribution g obtained from the previous one on the distribution g. k-1 (u|Q k-1 To estimate the next unknown distribution g from the sample. k (u|Q k The parameter Q k ={σ k ,p k};2) After obtaining Q k Next, how to generate g? k (u|Q k The sample size is used to ensure the stepwise estimation process can continue. The solutions to these two key problems will be detailed below. Furthermore, in the following text, the function symbol g... k (u|Q k (abbreviated as g) k (u) is used to make the expression more concise.

[0046] (1) Parameter Q k estimation method

[0047] g k (u) is rewritten in the following form:

[0048]

[0049] in,

[0050] Observation reveals that the intermediate smoothing joint IS-PDF sequence {g} constructed by this method k {(u), k = 0, 1, ..., K} have a smooth distribution sequence {h} k (x|B k The form of g is the same as that of k = 0, 1, ..., K, and g k (u) and h k (x|B k All of them contain the same unknown parameters, namely the smoothing parameter σk and the normalization parameter pk. Therefore, g k The parameter estimation method of (u) and h k (x|B k The process is the same as that of the previous one, and the specific process is as follows.

[0051] Normalized parameter p k With p k-1 The ratio is:

[0052]

[0053] ratio η k An approximate estimate is:

[0054]

[0055] Among them, u i From the previous known intermediate smooth joint IS-PDFg k-1 The i-th sample is drawn from (u), where N is the number of samples drawn.

[0056] The estimation accuracy can be achieved through The variance coefficient (Cov) is used to measure it:

[0057]

[0058] The approximate expression is:

[0059]

[0060] In obtaining after, It can be calculated using the following formula:

[0061]

[0062] As shown in formula (11), The estimation requires knowing the exact expression for π. k (u)=Φ(-S(u) / σ k ) / Φ(-S(u) / σ k-1 However, π k There exists an unknown smoothing parameter σ in (u). k Therefore, when using formula (11) to estimate Before that, σ must be obtained. k It is important to note that, in order to ensure The estimation accuracy, variance coefficient It should satisfy a pre-specified target value. As shown in equations (10)-(13), if given from the previous known intermediate smooth joint IS-PDF g k-1 Sample u drawn from (u) i If 1 ≤ i ≤ N, then It has only a single variable σ k The function. Considering this, the unknown parameter σ can be obtained by solving the following optimization model. k :

[0063]

[0064] in It is the target value, and ||·||2 represents the 2-norm.

[0065] In this method, to avoid significantly impacting the accuracy of subsequent estimation processes, it is necessary to simultaneously adjust the parameter σ. k and p k To make accurate estimates, this method is therefore effective. The value of is subject to relatively strict restrictions. It is worth noting that... The value of cannot be too large or too small. A larger value... It may lead to Inaccurate estimations, which in turn led to subsequent calculations of p k and g k (u) The estimate is inaccurate. However, the smaller... This may result in the acquisition of intermediate smooth joint IS-PDF g k (u) Combined with the previous intermediate smoothing IS-PDF g k-1 (u) is too similar, requiring a greater number of intermediate smoothed joint IS-PDF estimates to finally achieve g. k (u) for the theoretically optimal g op(u) gradually approximates, thus greatly increasing the computational cost. Typically, It can achieve a balance between accuracy and efficiency.

[0066] (2) For intermediate smoothing combined IS-PDF g k (u) sampling method

[0067] In addition to the initial intermediate joint IS-PDF g k (u) uses traditional sampling methods such as inverse transformation sampling and accept-rejection sampling to combine IS-PDFg from other intermediate samples. k Extracting samples from (u), k≥1, is extremely difficult; this is because it is impossible to obtain intermediate joint IS-PDF g containing the performance function S(u). k The analytical form of (u), k≥1. To address this issue, a sampling method based on sample resampling and MH sampling technique that avoids the curse of dimensionality is proposed.

[0068] Assuming that the previous intermediate smooth joint IS-PDF g has been obtained k-1 N iid samples {u1, u2, ..., u3} were obtained from (u). N}. Then sample u i The weight ω(u) i It can be defined by the following formula:

[0069]

[0070] If u i =u j ,but otherwise, use to represent u i The number of times it appears. It can be proven that as N approaches infinity, (N(u) i ) / N)ω(u i It will approach g k (u i Based on this, N can be obtained by performing a resampling operation. c <N obeying g k Seed sample of (u): with ω(u) i ), 1≤i≤N as u i The probability of occurrence of ,1≤i≤N, from the total sample {u1,u2,...,u N} Randomly select N with replacement c One sample.

[0071] After obtaining Nc seed samples {u [1] ,u [2] ,...,u[Nc] After that, MH sampling will be performed starting from these seed samples, while simulating N. c The length of the strip is L = N / N c A Markov chain. For each Markov chain, starting from the currently known sample u... i Obtain the target distribution g k The next random sample u in (u) i+1 The sampling process is as follows:

[0072] 1) From the proposed distribution function q(θ|u i Generate a candidate sample θ in the process.

[0073] 2) Calculate the ratio:

[0074]

[0075] 3) With P a =min{1,r(u i Accept θ: Generate a uniformly distributed random number r in the interval (0,1). If r is less than P a Then set u i+1 =θ; otherwise, set u i+1 =u i .

[0076] In general, the current sample u is selected. i The symmetric probability density function centered at θ, i.e., satisfying q(θ|u i )=q(u i The probability density function of |θ) is used as the proposal distribution function, such as the multivariate Gaussian probability density function. Therefore, the ratio in equation (17) simplifies to:

[0077]

[0078] As shown in formula (18), if the probability density function If the dimension n is high, then the ratio It may be very small, thus causing r(u) i The θ is also very small. For this reason, the probability of the MH algorithm obtaining duplicate samples is very high, leading to its low efficiency. To address this issue, the conditional MH sampling method is introduced, the details of which are described below.

[0079] In the conditional MH sampling method, the current state vector u is selected. i The conditional multivariate Gaussian probability density function is used as the proposal distribution function:

[0080]

[0081] Where I is the identity matrix, and γ∈[0,1] is a coefficient that can be adaptively adjusted.

[0082] Substituting formula (19) into formula (18), r(u i ,θ) can be further simplified to:

[0083]

[0084] As can be seen, the ratio no longer exists in formula (20). And r(u) i ,θ) becomes a ratio that depends only on a one-dimensional smooth function. r(u i ,θ) depends only on the current state u i The conditional MH sampling method is suitable for smoothing joint IS-PDFs from high-dimensional intermediate states, based on the similarity of the system performance function with candidate states θ, rather than the dimensionality of the random variable space. k Sampling is performed in (u).

[0085] It is worth noting that the choice of parameter γ affects the performance of the MH algorithm. If γ is chosen too small, many candidate samples will be rejected; while if γ is chosen too large, the obtained samples will have excessively high correlations. Therefore, an adaptive method is used to adjust γ:

[0086]

[0087] in, It is the newly generated N c The average acceptance probability of a sample is the ratio of the number of accepted candidate samples to the total number of candidate samples. a The selection of ζ should ensure that the variation in γ gradually decreases with each adjustment. For example, ζ can be set... a =(i a +1) -1 , where i a This represents the number of adjustments performed so far. γ max ≤1 is the upper limit of the value of γ. For example, γ can be set to... max =0.99. a * =0.44 is the ideal acceptance probability for the MH algorithm.

[0088] (3) Initial intermediate smoothing combined IS-PDF g0(u)

[0089] Starting from the known initial intermediate joint IS-PDF g0(u), the subsequently estimated intermediate smooth joint IS-PDF will gradually shift towards the failure domain and progressively approach the theoretical optimal IS-PDF g. op (u).

[0090] Clearly, the choice of the initial intermediate joint IS-PDFg0(u) has a certain impact on the performance of SSAM. Typically, one can choose... As g0(u), this is only effective when evaluating systems with low reliability. However, if the system being evaluated has high reliability, It will move away from the failure domain, making it possible to... Extracting fault samples becomes extremely difficult. Furthermore, to avoid excessive computational burden during estimation, the number of samples N generated from each intermediate smoothed joint IS-PDF is finite, thus... The N samples generated may not include fault samples. Because the initial stage of stepwise estimation lacks fault sample information, intermediate steps using IS-PDFg are necessary in subsequent estimation processes. k For (u), 1≤k≤K, moving towards the fault domain can become slow and difficult. Therefore, the final obtained g K Fault samples (u) may be concentrated primarily in fault regions where the performance function S(u) doesn't perform too poorly, leading to inaccurate reliability estimates. From this perspective, given a finite sample size N, an initial intermediate smoothed joint IS-PDF g0(u) that can easily generate fault samples is a better choice. Based on this idea, this example constructs a truncated normal distribution that ignores small load levels and uses it as g0(u) to more easily obtain fault samples. Specifically, in this example, the truncated normal probability density function... Used as the initial intermediate joint IS-PDF, where, This is the load distribution of the truncated system, and its specific form is shown below.

[0091]

[0092] in, It depends on the truncated load level (for example The threshold of ).

[0093] The parameters of g0(u) are Where F L (x L ) is the original cumulative distribution function of the system load. It is worth noting that the truncated load level... It can't be too big. If The value is very large, and in the initial stage of SSAM, a large number of fault samples with poor performance function S(u) may be extracted from g0(u). Although the subsequent intermediate joint IS-PDF estimate can quickly approximate the fault domain, the final obtained g... K (u) may be over-concentrated in the fault area where S(u) performs poorly, which may lead to a deviation in the estimated reliability index.

[0094] Therefore, in some embodiments of SSAM, the intermediate joint IS-PDF sequence {g1(u),…,g} will be processed sequentially. K Estimate (u)} until the variance coefficient (COV) of the reliability index obtained in the Kth estimation meets the preset accuracy, that is, the variance coefficient is less than or equal to β. max Based on the final obtained g K The estimated values ​​of the probability of load shedding (LOLP), expected power shortage (EENS), and COV of (u) are shown below:

[0095]

[0096] Among them, u i It is a joint IS-PDFg from the middle K The sample drawn from (u), N K L is the number of samples drawn. c (u i ) represents the system state u determined through optimal load shedding analysis. i The amount of load reduction, It is the likelihood ratio.

[0097] In summary, in some embodiments of the present invention, the SSAM algorithm can be implemented according to the following specific steps:

[0098] (1) Setting target value The maximum variance coefficient β of the reliability index max And let k = 0.

[0099] (2) From the initial intermediate joint IS-PDF N samples {u1,…,u} are drawn from the sample. N}, where σ0=∞, p0=0.5.

[0100] (3) Let k = k + 1. Using the previously known g k-1 N samples {u1,…,u} drawn from (u) N}, calculate the next intermediate joint IS-PDF using formula (15). k The smoothing parameter σ of (u) k .

[0101] (4) Estimate the ratio using formula (11)

[0102] (5) Use formula (14) to calculate g k Normalized parameters of (u)

[0103] (6) Use formula (16) for g k-1 The sample {u1,…,u} of (u) N} Calculate the weights {ω(u1),…,ω(u... N Resampling is performed to obtain samples that follow a distribution g. k (u) of N c Seed {u [1] ,u [2] ,…,u [Nc] Let j = 1, L = N / N c .

[0104] (7) Parallel execution of N c The Markov chain is simulated once, and each Markov chain is advanced one state forward using the conditional MH sampling algorithm to obtain N. c A new sample is generated. Then, using formula (21), the parameter γ is updated using these samples, and N is set to N. k =j·N c .

[0105] (8) Let K = k, then evaluate the reliability indices LOLP and EENS using formulas (23) and (25) respectively, and then obtain the results using formulas (24) and (26). and Cov(EENS). If Then set j = j + 1 and go to step (9); if The obtained intermediate joint IS-PDF g k (u) is the last distribution of the sequence, and the calculation result is output and the algorithm ends.

[0106] (9) If j≤L, then go to step (7); otherwise, obtain N=L·N c Each sample is from g k (u) sample and proceed to step (3).

[0107] To further improve the convergence speed of SSAM, some embodiments of the present invention also propose a comprehensive reliability index that can make full use of samples drawn from all intermediate joint IS-PDFs.

[0108] Assuming that, as SSAM proceeds, from the intermediate joint IS-PDF g k (u), 0≤k≤K, generated N from it. k Sample u i , (1≤i≤N) k Based on g k For a sample of (u), 0≤k≤K, the expected value and variance of various reliability indicators (such as LOLP, EENS, etc.) can be estimated using the following unified expression:

[0109]

[0110] in, It is the likelihood ratio, u i From g k The samples generated in (u) are F(u), which is the reliability index function in the u-space. If F(u) = I {S(u)＜0} , This refers to the LOLP metric. If F(u) = E c (u), This indicates the EENS index.

[0111] Let κ k for The weights, and satisfying The definition of the comprehensive reliability index is as follows:

[0112]

[0113] Definition makes The weight κ with the smallest variance k 0≤k≤K represents the optimal weight, and its specific form is shown below:

[0114]

[0115] In the optimal weight κ k Under (0≤k≤K), variance for:

[0116]

[0117] coefficient of variance The following formula can be used to estimate:

[0118]

[0119] An exemplary step of the SSAM algorithm using a comprehensive reliability metric is shown below:

[0120] Step (1) Set the following parameters:

[0121] ① The next estimated intermediate joint IS-PDF is numbered k = K = 1, the maximum number of samples N drawn from the intermediate joint IS-PDF each time is N (usually 10000-20000), and the number of seed samples is N. c (e.g., 10%-20% of N), and the Markov chain length L = N / N c ;

[0122] ② Maximum variance coefficient β max(usually 1%-2%) target value (usually 1%-2%)

[0123] ③ Cut-off load level (e.g., 0.8pu), the initial value of parameter γ (usually 0.5).

[0124] Step (2) transforms the x-space, which has different types of random variables, into the u-space, which has independent standard normal random variables.

[0125] Step (3) uses truncation As an intermediate joint IS-PDF, in which Draw N0 = N samples {u1,…,u} from g0(u). N}, and by configuring the system state x in the x space i =T(u i ) Perform power flow or optimal power flow analysis to determine S(u i (1≤i≤N) s ).

[0126] Step (4) Calculate g sequentially k The smoothing parameter σ of (u) k ,ratio and normalization parameters

[0127] Step (5) calculates the sample {u1,...,u} using formula (16). N}~g k-1 The weights of (u) are {ω(u1),...,ω(u)}. N Let ω(u1) be u i The probability of occurrence is determined, and a resampling operation is performed to obtain the probability from {u1,...,u Ns N are randomly selected with replacement. c Seeds Let j = 1.

[0128] Step (6) Run N in parallel c The Markov chain is simulated once, and each Markov chain is advanced one state forward using the conditional MH sampling algorithm to obtain N. c A new sample is drawn. During the sampling process, the candidate sample θ is drawn from the suggested PDF in formula (19), S(θ) is determined by the power flow or optimal power flow analysis of the system state x=T(θ) in x space, and the acceptance rate r(u i ,θ) is calculated from (20).

[0129] Step (7) uses formula (21) to utilize the newly generated N cEach sample updates the parameter γ, and N is set to... k =j·N c .

[0130] Step (8) calculates the expected values ​​of the LOLP and EENS indices using formulas (27) and (28) respectively. and variance

[0131] Step (9) calculates the expected values ​​of the LOLP and EENS indices using formula (30). weight κ k , (0≤k≤K).

[0132] Step (10) calculates the combined reliability index of LOLP and EENS using formulas (29), (30) and (31). variance and variance coefficient

[0133] Step (11) If the LOLP metric Greater than β max If yes, proceed to step (12); otherwise, proceed to step (13).

[0134] Step (12) Let j = j + 1. If j ≤ L, then go to step (7); otherwise, obtain N = L·N c Each sample is from g k Take a sample of (u) and set k = k + 1, K = k, then proceed to step (4).

[0135] Step (13) If the EENS index is Less than β max If the result is positive, output the obtained reliability index and stop the simulation; otherwise, set j = j + 1, and then go to step (6).

[0136] Calculation Example

[0137] In this example, reliability assessments were performed on the IEEE-RTS79 and IEEE-RTS96 test systems and their modified versions. The cut-in wind speed, rated wind speed, and cut-out wind speed of the wind turbine were 3 m / s, 13 m / s, and 25 m / s, respectively. This example comprehensively examines the efficiency and performance of the original Non-Sequential Monte Carlo Simulation (NCMCS), the SSAM method proposed in this invention, and the improved cross-entropy method (iCEM) in reliability assessment. In iCEM, normal, log-normal, and GMM were selected as the IS-PDF of the system load, and three different iCEM methods were constructed based on these. N iCEM L and iCEMG And set the multi-level number parameter ρ a =0.1, target value In SSAM, the target value, seed sample size, and cutoff load level are set to... N c =0.2N, Furthermore, unless otherwise specified, the sample size N is set to 20000 by default when updating parameters in iCEM for each iteration or when estimating parameters for each intermediate joint IS-PDF in SSAM. Additionally, the maximum variance coefficient β is set... max =0.02 is used as the convergence criterion for all methods.

[0138] (1) Reliability assessment results of the IEEE-RTS79 test system and its modified version

[0139] a. Original IEEE-RTS79 test system

[0140] To verify the effectiveness of the proposed SSAM method, a reliability assessment was first performed on the original IEEE-RTS79. Since the reliability of the original IEEE-RTS79 is not very good, the number of samples N required for each iteration was set to a specific value. s =10000 is more appropriate. The evaluation results are shown in Table 1, where K represents the total number of iterations in iCEM or the total number of consecutive estimates in SSAM, and N... K N represents the number of samples required for the iCEM main sampling phase or the final estimation phase of SSAM. t =K·N+N K This represents the total number of samples throughout the entire simulation process.

[0141] Table 1. Reliability assessment results of the IEEE-RTS79 system obtained by different methods.

[0142]

[0143] As shown in Table 1, all iCEM N iCEM L iCEM G Both the SSAM and NCM methods can obtain reliability metrics similar to NCMCS, but their total sample size N tBoth the CPU time T and the time required for IS-PDF are significantly less than those for NCMCS, because IS-PDF significantly improves their ability to capture rare failure events. The results also show that SSAM has higher efficiency than CEM. Furthermore, the proposed SSAM method can obtain a "true" joint IS-PDF with a similar functional form to the theoretical zero-variance IS-PDF, while the iCEM method can only obtain a "pseudo" joint IS-PDF that simply combines the IS-PDFs of all random variables. To examine the difference between the actually obtained IS-PDF and the theoretical zero-variance IS-PDF, this paper uses some commonly used distance metrics, such as Kullback-Leibler distance, Hellinger distance, and Pearson χ². 2 Metrics and α-divergence (parameters α = 0.5, 0.0, and -0.5) are used to measure the closeness of the IS-PDF obtained in the fault domain to the theoretically optimal IS-PDF.

[0144] As shown in Table 2, SSAM is significantly better than the three iCEM methods because all distance metrics of SSAM are smaller than the corresponding metrics of iCEM methods. This indicates that SSAM can obtain an IS-PDF that is closer to the theoretical optimal IS-PDF.

[0145] Table 2. Distance between IS-PDF obtained by different methods and the optimal IS-PDF.

[0146]

[0147] b. IEEE-RTS79 test system with reduced peak load

[0148] To demonstrate the superiority of SSAM in more reliable scenarios, the original IEEE-RTS79 was modified to reduce the load level to 80% of its original level, thereby significantly improving reliability. This modified IEEE-RTS79 is called MRTS79-1, and the evaluation results are shown in Table 3.

[0149] Table 3. Reliability assessment results of the MRTS79-1 system obtained by different methods.

[0150]

[0151] As shown in Table 3, the LOLP metric is on the order of 10. -5 This indicates that MRTS979-1 has very high reliability, therefore the total sample size N required to use the NCMCS method is low. tThe CPU time T also increases significantly. Compared to the NCMCS method, which suffers from heavy computational burden, the iCEM and SSAM methods still demonstrate significant performance advantages, even when the evaluated system has extremely high reliability. They greatly reduce the total sample size N required for evaluation while ensuring accurate reliability indicators. t And CPU time T. However, the results in Table 3 show that iCEM G It is less computationally efficient than other iCEM methods. The root cause lies in the fact that iCEM... G SSAM relies on a complex multi-step EM algorithm to update the parameters of the IS-PDF constructed based on GMM. Compared with the iCEM method, SSAM reduces the required sample size and computation time while obtaining accurate reliability metrics. Unlike the iCEM method, which requires K=3 parameter iterations for updates, SSAM only requires K=2 consecutive estimations of the intermediate joint IS-PDF.

[0152] Because the intermediate joint IS-PDF is a high-dimensional PDF, it is difficult to visualize its actual movement. Figure 2 The marginal PDF curves of the system load derived from the high-dimensional intermediate joint IS-PDF are plotted to further visualize the performance of the SSAM method.

[0153] exist Figure 1 middle, and g0(u L g1(u) represents the raw and truncated PDFs of the system load in u-space, respectively. L ) and g2(u L ) are the edge PDFs of the system load derived from the estimated intermediate joint IS-PDF, and f(x) are the edge PDFs of the system load derived from the intermediate joint IS-PDF. L ),g0(x L )~g2(x L ) represents the corresponding PDF in the x-space. For example... Figure 1 As shown, as the estimation progresses, the density peak of the system load edge PDF will gradually move to the risk region (high load region), which indicates that the SSAM method can effectively identify and highlight the important regions of variables.

[0154] c. IEEE-RTS79 test system containing multiple wind farms

[0155] To verify the performance of SSAM in scenarios with correlated random variables, the following modifications were made to IEEE-RTS79:

[0156] ①Increase the load level to 110% of the original level.

[0157] ② Six 200MW wind farms were added at nodes 13, 15, 16, 18, 21 and 23.

[0158] The revised version was named MRTS79-2, and the number of samples required for each iteration was set to N = 15000. The reliability evaluation results of MRTS79- are shown in Table 4.

[0159] Table 4. Reliability assessment results of the MRTS79-2 system obtained by different methods.

[0160]

[0161] As shown in Table 4, the SSAM method can obtain the same accurate reliability index as the CEM method, but requires a smaller sample size N. t Fewer steps are required, and the CPU time T used is also shorter. Furthermore, the number of consecutive estimations required in SSAM is K=2, while the number of parameter iterations required in iCEM is K=3, indicating that SSAM is more efficient than iCEM in finding IS-PDFs.

[0162] (2) Reliability assessment results of the IEEE-RTS96 test system and its modified version

[0163] a. Original IEEE-RTS96 test system

[0164] To further demonstrate the effectiveness of the SSAM method in reliability analysis of power systems with higher reliability and larger scale, this example further evaluates the reliability of the IEEE-RTS96 test system, and the evaluation results are shown in Table 5.

[0165] Table 5. Reliability assessment results of the IEEE-RTS96 system obtained by different methods.

[0166]

[0167] As shown in Table 5, when conducting reliability assessments for power systems with higher reliability and larger scale, the NCMCS method requires extracting a large number of system state samples, increasing the computational burden and resulting in a considerably long time required to complete the reliability assessment. However, thanks to the acceleration effect of the importance sampling method on reliability assessment, the iCEM method and the proposed SSAM method converge rapidly, using approximately 0.1% of the CPU time of NCMCS. Compared to the iCEM method, SSAM not only achieves a reliability index very close to that of NCMCS, but also uses the least CPU time T among all methods, thus exhibiting better performance and efficiency. Furthermore, iCEM… N and iCEM LThe obtained reliability metrics contain significant errors, especially the LOLP metric. This is because the normal and log-normal distributions are symmetrical and right-skewed, respectively, and can only construct IS-PDFs of fixed shapes, thus affecting iCEM. N and iCEM L The obtained IS-PDF cannot fully highlight the important areas of the system, which means that some samples with low load levels that contribute to the system reliability index cannot be obtained, thus resulting in a lower obtained reliability index.

[0168] To more intuitively demonstrate the differences between SSAM and iCEM, and to highlight the superiority of the SSAM method, in Figure 2 The IS-PDF curves of system load obtained by different methods are plotted.

[0169] Figure 2 (a) Shows the IS-PDF of the system load initially obtained by the iCEM and SSAM methods. To highlight load levels that significantly contribute to system failures, iCEM... N and iCEM L The IS-PDF curves obtained by this method have peak values ​​at relatively high load levels (approximately 0.95 pu), and a large portion of their right tail exceeds the annual peak load (1.00 pu). However, the corresponding system load samples cannot be extracted in practice. Therefore, iCEM... N and iCEM L The method actually uses a truncated IS-PDF in the master sampling phase, such as... Figure 3 As shown in (b). Since the normal and log-normal distributions are symmetrical and right-skewed distributions, respectively, they are only suitable for describing important regions with similar shapes. Therefore, iCEM N and iCEM L The left tail of the IS-PDF obtained by this method only extends to the region with a high load level (approximately 0.85 pu). This means that during the main sampling phase, iCEM... N and iCEM L These methods almost never obtain samples at low load levels, especially those low load levels (e.g., 0.80–0.85 pu) that have a non-negligible contribution to system reliability metrics, thus resulting in underestimating reliability metrics. In contrast, GMM can construct IS-PDFs of arbitrary shapes, enabling iCEM to... G The resulting IS-PDF ensures that the right tail does not exceed the annual peak load, while allowing the left tail to extend into areas with lower load levels (approximately 0.80 pu). Therefore, iCEM G The reliability index has a small error. However, the GMM parameter update is complex, involving multiple EM algorithm steps, which leads to the limitations of iCEM.G It is more time-consuming. Compared to iCEM, SSAM does not require selecting the IS-PDF type of the system load and has no complex parameter update process, and can also obtain an IS-PDF with a reasonable edge shape. Therefore, SSAM can obtain accurate reliability indicators in a shorter time.

[0170] b. IEEE-RTS96 test system with multiple wind farms

[0171] To further verify the performance of SSAM in scenarios with correlated random variables, the following modifications were made to IEEE-RTS96:

[0172] ①Increase the load level to 110% of the original level.

[0173] ② Six 200MW wind farms were added at nodes 13, 15, 16, 18, 21 and 23.

[0174] The revised version is named MRTS79-2, and its reliability assessment results are shown in Table 6.

[0175] Table 6. Reliability assessment results of the MRTS96 system obtained by different methods.

[0176]

[0177] As shown in Table 6, consistent with expectations, SSAM can derive accurate reliability metrics with minimal computation. However, because IS-PDF, constructed based on normal and log-normal distributions, cannot comprehensively describe the critical regions of the system, some samples at low load levels that contribute to the system's reliability metrics cannot be obtained. Therefore, iCEM... N and iCEM L The obtained reliability metrics have a large negative error.

[0178] (3) Analysis results of factors influencing SSAM

[0179] target value The value of σ affects the smoothing parameter σ of the intermediate joint IS-PDF. k and normalized parameter p k The estimation has a direct impact. Cut-off load level. The choice of sampling method plays a crucial role in constructing the initial intermediate joint IS-PDF. Furthermore, the choice of MH sampling method significantly impacts the effectiveness of sampling from the intermediate joint IS-PDF. The effects of all these factors will be discussed below.

[0180] a. Target value Impact analysis results

[0181] In order to analyze the target value The impact on SSAM performance, used in different The SSAM under the given values ​​is evaluated for IEEE-RTS96. As continuous estimation proceeds, the smoothing parameter σ... k and normalized parameter p k The change curve is as follows Figure 3 As shown in Table 7, the reliability assessment results are as follows.

[0182] Table 7 Differences Evaluation results of IEEE-RTS96

[0183]

[0184] like Figure 3 As shown, in each given Below, parameter σ k (Note σ0=+∞) and p k Both decrease as expected with the increase of the successive estimation number k, indicating that the estimated intermediate joint IS-PDF g decreases as the successive estimation proceeds. k (u)(0≤k≤K) will approach the theoretically optimal IS-PDF as closely as possible. However, as shown in Table 9, The choice of [specific value] does indeed have a significant impact on the performance of SSAM. Due to the target value... Essentially, it is about estimating the parameter σ. k and p k The required precision. Therefore, a very small target value, such as... This implies very high accuracy requirements. To meet such stringent requirements, the next estimation will use intermediate joint IS-PDFg. k (u) will not differ significantly from the previous gk-1(u), therefore SSAM needs to perform K=5 intermediate joint IS-PDF estimations to achieve an approximate estimate of the theoretically optimal IS-PDF, resulting in a significant increase in computational burden and a longer evaluation time. Conversely, a larger objective value, such as This means that the accuracy requirement is relatively low, so the next estimate of g... k (u) will be combined with the previous g k-1 The significant differences in (u) mean that SSAM only needs K=2 intermediate joint IS-PDF estimations to converge quickly, thus reducing the evaluation time. However, the lower accuracy requirement leads to inaccurate estimation of the intermediate joint IS-PDF, resulting in a large error in the obtained reliability index. Therefore, as shown in Table 7, to achieve a balance between accuracy and efficiency, The value should be kept between 1% and 3%.

[0185] (2) Cut-off load level Impact analysis results

[0186] In order to analyze The impact on SSAM performance, used in different The SSAM under the given values ​​was used to evaluate the IEEE-RTS96, and the reliability evaluation results are shown in Table 8. It is worth noting that the minimum system load is 0.34, therefore... This means choosing the original distribution. As the initial intermediate joint IS-PDF g0(u).

[0187] Table 8 Differences Evaluation results of IEEE-RTS96

[0188]

[0189] As shown in Table 8, when the truncated load level is very low, for example... or SSAM not only requires a long CPU time T, but also yields relatively large errors in the reliability metrics obtained. Smaller... This means that the initial intermediate joint IS-PDFg0(u) is far from the theoretically optimal IS-PDF. It should also be noted that IEEE-RTS96 is a highly reliable system because its LOLP is very low, reaching 10. -5 The order of magnitude, therefore, in a sample size of N s When the threshold is 20000, capturing fault events using g0(u) is difficult, causing the subsequently estimated intermediate joint IS-PDF to slowly approach the fault domain, resulting in a long estimation time. Furthermore, due to the lack of sufficient fault samples to estimate the subsequent intermediate joint IS-PDF, the final obtained intermediate joint IS-PDF g... K (u) might overemphasize fault areas with lower severity, i.e., fault areas where the performance function S(u) value is closer to 0, thus leading to an underestimation of reliability metrics. Conversely, if the truncated load level is very high, for example... The initial intermediate joint IS-PDF g0(u) will be very close to the theoretical optimal IS-PDF. Furthermore, due to the high load level, it is easy to obtain samples with very severe fault conditions from g0(u). Therefore, the subsequently estimated intermediate joint IS-PDF will quickly approach the fault domain, resulting in a shorter estimation time. However, because too many samples with very severe fault conditions are used to estimate the intermediate joint IS-PDF, the final obtained intermediate joint IS-PDF g... K(u) may overemphasize fault areas with very severe failures, i.e., fault areas with smaller performance function S(u) values, which will also lead to an underestimation of the reliability index. Therefore, as shown in Table 8, A value between 0.6 and 0.8 is very reasonable.

[0190] (3) Results of the impact analysis of the MH sampling method

[0191] To verify the advantages of conditional MH sampling, Table 9 compares the reliability assessment results of IEEE-RTS96 using SSAM based on the traditional MH sampling method and SSAM based on the conditional MH sampling method. For distinction, SSAM based on the traditional MH sampling method is referred to as SSAM. T .also, Figure 4 The study also compared the use of two MH sampling methods for intermediate joint IS-PDF g. k (u), 1≤k≤2, represents the change in acceptance rate during sampling. Since each time g... k (u), 1≤k≤2, N is retained during sampling. c =0.2N=4000 seed samples, therefore in Figure 4 The number of samples drawn by both conditional MH sampling and conventional MH sampling is the number after removing these seed samples.

[0192] Table 9 Differences Evaluation results of IEEE-RTS96

[0193]

[0194]

[0195] like Figure 4 As shown, due to the dimensionality curse, in SSAM based on traditional MH sampling... T In this method, the acceptance rate is significantly reduced, almost approaching zero. This phenomenon indicates that newly generated candidate samples are frequently rejected during traditional MH sampling, resulting in a large number of duplicate samples. This low acceptance rate may cause samples generated by traditional MH sampling methods to fail to effectively follow the target distribution g. k (u), 1≤k≤2, thus enabling SSAM T The obtained reliability metrics are significantly biased compared to those obtained using the NCMCS method. However, in SSAM employing conditional MH sampling, the sampling acceptance rate is maintained at a reasonable level, thus ensuring the accuracy of the obtained reliability metrics.

[0196] The above descriptions are merely embodiments of the present invention, and common knowledge such as specific technical solutions and / or characteristics are not described in detail here. It should be noted that those skilled in the art can make various modifications and improvements without departing from the technical solutions of the present invention, and these should also be considered within the scope of protection of the present invention. These modifications and improvements will not affect the effectiveness of the implementation of the present invention or the practicality of the patent. The scope of protection claimed by the present invention should be determined by the content of its claims, and the specific embodiments described in the specification can be used to interpret the content of the claims.

Claims

1. A power grid reliability assessment method based on successive smoothing estimation of the optimal importance sampling density function, characterized in that, Includes the following steps: Step 1: Using the normal transformation method, the mixed variable space of the original discrete and continuous variables is transformed into a multivariate standard normal space, and on this basis, a unified importance sampling framework for all random variables in the multivariate standard normal space is constructed. Step 2: Based on the important sampling framework, construct an intermediate smooth joint IS-PDF sequence, using a truncated Gaussian distribution as its initial distribution, and sequentially estimate the subsequent distributions in the intermediate smooth joint IS-PDF sequence. At the same time, use the estimated smooth joint IS-PDF for sampling and power grid reliability index evaluation. Until the variance coefficient of the reliability index obtained at the Kth estimation meets the preset accuracy, the current power grid reliability index is used as the final evaluation index. The process of constructing the intermediate smooth joint IS-PDF sequence and sequentially estimating the subsequent distributions in the intermediate smooth joint IS-PDF sequence includes: finding the optimal IS-PDF in the discrete x-space. Transformation into optimality in continuous u space : （6） in = Let n represent the n-dimensional standard normal probability density function, n = D + W + 1. , and Let D-dimensional, W-dimensional, and one-dimensional standard normal probability density functions be represented respectively. It is an indicator function: if ,but ,otherwise, , It is a system performance function defined in space u, which can be used to analyze the system state in space x. Determined by power flow or optimal power flow analysis; Combine formula (6) with the smoothing function Combining these, we obtain the optimal zero-variance IS-PDF for smoothing: （7） in, It is a sufficiently small smoothing parameter. It is a normalized parameter. ={ , }; Establish intermediate smooth joint IS-PDF sequence The initial distribution is a truncated Gaussian distribution. And estimate subsequent parameters in turn. ; parameter The estimation methods include: Will Rewrite it in the following form: （8） in, ; ratio The approximate estimate is: （11） in, From the previous known intermediate smooth joint IS-PDF The i-th sample is drawn from the sample, and N is the number of samples drawn; the unknown parameters are obtained by solving the following optimization model. , and thus obtain : (13) ， in It is the target value, ||·||2 represents the 2-norm, where yes The variance coefficient; In obtaining after, Calculated using the following formula: (14)。 2. The method according to claim 1, characterized in that: Step 1 includes: 1.1) System component state vector Each discrete variable in Equation (2) is converted into the corresponding standard normal variable. : （2） in, It is the threshold between the fault state and the normal state. It is the original unavailability rate of the j-th component. It is the inverse function of the standard normal cumulative distribution function. It is an indicator function: if ≤ ,but ,otherwise, ; 1.2) System load variables Transform it into the corresponding standard normal variable through equal probability transformation. : （3） in, It is the original cumulative distribution function The inverse function; 1.3) Correlated wind speed vectors The corresponding independent standard normal random vector is transformed by the Nataf transform. .

3. The method according to claim 2, characterized in that, Step 1.3) involves converting the vector into a corresponding independent standard normal random vector using the Nataf transform, which includes: First, through equal probability transformation Transform into a vector with zero mean and a correlation matrix The normally distributed random vector z, i.e. ,in It is the cumulative distribution function inverse function; set up yes The correlation matrix and They are and The element in the i-th row and j-th column, It can be obtained by solving the following formula: （4） in, and They are and The mean, and They are and standard deviation It has a correlation coefficient The bivariate standard normal probability density function; After obtaining Then, the results were calculated using Cholesky decomposition. lower triangular matrix ,Right now Then, the correlated normal random vector z is formed by independent standard normal random vectors. Indicates, that is ; Correlated wind speed vectors The corresponding independent standard normal random vector The conversion relationship between them can be expressed by equation (5): （5）。 4. The method according to claim 1, characterized in that, In step 2, for intermediate smoothing combined IS-PDF The sampling method employs MH sampling techniques based on sample resampling and avoiding the curse of dimensionality, including: The previous intermediate smooth joint IS-PDF has been completed. N iid samples were obtained. ,sample weight Defined by the following formula: （16） set up It is an indicator function such that: if ,but ,otherwise, ;use to indicate Number of times it appears; Obtained by performing a resampling operation Obedience Seed samples: As The probability of occurrence of , from the total sample Random selection with replacement One sample; In obtaining Seed samples Then, MH sampling is performed starting from these seed samples, while simultaneously simulating... The length of the strip is A Markov chain, for each of which, starting from the currently known samples Obtain the target distribution The next random sample The sampling process is as follows: 1) From the proposed distribution function Generate a candidate sample ; 2) Calculate the ratio: （17） 3) with accept Generate a uniformly distributed random number r within the interval (0,1). If r is less than P... a Then set Otherwise, set .

5. The method according to claim 4, characterized in that, The multivariate Gaussian probability density function is selected as the proposed distribution function, and the conditional MH sampling method is adopted. The specific details are as follows: Choose the current state vector The conditional multivariate Gaussian probability density function is used as the proposal distribution function: （19） Where I is the identity matrix, It is a coefficient that can be adaptively adjusted; Will Simplified to: （20） Use an adaptive method to adjust γ: （21） in, It is newly generated The average acceptance probability of a sample is the ratio of the number of accepted candidate samples to the total number of candidate samples. It is the upper limit of the value of γ. This is the ideal acceptance probability of the MH algorithm. This is an optional parameter. The selection of γ should ensure that the variation of γ gradually decreases with each adjustment.

6. The method according to claim 4, characterized in that, set up ,in This represents the number of adjustments made so far.

7. The method according to claim 1, characterized in that, In step 2, the truncated normal probability density function is used. Used as the initial intermediate joint IS-PDF, where, This is the load distribution of the truncated system, and its specific form is shown below; （22） in, It depends on the truncation load level The threshold; have, The parameters are ,in It is the original cumulative distribution function of the system load.

8. The method according to claim 5, characterized in that, In step 2, based on the obtained The power grid reliability assessment includes: probability of load shedding (LOLP), expected power shortage (EENS), and its variance coefficient (Cov), which are calculated as follows: With each estimation, a maximum sample size N is preset to obtain the maximum Markov chain length. ,set up K is the number of samples drawn, where K=k. If j=1, then: （23） （24） （25） （26） in, From the middle of IS-PDF The sample drawn from This represents the system state determined through optimal load shedding analysis. The amount of load reduction, It is the likelihood ratio; If the variance coefficient of the LOLP metric is greater than the maximum variance coefficient Let j = j + 1. If at this time... Then return to the seed sample extraction step and use the newly generated sample. Update parameter γ for each seed sample, and then re-sampling. The sample is used to calculate the index; if at this time Then the obtained Each sample is from The sample is used to continue the next round of estimation; If the variance coefficient of the LOLP index is less than or equal to the maximum variance coefficient The obtained intermediate joint IS-PDF will then be As the last distribution in the sequence, output the calculation results.

9. The method according to claim 5, characterized in that, The power grid reliability assessment index is a comprehensive reliability index that fully utilizes samples drawn from all intermediate joint IS-PDFs. Its specific calculation includes: Based on from For a sample of 0 ≤ k ≤ K, the expected values ​​and variances of various reliability metrics are estimated using the following unified expression: （27） （28） set up for The weights, and satisfying The definition of the comprehensive reliability index is as follows: （29） Definition makes The weight with the smallest variance 0≤k≤K represents the optimal weight, and its specific form is shown below: （30） In optimal weight Under (0≤k≤K), Variance Var( )for: （31） coefficient of variance The estimation is performed using the following formula: ；(32) in, It is a reliability index function in u-space, if , This refers to the overall reliability index of LOLP. , This represents the overall reliability index of EENS; If the variance coefficient of the overall reliability index of LOLP is greater than the maximum variance coefficient Let j = j + 1. If at this time... Then return to the seed sample extraction step and use the newly generated sample. Update parameter γ for each seed sample, and then re-sampling. The sample is used to calculate the index; if at this time Then the obtained Each sample is from The sample is used to continue the next round of estimation; If the variance coefficient of the LOLP comprehensive reliability index is less than or equal to the maximum variance coefficient Then determine whether the variance coefficient of the EENS index is less than If so, the obtained intermediate joint IS-PDF As the last distribution of the sequence, output the calculation result; otherwise, set... and rerun Sub-Markov chain simulation yields... A new sample is used to update the parameter γ, and then a new sample is drawn. The sample is used to calculate indicators.

Images