In recent years, China has provided policy support for the development and utilization of clean energy, which has promoted the development of renewable energy power generation (Zhang et al., 2022a; Zhang et al., 2022b; Song et al., 2022). The access to distributed photovoltaic energy in the distribution network brings the operation of the power system more uncertainty and puts forward higher requirements for distribution network planning (Zhang et al., 2019a). In order to optimize the operation of the distribution network, (Xiao et al. 2021) proposed a comprehensive control method for negative sequence current suppression and reactive power compensation. The optimal placement and sizing of distribution static compensator in radial distribution networks as multi-objective optimization with the objective of power loss reduction using whale optimization algorithm is discussed by (Noori et al. 2021). However, the compensation effects are rarely considered when evaluating distribution network planning options. In addition, the distribution network planning evaluation system based on a deterministic model and deterministic power flow will be difficult to be effectively applied to an active power distribution network (Meera and Hemamalini, 2021).

At present, a relatively standardized comprehensive evaluation system for distribution network planning (Song and Zhu, 2021; Dong et al., 2022) has been formed, covering flexibility, coordination, economy, and power supply reliability. In order to consider the characteristics of renewable energy, (Zeng et al. 2016) considered the environmental protection effects brought by a large number of distributed generators connected to the distribution network, and the index evaluation system of the new energy distribution network considering environmental protection and the comprehensive evaluation model of environmental benefits of intelligent distribution network were established, respectively. Based on the big data analysis method, (Visser et al. 2022) realized the comprehensive evaluation of the distribution network including “performance” and “benefit.” On this basis, (Li et al. 2018) constructed a four-tier index system including user energy efficiency and environmental pollution.

However, the existing research has not considered the volatility and intermittence brought by photovoltaic access to the power grid. Different from the traditional distribution network, for the fluctuation of the output power of distributed photovoltaic power stations and the correlation of the output of multiple photovoltaic electric fields (Murata et al., 2009), it is necessary to consider the power flow limit and voltage limit of the system caused by the aforementioned characteristics in the evaluation of distribution network (Fermín et al., 2022).

At present, the description of the fluctuation of photovoltaic output is usually based on the probability density function method (Liu et al., 2016; Zhang et al., 2019b; Chen et al., 2020). (Yang et al. 2017) and (Xiang et al. 2019) show the robust optimization model that effectively reflects the probability distribution characteristics of new energy. (Wang et al. 2022) proposed a spectrum analysis method based on an autoregressive model to analyze the fluctuation characteristics of photovoltaic output. (Yu et al. 2019) showed that the random characteristics of the actual output of photovoltaic power station approximately obey the beta distribution. (Wu et al. 2015) used the exponential distribution to simulate the variation characteristics of photovoltaic power and evaluate the short-term power fluctuation amplitude.

In general, the Pearson correlation coefficient method (Jia et al., 2021) and the copular function fitting method (Rajabalizadeh and Tafreshi, 2020) are extensively used to describe the output correlation of photovoltaic power stations. The Pearson correlation coefficient method can only characterize simple linear correlation (Singh et al., 2010). Compared with that, copular function has higher flexibility, but it is more complex to find the optimal parameters of copular function with actual data (Jin et al., 2021).

Existing methods can be categorized into two groups: 1) the refinement modeling of random characteristics for the output power of photovoltaic power station; 2) the refinement modeling of correlation characteristics for the output power of photovoltaic power station. However, it is hard to simultaneously depict them. Gaussian mixture model can effectively simulate the random characteristics of photovoltaic output and can describe the correlation of output power of multiple photovoltaic power stations based on their joint distribution (Jiang et al., 2015). In addition, the voltage and power flow limits will endanger the safe operation of the power grid. Therefore, when establishing the distribution network planning evaluation system of distributed photovoltaic access, it is necessary to consider the over-limit risk of voltage and power flow at the same time.

To address these important issues, this study proposes a comprehensive evaluation index system for distribution networks for distributed photovoltaic access. Relative to the state-of-the-art, the contributions of this study are threefold:

1) In the modeling of the randomness and correlation of output power of multiple photovoltaic power stations, the GMM based on the improved particle swarm optimization algorithm is proposed.

