1 Introduction

In pursuit of the economic and environmental benefits, an increasing popularity has been witnessed in distributed generators (DGs) such as photovoltaic (PV) panels. The PV owners aim to gain benefits by exporting surplus energy back to the grid [1], which causes a number of operational challenges for the distribution networks. These challenges include power fluctuations, renewable output curtailment and voltage stability risk [14]. In order to solve such problems, utility-scale battery energy storage system (BESS), is believed to be a promising alternative. BESS is a centralized controllable storage system with the capacity ranging from hundreds of kW to dozens of MW. It is able to improve the reliability, optimize voltage level, shift load profile, make monetary profits and defer infrastructure upgrade [511].

The optimal management of BESS in a distribution system has attracted significant attention in recent years. A number of studies on determining the storage size have been presented based on short-term scheduling [1214]. Long-term BESS planning, combined with charging/discharging state management, have been reported in [7, 15, 16]. Optimal allocation of BESS has been investigated under the scenario of high renewable penetration aiming for the mitigation of operational risk from price volatility [15]. In price-taker mechanism, fluctuating power price may influence the decision of decentralized entities, causing additional risk in daily operation. A novel planning method, aiming at optimal utilization of BESS for mitigating price-related risk was also presented in [15]. From the aspects of operation, maintenance and finance, a comprehensive cost-based formulation was proposed in [16] to determine the optimal size of BESS. For the long-term planning, the dynamically increasing penetration of DG is an essential factor to be taken into account. However, few studies have considered the sustainable energy growth pattern which describes the nature of dynamic development of DG.

Due to the long-range simulation of BESS planning, the number of scenarios that need to be adapted is the main computational burden for optimization algorithms. There exist a wide range of approaches in the literature of power system scenarios reduction. For instance, regarding the high uncertainty of wind generation and classical seasonal variation, Sedghi et al. [7] proposed a BESS allocation approach to optimize global benefits based on four typical operation modes to represent annual scenarios. In [17], a probability density function was applied to divide operational state clusters in order to reduce scenarios for planning verification. In most of the previous literature, the operational uncertainty is only associated with DG. However, in this paper, k-means is adopted to cluster annual aggregated profiles (demand and PV) along the time axis. This reduces the long-term continuous scenarios into a limited number of classical operation modes.

The primary challenge in BESS planning lies in the determination of the appropriate capacity and rating power due to limited space and cost. The conventional methods facilitating optimal sizing of utility-scale BESS usually ignore the effects related to the dynamic growth of dispersed PVs. To coordinate access of BESS and renewable energy, a multi-stage model is employed in this paper. This separates the BESS sizing task into several continuous decision-making phases and enhances adaptability of BESS planning strategy in terms of the trade-off between technical and economic benefits. There is a number of literature aiming at multi-stage allocation in distribution network with DG and BESS in [6, 18, 19]. A multi-stage distribution network expansion approach using DG and storage units, for both peak shaving and reliability improvement, is presented in [18]. A modified particle swarm optimization algorithm was proposed to solve the high dimensional nonlinear model [18]. In order to accommodate highly penetrated DG, an optimization approach in [19] was established to model multi-stage distribution expansion problem under uncertain conditions. However, the BESS was not taken into account in that research. The coordination of multiple-objectives is one of the main challenges in this field. For example, Tant et al. [6] investigated a particular model using voltage regulation, peak load curtailment and annual costs as characteristics to optimize the size of storage system with high PV penetration. However the dynamic evolution of PV in the long term was not discussed in that paper.

