Multi-Angle Reliability Evaluation of Grid-Connected Wind Farms with Energy Storage Based on Latin Hypercube Important Sampling

: Aiming to combat the problems of slow speed and poor accuracy of reliability evaluation of the power system in wind farms with energy storage, this paper proposes a method of reliability evaluation based on Latin hypercube important sampling (LHIS). Firstly, we aimed to establish the Latin hypercube important sampling evaluation model by combining the Latin hypercube sampling method with the important sampling method. Secondly, we aimed to optimize the sample probability distribution of the components and conduct hierarchical sampling of the system. Then, the comprehensive risk indicator (CRI) was proposed to evaluate the operational risk and the wind storage generation interrupted energy beneﬁt (WSGIEB) was proposed to evaluate the contribution of the reliability. Finally, simulation experiments were carried out through various power system operation scenarios. The simulation results show that the proposed method is 47% higher than the improving importance sampling method (IM-IS) in evaluation speed and 33% higher than the improving importance sampling method in calculation accuracy.


Introduction
Energy storage systems can improve the fluctuation of wind farm output; however, the operation conditions of the wind turbines and energy storage systems, the energy storage capacity and charging and discharging strategy, the wind power access capacity, and so on, will affect the reliability of the power system to varying degrees. Therefore, it is of great significance to study how to accurately evaluate the reliability of grid-connected wind farms with energy storage systems under different operating conditions for the stable operation of the grid under the background of the high proportion of wind power access [1]. The key to the accuracy of power grid reliability evaluation of wind farms and energy storage systems lies in selecting the power grid reliability evaluation method and constructing the evaluation index.
Aiming toward the selection of evaluation methods, early power system reliability evaluation mainly relied on the analytical method, whose computational complexity depended on the complexity of the power system solution space. However, as a large number of renewable energy sources were connected to the grid, the random variables in the system increased greatly, and their correlation structures became more and more complex. Therefore, the bottleneck of traditional analytical methods cannot be completely solved. Therefore, the simulation method of reliability evaluation based on Monte Carlo sampling (MCS) has received much attention [2,3]. Li. [4] used the sequential Monte Carlo method and the non-sequential Monte Carlo method to evaluate the power system and to evaluate the reliability of the power system. Finally, a variety of differentiated power grid operation scenarios are set up to verify the effectiveness of the proposed method. The results of one example show that this method can effectively evaluate the reliability of the grid after the wind power and storage systems are connected to the system and provide a useful reference for the design and dispatch strategy of the power system of wind farms and energy storage systems.
The main contributions of this paper are as follows: (1) To solve the issues of poor convergence and low accuracy of grid reliability evaluation of wind farms and energy storage, this paper proposes a power system reliability evaluation method based on LHIS. First, the sample space probability distribution of the sample object is optimized, and the optimized probability distribution is stratified sampled. Then, the correlation of the sampling matrix is reduced. (2) To solve the problem of the lack of risk analysis and the severity of the consequences in the case, this paper evaluates the power grid from two perspectives: the operational risk and the contribution to system reliability. Three operation scenarios are set for connecting wind farms with energy storage systems at different nodes, connecting wind farms of different capacities with energy storage systems, and connecting wind farms with different energy storage capacities. It has been verified that the method can evaluate the reliability of the grid more comprehensively and accurately.
The other parts of this paper are arranged as follows: Section 2 introduces the wind farm output model. Section 3 establishes the LHIS evaluation model and provides the specific evaluation process then defines the evaluation indicators. In Section 4, the LHIS method is compared with LHS method and IM-IS method in order to verify the efficiency and universality of this method. Simultaneously, the different power system connection modes based on the LHIS method are simulated and analyzed. Finally, Section 5 ends the paper by summarizing the main conclusions.

Wind Farm Output Model
The statistics on a large amount of wind speed data show that the Weibull distribution of two parameters can describe the variation law of wind speed. This distribution is represented by Formula (1) [18]: The distribution function is expressed by Formula (2): where λ is the scale parameter, representing the average wind speed of the wind farms, v is the wind speed value, and k is the shape parameter. The probability sampling value of wind speed can be obtained by changing the distribution function, and the distribution function is transformed into the following: Because the range of 1 − x and x in the formula is [0, 1], the distribution function can be rewritten as follows: Energies 2023, 16, 6427 4 of 18 The relationship between the fan output power and the wind speed can be expressed by formula [19], as follows: where P t (v) is the real-time output power of the fan; v ci and v co are the cut-in wind speed and the cut-out wind speed of the fan, respectively; and v r is the rated wind speed.

