A Stochastic Diffusion Process-based Multi-Level Monte Carlo for Predictive Reliability Assessment of Distribution System

Reliability assessment of electrical distribution systems is an important criterion to determine system performance in terms of i nterruptions. Probabilistic assessment methods are usually used in reliability analysis to deal with uncertainties. These techniques require a longer execution time in order to account for uncertainty. Multi-Level Monte Carlo (MLMC) is an advanced Monte Ca rlo Simulation (MCS) approach to improve accuracy and reduce the execution time. This paper provides a systematic approach to model the static and dynamic uncertainties of Time to Failure (TTF) and Time to Repair (TTR) of power distribution components using a Stochastic Diffusion Process. Further, the Stochastic Diffusion Process is integrated into MLMC to estimate the impacts of uncertainties on reliability indices. The Euler Maruyama path discretization applied to evaluate the solution of the Stochastic Diffusion Process. The proposed Stochastic Diffusion Process -based MLMC method is integrated into a systematic failure identification technique to evaluate the distribution system reliability. The proposed method is validated with analytical and Sequential MCS methods for IEEE Roy Bill inton Test Systems. Finally, the numerical results show the accuracy and fast convergence rates to handle uncertainties compared to Sequential MCS method. Author


Introduction
Power distribution system reliability can usually be evaluated through analytical or s imulation methods (Conti and Rizzo 2019;Brown 2008). Analytical methods can determine the average values of reliability indices, but they cannot handle uncertainties due to certain assumptions used in mathematical calculations. Meanwhile, simulation methods have the advantage of being able to handle uncertainties and are generally not restricted by the complexity of the power distribution network (IEEE Guide for Electric Power Distribution Reliability Indices 2004). There are two types of MCS methods: Non-sequential and Sequential MCS (Conti and Rizzo 2015). Non-Sequential MCS, also known as Random MCS, is a commonly used approach in the literature to assess the reliability of generation, transmission, and distribution systems. Further, Variance Reduction Techniques (VRT) has enhanced the MCS method convergence rate (Billinton and Jonnavithula 1996). With VRT methods, the expected value of a stochastic variable can be estimated for a specified accuracy level, increasing the computation speed.
Other methods for improving MCS have been explored, including the dagger sampling technique (Sun et al. 2010), accelerated state evaluation approach (Shu et al. 2014), and discrete event simulation (Aldhubaib and Kashef 2020). Among a different type of VRTs applied to power system reliability analysis, the beneficial methods are Quasi-Monte Carlo (Hou, Wang, and Guo 2017), Importance Sampling (Urgun 2019), Control Variates (Wang et al. 2018), Stratified Sampling (Wang et al. 2011, and Antithetic Variables (Benidris and Mitra 2014). A further meaningful way of variance reduction is by using Latin Hypercube Sampling (Jirutitijaroen and Singh 2008) and integrating these VRTs into MCS to enhance the process of estimating the reliability metrics. The MCS method's computation efficiency has been improved in recent literature using MLMC approach (Huda 2018;Tindemans and Strbac 2020). Notable literature is available on the application of the MLMC method for distribution system reliability analysis. The MLMC method constructs the TTF and TTR using Stochastic Differential Equation (SDE) (Huda 2018). The SDE of a stochastic variable is usually driven by its drift and diffusion coefficients. The reliability assessment of a simple distribution system has been performed using the MLMC method (Huda and Živanović 2017a). The Milstein discretization method is applied to evaluate stochastic differential equations (Huda 2018), and reliability indices are evaluated. The impact of two different discretization methods, i.e., Milstein and Euler-Maruyama, are examined on convergence criteria of MLMC in reliability assessment (Huda and Živanović 2017b). Furthermore, the impact of different power system components availability on reliability performance is estimated through the MLMC method (Huda and Živanović 2019b). The effect of time-varying failure rates, repair times, and load models on reliability indices; Expected Energy Not Supplied (EENS) (Huda and Živanović 2019a) and Expected Interruption Cost (ECOST) (Huda and Živanović 2019c) are estimated through MLMC method. In addition, the impact of adverse weather events like wind and lightning on reliability indices are assessed through the MLMC method (Huda and Živanović 2019d). From the above literature on the MLMC method, it is observed that the diffusion and drift coefficients of stochastic differential equations are assumed as constants tuned with reference to the analytical techniques for a specific level of accuracy. However, the coefficients of SDE indicate the nature and characteristics of the stochastic process. The above-specific assumption is not valid when we want to model the different uncertainties like aleatory and epistemic uncertainties (Teh 2018; Abunima and Teh 2020) associated with the failure rates and repair rates of power system elements. This paper proposes a new model that improves the MLMC method further to develop drift and diffusion coefficient expressions considering the various uncertainties associated with TTF and TTR.
The key contributions of the paper are listed below: a) Derivation of drift and diffusion coefficients of stochastic differential equations for the corresponding probability distribution function of stochastic variables (TTF and TTR) using the Stochastic Diffusion Process and Fokker-Planck equation. b) Integration of Multi-Level Monte Carlo simulation into a generalized approach for a computationally efficient and accurate evaluation of system reliability indices. c) Modeling of the stochastic variables (TTF and TTR) in two periods of the bathtub curve using the proposed method and case studies are performed on IEEE RBTS test systems. d) The efficiency of the improved method in handling various uncertainties over the other techniques are confirmed with the results. This paper is organized as follows. Section 2 presents the overview of MLMC simulation, Stochastic Diffusion Process, and generalized systematic procedure to evaluate the effect of branch failures on load points. Section 3 presents the proposed algorithm of Stochastic Diffusion Process-based MLMC for predictive reliability assessment. Section 4 presents the comparison of results between analytical, Sequential MCS (SMCS), and proposed method for IEEE RBTS test systems. Section 5 presents the conclusions and future scope of the research work.