In order to determine the optimal size of BESS to adapt the dynamic evolution of PV, a multi-stage sizing model for utility-scale BESS is presented in this paper. Based on the prediction of PV penetration, growth patterns are modelled using curve fitting, which enables cost-effective planning stage segmentation. With the goal of decreasing the computational burden, k-means clustering is applied to obtain the reduced annual operational scenarios without a significant compromise of generalization. Furthermore, a short-term operational model of BESS aiming at optimizing charge/discharge profit is established. This is a sub-problem of the long-term planning model which pursues maximization of comprehensive benefits in the long-run. The details of techno-economic benefits and chance-constrained security conditions are also provided. The contribution of this paper can be summarised as follows:

  1. 1)

    A comprehensive multi-stage storage planning model considering the increasing penetration of DG model is developed. The multi-stage planning can decrease the initial investment, and increase the liquidity of funds from the investor side.

  2. 2)

    The remaining life of the storage units is converted into the corresponding surplus value, making the multi-stage sizing approach for development of utility-scale BESS more applicable.

The remainder of this paper is organized as follows. Section 2 details the existing problems and provides the complete model of multi-stage sizing for BESS. Section 3 describes the methodologies applied to deal with short-term scheduling embedded in the sizing procedure. In order to illustrate the effectiveness of the proposed method, a case study is shown in Section 4. The conclusion is presented in Section 5.

2 Problem statement

2.1 Conditions and assumptions

In order to address the technical and economic issues of long-term planning for DGs integrated distribution, the following factors need to be taken into account in this research:

  1. 1)

    The location and size of DG are largely depended upon the decisions made by their owners. However, the planning of utility-scale storage units is carried out by distribution network company.

  2. 2)

    The planning horizon is segmented to the determined stages, enabling the enhancement of strategy adaptability of BESS accommodation. Particularly, the increase mode of distributed PVs connection can be investigated for each stage.

  3. 3)

    The daily operation scheduling of BESS can be seen as the typical strategy of state of charge (SOC).

  4. 4)

    In each stage, the interacting effects of long-term planning and short-term scheduling are modelled based on the integrated methodology, ensuring the conditions closer to the reality.

  5. 5)

    Limited by the charge/discharge cycles in the lifetime of BESS, the number of alternative state change is assumed to be once in a day.

For the sake of simplicity and generality, the following assumption are made:

  1. 1)

    The distribution network topology is considered to be constant during the BESS sizing procedure.

  2. 2)

    All the elements operate in the three-phase balanced distribution system. Therefore, one-line diagram is used to investigate the proposed method.

  3. 3)

    The growth pattern of self-invested DGs during the planning horizon can be estimated.

2.2 Objectives modelling

The proposed multi-stage sizing approach for the long-term development of BESS with the concern of short-term scheduling is shown in Fig. 1.

Fig. 1
figure 1

Sketch of multi-stage sizing approach

The primary goal in the proposed technique aims at maximizing the batteries profit in the complete planning horizon. For each planning stage i, the profit is evaluated via the difference value of BESS operational benefits and periodic investment cost. It is worth noting that due to the life-cycle constraints caused by charge state alternation limit, the surplus operational value B sur i modelled by remainder charge/discharge cycles ought to be included in the profit expression.

$$\hbox{max} \left\{ {\sum\limits_{i = 1}^{n} {(B_{i} - C_{i} + B_{i}^{sur} )\left( {\frac{{1 + r_{infr} }}{{1 + r_{intr} }}} \right)^{{t_{i} }} } } \right\}$$
(1)

where n is the number of the segmented stages; r infr and r intr respectively are the inflation rate and interest rate; t i is the beginning of stage i and appears according to the annual horizon.

The investment cost of stage i is the function of the storage capacity E i and rating power P i which are configured at the beginning of i th stage, shown as follow:

$$C_{i} = c_{E} E_{i} + c_{P} P_{i}$$
(2)

where c P is the unit cost of power rating for storage and c E represents the unit cost of battery investment.

Specifically, the batteries benefit is calculated based on four sub-factors: arbitrage income B arb i , environmental benefit B env i , reliability enhancement profit B rel i and the saving of upgrade deferring B del i .

$$B_{i} = B_{i}^{arb} + B_{i}^{\textit{env}} + B_{i}^{\textit{rel}} + B_{i}^{\textit{del}}$$
(3)