Charge and Discharge Model of Energy Storage System
According to the output and load of the wind farm, there are two charging and discharging strategies of the energy storage system. Strategy 1: When the output of the wind farm is less than P L (t)·η, the energy storage system begins to discharge. In order to limit the output power of the wind farm and the energy storage system, it is required that the sum of the output of the wind farm and the energy storage battery should not exceed P L (t)·η. The discharge power of the energy storage battery is expressed as follows [20,21]: where P max disch represents the maximum discharge power of the energy storage system, η is the wind power permeability, P W (t) is the power generated by the wind farm, and P L (t) is the maximum load of the system. When the output of the wind turbine is greater than P L (t)·η, the energy storage system begins to charge. The charging power of the energy storage system is expressed as follows: where P max ch represents the maximum charging power of the energy storage battery. Strategy 2: When the output of the wind farm is less than P L (t)·η, and the output of the conventional power plant cannot meet the demand of the residual load, the energy storage system starts to discharge. The discharge power of the energy storage system is expressed as follows [20,21]: When the output of the wind farm is greater than P L (t)·η, the energy storage system begins to charge. The charging power of the energy storage battery is expressed in Formula (7).

Latin Hypercube Important Sampling (LHIS)
Latin hypercube sampling (LHS) adopts a stratified sampling method, which can effectively reflect the overall distribution of the random variables through sampling values [7]; however, the sampling method of Latin hypercube sampling at the tail of probability density function will affect the accuracy of the algorithm. In this paper, the LHIS method is established by combining the LHS method with the IS method.

Sampling Principle
The LHIS method is carried out as follows: Firstly, the probability distribution of the sample space of the sampling object is optimized based on the IS method. Then, the optimized probability distribution is sampled by the LHS method. Finally, the correlation of the sampling matrix is reduced by the Cholesky decomposition method. The specific steps are as follows: (1) Optimize the spatial probability distribution of the sampling objects based on the IS method. Then, make the mathematical expectation E(X) of the original sample of the sampling object unchanged and replace the original sample space probability distribution P(X) with the new probability distribution P * (X). If the new variance is 0, then the P * op (X) is the optimal probability distribution at this time. The method is as follows [22]: where F(X) is the sampling function and E(F) is the sampling expectation value. Because it cannot be determined in the sampling process, the probability distribution of the optimized sample space can only be approximately optimal. From Formula (9), P * (X) is proportional to P(X), let m be the proportional coefficient, then Formula (9) is rewritten as: Let a component i in the power system have a sample space probability distribution of P(X i ), and take the optimal multiplier, γ, of the intermediate variable to construct its new probability distribution, P * (X i ), as shown in formula [23], as follows: Therefore, the problem of calculating the proportional coefficient, m, to obtain the approximate optimal probability distribution is transformed into the problem of finding the optimal multiplier, γ. The calculation formula of the proportional coefficient, m, is as follows [24]: where n 1 is the number of faulty components in the sampling process and n 0 is the number of running components. Based on the iterative method, γ is calculated. Firstly, it is given an initial value, and then iterated through small-scale sampling until the adjacent γ reaches the specified relative error accuracy. The calculation formula is as follows [24]: where A, B, and C are the calculation coefficients and λ the average failure rate of the components.
(2) Sample the optimized probability distribution based on the LHS method. Stratify the input random variables and select the random sampling points in each hypercube as the sampling points, according to the importance shown in the original probability density function, usually the boundary of the hypercube closest to the expected value. The calculation formula of X kn is as follows (assuming that the sampling number, n, is even) [11]: The LHIS method comprehensively considers the two ideas of "stratification" and "importance" in sampling, which makes the sampled points closer to the expected value and accelerates the convergence of LHIS.
(3) Reduce the correlation of the sampling matrix based on the Cholesky decomposition method. If K components are sampled, and the sampling times of each component is N, then all of the sampled values will form an initial sampling matrix, X KN , of order K × N. The accuracy of the evaluation results is greatly affected by the correlation of the sampled values in the sampling matrix, therefore, the correlation between the lines of the sampling matrix should be reduced.
Firstly, construct an arrangement matrix, L KN , and the element value of each row will represent the arrangement position of the corresponding row elements in each row in the X KN .
Then, construct the row correlation coefficient matrix ρ L of the L KN . ρ L is a matrix of K × K order, in which the elements are the correlation coefficients between each row, which can be expressed as follows [25]: Obtain the lower triangular matrix, D, by the Cholesky decomposition method, and then construct a matrix with less correlation, G KN , by combining the D matrix, which can be expressed as follows [25]: Finally, rearrange the elements in the X KN according to the G KN from the largest to the smallest, and repeat the above steps until the L KN is less than the required value.