MLMC Theory
The outcome of the MLMC method is to estimate the expected value of a random variable (reliability index) denoted by E [R]. In this section, the coming formulation follows the approach presented in Giles (2015). Let E[Ra] be an approximate function of E[R], where R is the reliability index to be estimated and Ra is a random variable that approximates R. The unbiased estimator of E[Ra] in the SMCS method is given by Giles (2015): (1) where ( ) is the i th sample of Ra and Ns is the samples on a fine level. The Root Mean Square Error (RMSE) is equal to 1/√ and the variance is equal to −1 V[Ra].
Mathematically, the concept of MLMC for different levels can be expressed by Equation (2) shown as follows (Giles 2015): Hence the outcome of the MLMC method for each performance index of the power distribution system is specified by Equation (3) (Giles 2015): where ℎ is the fine level timestep size and ℎ is coarse level timestep size.
The simulation convergence criteria in the MLMC method depends on the target Root Mean Square Error (RMSE) value specified by Equation (5) (Giles 2015) is: To meet the convergence of MSE ≤ 2 , it is needed to guarantee that both the variance and square root errors are distinctly less than 2 /2. The ideal number of samples is evaluated using Equation (6) (Giles 2015): where is the sample cost of level l. In the case of weak convergence error, the convergence criteria are specified as:

Stochastic diffusion process
In this section, the coming formulation follows the approach presented in Zárate-Miñano and Milano (2016). The stochastic differential equation usually defines the Stochastic Diffusion Process (SDP). A general form of single-dimensional SDE is as shown in Equation (8) where is the stochastic process (TTF and TTR) and ( ) is a Weiner process, also known as Brownian motion. The numerical solution of Equation (8) is known as a stochastic diffusion process, and the terms ( , ) and ( , ) are called drift and diffusion coefficients of the SDE. The Fokker-Planck equation generally describes the time-dependent expansion of the probability distribution function of a stochastic process. For the stochastic process defined by the SDE shown in Equation (8), the corresponding Fokker-Planck equation is shown in Equation (9) where ( , ) is the transient probability distribution function of the stochastic process . In general, the failure rates of the power system components are expressed by a two-state Markov model. Hence the Fokker-Planck Equation for a Markov process is given by: Therefore, if any of the coefficients ( , ) or ( , ) of a stochastic process is known, the other coefficient can be derived using Equation (10) for a specified probability distribution function ( , ).

