An Analytical Method for Delineating Feasible Region for PV Integration Capacities in Net-zero Distribution Systems Considering Battery Energy Storage System Flexibility

Shida Zhang; Shaoyun Ge; Hong Liu; Guocheng Hou; Chengshan Wang

网刊加载中。。。

使用Chrome浏览器效果最佳，继续浏览，你可能不会看到最佳的展示效果，

确定继续浏览么?

复制成功，请在其他浏览器进行阅读

An Analytical Method for Delineating Feasible Region for PV Integration Capacities in Net-zero Distribution Systems Considering Battery Energy Storage System Flexibility PDF

- ORCID：
Shida Zhang
✉
- ORCID：
Shaoyun Ge
✉
- ORCID：
Hong Liu
✉
- ORCID：
Guocheng Hou
✉
- ORCID：
Chengshan Wang
✉

the School of Electrical and Information Engineering, Tianjin University, Tianjin 300072, China

Updated：2024-03-26

DOI：10.35833/MPCE.2023.000633

OUTLINE

Abstract

To provide guidance for photovoltaic (PV) system integration in net-zero distribution systems (DSs), this paper proposes an analytical method for delineating the feasible region for PV integration capacities (PVICs), where the impact of battery energy storage system (BESS) flexibility is considered. First, we introduce distributionally robust chance constraints on network security and energy/carbon net-zero requirements, which form the upper and lower bounds of the feasible region. Then, the formulation and solution of the feasible region is proposed. The resulting analytical expression is a set of linear inequalities, illustrating that the feasible region is a polyhedron in a high-dimensional space. A procedure is designed to verify and adjust the feasible region, ensuring that it satisfies network loss constraints under alternating current (AC) power flow. Case studies on the 4-bus system, the IEEE 33-bus system, and the IEEE 123-bus system verify the effectiveness of the proposed method. It is demonstrated that the proposed method fully captures the spatio-temporal coupling relationship among PVs, loads, and BESSs, while also quantifying the impact of this relationship on the boundaries of the feasible region.

Keywords

Net-zero distribution system; photovoltaic (PV) integration capacity; feasible region; battery energy storage system (BESS).

I. Introduction

A net-zero energy system, where the total generation surpasses overall demand, signifies an ecologically sustainable paradigm that supports decarbonization procedure and potentially serves as a crucial form for future energy systems [

1]-[3]. The rapid development of photovoltaic (PV) technology in recent years has laid the foundation for the realization of net-zero distribution systems (DSs) [1]. The investigation into PV integration schemes provides invaluable guidance for establishing net-zero DSs.

Related to the mentioned topic, there are many outstanding research works on the PV hosting capability. Based on considerations of multiple PV locations, the current research works can be divided into two aspects: total PV integration capacity (PVIC) of all PV locations, and individual PVICs for different PV locations. In the first aspect, optimization models for total PVIC are employed to determine the maximum total PVIC, which have been explored in many literatures [

4], [5].

In contrast to quantifying the maximum total PVIC, the evaluation of individual PVICs places an emphasis on delineating the interaction relationships between PVICs at different PV locations, which can provide a more comprehensive information about PV hosting capability. However, the description and solution methods for individual PVICs are the main issues. References [

6] and [7] used the Monte Carlo method and affine arithmetic method, respectively, to find feasible PV integration schemes including locations and their capacities. Reference [8] extended the security region to the total quadrant by treating PVs as negative loads. In our prior work [9], the PV hosting capacities were characterized as a feasible region within a high-dimensional space, which provided explicit boundaries for the interdependent maximized PV capacities at different buses. It becomes evident that the region-based method is efficient for solving such problems.

In observation of the feasible region for PVICs, simulation techniques [

8] and multi-objective optimization [9] can be employed. However, they provide reliable solutions only at the Pareto points, and the accuracy of this region heavily depends on the accuracy of the high-dimensional fitting algorithm, which may be difficult to achieve in larger systems. Reference [10] defined the feasible operation region of a DS, which indicates the allowable power injections at different buses. Quadratic analytical expressions can effectively represent its boundaries. But the time-series PV outputs and temporal coupling flexible resources such as battery energy storage systems (BESSs) cannot be addressed in this method. Therefore, it is of significant value to employ analytical methods to accurately derive the expression of the feasible region for PVICs, which poses significant challenges.

For the constraining factors of PV integration, current PV hosting capacity studies primarily revolve around system security such as bus voltages [

9] and branch power flows [8]. However, a certain amount of PV integration is necessary to achieve net-zero. Specifically, for a net-zero system, the feasible PV integration capacity must not only adhere to the upper limits dictated by network security constraints but also satisfy a lower threshold. In this way, it is more appropriate to conduct the research on “feasibility of PV integration schemes” rather than “the maximum PV hosting capability” in net-zero systems. Furthermore, net-zero systems may exhibit varying requirements depending on their definitions of “net-zero” such as energy [2], [3] and carbon emission [1]. These requirements profoundly affect the configuration of the feasible region.

Furthermore, the importance of “flexibility” in the establishment of net-zero systems has been demonstrated [

1]. Flexibility is the capability of a system in mitigating the intermittency and volatility of renewable energy generation (REG) [11], thus ensuring secure and low-carbon system operation [12]. BESSs facilitate energy temporal transfer, substantially enhancing the PV hosting capability in DSs [5], and providing flexibility for net-zero DS construction. However, the interaction effects between viable PVICs and BESSs are complicated. Specifically, loads, PV outputs, and regulating power of BESSs are time-varying, but power balance constrains their mutual influences across temporal sequences; and the power of loads, PVs, and BESSs during a given time period is also limited by their locations and network security constraints. Consequently, quantifying the influence of BESSs on feasible PVICs relies on the formulation and solution of the spatio-temporal coupling among loads, PVs, and BESSs.

To address this issue, it is imperative to analyze the spatio-temporal coupling flexibility of a DS using a region-based method. References [

13] and [14] separately obtained the temporal coupling regulation power of aggregated flexible resources at the point of common coupling (PCC) by using external approximations and analytical methods, respectively. However, these methods may not be able to analyze the effect of BESSs on multiple PV locations. References [15] and [16] examined the capability of accommodating uncertainties for each REG buses, presenting the solution to spatial coupling. Nevertheless, this method may not consider temporal coupling constraints of BESSs (state-of-charge (SOC) constraints). Reference [17] analyzed the operation flexibility of mitigating uncertainties in a DS using basic feasible solution method, which is capable of addressing spatio-temporal coupling issues. Thus, it may provide an effective tool to formulate the feasible region for PVICs considering spatio-temporal coupling among loads, PVs, and BESSs in a DS.

In addition, the uncertainty of PV outputs, which mainly influence network security, is essential in the feasibility assessment of PV integration. Based on the data set of PV outputs, it is feasible to generate empirical probability distributions or assume certain parametric distributions to deploy stochastic programming [

18], but these distributions may differ from the true distribution when data samples are limited, potentially leading to an underestimation of the violation level and posing security risks to the safe operation of the system [19]. Therefore, we adopt the distributionally robust method to model the uncertainty. Currently, two main types of distributionally robust ambiguity sets are adopted: moment- and statistical-distance-based sets [20]. However, using two moments (the first and second moments) can capture only part of the true distribution information, sometimes leading to over-conservative decisions [12], [21]. To construct statistical-distance-based sets, the Kullback-Leibler (KL) divergence or the Wasserstein metric can be applied. Nevertheless, the KL divergence also has some limitations, and the main shortcoming is that it cannot be used for the models with heavy-tailed random functions [21], [22]. The Wasserstein metric has the following statistical guarantee property: a Wasserstein ball based on finite samples can contain the true probability distribution with a high probability level [23]. In this way, the traditional chance constraint can be reformulated as data-driven distributionally robust chance constraint (DRCC), and it is well-performed in the existing literature [21], [22].

Based on the above analysis, the contributions of this paper are presented as follows.

1) This paper introduces a feasible region model for PVICs under DRCC-based network security constraints, the network loss constraint, and energy or carbon-emission net-zero constraints, which define the upper and lower boundaries, respectively. This model fully clarifies the requirements for multiple PV locations, which are necessary for establishing a well-performed net-zero system. Moreover, the proposed model also presents the influence of BESS flexibility on enhancing the PV hosting capability and facilitating the construction of net-zero systems.

2) A novel analytical method for delineating the feasible region is presented. Utilizing the basic feasible solution concept, it deduces an analytical expression of the region from the DRCC-based linearized model, which is not reliant on existing optimization methods. A procedure is designed to verify and adjust the feasible region, ensuring that it satisfies network loss constraints under alternating current (AC) power flow. The systematic solution method captures the spatio-temporal coupling among PVs, loads, and BESSs, while also quantifying the impact of this relationship on the boundaries of the feasible region.

The significant differences between this paper and our prior work [

9] are elaborated as follows.

1) The most crucial difference lies in the distinct methods and their outputs. In terms of methodologies, [

9] employed a multi-objective optimization framework, treating the Pareto frontier as the boundary of the feasible region. In contrast, this paper proposes a novel analytical method, enabling us to derive a mathematical expression for the feasible region. The proposed method in this paper is entirely independent of optimization methods. In terms of outputs, the method in [9] produced the Pareto points. While it can generate boundary expressions via polynomial fitting or convex hull techniques, achieving higher precision necessitates more densely packed Pareto points. This leads to exponentially increased solution time and complexity of expressions, rendering them impractical. The proposed method in this paper directly yields an explicit and relatively simple mathematical expression for the boundaries, which is free from these issues.