Sampling the Power Systems with Wind Farms and Storage Systems
Sample the distribution function (2) based on the LHIS method and obtain the n-th wind speed sampling value as follows: Considering the influence of the wake effect, construct the sampling matrix sequence of the wind speed by completing the set N sampling times, as follows [26]: Based on the failure rate and repair rate of the fan, obtain the distribution function of the fan state as follows: By sampling Equation (20) N times by LHIS method, obtain the sampling matrix sequence of the wind turbine state as follows: Then, obtain the initial sampling matrix, X KN , after sampling the power system with the wind farm and energy storage system N times, which can be expressed as follows: In the formula, X WTi , X Gg , X Tk , X ESp , and X Lm are, respectively, represented as sampling vector groups that sample wind turbines, conventional generators, transformers, lines, and energy storage systems N times, where i + g + k + m + p + 1 = K.
The correlation between the sampling values of each row of X KN is large, therefore, it is necessary to rearrange X KN according to L KN , as proposed in Section 3.2.1, to reduce the correlation of X KN .

Evaluation Indicators
To analyze the influence of wind storage systems on the power grid comprehensively, comprehensive evaluation indexes are needed in order to reflect the situation of the power grid. In this paper, a comprehensive risk index (R CRI ) [27,28] and the wind storage generation interrupted energy benefit (WSGIEB) are used as evaluation indicators.

Comprehensive Risk Indicator (R CRI )
Considering the evaluation results of the risk of load cutting, R CUT , the risk of line active power exceedance, R LI NE , and the risk of voltage exceedance, R V , the comprehensive risk index, R CRI , is proposed, which can be expressed as follows [29]: In the formula, ω 1 , ω 2 , and ω 3 are the weights of the above three risk indicators, which can be obtained by the analytic hierarchy process. The R CUT represents the risk of load cutting under certain operating conditions of the system. The lower the risk number, the more reliable the system. The R LI NE represents the risk of the line active power exceeding the maximum allowable active power under a certain operating condition. The R V represents the risk of voltage exceedance of the upper or lower limit under a certain operating condition. The calculation formula of R CUT , R LI NE , and R V is as follows [29]: where p wt (E k ), p g (E k ), p ling (E k ), and p es (E k ) are the real-time state probabilities of the wind turbine, conventional generator, transmission line, and energy storage system when the system state is E k . K is the number of simulations to run the risk assessment. S ev−CUT (E k ) is the severity of the consequences when load cutting occurs, S ev−LI NE (E k ) is the severity of the consequences when the line active power exceeds the limit, and S ev−V (E k ) is the severity of the consequences when the voltage exceeds the limit. The calculation formula is as follows [29]: where C c (E k ) is the unit value of load cutting; P load is the unit value of the system load; M is the number of lines; p line (m) is the unit value of the active power of the m-th line; p line,max (m) is the unit value of the maximum active power allowed to flow through the mth line; N is the number of power system nodes; and U i is the unit value of the node voltage.
Based on the principle of as low as reasonably practicable (ALARP) and the requirements of the State Grid Corporation's Safety Accident Investigation Regulations, this paper classifies the risk levels of power system operation. The ALARP standard divides the risk into three zones, namely, the unacceptable zone, the reasonably feasible minimum zone, and the negligible zone. The ALARP criterion consists of two risk boundaries; one is an unacceptable risk level and the other is a negligible-risk level [30,31]. In this paper, the power system is divided into three levels and contains two risk lines, as follows: high-risk level, medium-risk level, low-risk level, the unacceptable high-risk line and the negligible low-risk horizontal line.
The negligible low-risk horizontal line is defined as the product of the lowest level of failure probability and the lowest level of failure consequence [30,31]. The lowest level of failure probability is the average value of national line failure over the past years. The reliability data released by China Electricity Council over the years are shown in Table 1. From Table 1, it can be seen that the average forced outage rate, f, is 0.192 times/100 km·year, and the average MTTR is 2.73 h/100 km·year. According to the following formula, the lowest level of failure probability is 5.98 × 10 −5 : According to the safety accident investigation procedures of the State Grid Corporation, a reduction of 4% or less of the power supply load does not constitute a level 4 or above power grid event; therefore, the lowest level of failure consequence is defined as a 4% reduction of the power supply load. Therefore, the negligible low-risk horizontal line is: The unacceptable high-risk horizontal line is defined as the product of the average failure probability and the unacceptable consequence [30]. According to the safety accident investigation procedures of the State Grid Corporation, the unacceptable consequence Energies 2023, 16, 6427 9 of 18 is 30% load reduction. Therefore, the unacceptable high-risk horizontal line is: R max = 5.98 × 10 −5 × 30% = 1.79 × 10 −5 .
Based on the ALARP principle and the requirements of the State Grid Company's safety accident investigation procedures, the standards for the classification of power system operation risk are shown in Table 2.