The first factor is composed of charging cost and discharging earning. Both energy of charge p cha i and energy of discharge p dis i depend upon the rating power of storage units existing in i th stage. The charge power is purchased from distribution network for energy storing while the discharge power is sold back for peak load shaving. Therefore, the arbitrage income can be evaluated as follow:

$$B_{i}^{arb} = \sum\limits_{y = 1}^{{Y_{i} }} {\sum\limits_{c = 1}^{{N_{y} }} {\sum\limits_{t = 1}^{24} {(c_{t} p_{i}^{dis} - c_{t} p_{i}^{cha} )} } \omega_{c,y} }$$
(4)

where c t is the active energy price at time t; N y is the number of cluster scenarios in y th year and ω c,y is the weight of scenarios c in y th year; Y i is the number of year at stage i.

From the environmental point of view, the utilization of storage units is able to reduce pollutant emission, as reported in [20]. Therefore, the environmental benefit is taken into account in this research using the following formula:

$$B_{i}^{\textit{env}} = \sum\limits_{y = 1}^{{Y_{i} }} {\sum\limits_{c = 1}^{{N_{y} }} {\sum\limits_{t = 1}^{24} {(C(P_{t,c,y}^{bef} ) - C(P_{t,c,y}^{aft} ))} } \omega_{c,y} }$$
(5)
$$C(P_{t,c,y}^{R} ) = a + bP_{t,c,y}^{R} + c(P_{t,c,y}^{R} )^{2} \quad R \in (\textit{bef},\textit{aft})$$
(6)

where C(  ) is the unit production cost of the conventional power; P bef t,c,y and P aft t,c,y are the power purchase from substation before and after the installation of BESS and PV; a, b, c are the fuel cost coefficients of the conventional power units.

In this work, the enhanced reliability profit is quantified as the savings from reduced probability of electricity outage. Additionally, the PV power is not used to compute the reliability enhancement in this paper, therefore the potential reliability-improvement benefits are not involved in the proposed model. The interrupted energy assessment rate can be applied in this paper to calculate the profits BESS brought in distribution network [21]:

$$B_{i}^{\textit{rel}} = \sum\limits_{y}^{{t_{i} }} {\sum\limits_{k}^{{N_{\textit{ev}} }} {\hbox{min} (P_{i} ,P_{k,t}^{\textit{out}} )c\;\hbox{min} (h_{k} ,h_{\textit{ESS,k}} )\lambda_{k} } }$$
(7)
$$h_{ESS,k} = \frac{{E_{t} \eta_{\textit{dis}}}}{{\hbox{min} (P_{i} ,P_{k}^{\textit{out}} )}}$$
(8)

where N ev is the number of event; λ k is the average failure rate; P out k,t is the power not supplied at time t in node k; c is interrupted energy assessment rate (in ¥/kWh); h k is the time to repair/replacement; h ESS,k is the time BESS can supply; E t is the energy storage had at time t.

The approach to quantify the distribution network deferral value is shown in [22], which can be adopted to the system with BESS. The deferral saving is associated with the distribution network construction cost C e and deferral time τ, related with the rate of peak shaving α and growth rate of load ε. In order to avoid repeated calculation, the deferral saving, will be counted only at the final stage.