The first section of the article is the introduction. Section 2 provides the uncertainty model considering the randomness and correlation of photovoltaic power stations; Section 3 discusses the analytical probabilistic power flow analysis method based on a linearized power flow equation. In Section 4, the evaluation system of distribution network planning considering power flow uncertainty are shown. Simulation analysis is presented in Section 5. The conclusion and limitations are provided at the end of this article.

Gaussian mixture model (GMM) can accurately model non-Gaussian random variables and can approximate any probability density distribution with the help of a linear combination of a certain amount of Gaussian density function (Angelim and Affonso, 2020), shown as follows: f X ˜ ( X ) = ∑ m = 1 M ω m N m ( X | μ m , Σ m ) , ∑ m=1 M ω m = 1 , (1) N m ( X | μ m , Σ m ) = exp ( − 1 2 ( X − Σ m ) T Σ i − 1 ( X − μ m ) ) ( 2 π ) W / 2 det ( Σ m ) W / 2 , (2) where X is the input variable of the expression and represents the output power vector of the photovoltaic power station; ω m represents the weight coefficient of the mth Gaussian distribution sub-component; μ m and Σ m represent the expectation vector and covariance matrix of the mth Gaussian sub-component, respectively, in which the non-diagonal elements of the covariance matrix are used to describe the correlation between different photovoltaic power stations. The parameters of GMM that needs to be further estimated are ω m , μ m , and Σ m .

In parameter estimation, the samples X = { X 1 , ⋯ , X N } based on the actual collected output active power of the photovoltaic power station are formed first, and the likelihood function is built then, as shown in Eqs. 3, 4.the following equations: L ( X 1 , ⋯ X N ) = ∏ j = 1 N f ( X j ) , (3) ln L ( X 1 , ⋯ X N ) = ∑ j =1 N ln ∑ m = 1 M ω m N m ( X j ) . (4)

Then, iterative optimization is carried out with the maximum likelihood function as the optimization objective, and the optimization result parameters are the parameters that GMM needs to be estimated.

In the traditional GMM parameter solving process, the expectation maximization (EM) algorithm is extensively used. Considering that the accuracy of the algorithm depends on the selection of the initial value, the improper initial value may lead to the local optimal solution during the iterative optimization, which will have a negative impact on the next settlement results. In contrast, the improved optimal particle swarm optimization algorithm can greatly reduce the influence of initial value and realize the balance of global optimization and local optimization ability, so it is easier to obtain the global optimal solution (Peng et al., 2017). Therefore, this section adopts the improved particle swarm optimization algorithm.

1) First, the particle swarm x = [ ω m , μ m , Σ m ] to be determined is formed to initialize the particle position and velocity.

The structural configuration of the traditional distribution network is that there is only one generator as the swing bus, and the rest are all PQ nodes, excluding PV nodes. The research on the linearized power flow equation of the distribution network has been relatively mature. However, since the distribution network is connected to the distributed photovoltaic power station, sometimes it needs to be processed as a PV node in the power flow calculation, and the ZIP load is getting more and more attention in the analysis of the distribution network. Thus, it is necessary to establish a more applicable power flow calculation method for the distribution network.

In this section, the DG access node is regarded as a PV node, the load is set as the ZIP model, the node injection power equation is used as the original power flow equation, the power flow equation is linearized, and the linear power flow equation between node voltage amplitude and line power and injection power is derived. The final linearized power flow equation is as follows: V = A 1 X + B 1 , S = A 2 X + B 2 , (8) where V and S are the output variables and represent the node voltage amplitude and branch power flow, respectively. X is the input variable and represents the photovoltaic active output. A ₁ , A ₂ , B ₁ , and B ₂ are linearization expression coefficients. In order to simplify the expression, all linearized expressions are uniformly expressed as shown in Eq. 9 the following equation: W = B X + C , (9) where W is the output variable, representing the node voltage amplitude and branch power; X is the input variable and represents the active output of new energy. B and C are the coefficients of linearization expression.