Wind Storage Generation Interrupted Energy Benefit (WSGIEB)
This indicator is the ratio of the change degree of the E EENS to the sum of C w and C s . This indicator can directly reflect the contribution of the power system to the reliability. The expression is as follows: In the formula, C w and C s represent the wind farm capacity and the energy storage system capacity, respectivley. E EENS0 and E EENS1 represent the expected power shortage of the wind farm before and after the energy storage system is connected, respectivley. The calculation formula of is as follows [5]: where S represents the set of system states with load cutting conditions, P i represents the probability that the system is in state i, and C i represents the load cutting power when the system is in state i.

Evaluation Process
Based on LHIS method, the state of the power system is sampled to form a time series of each component state, the system state is extracted, and the probability power flow is calculated to obtain the evaluation index. A flow chart is shown in Figure 1 and the specific evaluation process is as follows: (1) Input the original data, calculation accuracy, and evaluation time of the power system. Sample the power system based on the LHIS method and obtain the initial sampling matrix X KN . total output of the conventional power plants, wind farms, and energy storage systems. (4) Calculate the probability power flow of the system, the load cutting, the line power exceeding the limit, the node voltage exceeding the limit, and the power shortage. Calculate the evaluation indicator and judge whether the calculation results converge [32]. If the results do not converge, increase the time by one and return to Step 3.

Simulation Analysis
This paper analyzes the reliability test system of IEEE-RTS79 based on MATLAB R2022a software. The system consists of 32 generator sets, including 24 buses, a strip camera, 33 transmission lines, and 5 transformers. The total installed capacity is 3405 MW, and the peak load of the system is 2850 MW.

Verification of Evaluation Methods
In this paper, the accuracy of the evaluation method is measured by the relative error between the expected value and the standard deviation of the results of the probabilistic power flow calculation, and the calculation method is as follows [30]: where µ h,acc is the expected value, σ h,acc is the standard deviation of the calculated results of the probabilistic power flow under a certain sampling scale, µ h,sim is the expected value, and σ h,sim is the standard deviation of the accurate value. In this paper, the result of simple random sampling (SRS) 10,000 times is adopted as the accurate value. h is the physical quantity of the result of the probabilistic power flow calculation and voltage amplitude, u, and line active power, p, are selected as the verification parameters. Connect a 350 MW wind farm and a 300 MW energy storage system to 15 nodes and use the Latin hypercube important sampling method (LHS) proposed by Huang et al. [7], the improving importance sampling method (IM-IS) proposed by Tómasson et al. [8], and the LHIS method proposed in this paper to perform the probabilistic power flow calculations on the power system. Then, calculate the relative errors of the node voltage, the branch active power expectation, and standard deviation, ε u,µ , ε u,σ , ε p,µ , ε p,σ , according to Equations (29) and (30). The calculation results are shown in Figure 2. When the relative error of ε u,µ and ε p,µ reaches 0.2%, and the accuracy of ε u,σ and ε p,σ reaches 3% of the three calculation methods, the calculation time is calculated, as shown in Table 3.      The active power of branch 6-8 is selected as the research object. The probability density function and cumulative distribution function of the calculated results of the LHIS method, the LHS method, and the SRS method are shown in Figure 3.  The active power of branch 6-8 is selected as the research object. The probability density function and cumulative distribution function of the calculated results of the LHIS method, the LHS method, and the SRS method are shown in Figure 3. The active power of branch 6-8 is selected as the research object. The probability density function and cumulative distribution function of the calculated results of the LHIS method, the LHS method, and the SRS method are shown in Figure 3. It can be seen from Figures 2 and 3 that, after the calculation results converge, the relative error of expectation of the LHS and IM-IS calculation results reaches 0.44%, and the relative error of the standard deviation of the LHS and IM-IS calculation results exceeds It can be seen from Figures 2 and 3 that, after the calculation results converge, the relative error of expectation of the LHS and IM-IS calculation results reaches 0.44%, and the relative error of the standard deviation of the LHS and IM-IS calculation results exceeds 4.05%. It can be seen from Table 3 that, when the relative error of the expected value and the relative error of the standard deviation of the calculation result reach 0.2 and 3, respectively, the calculation time of the LHS method and the IM-IS method should be more than 9 s. According to the analysis of Figures 2 and 3 and Table 3, it can be concluded that the LHS method and the IM-IS method still have the problem of poor accuracy of calculation results and a long calculation time.