$$B_{i}^{\textit{del}} = \left\{ \begin{array}{cl} \!\!\!\!\!\! \!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!0 &\quad \quad i = 1,2, \cdots,n - 1 \hfill \\ \!C_{e} (1 - \frac{1}{{e^{{r_{\textit{intr}} \tau}} }}) & \quad \quad i = n \end{array} \right.$$
(9)
$$\tau = \frac{\lg (1 + \alpha )}{\lg (1 + \varepsilon )}$$
(10)
$$\alpha = \frac{{P_{\textit{max}}^{\textit{bef}} - P_{\textit{max}}^{\textit{aft}} }}{{P_{\textit{max}}^{\textit{bef}} }}$$
(11)

where \(P_{\textit{max}}^{\textit{bef}}\) and \(P_{\textit{max}}^{\textit{aft}}\) are the maximum power purchasing form the substation before and after configured BESS respectively.

In this paper, the benefit of BESS remaining operation B sur i is quantified by the charge/discharge value out of the planning cycle. The surplus value ought to be taken into account due to the fact that the value in use of BESS remained after planning period is still available with assumption that the storage ability is evaluated by the number of charge/discharge actions. Therefore, the surplus value can be written as:

$$B_{i}^{\textit{sur}} = C_{i} \left(1 - \frac{{n_{\textit{opi}} }}{{\chi_{\textit{bat}} }}\right)$$
(12)

where χ bat is the number of charge/discharge cycles in storage life; n opi is the cycles that the storage installed in stage i have operated.

To improve the performance of the proposed model, three new factors are taken into account: i) The balance of active/reactive power shown in (13); ii) limitation of equipment output (14) are obliged to satisfy; iii) the chance constraint is introduced to deal with voltage constraint present in (15):

$$\sum\limits_{i = 1}^{{n_{LD} }} {S_{\textit{LD},t} + S_{\textit{BESS,t}} + S_{\textit{LOSS,t}} = \sum\limits_{d = 1}^{{n_{\textit{DG}} }} {S_{\textit{DG,t}}}}+ S_{\textit{SUB}}$$
(13)
$$P^{R} \le P_{\textit{max}}^{R} \quad R \in \{ {\text{DG}},{\text{BESS}},{\text{Grid}}\}$$
(14)
$$p\{ V_{\textit{min}} \le V \le V_{\textit{max}}\} > \beta$$
(15)

where S LD,t , S BESS,t , S DG,t and S SUB respectively are the electric load, charge/discharge load of BESS, DG output and power exchange of balancing bus at time t; S LOSS,t is the power loss of network; V min and V max are the permitted minimum and maximum voltage; β is the chance limit factor; P R means the power limits of the specific devices including DG and BESS.

3 Methodology

The implementation framework shown in Fig. 2 can achieve the proposed multi-stage sizing approach for utility-scale BESS in the long-time scale through the following steps.

Fig. 2
figure 2

Implementation framework of multi-stage BESS configuration

  1. 1)

    Segment the whole planning period into a number of discrete stages according to the estimated growth pattern of PVs connection.

  2. 2)

    Due to the similarity of operating features associated with power profiles in different days, k-means clustering which is able to assess the optimum number of clusters is employed to reduce the short-term operation scenarios [9]. The clustering centres which represent the specific groups of operation points are obtained to generate classical scenarios for short-term scheduling in the years. In addition, the weight of centre is defined as the appearance probability of similar scenarios, which can be determined by the number of scenarios in the same group. Applying this method is able to enhance the generalization of scenarios segmentation by using the scalar variables ω c,y in (4) and (5).

  3. 3)

    In order to gain the solution with higher efficiency, the short-term scheduling, a day-ahead operation optimization, is performed. The strategy is designed for daily operation: charge when the charging price is in valley and vice versa. Aiming at maximizing the BESS utilization, the charge/discharge power is assumed to be operated under the rating power condition. The short-term scheduling is detailed in Section 3.3.

  4. 4)

    An adaptive particle swarm optimization (APSO) algorithm [23] is adopted to solve the proposed model.

3.1 Growth pattern of DG and multi-stage definition

The growth pattern of general products is modelled by innovative and imitative behaviour in [24]. It is worth noting that, distributed sustainable energy, as the newly developing object, is commonly assumed to conform the likewise growth law. Based on the related theory, three growth patterns of PV forecasted in a long range in the distribution network is shown in Fig. 3. Pattern 1 is the PV grows at a constant velocity. In pattern 1, PV panels are allocated with the constant growth rate. Pattern 2 describes the situation that PV grows very fast at the very beginning of the initial date while growth speed gradually decreases. Pattern 3 provides the tendency that a large amount PV panels are connected in the middle of the forecasting interval.