MLMC method for Distribution System Reliability
In this section, the coming formulation follows the approach presented in Huda and Živanović (2019a). In the distribution system reliability evaluation, a two-state Markov model of a power system component is used. The two states are represented as up state and down state. The average duration during which a component is in the up state is called TTF, and the average duration during which a component is in the down state is called TTR, respectively (Huda and Živanović 2019a).
where represents failure rate of k th branch µ represents repair rate of k th branch In conventional SMCS, the up time (Tup) and down time (Tdn) of a component k are usually sampled using failure and repair rates of components using Equation (13) and Equation (14) (Huda and Živanović 2019a) is: where U1 and U2 are uniform random numbers sampled using uniform probability distribution between [0, 1]. In conventional MLMC method, the up time (Tup) and down time (Tdn) of a component k are usually sampled using failure and repair rates of components using Equation (15) (Huda 2018) is: where is evaluated using SDP represented by Equation (8) as shown below (Huda 2018): where and are the drift and diffusion coefficients assumed as constants (Huda 2018) for modeling of TTF and TTR tuned with reference to the analytical method for a specific level of accuracy. However, the coefficients of SDP indicate the nature and characteristics of the stochastic process. The above-specific assumption is not valid when we want to model the different uncertainties associated with the failure rates and repair rates of power system elements. As the main contribution, this paper develops the expressions for drift and diffusion coefficients considering the various uncertainties.

Methodical approach for cases identification
A generalized methodical approach (Conti, Nicolosi, and Rizzo 2012) is integrated into the MLMC method to evaluate the impact of branch failures on load point indices. The effect of branch failure on customers is categorized into different failure cases depending on load point location, failure branch location, protection device type, and installation location. Table 1

Distribution system reliability indices
The reliability performance of an electrical distribution system is examined thanks to two important metrics that are the System Average Interruption Frequency Index (SAIFI), and System Average Interruption Duration Index (SAIDI) (IEEE Guide for Electric Power Distribution Reliability Indices 2004). These indices are calculated from the average load point outage frequency ( ) and outage duration ( ) of the distribution system and are expressed as follows: where represents the number of load points (LPs), , represents the number of customers at LP-i, represents the LP-i outage frequency and represents the LP-i outage duration.

Stochastic Diffusion Process Based MLMC Method
In the MLMC method, the stochastic uncertainty of TTF of component k is modeled using an exponential distribution. The Weiner process drives SDE, Wt on the time interval [0, T]. Then, the SDE based TTF model with specified values of drift and diffusion coefficients and an initial TTF [ 0 = 1⁄ ] is expressed using Equation (19): where ( ) is the TTF value at time t. The numerical solution of SDE-based TTF model in Equation (19) is determined by applying Euler Maruyama discretization, as shown in Equation (20) Then, the SDE based models of TTF and TTR are represented using Equation (21) and Equation (22), as shown below: = − µ ( +1 ) * ( 2) As part of the main objective of this paper, the drift and diffusion coefficients of stochastic variables representing TTF and TTR are derived using the Stochastic Diffusion Process and Fokker-Planck Equation explained as follows.

Estimation of drift and diffusion coefficients of TTF using stochastic weibull diffusion process
The failure rates of components are usually represented by the Weibull distribution function in three periods of the bathtub curve, i.e., the infant mortality period, the useful life period, and the aging period. The probability distribution function of Weibull distribution is given by Equation (23) (Huda 2018) is: where is the scale parameter and is the shape parameter. The shape parameter for three periods of the bathtub curve is: The Stochastic Diffusion Process can be generated in different ways by considering the probability distribution function of a random variable. A Stochastic Diffusion Process representing the Weibull distribution has been developed to model the uncertainty of a random variable. The drift coefficient is considered as proportional to the probability distribution function of Weibull distribution (Nafidi et al. 2019). Therefore, the drift coefficient representing the random variable, TTF is expressed as Equation (24) shown below: The expression for the diffusion coefficient ( , ) is determined by substituting Equation (24) into Equation (10) for Weibull distribution. Equation (25) Thus, the expression to sample TTF in the proposed method is determined by substituting Equation (27) into Equation (21) and is given by Equation (28)