2) Unlike [

9], which only considered system security constraints, this paper takes into account the support of PV in building net-zero systems. The feasible region in this paper not only includes an upper limit for PV integration but also a lower one.

3) For addressing uncertainty, this paper uses a linearization transformation method for Wasserstein-metric-based distributionally robust chance constraints (WDRCCs) to derive analytical expressions. In contrast, [

9] used an iterative method to solve distributionally robust chance-constrained programming.

The rest of this paper is organized as follows. Section II describes the specific problem (feasible region for PVICs) that we try to address in this paper. Operation constraints and net-zero requirements for DSs are presented in Section III. Section IV presents the formulation and solution of the feasible region for PVICs. Section V presents case studies to verify the effectiveness of the proposed method on the 4-bus system, the IEEE 33-bus system, and the IEEE 123-bus system. Finally, Section VI concludes this paper.

II. Problem Description

To ensure generality, the research object of this paper is the set of feasible PVICs at multiple PV locations. Due to the temporal coupling BESSs, scattering PV locations, network security constraints, and net-zero operation requirements, the feasible PVICs at different locations are coupled. To fully represent this relationship, the feasible region of all PVICs forms a high-dimensional geometric as follows.

Ω^{P V} = \{{(S^{P V})}^{T} \in R^{N^{P V}} |F (S^{P V}) \leq 0\}

(1)

S^{P V} = {[\begin{matrix} S_{1}^{P V} & S_{2}^{P V} & \dots & S_{N^{P V}}^{P V} \end{matrix}]}_{1 \times N^{P V}}

(2)

where $F (S^{P V}) \leq 0$ indicates all constraints that only involve $S^{P V}$ .

A schematic diagram shown in Fig. 1 is used to clarify the feasible region for PVICs introduced in this paper, which shows the PV hosting capability for the three integrated locations, with each dimension indicating the PVIC for the corresponding location. Figure 1 demonstrates that the feasible region for PVICs provides sufficient information on the PV hosting capability in a DS.

Fig. 1 Schematic diagram for feasible region for PVICs.

Note that in this paper, we do not impose any restrictions on the potential PV locations. In practice, the PV locations are influenced by geographical conditions and investor’s preferences. Therefore, an appropriate method should be adaptable to various user inputs (i.e., different sets of PV integration locations) without preset conditions. The proposed method in this paper accommodates any subset of system buses, including the set comprising all buses. While users ensure the rationality of chosen PV locations, the proposed method guarantees adaptability to these inputs.

III. Operation Constraints and Net-zero Requirements for DSs

In this section, we introduce the operation constraints and net-zero requirements that feasible PVICs need to satisfy in DS operation.

A. Linearized Power Flow Model

Given that our study considers the impact of temporal coupled BESSs, the calculation of the feasible region must be based on an accurate measure of flexibility. However, embedding the nonlinear AC power flow equations and the temporally coupled constraints in the flexibility analysis problem is computationally intractable [

14]. Consequently, existing literature considering nonlinear power flow does not take into account the temporal coupling in flexibility [24], [25]. Many research works in the related field use linear power flow models, indicating that linear power flow can achieve high accuracy and computational efficacy [26]-[29]. Therefore, researchers have highlighted that the error of whether to precisely consider temporal coupled flexibility or not is more significant than that of the choice of linear or nonlinear power flow equations [14]. To accurately measure the impact of BESSs on PV hosting capability, we opt for linearized flow equations in this paper, which allows for obtaining analytical expressions for the feasible region with minimal error. The linearized AC power flow model [30] is given as:

[\begin{matrix} P_{t} \\ Q_{t} \end{matrix}] - [\begin{matrix} B_{2}^{c o l} \\ - B_{1}^{c o l} \end{matrix}] θ_{1, t} - [\begin{matrix} B_{1}^{c o l} \\ B_{2}^{c o l} \end{matrix}] V_{1, t} = [\begin{matrix} B_{2}^{s u b} & B_{1}^{s u b} \\ - B_{1}^{s u b} & B_{2}^{s u b} \end{matrix}] [\begin{matrix} θ_{t} \\ V_{t} \end{matrix}] \forall t

(3)

where $P_{t}$ , $Q_{t}$ , $θ_{t}$ , and $V_{t}$ indicate the vectors of $P_{i, t}$ , $Q_{i, t}$ , $θ_{i, t}$ , and $V_{i, t}$ , respectively, at all buses. For brevity, the detailed expressions are available in [

30]. This equation is not limited to radial DSs but applicable to meshed DSs as well. Due to its assumption, the bus voltages of this equation are more optimistic than AC power flow equations. Next, the equations of bus power injections are given as:

\{\begin{array}{l} P_{i, t}^{} = - P_{i, t}^{l o a d} + (1 - μ_{t}) κ_{t}^{m a x} S_{i}^{P V} (1 + ξ) - P_{i, t}^{B E S S} \\ Q_{i, t}^{} = - Q_{i, t}^{l o a d} \end{array} \forall i \in Ω^{P V L}, \forall t

(4)

\{\begin{array}{l} P_{i, t}^{} = - P_{i, t}^{l o a d} - P_{i, t}^{B E S S} \\ Q_{i, t}^{} = - Q_{i, t}^{l o a d} \end{array} \forall i \in Ω^{b u s} \ Ω^{P V L}, \forall t

(5)

The output power of PV panel has losses, including dust and dirt, inverter losses, cable losses, etc. [

31]. The value may vary under actual conditions.

Equation (4) only needs to consider $κ_{t}^{m a x}$ since violations of network security constraints are most likely to happen in the season corresponding to $κ_{t}^{m a x}$ . PV integration leads to an increased voltage, presenting a risk of exceeding network security constraints. Therefore, for a certain PV integration scheme, the season with the maximum PV output $κ_{t}^{m a x}$ is most likely to exceed limits, as the voltage will be higher and the reverse power flow will be greater compared with other seasons. Thus, in the power flow model, it is enough to consider the profile $κ_{t}^{m a x}$ .

For the above power flow equations, the linearized network loss calculation formula based on the operating point $(P_{t}^{*}, Q_{t}^{*})$ is expressed as:

\begin{array}{l} P_{t}^{l o s s} (P_{t}, Q_{t}) \approx P_{t}^{l o s s} (P_{t}^{*}, Q_{t}^{*}) + \sum_{i = 1}^{N} \frac{\partial P_{t}^{l o s s}}{\partial P_{i, t}^{}} (P_{i, t}^{} - P_{i, t}^{*}) + \\ \sum_{i = 1}^{N} \frac{\partial P_{t}^{l o s s}}{\partial Q_{i, t}^{}} (Q_{i, t}^{} - Q_{i, t}^{*}) \forall t \end{array}

(6)

The detailed expressions are available in [

30]. In this paper, the default operating point

(P_{t}^{*}, Q_{t}^{*})

is the net load without PV integration.

B. PV Output Uncertainty Modeling

Fluctuations in PV outputs have an impact on the maximum PVIC. In this paper, a normalized random error $ξ$ is added in (4) to represent the actual PV output. Assuming that the historical sample size for the error $ξ$ is S, the data set $\hat{ξ}$ can be expressed as:

\hat{ξ} = \{{\hat{ξ}}_{1}^{}, {\hat{ξ}}_{2}^{}, \dots, {\hat{ξ}}_{s}^{}, \dots, {\hat{ξ}}_{S}^{}\}

(7)

The empirical distribution $\hat{P}$ from the historical data can be defined as:

\hat{P} = \frac{1}{S} \sum_{s = 1}^{S} δ (ξ - {\hat{ξ}}_{s}^{})

(8)

Assuming that the true distribution $P$ is not too far away from the empirical distribution $\hat{P}$ [

21], an ambiguity set

P

based on the Wasserstein metric can be constructed to express the potentials of true distribution. The Wasserstein distance

d^{W}

is used to measure the distance between

\hat{P}

and

P

. The details of definition of

d^{W} (\hat{P}, P)

can be found in [23]. Therefore, the ambiguity set

𝒫

can be defined by the Wasserstein ball of radius r centered at the empirical distribution

\hat{P}

𝒫 = \{P : d^{W} (\hat{P}, P) \leq r\}

(9)

The probability that $P$ is in the Wasserstein ball $𝒫$ is greater than confidence level $1 - ε^{W}$ .

P r \{d^{W} (\hat{P}, P) \leq r\} \geq 1 - ε^{W}

(10)

C. Operation Constraints of BESSs

- P_{i}^{d i s, m a x} \leq P_{i, t}^{B E S S} \leq P_{i}^{c h, m a x} \forall i, \forall t

(11)

E_{i}^{m i n} \leq E_{i}^{0} + \sum_{t_{1} = 1}^{t} P_{i, t_{1}}^{B E S S} Δ t \leq E_{i}^{m a x} \forall i, \forall t

(12)

\sum_{t = 1}^{T} P_{i, t}^{B E S S} Δ t = 0 \forall i

(13)

Constraint (11) limits the charging and discharging power of the BESS. Constraint (12) ensures that the residual energy of BESS meets the SOC constraints at every period throughout the day. Constraint (13) ensures that the initial residual energy of BESS for a day is equal to its ending residual energy, thereby ensuring the feasibility of BESS utilization in the following day. The above constraints delineate the feasible operation range of BESSs.

D. DRCC-based Network Security Constraints and Their Transformation

In this paper, considering the uncertainty of PV outputs based on Wasserstein ambiguity set $𝒫$ , the network security constraints are formulated as DRCCs to accommodate the worst-case distribution of the normalized random error $ξ$ . The security constraints indicate the limitation of bus voltages and branch flows [

14] as follows.

\underset{P \in 𝒫}{i n f} P r \{V_{i}^{m i n} \leq V_{i, t}^{} \leq V_{i}^{m a x}\} \geq 1 - ε^{C} \forall i, \forall t

(14)

\underset{P \in 𝒫}{i n f} P r \{- P_{i j}^{m a x} \leq P_{i j, t}^{} \leq P_{i j}^{m a x}\} \geq 1 - ε^{C} \forall i j, \forall t

(15)

The constraints for $Q_{i j, t}^{}$ are simplified since the variation of PVICs may not affect $Q_{i j, t}^{}$ . Based on the approximation method proposed in [

21], the DRCCs are formulated as:

\underset{P \in 𝒫}{i n f} P r \{a {(x)}^{T} ξ \leq b (x)\} \geq 1 - ε^{C}

(16)

It can be transformed into tractable linear constraints as:

\{\begin{array}{l} r v - ε^{C} γ \leq \frac{1}{S} \sum_{s = 1}^{S} \{b (x) - a (x) {\hat{ξ}}_{s} - γ\} \\ {‖a (x)‖}_{*} \leq v \\ v \geq 0, γ \geq 0 \end{array}

(17)

In this paper, Wasserstein distance defines on 1-norm.

Thus, bus voltage $V_{i, t}^{}$ and branch flow $P_{i j, t}^{}$ can be linearly expressed by $κ_{t}^{m a x} S_{i}^{P V} (1 + ξ)$ and $P_{i, t}^{B E S S}$ , respectively, and power flow equation (3) and DRCCs (14), (15) can be transformed into tractable linear constraints. The details are presented in Appendix A.

E. Operation Requirements

This subsection includes net-zero requirements and the network loss constraint. For net-zero requirements, given the significant variation in PV outputs across seasons, it is imperative to assess the net-zero state of the DS from an annual perspective, which needs to consider the annual energy generation of a PV integration scheme. Based on the existing literature on net-zero investigations, this paper models net-zero requirements from two distinct dimensions as follows.

1) Net-zero Energy

Some research works such as [

3] stated that “net-zero energy” represents a balance between the energy produced and consumed during a single year. The corresponding requirement is expressed as:

\sum_{s e} \sum_{i} \sum_{t} (P_{i, t}^{l o a d} - κ_{t, s e}^{} S_{i}^{P V}) \leq 0

(18)

2) Net-zero Carbon Emission

Similarly, for net-zero carbon emission DSs [

1], the requirement is expressed as:

\sum_{s e} \sum_{i} \sum_{t} π_{t, s e}^{c a r b o n} (P_{i, t}^{l o a d} - κ_{t, s e}^{} S_{i}^{P V}) \leq 0

(19)

For the network loss constraint, as PV integration capacity increases, the reverse power flow in the DS also increases, potentially leading to significant energy losses. To avoid excessive network losses that may cause inefficient operation of the DS, it is necessary to add network loss constraints is expressed as:

\sum_{t} P_{t}^{l o s s} Δ t \leq W_{}^{l o s s, r e q}

(20)

IV. Formulation and Solution of Feasible Region for PVICs

In this section, we summarize the equations and inequalities based on the DS operation model. Then, we use the concept of basic feasible solution to obtain the affine relationships between PVICs and constrained variables, which indicate that the variables with known boundaries (such as voltages) can be expressed by the linear combinations of PVICs. Therefore, the feasible region can be formulated from these linear combinations and their ranges. Finally, the solution procedure can be concluded.

A. Obtaining Affine Relations for PVICs Using Basic Feasible Solution Concept

The affine relations between PVICs $S_{}^{P V}$ and constrained variables (such as branch power $P_{t}^{}$ , voltage V_t, and BESS power $P_{t}^{B E S S}$ ) play a crucial role in determining the feasible region. Using these identity relations, we can derive the constraints that indicate the affine mappings of $S_{}^{P V}$ are constrained by the known boundaries of $P_{t}^{}$ , $V_{t}^{}$ , and $P_{t}^{B E S S}$ . Then, the feasible region of $S_{}^{P V}$ can be determined. Therefore, we need to obtain the affine relations.

For a given time period t, the original constraints in Section III can be summarized as follows.

1) For BESSs: (11). However, the bounds of (11) may be temporally coupled during different time periods, due to constraints (12) and (13). The treatment of this will be elaborated in Section V-B.

2) For network security: (3), (14), and (15). Note that (14) and (15) can be further transformed into tractable linear constraints and $ξ$ in (4) can be removed. The details are presented in Appendix A.

3) For net-zero operation requirements: constraints (18) and (19) can be rewritten as (21), (22) and (23), (24), respectively.

F_{t}^{e n e r g y} = \sum_{s e} \sum_{i} (P_{i, t, s e}^{l o a d} - κ_{s e, t}^{} S_{i}^{P V})

(21)

\sum_{t} F_{t}^{e n e r g y} \leq 0

(22)

F_{t}^{c a r b o n} = \sum_{s e} \sum_{i} π_{t, s e}^{c a r b o n} (P_{i, t, s e}^{l o a d} - κ_{s e, t}^{} S_{i}^{P V})

(23)

\sum_{t} F_{t}^{c a r b o n} \leq 0

(24)

Thus, equations for net-zero operation requirements are (21) and (23); and inequalities are given in (22) and (24).

4) For network loss: constraint (20) can be rewritten as (25) and (26).

F_{t}^{l o s s} = \sum_{t} P_{t}^{l o s s} Δ t

(25)

F_{t}^{l o s s} \leq W_{}^{l o s s, r e q}

(26)

Thus, equations for network loss constraint are (6) and (25); and inequality is (26).

It is assumed that $F_{t}^{e n e r g y}$ , $F_{t}^{c a r b o n}$ , and $F_{t}^{l o s s}$ during time period t have bounds ${\bar{F}}_{t}^{e n e r g y}$ , ${\bar{F}}_{t}^{e n e r g y}$ , and ${\bar{F}}_{t}^{l o s s}$ , respectively. Strictly speaking, these bounds are temporally coupled during different time periods according to (22), (24), and (26), which will be elaborated in Section V-B.