If the random variable obeys Gaussian distribution, it still obeys Gaussian distribution after linear transformation. Therefore, the output power X of a photovoltaic power station is modeled based on multivariable Gaussian distribution Nm (x) first, and after linear transformation in Eq. 9, the node voltage and line power flow also obey Gaussian distribution. Here, the node voltage is selected for derivation. The amplitude of the node voltage follows Gaussian distribution, in which the expectation vector is Aμ _m + B _1, the covariance matrix is AΣ _m A ^T , so the joint probability density function (PDF) of node voltage can be expressed as N m ( W| μ m , Σ m ) = 1 ( 2 π ) K / 2 det ( B Σ m B T ) 1 / 2 × e ( − 1 2 ( W - B μ m - C ) T ∑ i − 1 W -B μ m − C ) , (10) where K represents the number of nodes. By multiple integrations of Eq. 10, the joint cumulative distribution function (CDF) of node voltage V can be obtained, as shown in the following equation: Eq. 11, F m ( W ) = ∫ ⋯ ∫ − ∞ W N m ( W | μ m , Σ m ) d W 1 ⋯ d W k . (11)

When photovoltaic power X is modeled by GMM with multiple weights and variables, the joint CDF of node voltage can be expressed as F ( W ) = ∫ ⋯ ∫ − ∞ W ∑ m = 1 M ω m N m ( W | μ m , Σ m ) d W 1 ⋯ d W K , = ∑ m = 1 M ω m [ ∫ ⋯ ∫ − ∞ W ∑ m = 1 M N m ( W | μ m , Σ m ) d W 1 ⋯ d W K ] , = ∑ m = 1 M ω m F m ( W ) . (12)

Accordingly, the joint probability density function of node voltage V can be obtained by differentiating Eq. 12. f ( W ) = ∂ K ∂ W 1 ⋯ ∂ W K ∑ m = 1 M ω m F m ( W ) , = ∑ m = 1 M ω m ∂ K ∂ W 1 ⋯ ∂ W K F m ( W ) , = ∑ m = 1 M ω m N m ( W ) . (13)

Eqs 12, 13 are the node voltage joint CDF and joint PDF when the photovoltaic output power is modeled by GMM.

Under the influence of light intensity, the output power of the photovoltaic power station meets a certain probability distribution, which will lead to fluctuations in the voltage of key nodes in the power grid and the transmission power of key power flow section. Distributed photovoltaic grid connection may also cause the bus voltage of the power grid to deviate to varying degrees at different times. The degree of system voltage offset is an important factor in evaluating voltage quality. Excessive voltage offset may lead to “voltage collapse” and eventually large-scale power outage (Zhang et al., 2021). In addition, because the random fluctuation of line power flow is not considered, the distribution network planning evaluation based on traditional power flow will not accurately measure the risk. Therefore, in order to form a distribution network planning evaluation system that is suitable for distributed photovoltaic access, the aforementioned factors are added here, and the following evaluation indexes are defined:

1) The risk index of voltage over-limit Risk (V). In order to measure the voltage quality of the distribution network at the key bus under the condition of random fluctuating power output, based on the voltage joint cumulative distribution function obtained in Section 3, the voltage out-of-limit risk index is defined as follows:

The voltage out-of-limit risk index, voltage deviation index, and line overload risk index are incorporated into the existing comprehensive evaluation index system of distribution network planning, that is, considering reliability, economy, safety, and adaptability at the same time, a three-tier distribution network index evaluation system considering distributed photovoltaic access is formed, as shown in Figure 1.

In addition to the three indicators defined in Section 4.1, the definitions of relevant traditional indicators in the figure are shown in the study by (Guo et al. 2021), and the weight of indicators shown in Figure 1 is calculated by the analytic hierarchy process.

1) Establish the importance matrix of the benchmark layer and the index layer under each benchmark layer. When constructing the importance matrix, the relative importance of the elements of each layer is quantified by numbers 1∼9, and finally the importance matrix P = ( p i j ) n × n is obtained. According to the Figure 1, there are four elements of the benchmark layer in the study. Therefore, four important matrices needed to be formed in this study, that is the attribute of the criterion layer to the target layer M = ( m i j ) 4 × 4 , the power supply reliability A = ( a i j ) 4 × 4 , voltage quality B = ( b i j ) 3 × 3 , adaptability C = ( c i j ) 3 × 3 , and economy D = ( d i j ) 2 × 2 .

In summary, this section forms the following flow chart of comprehensive evaluation score of active distribution network planning based on probabilistic power flow analysis, as shown in Figure 2.

In order to verify the rationality of the active distribution network planning evaluation method based on probabilistic power flow analysis proposed in this study, two areas, M and N, of an urban distribution network in China are used as simulation examples in this section. The distribution network in this area is connected with a large number of distributed photovoltaics. The system structure diagram is shown in Figure 3, and the dotted line in the figure indicates that it is not directly connected.