Fig. 3
figure 3

Three PV growth patterns

On the account of the dynamically increasing growth of PV, we propose a method to segment the whole horizon into several stages using the annual forecast data of PV. The Fitting Curve Tool of MATLAB is applied to obtain the fitting curve which contains the continuous information of PV growth pattern.

The segmentation of multi-stages is determined by the following steps. First, the derivation of the fitted growth curve is calculated in order to get the growth rate of PV. Second, confirm the starting and ending point of the most rapidly growing period of PV. Third, divide the focused period into several segments equally, which enables the stage segmentations. The detailed example is illustrated in Section 4.1.

3.2 Scenarios clustering

The long-term PV generation and load simulation approach are needed to execute the planning in a long time.

Tina et al. [25] proposed a model in terms of solar radiation and temperature to simulate the long-term generation of PV:

$$I_{r,t} = \frac{1}{2}\sin \alpha_{t} k$$
(16)
$$P_{t} = P_{s} \frac{{I_{r,t} }}{{I_{s} }}[1 + \alpha_{T} (T_{t} - T_{s} )]$$
(17)

where α t is the solar altitude; k is the random factor; I r,t is the radiation intensity at time t; T t is the temperature of photovoltaic panel at time t; α T is the power temperature coefficient of photovoltaic panels; P s is the PV power at the standard condition (the solar radiation intensity I s is 1000 W/M2, temperature T s is 25°C).

As described in [26], the load is modelled by the nonlinear equation consisted of several factors. In the current paper, we introduce monthly and daily factors based on load type to replace the demand level:

$$S_{i,m,t} = S_{i,m,b} F_{m,t} F_{D,t} (1 + \varepsilon )^{T} R$$
(18)

where S i,m,b is the base power of node i at m th month; F m,t is the monthly factor whiles F D,t is the daily factor; ε is the annual increment of load and R is the random factor .

In order to reduce the computational burden of a full time-series analysis, in this paper, k-means is adopted to deal with photovoltaic availability and demand of one year into a manageable number of photovoltaic/demand scenarios based on their joint probability of occurrence. In addition, the mean index adequacy (MIA) [27] defined as the average of the distances between each input vector, assigned to the cluster and its centre, is employed to determine the number of scenarios clustering:

$$r_{MIA} = \sqrt {\frac{1}{S}\sum\limits_{sc = 1}^{S} {d^{2} (\vec{\omega }_{\textit{sc}} ,\vec{\varOmega }_{\textit{sc}} )} }$$
(19)
$$d(\vec{\omega }_{sc} ,\vec{\varOmega }_{sc} ) = \sqrt {\frac{1}{{N_{sc} }}\sum\nolimits_{{\vec{x}_{l} \in \varOmega_{sc} }} {\frac{1}{K}\sum\limits_{k = 1}^{K} {(x_{{l_{1} i}} - x_{{l_{2} i}} )^{2} } } }$$
(20)

where S is the corresponding set of the classes; d 2(ω,Ω) is the distance between the representative vector ω sc of sc th cluster and the subset Ω, calculated as geometric mean of the Euclidean distances between ω and each member of Ω.

The center scenario of each group is used to represent this group of scenarios. Therefore, the number of members, namely ω c,y , of this specific group is defined as the weight of scenarios in (4) and (5).

3.3 Short-term scheduling

Several typical photovoltaic/demand scenarios can be obtained to simulate the short-term scheduling based on Section 3.2.

With the introduction of the DGs and BESS, the power flow may be reversed in the radial distribution network. Therefore, the forward-backward algorithm that based on basic circuit theory is adopted to calculate power flow and the details can be found in reference [28].