Estimation of drift and diffusion coefficients of TTR using stochastic exponential diffusion process
The Exponential distribution function is usually used to represent the failure rates of components in the useful life period of the bathtub curve and repair times of the power system components. The probability density function of Exponential distribution is given by Equation (29) The drift coefficient is considered proportional to the probability distribution function of Exponential distribution. Therefore, the drift coefficient representing the random variable TTR is expressed as Equation (30): The expression for the diffusion coefficient ( , ) is determined by substituting Equation (30) into Equation (10) for exponential distribution. Equation (31) represents the derived expression for the diffusion coefficient ( , ): Hence, the Stochastic Exponential Diffusion Process (SEDP) for TTR modeling is given by: The numerical solution of Equation (32) is determined using an Euler Maruyama discretization scheme represented by Equation (20) is: Thus, the expression to sample TTR in MLMC is determined by substituting Equation (33) into Equation (22) and is given by Equation (34):

Integrated stochastic diffusion process based MLMC and generalized methodical approach
The algorithm and flowchart of the integrated SDP based MLMC (SDP_MLMC) and generalized methodical approach for predictive reliability assessment is summarized in this section:

Algorithm:
1. Initialization of the component reliability data, network data and MLMC simulation parameters: Number of samples (N) Maximum levels (Lm) Desired accuracy (ε) and Initial samples (Ns) 2. If L ≤ Lm, then set l = 0 and go to step (3) or else set Lm = Lm + 1 and repeat step (2). 3. If l ≤ L, then set n = 1, N = Ns and go to step (4) or else go to step (10). 4. If n < N, then set k = 1 and go to step (5) or else go to step (9). 5. If k ≤ NB, then set Tt = 0 and go to step (6) or else calculate , and reliability index and set n = n + 1 and go to step (4). 6. If Tt < T, then go to step (7) or else set k = k + 1 and go to step (5). 7. Set i = 1 and for branch k, evaluate TTF using SWDP given by Equation (28) and TTR using SEDP given by Equation ( 8. If i ≤ NLP, then determine the failure cases from Table 1 and evaluate , , , and set i = i + 1 and repeat step (8). Else update Tt = Tt + TTFk + TTRk and go to step (6). 9. If RMSE ≤ √2 ⁄ , then set l = l + 1 and go to step (3) or else set N = Nopt calculated using (6) and go to step (4). 10. If WE ≤ √2 ⁄ , then stop the algorithm or else L = L + 1 and go to step (2).

1) Flowchart
The flowchart of the integrated SDP based MLMC and generalized methodical approach is presented in Figure 1.

Case Studies and Results
The SDP based MLMC method is validated using three distribution test systems (Wang 1995): IEEE RBTS Bus-2, IEEE RBTS Bus-5, and IEEE RBTS Bus-6. The proposed SDP_MLMC method is implemented in Visual C++, and description of the test systems along with numerical results are discussed in this section:

Test systems
A complex 11 kV IEEE RBTS Bus-2 (Wang 1995) distribution system is considered to validate the proposed method, as shown in Figure 2.
The test system consists of 4 main feeders, 56 nodes, 36 sections, 20 transformers, and 22 load points. The reliability data, line lengths, average load data, and consumers' data are in (Wang 1995). Figure 3 presents the IEEE RBTS Bus-5 (Wang 1995) urban distribution system. The test system consists of 4 feeders, 69 nodes, 43 lines, 26 transformers, and 26 load points. Figure 4 presents the IEEE RBTS Bus-6 (Wang 1995) rural distribution system. The test system consists of 4 feeders, 102 nodes, 64 lines, 38 transformers, and 40 load points.