\{\begin{array}{l} F_{t}^{e n e r g y} \leq {\bar{F}}_{t}^{e n e r g y} \\ F_{t}^{c a r b o n} \leq {\bar{F}}_{t}^{c a r b o n} \\ F_{t}^{l o s s} \leq {\bar{F}}_{t}^{l o s s} \end{array}

(27)

The aforementioned model can be summarized as:

Ψ_{t} ψ_{t} = Θ_{t} S_{}^{P V} + τ_{t}

(28)

{\underset{̲}{ψ}}_{t} \leq ψ_{t} \leq {\bar{ψ}}_{t}

(29)

where $ψ_{t} = [P_{t}^{B E S S}, V_{t}^{}, P_{t}^{}, F_{t}^{e n e r g y}, F_{t}^{c a r b o n} o r F_{t}^{l o s s}]^{T}$ is the vector of variables with known boundaries, and $V_{t}$ and $P_{t}^{}$ are actually replaced by are auxiliary variables mentioned in Appendix A.

The number of (28) is usually smaller than the number of constrained variables in $ψ_{t}$ . Also, there may be collinear data in the calculated matrices $Ψ_{t}$ . Thus, the length of $ψ_{t}$ is usually greater than the rank of $Ψ_{t}$ , causing difficulty in obtaining affine relations.

To address these issues, inspired by [

17], we employ the concept of basic feasible solution (BFS) from linear programming theory [32] to obtain the affine relations between

ψ_{t}

and

S_{}^{P V}

. The affine relationship based on a BFS of the linear programming problem constrained by (28) and (29) can be expressed as:

ψ_{l, t}^{a *} = {(Ψ_{l, t}^{a})}^{- 1} (Θ_{l, t}^{a} S_{}^{P V} + τ_{l, t}^{a})

(30)

where $l \in Ω^{b}$ is the index of BFSs that satisfy the condition presented in Appendix B. The elements in $ψ_{l, t}^{a *}$ indicate a part of elements in $ψ_{t}$ . D-dimensional matrix $Ψ_{l, t}^{a}$ is a full-rank submatrix of matrix $Ψ_{t}$ . $D \times N^{P V}$ matrix $Θ_{l, t}^{a}$ and $D \times 1$ matrix $τ_{l, t}^{a}$ denote submatrices of $Θ_{t}$ and $τ_{t}$ , respectively. The detailed explanation of (30) can be found in Appendix B.

Based on (30) of a BFS, the affine relations between $ψ_{l, t}^{a *}$ and $S_{}^{P V}$ are derived. Therefore, we obtain the whole affine relations from all BFSs $Ω^{b}$ .

B. Formulation of Feasible Region for PVICs

Due to the temporal coupling of the variables in $ψ_{t}$ (such as $P_{t}^{B E S S}, F_{t}^{e n e r g y},$ $F_{t}^{c a r b o n}, F_{t}^{l o s s}$ ), the operation constraints for each $ψ_{c, t}$ in $ψ_{t} = [ψ_{1, t} ψ_{2, t} . . . ψ_{c, t} . . . ψ_{C, t}]$ need to be further reformulated as:

{[\begin{matrix} {\underset{̲}{ψ}}_{c, 1}^{c o u p} & {\underset{̲}{ψ}}_{c, 2}^{c o u p} & \dots & {\underset{̲}{ψ}}_{c, T}^{c o u p} \end{matrix}]}^{T} \leq L {[\begin{matrix} ψ_{c, 1}^{} & ψ_{c, 2}^{} & \dots & ψ_{c, T}^{} \end{matrix}]}^{T} \leq {[\begin{matrix} {\bar{ψ}}_{c, 1}^{c o u p} & {\bar{ψ}}_{c, 2}^{c o u p} & \dots & {\bar{ψ}}_{c, T}^{c o u p} \end{matrix}]}^{T}

(31)

The specific treatments of each kind of variables are discussed as follows.

For the temporal decoupling variables $ψ_{c, t}$ in $ψ_{t}$ (corresponding to $V_{i, t}^{}$ and $P_{i j, t}^{}$ ), the bounds in (31) are shown as follows. If there are only lower or upper bounds, virtual bounds based on the Big-M method can be useful.

\{\begin{array}{l} {[\begin{matrix} {\underset{̲}{ψ}}_{c, 1}^{c o u p} & {\underset{̲}{ψ}}_{c, 2}^{c o u p} & \dots & {\underset{̲}{ψ}}_{c, T}^{c o u p} \end{matrix}]}^{T} = L {[\begin{matrix} {\underset{̲}{ψ}}_{c, 1}^{1} & {\underset{̲}{ψ}}_{c, 2}^{2} & \dots & {\underset{̲}{ψ}}_{c, T}^{T} \end{matrix}]}^{T} \\ {[\begin{matrix} {\bar{ψ}}_{c, 1}^{c o u p} & {\bar{ψ}}_{c, 2}^{c o u p} & \dots & {\bar{ψ}}_{c, T}^{c o u p} \end{matrix}]}^{T} = L {[\begin{matrix} {\bar{ψ}}_{c, 1}^{1} & {\bar{ψ}}_{c, 2}^{2} & \dots & {\bar{ψ}}_{c, T}^{T} \end{matrix}]}^{T} \end{array}

(32)

Due to the temporal coupling characteristics of BESSs in (12) and (13), the bounds of variable $ψ_{c, t}$ (corresponding to BESS power $P_{i, t}^{B E S S}$ ) in (28) are as given by:

\{\begin{array}{l} {[\begin{matrix} {\underset{̲}{ψ}}_{c, 1}^{c o u p} & {\underset{̲}{ψ}}_{c, 2}^{c o u p} & \dots & {\underset{̲}{ψ}}_{c, T}^{c o u p} \end{matrix}]}^{T} = {[\begin{matrix} e & 0 \end{matrix}]}^{T} (E_{i}^{m i n} - E_{i}^{0}) / Δ t \\ {[\begin{matrix} {\bar{ψ}}_{c, 1}^{c o u p} & {\bar{ψ}}_{c, 2}^{c o u p} & \dots & {\bar{ψ}}_{c, T}^{c o u p} \end{matrix}]}^{T} = {[\begin{matrix} e & 0 \end{matrix}]}^{T} (E_{i}^{m a x} - E_{i}^{0}) / Δ t \end{array}

(33)

where $e = {[\begin{matrix} 1 & 1 & \dots & 1 \end{matrix}]}_{1 \times (T - 1)}$ .

Considering temporal coupling of (22), (24), and (26), the bounds of variable $ψ_{c, t}$ (corresponding to $F_{t}^{e n e r g y}$ , $F_{t}^{c a r b o n}$ , $F_{t}^{l o s s}$ ) in (31) are given by:

\{\begin{array}{l} {[\begin{matrix} {\underset{̲}{ψ}}_{c, 1}^{c o u p} & {\underset{̲}{ψ}}_{c, 2}^{c o u p} & \dots & {\underset{̲}{ψ}}_{c, T}^{c o u p} \end{matrix}]}^{T} = {[\begin{matrix} e & 1 \end{matrix}]}^{T} (- M) \\ {[\begin{matrix} {\bar{ψ}}_{c, 1}^{c o u p} & {\bar{ψ}}_{c, 2}^{c o u p} & \dots & {\bar{ψ}}_{c, T}^{c o u p} \end{matrix}]}^{T} = {[\begin{matrix} e & 0 \end{matrix}]}^{T} M \end{array}

(34)

Thus, for all time periods, the bounds for all $ψ_{c, t}$ during all time periods, named $ψ$ , can be summarized as:

{[\begin{matrix} \underset{̲}{ψ} & {\underset{̲}{ψ}}^{c o u p} \end{matrix}]}^{T} \leq Z ψ^{T} \leq {[\begin{matrix} \bar{ψ} & {\bar{ψ}}^{c o u p} \end{matrix}]}^{T}

(35)

where $ψ = [\begin{matrix} ψ_{1} & ψ_{2} & \dots & ψ_{T} \end{matrix}]$ , $Z = {[\begin{matrix} I^{T} & L^{T} \end{matrix}]}^{T}$ , $\underset{̲}{ψ} = [\begin{matrix} {\underset{̲}{ψ}}_{1} & {\underset{̲}{ψ}}_{2} & \dots & {\underset{̲}{ψ}}_{T} \end{matrix}]$ , $\bar{ψ} = [\begin{matrix} {\bar{ψ}}_{1}^{} & {\bar{ψ}}_{2}^{} & \dots & {\bar{ψ}}_{T}^{} \end{matrix}]$ , ${\underset{̲}{ψ}}_{c o u p} = [\begin{matrix} {\underset{̲}{ψ}}_{1}^{c o u p} & {\underset{̲}{ψ}}_{2}^{c o u p} & \dots & {\underset{̲}{ψ}}_{T}^{c o u p} \end{matrix}]$ , and ${\bar{ψ}}^{c o u p} = [\begin{matrix} {\bar{ψ}}_{1}^{c o u p} & {\bar{ψ}}_{2}^{c o u p} & \dots & {\bar{ψ}}_{T}^{c o u p} \end{matrix}]$ .

For a BFS $ψ_{l, t}^{*}$ mentioned in Section V-A, considering all time periods, $ψ_{l}^{a *} = {[ψ_{l, 1}^{a *} ψ_{l, 2}^{a *} \dots ψ_{l, T}^{a *}]}_{S \times T}$ corresponds to part of $ψ$ and meets the following constraint.

{[\begin{matrix} {\underset{̲}{ψ}}_{l}^{a} & {\underset{̲}{ψ}}_{l}^{c o u p, a} \end{matrix}]}^{T} \leq Z (ψ_{l}^{a *})^{T} \leq {[\begin{matrix} {\bar{ψ}}_{l}^{a} & {\bar{ψ}}_{l}^{c o u p, a} \end{matrix}]}^{T}

(36)

Therefore, the bounds for $(Ψ_{l, t}^{a})^{- 1} (Θ_{l, t}^{a} S_{}^{P V} + τ_{l, t}^{a})$ during all time periods are explicitly settled based on (36). The feasible region of PVICs derived from BFS l is described as:

Ω_{l} = \{S_{P V}^{T} |{\underset{̲}{γ}}_{l} \leq Z Y_{l} \leq {\bar{γ}}_{l}\}

(37)

where ${\underset{̲}{γ}}_{l} = {[\begin{matrix} {\underset{̲}{ψ}}_{l}^{a} & {\underset{̲}{ψ}}_{l}^{c o u p, a} \end{matrix}]}^{T}$ and ${\bar{γ}}_{l} = {[\begin{matrix} {\bar{ψ}}_{l}^{a} & {\bar{ψ}}_{l}^{c o u p, a} \end{matrix}]}^{T}$ are of order $2 T \times D$ , and $Y_{l} = [(Ψ_{l, 1}^{a})^{- 1} (Θ_{l, 1}^{a} S_{}^{P V} + τ_{l, 1}^{a}) (Ψ_{l, 2}^{a})^{- 1} (Θ_{l, 2}^{a} S_{}^{P V} + τ_{l, 2}^{a}) \dots (Ψ_{l, T}^{a})^{- 1} (Θ_{l, T}^{a} S_{}^{P V} + τ_{l, T}^{a} {)]}^{T}$ is of order $T \times D$ .

As a result, the obtained $Ω_{l}$ based on affine relation (30) is accurate and free from redundancy. Hence, the comprehensive feasible region for PVICs can be formulated as:

Ω^{P V} = \underset{l \in Ω^{b}}{\cup} Ω_{l} = \underset{l \in Ω^{b}}{\cup} \{{(S^{P V})}^{T} |{\underset{̲}{γ}}_{l} \leq Z Y_{l} \leq {\bar{γ}}_{l}\}

(38)

Note that the constraint ensuring the positivity of each $S_{P V, i}$ should be considered within the feasible region for PVICs.

Ω^{P V'} = Ω^{P V} ⋃ \{{(S^{P V})}^{T} |{(S^{P V})}^{T} \geq 0\}

(39)

C. Solution Procedure of Feasible Region for PVICs

The essence of calculating the feasible region for PVICs is obtaining all the affine relations (all BFSs). The specific solution procedure is as follows.

Step 1: input the DS data set and determine the operation requirements. Input the potential PV locations that need to be investigated.

Step 2: reformulate the model as (28) and (29).

Step 3: solve all BFSs from (28) and (29), obtaining all affine relationships between $ψ_{t}^{}$ and $S_{}^{P V}$ .

Step 4: obtain the upper/lower bounds of $ψ_{t}^{}$ during all time periods based on (35).

Step 5: construct the feasible region $Ω_{l}$ based on affine relations and bounds of the $l^{t h}$ BFS. Unite all $Ω_{l}, \forall l \in Ω^{b}$ , and obtain the complete feasible region for PVICs.

D. Verification and Adjustment Procedure to Satisfy Network Loss Constraint Under AC Power Flow

To ensure that the obtained PVIC feasible region satisfies network loss constraints under AC power flow constraints, we design the following verification and adjustment procedure.

Step 1: set iteration number $m = 1$ , and

W_{1}^{l o s s, r e q} = W_{}^{l o s s, r e q}

(40)

Step 2: set the net load point as the operating point, then obtain the feasible region under the linearized AC power flow model, considering the network loss constraint (41) and calculation formula (6) (Steps 1-5 in Section IV-C). Obtain a convex polyhedron, and determine all its vertices $\forall v_{m} \in Ω_{m}^{v e r t e x}$ .

\sum_{t} P_{t}^{l o s s} Δ t \leq W_{m}^{l o s s, r e q}

(41)

Step 3: for each vertex $v_{m}$ , with the corresponding PV integration scheme as the boundary condition, keep other system constraints unchanged, and solve the distributionally robust chance-constrained AC optimal power flow [

33] with network loss as the objective function. Obtain the accurate network loss value for each vertex

v_{m}

(PV integration scheme).

Step 4: verify whether the network loss value $W_{}^{l o s s, r e q} (v_{m})$ of each vertex $v_{m}$ satisfies the network loss constraint (41). If all vertices satisfy (41), the procedure ends. If there are vertices $v_{m}^{e x} \in Ω_{m}^{v e r t e x, e x}$ with network loss exceeding the limit and $m \leq m^{m a x}$ , calculate the adjustment amount of network loss constraint $Δ W_{m}^{l o s s, r e q}$ as:

Δ W_{}^{l o s s, r e q} (v_{m}) = W_{}^{l o s s, r e q} (v_{m}) - W_{m}^{l o s s, r e q} \forall v_{m}^{e x} \in Ω_{m}^{v e r t e x, e x}

(42)

Δ W_{m}^{l o s s, r e q} = m a x \{Δ W_{}^{l o s s, r e q} (v_{m}), \forall v_{m}^{e x}\}

(43)

Then, update the network loss constraint (41) as:

W_{m + 1}^{l o s s, r e q} = W_{m}^{l o s s, r e q} - Δ W_{m}^{l o s s, r e q}

(44)

Increase the iteration number $m = m + 1$ , and return to Step 2.

The verification and adjustment procedure described above is further shown in Fig. 2. This procedure may make the obtained feasible region slightly conservative, but it ensures that the feasible region satisfies all constraints under AC power flow.

Fig. 2 Verification and adjustment procedure to satisfy network loss constraints.

V. Case Study

In this section, the 4-bus system, the IEEE 33-bus system, and the IEEE 123-bus system are used to test the formulation and solution method of the feasible region for PVICs.

A. 4-bus System

The 4-bus system in Fig. 3, which is part of the IEEE 33-bus system, is used to validate the effectiveness of the proposed method for the feasible region for PVICs. There is no particular demand for integration locations and PV integrated at all locations in this system. Two BESSs (500 kW/1000 kWh for each BESS) are integrated at Bus 2 and Bus 3, respectively. The seasonal characteristics of PV outputs are from [

34]. The intraday carbon emission intensity profile is from [12], assuming that it does not change in different seasons.

Fig. 3 Configuration of 4-bus system.

The energy net-zero constraint is considered in the case study. We do not consider the network loss constraint at basic analysis section. For uncertainty consideration, to guarantee the reproducibility, we assume that the prediction error $ξ$ follows a normal distribution [

21], and the mean value is 0 and the standard deviation is set to be 10%. A data set with 100 samples is generated. The confidence level for the Wasserstein metric constraints

1 - ε^{W}

is set to be 90%, and the confidence level for the DRCCs

1 - ε^{C}

is set to be 95%. The PV panel loss percentage is set to be 6.5% based on [31].

All experiments are conducted on an Intel-i5 computer with 16 GB RAM and a basic frequency of 4.1 GHz using the MATLAB platform.

1) Basic Analysis

Figure 4 shows the analytical feasible region for PVICs proposed in this paper. It can be observed that, for the three PV integration locations, the analytical feasible region for PVICs forms a polyhedron in three-dimensional space, delineated by the red boundaries representing security constraints and the blue boundaries representing net-zero requirements. Specifically, the maximum feasible PVICs for the location PV1, PV2, and PV3 are 2360 kW, 1760 kW, and 957 kW, respectively. The net-zero constraint is satisfied when the total PV capacity of the three locations exceeds 827 kW.

Fig. 4 Analytical feasible region for PVICs in 4-bus system.

2) Comparison with Existing Method

Note that few existing methods can obtain the lower boundary of feasible region for PVICs due to the net-zero operation requirements.

For security constraints, the existing method [

9] is capable of determining the feasible region for PVICs, as depicted by the blue grid lines in Fig. 5. This method employs multi-objective optimization to obtain feasible region (the Pareto frontier). In the case study, the existing method yields 56 Pareto points. The resulting boundaries in Fig. 5 are obtained by a three-dimensional fitting algorithm. This method provides reliable solutions only at the Pareto points, and the accuracy of the region heavily depends on the accuracy of the fitting algorithm, which may be difficult to achieve in larger systems.

Fig. 5 Comparison with existing method.

In this paper, we propose a novel method based on the linear power flow to derive analytical PVIC boundaries under network security constraints, as depicted by the red lines in Fig. 5. The basic contours of the boundaries of the two feasible regions closely resemble each other, thus verifying the effectiveness of the proposed method.

In contrast to the existing method, the proposed method provides a relatively simple mathematical expression for the boundary and does not require the multi-objective optimization and high-dimensional fitting algorithm. It can be observed that the feasible region based on the proposed method is smaller than that based on multi-objective optimization. It may be due to the conservative strategy for network loss constraint addressing.

3) Sensitivity Analysis: BESS Allocation

Figure 6 provides a sensitivity analysis for BESS allocation by comparing between the scenario without BESSs and that with double BESSs. The BESS has a significant effect on the shape of the feasible region. However, the effect of BESSs on the maximum PVIC is non-linear due to the location of these resources and network constraints. In the case without BESSs, the security boundary intersects with the net-zero boundary, resulting in a substantial reduction in the size of the feasible region.

Fig. 6 Sensitivity analysis for BESS allocation.

It is worth noting that the proposed method is capable of capturing the boundary (blue line in Fig. 6) at the intersection of security boundary and net-zero boundary. The extended orange facet in Fig. 6, corresponding to the scenario without BESSs, is used to enhance the clarity of explanation.

4) Sensitivity Analysis: Confidence Level of DRCCs

Figure 7 provides a sensitivity analysis for the confidence level of DRCCs by comparing the feasible region at different confidence levels, corresponding to network security. A lower confidence level implies a relaxation of the security constraint requirements, resulting in an expanded range of the feasible region.

Fig. 7 Sensitivity analysis at different confidence levels of DRCCs.

5) Sensitivity Analysis: Net-zero Requirement