In the current work, some special constraints applied to the storage units are listed as follows.

  1. 1)

    On account of lifetime limitations and replacement cost, the storage units should be charged and discharged only once in each day.

  2. 2)

    At the status of charge, the SOC of batteries is limited by the maximum SOC. Meanwhile, to protect the storage units, BESS should not be allowed to discharge when the SOC is less than the minimum SOC.

  3. 3)

    Considering the storages installed are not able to change power quickly, the charge and discharge power in a daily operation remains unchangeable.

To calculate the arbitrage income and environmental benefit more acceptable, the SOC should be updated in each hour by the following equation [7]:

$$SOC_{t + 1} = \left\{ \begin{array}{l} SOC_{t} + \frac{{p_{t,j}^{dis} }}{{\eta_{dis} S}} \hfill \\ SOC_{t} + \frac{{\eta_{ch} p_{t,j}^{cha} }}{S} \hfill \\ \end{array} \right.$$
(21)

In the realistic application, the shape of electricity price in a day is quite similar with the variation of load profile. For example, the heaviest loaded period of power systems usually correspond to the highest electricity price. Therefore, from the supply theory point of view [29], it is reasonable to assume that the charge action is response to the electricity price instead of load profile. The flowchart of short-term operation strategy aiming at maximizing the arbitrage income and environmental benefits while computing by PSO is shown in Fig. 4. The particles are the time of start charge T cs and time of start discharge T ds in daily operation. The flowchart is explained in detail in the following:

  1. 1)

    Initialization. This section is tagged using symbol ① in Fig. 4. A number of factors including rating power P i and capacity E i of BESS in i th stage need to be confirmed in this block. For example, in this study, we set the current time point T = 1.

  2. 2)

    Charge strategy confirmation. This section is tagged using symbol ②. This section of flowchart describes the charging behavior of BESS in the process of short-term scheduling. More specifically, when the condition T ≥ T cs , is satisfied, the actual charge power P cha is determined by SOC. If SOC < SOC max , set P cha  = P i , else let P cha  = 0. Compute the comprehensive benefits including arbitrage income, environmental benefits and auxiliary information.

    Here, the auxiliary information refers to the other useful data involving the times of abnormal voltage occurrence and maximum power injection from the upper substation.

  3. 3)

    Discharge strategy confirmation: This section is tagged using symbol ③. This section describes the discharging behavior of BESS, which is the inverse procedure of charge strategy confirmation.

  4. 4)

    Fitness update: This part is tagged using symbol ④ This section is used to update the fitness of particles and related auxiliary information.

  5. 5)

    Optimal results: This section is tagged using symbol ⑤. The optimal T cs and T ds are determined. Meanwhile, the maximum comprehensive and related auxiliary information are computed.

The short-term scheduling strategy investigated in this paper aims at maximizing the comprehensive benefits integrated with arbitrage incomes. Particularly, it is widely accepted that the arbitrage incomes have a highly close relationship with charge/discharge price rather than environmental benefits. Therefore, the operation strategy of BESS is mainly affected by the given electricity price. The charge and discharge price are assumed to be fixed, hence resulting in the constant T ds and T cs .

Fig. 4
figure 4

Flow chart of short-term scheduling

3.4 APSO algorithm

The planning of BESS in distribution network considering the growth of DG is a mixed integer nonlinear optimization which cannot be easily solved by conventional optimization tools. Therefore, a heuristic algorithm is needed to compute both in long-term planning and short-time scheduling. Compared with traditional PSO algorithm, the APSO algorithm is more efficient in searching global optimization solution and the details of APSO algorithm is presented in literature [23].

4 Case study

The test system used for the case study is the modified IEEE 69-bus distribution system shown in Fig. 5. The available locations which are feasible for PV’s growth are assumed at bus 9, 17, 20, 33, 36, 48 and the BESS is located at bus 7. The system peak demand is 5.1746 + 2.7054 MVA at base year, and the annual increment rate of load ε is a constant 5%. In this paper, the planning period is set to be 20 years. The interest rate r intr and inflation rate r infr are further assumed to be 7% and 10% respectively. The interrupted energy assessment rate c is 29.84 ¥/kWh [21]. The price profiles of distribution network purchasing from upstream network and selling to consumers, adopted from [4], are shown in Fig. 6. The price of PV output is set to be 0.6 ¥/kW. Fuel cost coefficients of the conventional power units a, b, c are 2933.4, 300.0 and 0.6877 respectively [20]. The chance constraint of voltage magnitude is set to be 0.9. The parameter’s configurations of BESS performance as well as reliability indices of test system are listed in Table 1[7].