LHIS (s) LHS (s) IM-IS (s)
However, after the calculation results converge, the relative error of expectation of the LHIS method is less than 0.2%, and the relative error of standard deviation of the LHIS method is less than 3%. The calculation result of the LHIS method 150 times is similar to that of the SRS method 10,000 times, while the LHS requires 400 times and the IM-IS requires 500 times. When the evaluation index reaches the allowable accuracy of error, the time that the LHIS method takes is about 5 s, which is nearly 45% shorter than that of the LHS method, and nearly 50% shorter than that of the IM-IS method. It can be concluded that the calculation efficiency of the LHIS method is 40% higher than that of the LHS method and 47% higher than that of the IM-IS method. The calculation accuracy of the LHIS method is 20% higher than that of the LHS method and 33% higher than that of the IM-IS method.
In order to analyze the adaptability of the evaluation methods under different permeability rates, Figure 4 and Table 4 show the error indicators of the three methods under different permeability conditions. The calculation time of the three methods is shown in Table 5. LHIS method is 20% higher than that of the LHS method and 33% higher than that of the IM-IS method. In order to analyze the adaptability of the evaluation methods under different permeability rates, Figure 4 and Table 4 show the error indicators of the three methods under different permeability conditions. The calculation time of the three methods is shown in Table 5.     . ε u,µ , ε u,σ , ε p,µ , and ε p,σ relative error curves under different permeability conditions.
As can be seen from Figure 4 and Table 4, the relative error of expectation of the LHS and IM-IS calculation results exceeds 0.3% with the increase in permeability, and even the expected relative error of the IM-IS method exceeds 0.4% when the permeability reaches 60%. The relative error of standard deviation of the LHS and IM-IS calculation results exceeds 4% with the increase in permeability. As can be seen from Table 5, when the permeability reaches 60%, the calculation time of LHS and IM-IS exceeds 10 s and 12 s, respectively. According to the analysis of Figure 4 and Tables 4 and 5, when the power system has high permeability, the LHS and IM-IS methods have obvious lag, the calculation accuracy is poor, and the calculation period is long.
However, with the increase in permeability, the relative error of expectation of the LHIS calculation results does not exceed 0.3%, and the relative error of standard deviation of the LHIS calculation results does not exceed 3%. However, the algorithm still maintains a high precision. The calculation time of the LHIS method also fluctuates in a small range with the increase in permeability, and the calculation time does not exceed 6 s when the permeability reaches 60%. This shows that the LHIS method is universal under the high permeability of the power system and ensures the high precision of the calculation results and a fast calculation time.