Considering energy net-zero and carbon net-zero, Fig. 8 provides a sensitivity analysis for net-zero requirement by comparing the feasible region boundaries corresponding to different net-zero constraints.

Fig. 8 Sensitivity analysis for net-zero requirement.

For the net-zero requirement in energy, the presence or absence of BESSs does not alter the respective boundaries, since BESSs do not modify the energy balance of a year.

If the net-zero energy constraint is imposed during each time period, it still does not change the boundaries, but the demand for BESSs significantly increases. By testing, net-zero energy is achieved for each time period when each BESS is configured at 680 kW/1360 kWh.

For the net-zero requirement in carbon emissions, the influence of the presence or absence of BESS on the boundary is significant due to the intraday variability in the carbon emission intensity. During the periods of high PV generation around midday, the carbon intensity is lower. BESS can facilitate the temporal transfer of energy, thereby reducing the required PV integration capacity for achieving net-zero carbon emissions.

6) Validation of Network Loss Constraint

This part examines the impact of the network loss constraint (20). In fact, as shown in Fig. 4, the vertex with the maximum network loss under AC power flow is Point A in Fig. 9, with a network loss of 27.7 kWh.

Fig. 9 Validation of network loss constraint.

We set the network loss constraint for a typical day to be 25 kWh (about 90% of 27.7 kWh) to test the effectiveness of the network loss constraint (20). After two iterations, the results are depicted in Fig. 9. Since the network loss constraint value only affects the upper boundary of the feasible region, the upper boundary is illustrated in Fig. 9.