Fig. 5
figure 5

Modified IEEE 69-bus distribution system

Fig. 6
figure 6

Price profiles of charge and discharge

Table 1 Parameters of test system

4.1 Data pre-processing

According to the growth pattern model introduced in Section 3.1, there are three classical patterns of PV growth data used in the study, shown in Fig. 7. The commonly used Gaussian function is adopted to fit the forecast data. The root mean square errors (RMES) of three fitting results are 0.9999, 0.9967 and 0.9993 respectively, indicating the approach has well fitness for PV growth tendency. The curve fitting results and theirs differential coefficients illustrating the absolute growth rate are also shown in Fig. 7.

Fig. 7
figure 7

Three PV growth pattern segment

There is a strong relationship between increment of PV and the sizing of BESS. Therefore, the derivative curves of the fitting data are employed to segment the whole planning period into several stages. For instance, the PV panels with pattern B, which grow quickly in the early period the planning stages, can be segmented into 3 stages: 1st–4th year, 5th–7th year, 8th–20th year. The segment results for other two patterns are shown in Fig. 7.

Based on the long-range forecast model of PV and load, PV output and load for a day in a year can be generated. Taking the annual PV amount into consideration, the clustering would be executed 20 times in the planning. For example, the 5th year results with clustering number 30 under growth pattern A are presented in Fig. 8. Comparing trends between before and after clustering, it is worth noting that the trends of clustering are consistent with the trends of time-series scenarios. Thus it indirectly validates the choice of k-means to deal with the uncertainty of PV and load.

Fig. 8
figure 8

Time-series scenarios and clustering results

4.2 Results analysis

The proposed APSO algorithm has been computed several times and the optimal BESS sizing results in the distribution network with the increasing number of PV panels are shown in Table 2. Under the growth of pattern A, the power and capacity installed in each stage are all 500 kW. The allocation of rating power and capacity are different in each stage owing to the difference of PV increment in each stage of the pattern B and pattern C. The proportions of the rating power and capacity in each stage are approximately identical under the same PV growth pattern.

Table 2 Multi-stage BESS sizing results

In Fig. 9, the investment cost, surplus operation value, arbitrage income and environmental benefit, reliability enhancement profit and saving of upgrade deferring under three growth patterns of PV are given. The results indicate that the storage units are economical and helpful. Under the short-term operation strategy proposed in this paper, the arbitrage income plays the dominating role. And there is a positive correlation between the BESS installed factors (time and scale) and the three operation benefits (arbitrage income, environmental benefits and reliability enhancement profits). Comparing the saving of deferring upgrade under different growth patterns, it can be seen that the deferral saving is almost the same because the model proposed for calculating savings is based on the BESS had

Fig. 9
figure 9

Investment cost and benefits of BESS in the planning

been installed in last stage. It is worth to note that there are little difference among the total rating power and capacity configuration results during the planning under the three growth patterns, but the objective benefits are somewhat different. For instance, the total rating power (1500 kW) and capacity (4567.8 MWh) of BESS under pattern A are approximately to power and capacity under growth pattern B whose power and capacity are 1509.3 kW and 4609 MWh respectively, while the objective values are ¥2.3156 million and ¥3.043 million separately. The reason is that the time and amount of BESS allocated in each stage under three patterns are different. According to Fig. 9, the investors can gain more value under pattern B. The reason is that the multi-stage size of BESS under pattern B gains more operation benefits and less investment cost.