Results and discussions
This section discusses the proposed method applied for two scenarios, i.e., useful life and aging periods of failure rate model. The target RMSE values (ε) of 0.01% and 0.1% are considered as stopping criteria for estimation of SAIFI and SAIDI respectively. Table 2 presents the SAIFI values evaluated using analytical, SMCS, MLMC, and SDP_MLMC methods for the useful life period. The objective of the comparison in Table 2 is to assess the performance of the proposed method. As it is observed from the structure of the methods proposed, the SAIFI values evaluated with the analytical method does not consider the uncertainties whereas the SMCS method only considers randomness and constant failure rates and repair times for sampling TTF and TTR using exponential distribution. The conventional MLMC method modeled the TTF and TTR using SDEs with constant drift and diffusion coefficients tuned with reference to the analytical method for all three test systems. In general, the drift and diffusion coefficients represent the nature of a stochastic process and considering them constant does not take into account the uncertainty of failure rates and repair times. This limitation has been addressed with proposed SDP_MLMC method. The proposed method models the continuous and discrete uncertainties affecting the failure rates and repair times effectively. This method samples TTF and TTR using exponential probability distribution, SDP and Fokker-Planck equation, and calculates the drift and diffusion coefficients. The results of the proposed method are validated with the published results of analytical, SMCS and MLMC methods for the three test systems. It is concluded that the SAIFI values evaluated with the proposed SDP_MLMC method and outcomes are near to the published methods.     Table 4 presents the comparison of SAIFI and SAIDI values calculated using SMCS and SDP_MLMC methods for the aging period for the three test systems. It is evident that in the aging period, the failure rates are increased exponentially and are modeled with Weibul l distribution in SMCS method. Since the MLMC method is considering constant drift and diffusion terms and hence, it is limited to the useful life period only and not applicable for aging period. However, the advantage of the MLMC method is that it results in the improvement of execution time (Huda 2018) due to reduction in number of simulations for convergence compared to SMCS method. The proposed SDP_MLMC holds good for both aging period and useful life period as drift and diffusion coefficients are derived by representing TTF as Weibull distribution and exponential distribution for TTR. With respect to the computation performance of the proposed SDP_MLMC method, the number of simulations for convergence of SAIFI value is 18152, 111049, and 135850 for B2, B5, and B6 respectively for SDP_MLMC approach compared to 10 6 simulations for the SMCS method in the aging period. The main advantage of the proposed SDP_MLMC method is that it holds good for handling aging period same as SMCS method and the faster convergence rate as of MLMC method. From Table 4, it is concluded that the SAIFI and SAIDI values estimated for aging period with the proposed SDP_MLMC method are comparable with SMCS with less number of iterations and provides the solution where conventional MLMC fails.    From the analysis, it is concluded that the proposed method after including the uncertainty of TTF and TTR results in the indices values that are in acceptable accuracy levels with faster convergence rates. Further, it is observed that for two different scenarios, there is an increase in the value of SAIFI and SAIDI in the aging period compared to the useful life period. Results demonstrate the effectiveness of the proposed method to handle uncertainty in all periods of the bathtub curve compared to the MLMC method, which is limited to the useful life period.

Conclusions and Future Scope
This research illustrates the Stochastic Diffusion Process application for modeling uncertainties of TTF and TTR in the MLMC method for predictive reliability assessment. The SDP coefficients representing the stochastic variables TTF and TTR of power system elements are determined using SWDP and SEDP. The SDP-based MLMC method has been integrated into a generalized failure detection procedure to estimate the effect of branch failures on load points. The efficacy of the proposed method in terms of handling uncertainty and computation speed is validated by IEEE RBTS Bus-2, Bus-5, and Bus-6 test systems. Two different case studies have been performed to investigate the TTF and TTR modeling in useful life and aging periods of components. Numerical results demonstrate that the proposed method can model any uncertainty associated with TTF and TTR. The conventional MLMC method limitation of considering the constant drift and diffusion terms with reference to analytical method is overcome with the proposed SDP_MLMC approach. Also, the proposed method is within acceptable levels of accuracy with the improved convergence rate. Further, the proposed SDP_MLMC method is better suited to model various uncertainties associated with the DERs and perform the reliability analysis more effectively with faster execution times. As part of the future scope, the Stochastic Diffusion Process can be used to incorporate discrete random events and model the various uncertainties of Distributed Energy Resources (DER) in reliability assessment with less computational burden.