It can be observed that the network loss constraint (20) results in a plane within the feasible region. The proposed method guarantees that for any PV integration scheme within the derived feasible region, the network losses under AC power flow conditions remain within the prescribed limits.

B. IEEE 33-bus System

To validate the scalability of the proposed method, testing is conducted on the IEEE 33-bus system, as shown in Fig. 10. To enhance the significance of the results, the loads at each bus are reduced by half. Four candidate PV integration locations are considered, along with the five BESSs rated at 2000 kW/4000 kWh. The network loss constraint for the typical day is 500 kWh. Specifically, two sets of candidate locations are selected: ① Set 1: Buses 10, 24, and 32; and ② Set 2: Buses 17, 24, and 32. The feasible regions for these two cases are plotted in Fig. 11(a) and (b), respectively.

Fig. 10 Configuration of IEEE 33-bus system.

Fig. 11 Feasible region for PVICs in IEEE 33-bus system. (a) Set 1. (b) Set 2.

According to Fig. 11, with the increase in the number of buses and constraints, the configuration of the feasible region becomes more complex. Comparing the two cases, it can be observed that the PV integration potential at Bus 17 is lower than that at Bus 10. Consequently, for Set 2, this results in the intersection of the security boundary and the net-zero boundary, leading to a significant reduction in the size of the feasible region compared with that of Set 1.

C. IEEE 123-bus System

For a more comprehensive validation, we use the IEEE 123-bus system [

35], as shown in Fig. 12, which operates at a nominal voltage of 4.16 kV. Buses 83, 96, 250, and 450 are chosen for PV integration and three BESSs rated at 1000 kW/2000 kWh are installed at Buses 13, 25, and 67, respectively. The network loss constraint for the typical day is set to be 1000 kWh.

Fig. 12 Configuration of IEEE 123-bus system.

Given that the feasible region for the four PV integration locations cannot be depicted in a three-dimensional space, we perform dimensionality reduction for ease of presentation. When 4000 kW PV is integrated at Bus 250, the feasible region for Buses 83, 96, and 450 is shown in Fig. 13.

Fig. 13 Feasible region for PVICs in IEEE 123-bus system.

It takes 4945 s to complete the entire process including solution, verification, and iteration on an Intel-i5 computer with 16 GB RAM. This time cost is acceptable for planning purposes.

VI. Conclusion

This paper proposes an analytical method to delineate the feasible region for PVICs in net-zero DSs. First, operation constraints and net-zero requirements that feasible PVICs need to satisfy are introduced. Then, the formulation and solution of the feasible region using BFS concept are proposed.

The case study validates the effectiveness of the proposed method. The DRCC-based network security constraints set the upper boundaries of the feasible region, while the energy and carbon-emission net-zero constraints establish the lower boundaries. These boundaries may intersect under certain parameters, and the configuration of the feasible region may become more complicated when the system is larger.

Sensitivity analysis for BESS allocation, confidence levels of DRCCs, and net-zero requirements are conducted, and the network loss constraint is validated. It illustrates that the proposed method fully captures the spatio-temporal coupling relationship among PVs, loads, and BESSs, while also quantifying the impact of this relationship on the boundaries of the feasible region for PVICs.

NOMENCLATURE

Symbol	——	Definition
A.	——	Indices
c	——	Index of variables in $ψ_{t}$
i, j, k	——	Indices of buses for distribution system
ij	——	Index of branches
l	——	Index of basic feasible solutions
m	——	Iteration number of loss constraint adjustment
se	——	Index of seasons
s	——	Index of samples
t	——	Index of time periods
B.	——	Sets
$Ω^{b}$	——	Set of basic feasible solutions
$Ω^{b u s}$	——	Set of buses of distribution system
$Ω_{l}$	——	Set of partial feasible regions based on basic feasible solution l
$Ω^{P V}$	——	Set of feasible regions for photovoltaic (PV) integration capacities
$Ω^{P V L}$	——	Set of locations of PV integration
$Ω_{m}^{v e r t e x}$	——	Set of vertices of a feasible region
$Ω_{m}^{v e r t e x, e x}$	——	Set of vertices that exceed loss constraint
$S$	——	Set of samples of uncertainties
$𝒫$	——	Ambiguity set of distributions
C.	——	Parameters and Functions
$δ$	——	Dirac delta function
$1 - ε^{C}$	——	Confidence level of distributionally robust chance constraint (DRCC)
$1 - ε^{W}$	——	Confidence level of Wasserstein ball
$κ_{t}^{m a x}$	——	Profile of season with the highest PV outputs
$κ_{t, s e}^{}$	——	Profile of PV outputs in different seasons
$μ_{t}$	——	PV panel loss percentage
${\hat{ξ}}_{s}^{}$	——	Historical sample of error of PV output
$π_{t, s e}^{c a r b o n}$	——	Locational carbon intensity
$Ψ_{t}$ , $Θ_{t}$ , $τ_{t}$	——	Parameter matrices and vectors
a(x), b(x)	——	Affine mappings of x
$B_{1}^{s u b}$ , $B_{2}^{s u b}$ , $B_{1}^{c o l}$ , $B_{2}^{c o l}$	——	Matrices of distribution system impedances
C	——	Number of variables in $ψ_{t}$
D	——	Dimension of full-rank submatrix of $Ψ_{t}$
$d^{W}$	——	Wasserstein distant metric
$E_{i}^{m i n}$ , $E_{i}^{m a x}$	——	The allowable lowest and highest states-of-charge (SOCs) of an energy storage system (ESS)
$E_{i}^{0}$	——	Initial SOC of an ESS
${\bar{F}}_{t}^{e n e r g y}$ , ${\bar{F}}_{t}^{c a r b o n}$ , ${\bar{F}}_{t}^{l o s s}$	——	Virtual boundaries of auxiliary variables
L	——	Lower triangular matrix
M	——	A number large enough
m^max	——	The maximum iteration number
$N$ , $N^{P V}$	——	Number of buses and PV locations
$P_{i}^{c h, m a x}$ , $P_{i}^{d i s, m a x}$	——	The maximum charging and discharging power of ESS
$P_{i, t}^{l o a d}, Q_{i, t}^{l o a d}$	——	Active and reactive load demands
$P_{i j}^{m a x}$	——	Limitation of branch power flows
$S$	——	Number of samples
T	——	Number of time periods
$Δ t$	——	Duration of a time period
$V^{m i n}, V^{m a x}$	——	Limitations of bus voltage
$W_{}^{l o s s, r e q}$	——	Network loss limitation of a typical day
${‖\cdot‖}_{*}$	——	Norm dual to that in definition of Wasserstein distance
D.	——	Variables
$v, γ, λ_{1}, λ_{2}, λ_{3}$	——	Auxiliary variables in transformed DRCC
$ξ$	——	Normalized random error of PV output
$ψ_{t}^{}$	——	Variables with known boundaries
$F_{t}^{e n e r g y}$ , $F_{t}^{c a r b o n}$ , $F_{t}^{l o s s}$	——	Auxiliary variables in transformed operation requirements
$P$ , $\hat{P}$	——	True and empirical distributions
$P_{i, t}^{}$ , $Q_{i, t}^{}$	——	Active and reactive power of bus injection
$P_{i j, t}^{}$ , $Q_{i j, t}^{}$	——	Active and reactive power of branch
$P_{i, t}^{B E S S}$	——	Power of battery energy storage system (BESS)
$P_{t}^{l o s s}$	——	Network loss during a time period corresponding to $κ_{t}^{m a x}$
$r$	——	Radius of Wasserstein ball
$S_{i}^{P V}$	——	PV integration capacity
$V_{i, t}^{}$ , $θ_{i, t}$	——	Amplitude and angle of bus voltage
$v_{m}$	——	Vertex of a feasible region
$v_{m}^{e x}$	——	Vertex that exceeds network loss constraint
$W_{m}^{l o s s, r e q}$	——	Network loss constraint in an iteration
$Δ W_{m}^{l o s s, r e q}$	——	Adjustment amount of network loss constraint in an iteration
x	——	Decision variable
$z_{t}^{1}, z_{t}^{2}, z_{t},$	——	Auxiliary variables in application of basic feasible solution (BFS) theory
Jt,Kt

Appendix

Appendix A

Appendix A shows the transformation of DRCCs. For a time period t, power flow equation (3) can be reformulated as matrix forms.

[\begin{matrix} θ_{t} \\ V_{t} \end{matrix}] = {[\begin{matrix} Β_{2}^{s u b} & Β_{1}^{s u b} \\ - Β_{1}^{s u b} & Β_{2}^{s u b} \end{matrix}]}^{- 1} ([\begin{matrix} P_{t} \\ Q_{t} \end{matrix}] - [\begin{matrix} Β_{2}^{c o l} \\ - Β_{1}^{c o l} \end{matrix}] θ_{1, t} - [\begin{matrix} Β_{1}^{c o l} \\ Β_{2}^{c o l} \end{matrix}] V_{1, t}) \forall t