Connecting Wind Farms with Energy Storage Systems at Different Nodes
In order to fully evaluate the influence of wind farms with energy storage systems on the power system at different nodes, this paper chooses to connect a 350 MW wind farm and a 300 MW energy storage system at nodes 1, 2, 13, 15, 16, 18, 21, 22, and 23, and, at the same time, lets the energy storage systems run in charge-discharge strategy 1 and charge-discharge strategy 2, respectively, in order to calculate the R CRI and WSGIEB. The calculation results are shown in Figure 5. the power system at different nodes, this paper chooses to connect a 350 MW wind farm and a 300 MW energy storage system at nodes 1, 2, 13, 15, 16, 18, 21, 22, and 23, and, at the same time, lets the energy storage systems run in charge-discharge strategy 1 and charge-discharge strategy 2, respectively, in order to calculate the CRI R and WSGIEB The calculation results are shown in Figure 5. It can be concluded from Figure 5 that the R CRI of Node 1, Node 2, and Node 18 is very high after they are connected to wind farms with energy storage systems. Compared with Table 2, the R CRI exceeds the unacceptable high-risk line and the system reaches the high-risk level after connecting to the wind farms with energy storage systems at Node 2 and Node 18. At the same time, the WSGIEB is negative, therefore, the reliability of the power system is reduced after these three nodes are connected to wind farms with energy storage systems. However, the R CRI of Node 15 is very low and is closest to the negligible low-risk horizontal line, the system operation is in the negligible-risk level, and the WSGIEB is relatively high. Therefore, the reliability of the power system is high after the node is connected to the wind farm with energy storage. Based on the above situation, we can see that it is not suitable to connect wind power at Node 1, Node 2, or Node 18, and the reliability of connecting the wind farm at Node 15 is high.
The reliability of the energy storage systems under charge and discharge strategy 2 is higher than that under charge and discharge strategy 1. This is because, in charging and discharging strategy 2, when the output of the wind turbines cannot reach P L (t)·η, the gap will be filled by the output of the conventional power plants first, and the energy storage will be discharged if the sum of the output of the conventional power plants and the wind turbines cannot meet the current power. This discharge strategy can put the energy storage battery in a state of full capacity, or more capacity, for a long time, and, once the power supply is insufficient, there is more energy in the energy storage system to make up for the deficiency. Under charge and discharge strategy 1, the energy storage battery will be discharged for a long time, with relatively low capacity and poor reliability.

Connecting Wind Farms of Different Capacities with Energy Storage Systems
According to the analysis results shown in the previous section, the wind farms with energy storage systems are selected to be connected at Node 15, and the energy storage systems are operated under charging and discharging strategy 2, so that the power grid has a higher reliability. The capacity of the energy storage systems is kept constant at 350 MW, and wind power capacities of 100 MW, 200 MW, 300 MW, 400 MW, 500 MW, and 600 MW are connected, respectively. The calculation and evaluation results are shown in Figure 6. With the increase in the connected wind power capacity, the WSGIEB gradually increases, indicating that the increased connected wind power capacity contributes to the reliability of the power system in a certain range; however, when the connected wind power capacity exceeds 500 MW, the WSGIEB decreases, indicating that the newly added wind turbines contribute less and less to the reliability of the power system, resulting in the saturation phenomenon, therefore, further adding wind power capacity will cause a waste of resources.

Connecting Wind Farms with Different Energy Storage Capacities
The capacity of wind farms is kept unchanged at 450 MW, and the energy storag systems capacities of 100 MW, 200 MW, 300 MW, 400 MW, and 500 MW are connected to Node 15, respectively. The calculation results of the evaluation indicators are shown in Figure 8. When the capacity of 15-node-connected wind farms reaches 500 MW, the R CRI increases rapidly. The settlement results of indicators such as R CUT , R LI NE , and R V are further analyzed, as shown in Figure 7. When the capacity exceeds 500 MW, the R V and the R LI NE increase sharply, which exceeds the unacceptable high-risk line. Through the analysis of Figures 6 and 7, it can be concluded that, when the storage system capacity is 350 MW, it is more appropriate to connect 450 MW of wind power capacity at 15 nodes.

Connecting Wind Farms with Different Energy Storage Capacities
The capacity of wind farms is kept unchanged at 450 MW, and the energy storage systems capacities of 100 MW, 200 MW, 300 MW, 400 MW, and 500 MW are connected to Node 15, respectively. The calculation results of the evaluation indicators are shown in Figure 8.

Connecting Wind Farms with Different Energy Storage Capacities
The capacity of wind farms is kept unchanged at 450 MW, and the energy storage systems capacities of 100 MW, 200 MW, 300 MW, 400 MW, and 500 MW are connected to Node 15, respectively. The calculation results of the evaluation indicators are shown in Figure 8. With the increase in the energy storage system capacity, the risk of power system failure decreases. In a certain range, the energy storage system capacity increases and the reliability improves. However, when the capacity reaches 400 MW, the value of CRI R is basically unchanged, and the value of WSGIEB will be slightly reduced, mainly because, when the energy storage capacity is too large, the excess output of the wind turbines is completely stored in the energy storage system. Moreover, when the capacity of the energy storage system is increased, the capacity response of the energy storage system reaches saturation, which will not significantly improve the system reliability. When the output of the wind power capacity is kept at 450 MW, the energy storage system capacity of 400 MW is selected in order to optimize the system economy and reliability.