The maximum load-shaving, the average load shift and the rate of peak-shaving with three growth patterns of PV in different years are listed in Table 3, illustrating the distinct effects of BESS allocation in different stages. It can be learned from Table 3 that the stages segmentation is related to the growth pattern of PV. The segmentation point is assumed to be located at the starting, middle and ending years of rapid PV growing period in this paper. In addition, along with the increased sizing of utility-scale BESS, the global load shift and shaving performance is promoted gradually.

Table 3 Peak-shaving index

Figure 10 presents the active power curves that distribution network company purchases from the substation with and without PV and BESS in three scenarios under three growth patterns. The plot reveals that the storage units have certain influence on peak-shaving, even if the penetration of PV is low. The reason is that the BESS charge and discharge earns benefits meanwhile it also has a good performance in shifting load. However, with the allocation of BESS and PV, the peak-shaving rate rises considerably. For instance, under the growth pattern B, the rates of peak-shaving in the 2nd year, 5th year and 9th year are 6.01%, 10.7%, 15.82% respectively. In a word, with the cooperation of PV and BESS, the performance of peak-shaving is substantially improved.

Fig. 10
figure 10

Exchange power of feeder root node curves

In order to detail the operation of BESS with three growth patterns, the daily storage capacity profiles under the scenarios of Fig. 10 are depicted in Fig. 11. The storage capacity curves reflect that the BESS are efficiently utilized. The energy price profile in this study is assumed to keep constant in the whole planning horizon, resulting in the same operation strategies of BESS in a particular stage. In addition, under the specific charging/discharging price, the self-invested PV power is not used to charge the centralized storage units in this paper.

Fig. 11
figure 11

Storage capacity curves in short-term scheduling

In IEEE 69-bus test system, the branch from bus 43 to bus 55 is heavily loaded. Without installing the PV and storage, the number of voltage violation is almost 28000 times in the coming 20 years and most of them are violated towards the lower limitation. While with the placement of the PV and BESS, the simulated results see a sharp decrease in the occurrence of voltage violation at only 4672 times in 20 years. The reason for such voltage quality improvement is that PV installed in bus 48 and BESS placed at bus 7 can provide a part of the active and reactive power to compensate the power in this branch. It can be concluded that the cooperation of BESS and PV can improve the voltage profile. Moreover, the BESS can work as backup power to enhance the reliability of distribution system when power is not supplied as shown in Fig. 9.

4.3 Discussion

According to the Fig. 9, it can be learned that the investment cost, surplus operation value, arbitrage income play crucial roles in computing the optimal results. In other words, the spot price of active power and the units cost of batteries have great influence on optimal BESS sizing results. When the approach is used to calculate the sizing of BESS in one-stage planning, the uncertainties of batteries cost and energy price rise the economic risk and the configuration results may not suit the future distribution network. However, the model and optimal result of multi-stage planning can be modified along with the time.

The allocation results under different patterns of PV increase are different. Therefore, it is critical to take the growth patterns into consideration due to the uncertainty of PV amount in the coming years. As a result, the optimal allocation results are evaluated considering different PV growth patterns. For instance, based on the similarity of allocation under different growth patterns, the most suitable allocation for the first stage can be chosen.

5 Conclusion

In this paper, the multi-stage planning approach for BESS applied in distribution network is presented, enabling the best trade-off between storage system cost and benefits including maximizing the PV connections. The short-term operation behaviour and long-term planning process are modelled in the proposed approach in order to investigate a number of interaction effects. For the sake of computational simplicity, the scenarios based on PV-load joint probability are diminished by applying k-means clustering. The results of numerical study illustrate that the proposed methodology is able to determine the optimal multi-stage strategy for utility-scale BESS sizing, adapting to the diverse growth pattern of disperse PVs. Especially, peak load shaving is facilitated due to the collaborative consideration of time-serial charge action and probabilistic PVs generation. In the future, additional factors including the volatility of charge/discharge tariff and detailed life cycle cost of utility-scale BESS will be taken into account.