(A1)

It can be further formulated as:

\begin{array}{l} [\begin{matrix} θ_{t} \\ V_{t} \end{matrix}] = [\begin{matrix} D_{11} & D_{12} \\ D_{21} & D_{22} \end{matrix}] \cdot \\ ([\begin{matrix} P_{t}^{l o a d} + (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} (1 + ξ) + P_{t}^{B E S S} \\ Q_{t}^{l o a d} \end{matrix}] - [\begin{matrix} C_{1} \\ C_{2} \end{matrix}]) \end{array}

(A2)

For a bus,

\begin{array}{l} V_{i, t}^{} = D_{21} (i) [P_{t}^{l o a d} + (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} (1 + ξ) + P_{t}^{B E S S}] + \\ D_{22} (i) Q_{t}^{l o a d} - C_{2} (i) \end{array}

(A3)

For a branch [

29],

\begin{array}{l} P_{i j, t}^{} = \frac{r_{i j} x_{i j}}{r_{i j}^{2} + x_{i j}^{2}} \frac{V_{i, t} - V_{j, t}}{x_{i j}} + \frac{x_{i j}^{2}}{r_{i j}^{2} + x_{i j}^{2}} \frac{θ_{i, t} - θ_{j, t}}{x_{i j}} = A_{1} (i j) V_{t} + A_{2} (i j) θ_{t} = \\ (A_{1} (i j) D_{21} + A_{2} (i j) D_{11}) [P_{t}^{l o a d} + (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} (1 + ξ) + P_{t}^{B E S S}] + (A_{1} (i j) D_{22} + A_{2} (i j) D_{12}) Q_{t}^{l o a d} - (A_{1} (i j) C_{2} + A_{2} (i j) C_{1}) = E_{1} (i j) [P_{t}^{l o a d} + (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} (1 + ξ) + P_{t}^{B E S S}] + \\ E_{2} (i j) \end{array}

(A4)

where $A_{1} (i j)$ , $A_{2} (i j)$ , and $E_{1} (i j)$ are row vectors； and $E_{2} (i j)$ is the element of vector.

Then, for bus i, the corresponding DRCC (14) can be reformulated as (A5). Following the formulation of (16), (A5) can be concluded as (A6). Then, for branch ij, the corresponding DRCC (15) can be reformulated as (A7).

\begin{array}{l} \underset{P \in 𝒫}{i n f} P r \{D_{21} (i) [P_{t}^{l o a d} + (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} (1 + ξ) + P_{t}^{B E S S}] + \\ D_{22} (i) Q_{t}^{l o a d} - C_{2} (i) \leq V_{i}^{m a x}, - [D_{21} (i) \cdot \\ (P_{t}^{l o a d} + (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} (1 + ξ) + P_{t}^{B E S S}) + \\ D_{22} (i) Q_{t}^{l o a d} - C_{2} (i)] \leq - V_{i}^{m i n}\} \geq 1 - ε^{C} \forall i, \forall t \end{array}

(A5)

where $D_{21} (i)$ and $D_{22} (i)$ are row vectors from $D_{21}$ and $D_{22}$ , respectively.