Conclusions
In this paper, the LHIS reliability evaluation method was proposed by combining the LHS method and the IS method. Then, the reliability of the power system with wind farms and energy storage systems was evaluated from two aspects of the power system risk and the contribution coefficient to the power system reliability. Through the analysis of the IEEE-RTS79 example, the following conclusions have been obtained: (1) Through comparative experiments with the LHIS method, the LHS method, and the IM-IS method, it can be concluded that the calculation efficiency of the LHIS method is 40% higher than that of the LHS method and 47% higher than that of the IM-IS method. The calculation accuracy of the LHIS method is 20% higher than that of the LHS method and 33% higher than that of the IM-IS method. (2) The LHIS method is universal to the reliability analysis of power systems under different permeability conditions. When the permeability reaches 60%, the relative error of expectation does not exceed 0.3%, and the relative error of standard deviation does not exceed 3%. Simultaneously, the calculation time does not exceed 6 s. The LHIS method ensures the high precision of the results and a fast calculation time under high permeability conditions. (3) This paper proposed a comprehensive risk indicator to evaluate the operational risk and the wind storage generation interrupted energy benefit (WSGIEB) to evaluate the contribution of the reliability, which solved the problem of a lack of risk analysis of operation state in traditional power system reliability evaluation. Combined with the operation risk division table, the operation of the power system in this state can be more intuitively obtained. With the increase in the energy storage system capacity, the risk of power system failure decreases. In a certain range, the energy storage system capacity increases and the reliability improves. However, when the capacity reaches 400 MW, the value of R CRI is basically unchanged, and the value of WSGIEB will be slightly reduced, mainly because, when the energy storage capacity is too large, the excess output of the wind turbines is completely stored in the energy storage system. Moreover, when the capacity of the energy storage system is increased, the capacity response of the energy storage system reaches saturation, which will not significantly improve the system reliability. When the output of the wind power capacity is kept at 450 MW, the energy storage system capacity of 400 MW is selected in order to optimize the system economy and reliability.

Conclusions
In this paper, the LHIS reliability evaluation method was proposed by combining the LHS method and the IS method. Then, the reliability of the power system with wind farms and energy storage systems was evaluated from two aspects of the power system risk and the contribution coefficient to the power system reliability. Through the analysis of the IEEE-RTS79 example, the following conclusions have been obtained: (1) Through comparative experiments with the LHIS method, the LHS method, and the IM-IS method, it can be concluded that the calculation efficiency of the LHIS method is 40% higher than that of the LHS method and 47% higher than that of the IM-IS method. The calculation accuracy of the LHIS method is 20% higher than that of the LHS method and 33% higher than that of the IM-IS method. (2) The LHIS method is universal to the reliability analysis of power systems under different permeability conditions. When the permeability reaches 60%, the relative error of expectation does not exceed 0.3%, and the relative error of standard deviation does not exceed 3%. Simultaneously, the calculation time does not exceed 6 s. The LHIS method ensures the high precision of the results and a fast calculation time under high permeability conditions. (3) This paper proposed a comprehensive risk indicator to evaluate the operational risk and the wind storage generation interrupted energy benefit (WSGIEB) to evaluate the contribution of the reliability, which solved the problem of a lack of risk analysis of operation state in traditional power system reliability evaluation. Combined with the operation risk division table, the operation of the power system in this state can be more intuitively obtained.
(4) By analyzing different wind farms with different energy storage system connection strategies, it can be confirmed that the LHIS method can accurately evaluate the reliability of power systems and provide useful references for power system design. Funding: This work was supported by the scientific and technological project of North China Electric Power Research Institute Co., Ltd. Technology for Improving the Operation Reliability of Wind Turbines Considering the Influence of Regional Extreme Weather (project number Z2022KJ041).

Data Availability Statement:
The data are available from the corresponding author upon reasonable requests.

Conflicts of Interest:
The authors declare no conflict of interest.