\{\begin{array}{l} a (x) = D_{21} (i) (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} \\ b_{1} (x) = V_{i}^{m a x} - D_{21} (i) [P_{t}^{l o a d} + (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} + P_{t}^{B E S S}] - \\ D_{22} (i) Q_{t}^{l o a d} + C_{2} (i) \\ b_{2} (x) = - V_{i}^{m i n} + D_{21} (i) P [P_{t}^{l o a d} + (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} + P_{t}^{B E S S}] + \\ D_{22} (i) Q_{t}^{l o a d} - C_{2} (i) \end{array}

(A6)

\begin{array}{l} \underset{P \in 𝒫}{i n f} P r \{E_{1} (i j) [P_{t}^{l o a d} + (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} (1 + ξ) + P_{t}^{B E S S}] + \\ E_{2} (i j) \leq P_{i j}^{m a x}, - [E_{1} (i j) (P_{t}^{l o a d} + (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} (1 + ξ) + \\ P_{t}^{B E S S}) + E_{2} (i j)] \leq P_{i j}^{m a x}\} \geq 1 - ε^{C} \forall i j, \forall t \end{array}

(A7)

Following the formulation of (16), (A7) can be concluded:

\{\begin{array}{l} a (x) = E_{1} (i j) (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} \\ b_{1} (x) = P_{i j}^{m a x} - \{E_{1} (i j) [P_{t}^{l o a d} + (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} + P_{t}^{B E S S}] + \\ \overset{}{E_{2} (i j)}\} \\ b_{2} (x) = P_{i j}^{m a x} + \{E_{1} (i j) [P_{t}^{l o a d} + (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V} + P_{t}^{B E S S}] + \\ \overset{}{E_{2} (i j)}\} \end{array}

(A8)

Based on the conclusion (17) in [

21], the reformulated DRCCs (A5) and (A7) can be transformed as:

\{\begin{array}{l} r v - ε^{C} γ = λ_{1} + \frac{1}{S} \sum_{s = 1}^{S} (b_{1} (x) - a (x) {\hat{ξ}}_{s}^{} - γ) \\ r v - ε^{C} γ = λ_{2} + \frac{1}{S} \sum_{s = 1}^{S} (b_{2} (x) - a (x) {\hat{ξ}}_{s}^{} - γ) \\ {‖a (x)‖}_{*} = v + λ_{3} \\ v, γ, λ_{1}, λ_{2}, λ_{3} \geq 0 \end{array}

(A9)

Note that ${‖a (x)‖}_{*}$ in (17) denotes $\underset{i}{m a x} (|a {(x)}_{i}|)$ for the 1-norm DRCC, and the following equations hold in this paper.

For (A5),

{‖a (x)‖}_{*} = D_{21} (i) (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V}

(A10)

For (A7),

{‖a (x)‖}_{*} = E_{1} (i j) (1 - μ_{t}) κ_{t}^{m a x} S_{}^{P V}

(A11)

In conclusion, the original DRCC formulation is transformed into a set of linear constraints. Specifically, the bus voltage, injection power, and branch power variables are replaced by auxiliary variables. The original DRCCs (14), (15) and power flow equation (3) are transformed into linear equations and inequalities (A9) involving the auxiliary variables.

Appendix B

Appendix B shows the application of basic feasible solution concept. Firstly, introduce auxiliary variables as:

\{\begin{array}{l} z_{t}^{1} = ψ_{t} - {\underset{̲}{ψ}}_{t} \\ z_{t}^{2} = {\bar{ψ}}_{t} - ψ_{t} \end{array}

(B1)

Then, (28) and (29) can be reformulated as:

\{\begin{array}{l} J_{t} z_{t} = K_{t} \\ z_{t} \geq 0 \end{array}

(B2)

\{\begin{array}{l} J_{t} = [\begin{matrix} Ψ_{t} & 0 \\ I & I \end{matrix}] \\ z_{t} = [\begin{matrix} z_{t}^{1} \\ z_{t}^{2} \end{matrix}] \\ K_{t} = [\begin{matrix} Θ_{t} S_{}^{P V} + τ_{t} - Ψ_{t} {\underset{̲}{ψ}}_{t} \\ {\bar{ψ}}_{t} - {\underset{̲}{ψ}}_{t} \end{matrix}] \end{array}

(B3)

A BFS l of the linear programming problem, constrained by (B2), can be expressed as:

\{\begin{array}{l} z_{l, t}^{a *} = {(J_{l, t}^{a})}^{- 1} K_{l, t}^{a} \\ z_{l, t}^{b *} = 0 \end{array}

(B4)

The elements in $z_{l, t}^{a *}$ and $z_{l, t}^{b *}$ comprise all elements in $z^{t}$ , indicating a BFS $z_{l, t}^{*}$ . If $z_{l, t}^{a *}$ is all from $z_{t}^{1}$ , the first equation in (B4) can be further expressed as:

\begin{array}{l} ψ_{l, t}^{a *} - {\underset{̲}{ψ}}_{l, t}^{a} = {(Ψ_{l, t}^{a})}^{- 1} (Θ_{l, t}^{a} S_{}^{P V} + τ_{l, t}^{a} - Ψ_{l, t}^{a} {\underset{̲}{ψ}}_{l, t}^{a}) = \\ {(Ψ_{l, t}^{a})}^{- 1} (Θ_{l, t}^{a} S_{}^{P V} + τ_{l, t}^{a}) - {\underset{̲}{ψ}}_{l, t}^{a} \end{array}

(B5)

If $z_{l, t}^{a *}$ is not entirely from $z_{t}^{1}$ , no significant conclusions can be derived.

Thus, a group of affine relationships (B5) for BFS l is obtained. While the affine relationships combine the ranges in (29) under the other equations in (B4), the feasible region for BFS l can be constructed.

References

J. Wu, Y. Zhou, and W. Gan, “Smart local energy systems towards net zero: practice and implications from the UK,” CSEE Journal of Power and Energy Systems, vol. 9, no. 2, pp. 411-419, Mar. 2023. [Baidu Scholar]

H. Mehrjerdi, R. Hemmati, M. Shafie-Khah et al., “Zero energy building by multicarrier energy systems including hydro, wind, solar, and hydrogen,” IEEE Transactions on Industrial Informatics, vol. 17, no. 8, pp. 5474-5484, Aug. 2021. [Baidu Scholar]

Y. Yamaguchi, Y. Shoda, S. Yoshizawa et al., “Feasibility assessment of net zero-energy transformation of building stock using integrated synthetic population, building stock, and power distribution network framework,” Applied Energy, vol. 333, p. 120568, Mar. 2023. [Baidu Scholar]

H. Wu, Y. Yuan, J. Zhu et al., “Assessment model for distributed wind generation hosting capacity considering complex spatial correlations,” Journal of Modern Power Systems and Clean Energy, vol. 10, no. 5, pp. 1194-1206, Sept. 2022. [Baidu Scholar]

S. Zhang, Y. Fang, H. Zhang et al., “Maximum hosting capacity of photovoltaic generation in SOP-based power distribution network integrated with electric vehicles,” IEEE Transactions on Industrial Informatics, vol. 18, no. 11, pp. 8213-8224, Nov. 2022. [Baidu Scholar]

M. J. Chihota, B. Bekker, and T. Gaunt, “A stochastic analytic-probabilistic approach to distributed generation hosting capacity evaluation of active feeders,” International Journal of Electrical Power & Energy Systems, vol. 136, p. 107598, Mar. 2022. [Baidu Scholar]

S. Wang, Y. Dong, Q. Zhao et al., “Bi-level multi-objective joint planning of distribution networks considering uncertainties,” Journal of Modern Power Systems and Clean Energy, vol. 10, no. 6, pp. 1599-1613, Nov. 2022. [Baidu Scholar]

J. Xiao, G. Zu, H. Zhou et al., “Total quadrant security region for active distribution network with high penetration of distributed generation,” Journal of Modern Power Systems and Clean Energy, vol. 9, no. 1, pp. 128-137, Jan. 2021. [Baidu Scholar]

S. Zhang, S. Ge, H. Liu et al, “Model and observation of the feasible region for PV integration capacity considering Wasserstein-distance-based distributionally robust chance constraints,” Applied Energy, vol. 347, p. 121312, Oct. 2023. [Baidu Scholar]

X. Jiang, Y. Zhou, W. Ming et al., “Feasible operation region of an electricity distribution network,” Applied Energy, vol. 331, p. 120419, Feb. 2023. [Baidu Scholar]

H. Li, N. Zhang, Y. Fan et al., “Decomposed modeling of controllable and uncontrollable components in power systems with high penetration of renewable energies,” Journal of Modern Power Systems and Clean Energy, vol. 10, no. 5, pp. 1164-1173, Sept. 2022. [Baidu Scholar]

J. Li, Z. Xu, H. Liu et al., “A Wasserstein distributionally robust planning model for renewable sources and energy storage systems under multiple uncertainties,” IEEE Transactions on Sustainable Energy, vol. 14, no. 3, pp. 1346-1356, Jul. 2023. [Baidu Scholar]

S. Wang and W. Wu, “Aggregate flexibility of virtual power plants with temporal coupling constraints,” IEEE Transactions on Smart Grid, vol. 12, no. 6, pp. 5043-5051, Nov. 2021. [Baidu Scholar]

Y. Wen, Z. Hu, and L. Liu, “Aggregate temporally coupled power flexibility of DERs considering distribution system security constraints,” IEEE Transactions on Power Systems, vol. 38, no. 4, pp. 3884-3896, Jul. 2023. [Baidu Scholar]

W. Wei, F. Liu, and S. Mei, “Real-time dispatchability of bulk power systems with volatile renewable generations,” IEEE Transactions on Sustainable Energy, vol. 6, no. 3, pp. 738-747, Jul. 2015. [Baidu Scholar]

W. Wei, Y. Chen, C. Wang et al, “Nodal flexibility requirements for tackling renewable power fluctuations,” IEEE Transactions on Power Systems, vol. 36, no. 4, pp. 3227-3237, Jul. 2021. [Baidu Scholar]

S. Zhang, S. Ge, H. Liu et al., “Region-based flexibility quantification in distribution systems: an analytical approach considering spatio-temporal coupling,” Applied Energy, vol. 355, p. 122175, Feb. 2024. [Baidu Scholar]

R. Leng, Z. Li, and Y. Xu, “Two-stage stochastic programming for coordinated operation of distributed energy resources in unbalanced active distribution networks with diverse correlated uncertainties,” Journal of Modern Power Systems and Clean Energy, vol. 11, no. 1, pp. 120-131, Jan. 2023. [Baidu Scholar]

J. Wang and Y. Song, “Distributionally robust OPF in distribution network considering CVaR-averse voltage security,” International Journal of Electrical Power & Energy Systems, vol. 145, p. 108624, Feb. 2023. [Baidu Scholar]

Y. Zhou, Z. Wei, M. Shahidehpour et al., “Distributionally robust resilient operation of integrated energy systems using moment and Wasserstein metric for contingencies,” IEEE Transactions on Power Systems, vol. 36, no. 4, pp. 3574-3584, Jul. 2021. [Baidu Scholar]

A. Zhou, M. Yang, M. Wang et al., “A linear programming approximation of distributionally robust chance-constrained dispatch with Wasserstein distance,” IEEE Transactions on Power Systems, vol. 35, no. 5, pp. 3366-3377, Sept. 2020. [Baidu Scholar]

P. Li, M. Yang, and Q. Wu, “Confidence interval based distributionally robust real-time economic dispatch approach considering wind power accommodation risk,” IEEE Transactions on Sustainable Energy, vol. 12, no. 1, pp. 58-69, Jan. 2021. [Baidu Scholar]

W. Xie, “On distributionally robust chance constrained programs with Wasserstein distance,” Mathematical Programming, vol. 186, pp. 115-155, Nov. 2019. [Baidu Scholar]

D. Lee, K. Turitsyn, D. K. Molzahn et al., “Robust AC optimal power flow with robust convex restriction,” IEEE Transactions on Power Systems, vol. 36, no. 6, pp. 4953-4966, Nov. 2021. [Baidu Scholar]

N. Nazir and M. Almassalkhi, “Grid-aware aggregation and realtime disaggregation of distributed energy resources in radial networks,” IEEE Transactions on Power Systems, vol. 37, no. 3, pp. 1706-1717, May 2022. [Baidu Scholar]

Y. Liu, L. Wu, Y. Chen et al., “On accurate and compact model of high DER-penetrated sub-transmission/primary distribution systems in ISO energy market,” IEEE Transactions on Sustainable Energy, vol. 12, no. 2, pp. 810-820, Apr. 2021. [Baidu Scholar]

Y. Liu, L. Wu, Y. Chen et al., “Integrating high DER-penetrated distribution systems into ISO energy market clearing: a feasible region projection approach,” IEEE Transactions on Power Systems, vol. 36, no. 3, pp. 2262-2272, May 2021. [Baidu Scholar]

X. Chen, E. D. Anese, C. Zhao et al., “Aggregate power flexibility in unbalanced distribution systems,” IEEE Transactions on Smart Grid, vol. 11, no. 1, pp. 258-269, Jan. 2020. [Baidu Scholar]

X. Chen and N. Li, “Leveraging two-stage adaptive robust optimization for power flexibility aggregation,” IEEE Transactions on Smart Grid, vol. 12, no. 5, pp. 3954-3965, Sept. 2021. [Baidu Scholar]

H. Yuan, F. Li, Y. Wei et al., “Novel linearized power flow and linearized opf models for active distribution networks with application in distribution LMP,” IEEE Transactions on Smart Grid, vol. 9, no. 1, pp. 438-448, Jan. 2018. [Baidu Scholar]

Solar Empower. (2023, Jan.). Solar PV system losses – how to calculate solar panel efficiency [Online]. Available: https://www.solarempower.com/blog/10-solar-pv-system-losses-their-impact-on-solar-panel-output/ [Baidu Scholar]

G. Bernd and M. Jiří, Understanding and Using Linear Programming, Berlin: Springer, 2006. [Baidu Scholar]

A. Arab and J. E. Tate, “Distributionally robust optimal power flow via ellipsoidal approximation,” IEEE Transactions on Power Systems, vol. 38, no. 5, pp. 4826-4839, Sept. 2023. [Baidu Scholar]

H. M. A. Ahmed, A. B. Eltantawy, and M. M. A. Salama, “A planning approach for the network configuration of AC-DC hybrid distribution systems,” IEEE Transactions on Smart Grid, vol. 9, no. 3, pp. 2203-2213, May 2018. [Baidu Scholar]

J. Zhao, M. Zhang, H. Yu et al., “An islanding partition method of active distribution networks based on chance-constrained programming,” Applied Energy, vol. 242, pp. 78-91, May 2019. [Baidu Scholar]

Address:No.19 Chengxin Avenue, Jiangning District, Nanjing 211106, China

E-mail: mpce@alljournals.cn

Tel:86-25-81093060

Fax:86-25-81093040

Home

Introduction

Editorial Board

For Author

Call For Papers

APC

Sponsor & Publisher