Distributionally Robust Decentralized Scheduling Between the Transmission Market and Local Energy Hubs

The cross impacts between transmission and distribution systems have drawn extensive attention, where multi-energy carriers become increasingly dominant in local market operations. Local energy hubs, integrated with multi-energy carriers, e.g., electricity, gas, heat, present great potentials to provide adequate power and reserve support for the transmission system, which could mitigate issues such as boundary mismatches in the upstream market. However, energy hubs are typically equipped with abundant distributed renewable resources. Such local uncertainties can adversely affect the well-being of the coordinated multi-energy market hierarchy. This paper presents a novel multi-period scheduling framework that considers the coordination between the transmission network and distribution energy hubs. Each agent performs local scheduling operations capturing independent uncertainties via a distributionally robust formulation. We then apply a tailored accelerated augmented Lagrangian algorithm embedded with the column-and-constraint method to decentralize the overall operation with agents' privacy preserved. The fast and convergent feature of the algorithm ensures its scalability and reliability in real-world applications of the energy hub design with transmission coordination. Numerical experiments confirm the efficacy of the proposed method.

Mixture temperature at return node.

I. INTRODUCTION
E LECTRICITY markets have been in a transitional stage that follows from a top-down hierarchy to a bottom-up one. As the local distribution networks and microgrids become increasingly active, they attain a growing capability to provide both power and ancillary service support to the upstream operator. While the transmission system could not overlook the local impacts, distribution system operators (DSOs) are emerging to seek cooperation with the conventional independent system operator (ISOs) for better managing a reliable and efficient power grid [1]. The variable distributed energy resource (DER) plays a critical role in this transition with governmental incentives and regulations, e.g., FERC order 2222 [2]. Nevertheless, substantial operational uncertainties in DERs such as wind and solar variations shape a bottleneck when considering them in the market operation, especially with the transmission and distribution (T&D) coordination. The current research has conducted thorough investigations in the impacts of variable DERs in the power sector, e.g., [3]. However, another non-negligible transition, i.e., the proliferation of multi-energy carriers in the local energy network, is turning the power distribution system into a multi-energy hub. Its coordination with the T&D market hierarchy has not been sufficiently studied.
Existing literature has devoted considerable efforts to studying T&D coordination. Designing the future DSO frameworks gains prevalence since the ISO operation has been well structured and developed. As summarized in [4], there are three stages of the DSO design philosophy and its coordination scheme with the ISO, i.e., centralized, decentralized, and transactive models. From first to last, the system of systems becomes more flexible and distributed. Due to the natural compatibility with the existing market framework, based on the distribution utility model (Model 3) in [4], our previous work [5] focuses on an ISO-leading unit commitment (UC) problem capturing the individual system uncertainties. The ISO still dominates the T&D market in this model, while the DSO performs as a utility managing the local assets. In this paper, we extend that work to investigate the case that both the ISO and DSO act as independent entities based on Model 4 in [4], where the ISO cannot easily dominate the coordinated market with no vision and control on DSOs. Entities only exchange boundary information in the physical coupling node and operate independent market operations. There will no longer be an integrated UC problem as there is no central operator, but a multi-objective and multi-agent optimization problem, which could be formulated as a distributed economic dispatch (ED). We refer to [6] for deeper discussions on decentralized T&D cooperation. The local energy system in the distribution level, however, does not only serve electric power. The proliferating energy hub in smart distribution systems houses multi-energy carriers, within which the electric power, natural gas, and heat energies are the most common and substantial [7]. It is natural that the DSO takes the multi-energy carriers' operations into account when coordinating with the power transmission system since other energy networks exert impacts on both system dynamics and techno-economic values of the power distribution system. However, the characteristics of each energy carrier are different and complex. For instance, the electric power flow in the power distribution system and the gas pipeline flow in the gas distribution system inflict nonlinearity and intractability to the coordinated market model. In our paper, to retain convexity, we propose using the second-order cone programming (SOCP) to model the nonlinear AC power flow and the gas flow equations while keeping the heat network linear via the variable-temperature constant-flow (VTCF) model. Note that it is normal to consider multi-energy carriers in the distribution system to be operated by a single entity, i.e., the DSO. Many utilities manage and operate multi-energy assets, including power and gas, such as PG&E in California and ComEd in Illinois.
Operational uncertainties, e.g., variable DER outputs, pose another hurdle in the T&D coordination. Though many works focus on the deterministic T&D operation due to the relatively small variation in the distribution network [8], [20], [21], it has been confirmed that a stochastic optimization-based model yield greater economic values and better capture the operation changes [10], [13]. Various stochastic optimization techniques have been widely applied in T&D coordination, among which the most common ones are scenario-based stochastic programming (SP) [5] and robust optimization (RO) [9]. The issue of the curse of dimensionality arises when the number of scenarios increases to capture enough uncertain patterns in the SP, and the over-conservatism of the RO beclouds its applicability in normal system operations. In this context, distributionally robust optimization (DRO) attracts extensive attention to achieve a  I  COMPARISON BETWEEN THIS PAPER AND STATE-OF-THE-ART STUDIES tradeoff between conservativeness and uncertainty realization, which outperforms SP and RO as a general-use stochastic optimization framework [3]. It models the uncertainty by assuming an ambiguity set that approximates its true probability distribution, which enjoys broad applications in the literature. C. He et al. [22] adopted DRO modeling a power-gas energy hub operation considering electric and gas load uncertainties. D. Mohammadreza et al. [23] proposed a two-stage transactive energy framework using DRO to capture renewable DER uncertainties. P. Li et al. [10] built a DRO chance-constrained distributed T&D coordination problem with wind uncertainties.
The resulting T&D problem with energy hubs is highly intractable due to the multi-scale T&D system characteristics and the stochastic nature of DERs. Meanwhile, when the ISO and DSO act as independent entities, they should execute individual scheduling problems with the privacy of system information preserved during the operation. Hence, centralized decomposition techniques such as Benders decomposition [15] are not applicable. Instead, ISOs and DSOs require distributed algorithms to decentralize communications such as the well-known augmented Lagrangian relaxation (LR) method. It helps each agent formulate the individual subproblem and only updates Lagrangian multipliers, i.e., gradient information, in each iteration, protecting agents' confidentiality. A recent enhancement, i.e., the accelerated augmented Lagrangian (AAL) [24], becomes a more efficient alternative than conventional LR-based methods such as the alternating direction of multipliers method (ADMM). Though originally devised for solving nonconvex problems, the efficacy of the AAL has also been validated in distributed convex applications with better performance than ADMM [25]. To solve the DRO-based subproblems for the ISO and DSOs, we propose using the widely-adopted column-and-constrained generation (C&CG) method [26], which is embedded in the AAL procedure without jeopardizing the convergence.
Though the existing literature has conducted thorough research in the T&D operation and distribution-level energy hubs, to the best of the authors' knowledge, this is the first work investigating the impacts of multi-energy hubs in the T&D electricity market. We also extend the application of AAL to the DRO, which yields better performance and economic outcomes than deterministic [25] and RO [13] applications. We conduct a comprehensive comparison between this study and other recent works in Table I. Furthermore, we summarize the threefold contributions of this paper as follows.
r We propose a novel multi-period scheduling framework for coordination between the transmission system and local energy hubs. The energy hub considers electric, gas, and heat energies with explicit network and resource modeling. The quantitative analyses demonstrate the cost-effectiveness of considering multi-energy hubs in T&D coordination.
r We implement a T&D market paradigm based on decentralized operations, which serves as a potential future T&D market structure as discussed in [4]. Each agent considers individualized uncertainties such as renewable generation and multi-energy demand without impacting others.
r We tailor the AAL algorithm with the C&CG method and validates its applicability to a DRO-based multi-agent distributed operation. The proposed method meets the requirements for solving a stochastic T&D coordination problem, e.g., computational efficiency, convergence guarantee, and privacy protection. We organize the rest of the paper as follows. Section II provides a detailed problem formulation for the transmission system operation and the local energy hub scheduling while also introducing the ambiguity set construction for each system's uncertainties; Section III describes the decentralization between agents and the overall solution procedure of the distributionally robust T&D problem based on the AAL algorithm; Section IV demonstrates the effectiveness of the proposed framework using numerical experiments; Section V concludes the paper with several technical remarks.

II. PROBLEM FORMULATION OF THE DRO-BASED MULTI-ENERGY T&D COORDINATION
In this section, we detail the mathematical formulation of the transmission system and local energy hubs, with uncertainties involved in the variable renewable generation. A two-agent co-optimization describes the T&D coordination with different agent-specific modeling details. We consider the upstream agent problem as a conventional economic dispatch carried out for the transmission network with industry-calibrated DC power flow. The problem for downstream agents is an SOCP-based multi-energy hub optimization considering nonlinear network equations for more accurate AC power and gas pipeline flows. The heat network brings a linearized VTCF model. Both agents follow a two-stage distributionally robust optimization formulation with individualized ambiguity sets for uncertainties.

A. Upstream Agent Problem: Transmission Market
The upstream agent solves a transmission ED problem capturing utility-level wind uncertainties. We write a canonical form of the formulation as shown in (1). (1) follows a typical day-ahead ED model, including the power balance equation (1b), DC power flow equation (1c), flow capacity constraint (1d), ramping constraints (1e)-(1f), power capacity constraint (1g), upward/downward reserve limit constraint (1h), reserve requirements (1i), and renewable curtailment equation (1j). We slightly abuse the notation of subscript T in sets and variables for readability. We also defer the ESS model in the transmission system to Section III. C with a generalized formulation. Note that we adopt the DC power flow to increase the proposed model's compatibility with the current ISO's practice [27]. In the objective function (1a), the first-stage cost contains the power generation cost and energy storage system (ESS) degradation cost, and the renewable curtailment cost and load shedding cost remain in the second stage. It is natural to regard the first stage as a normal market clearing and the second stage as an ex post adjustment during the operation. Though we utilize linear cost coefficients in the objective for showcase, it could adapt to more practical designs such as the quadratic cost curve and more sophisticated ESS degradation processes, which is out of the scope of this study and requires unique treatments in the decentralization process. In this case, the upstream agent features two-stage stochastic linear programming.

B. Downstream Agent Problem: Multi-Energy Hub
The downstream agent solves a distribution ED problem capturing DER-level uncertainties. However, the synergy between electric power, gas, and heat energies far complicates the problem than the upstream agent.
denote the production costs (including power generation cost, gas procurement cost, and combined heat and power (CHP) generation cost), multi-energy ESS degradation costs, DER curtailment costs, and multi-energy load shedding costs, respectively. Note that the first stage and second stage follow the same structure as in the upstream agent with a similar ambiguity set construction, which will be detailed in Section III.D. Based on different energy systems' characteristics, we sectionalize the constraint space for electric power, natural gas, and heat systems as follows. Besides, note that the convexification of multi-energy network flow models is necessary for finding a convergent solution, yet retains exactness, which has been corroborated in various studies [5], [28], [29].

1) Electric Power Distribution Network:
We cast the SOCPbased AC branch flow model [5] in the power distribution system with generation constraints.
2) Natural Gas Distribution Network: We adopt the SOCPbased Weymouth equation in the gas network model.

3) Heat Distribution Network:
We consider a linear VTCF heat network model with guaranteed accuracy [29] as follows.
with the heat pipeline loss defined as . The heat network model follows a typical VTCF formulation constraining the thermal energy conversion with temperature (5a)-(5b), confluence nodal temperatures (5c)-(5d), mass flow balance (5e)-(5f), and heat inlet/outlet flow-temperature (5g)-(5h). We perform a heuristic two-step hydraulic-thermal decomposition to fix the mass flow rate M and hence reduce the heat network model to a linear formulation, which reportedly maintains accuracy. We refer to [29] for detailed descriptions of the heuristic and the derivation of (5). Note that, unlike the power and gas system, formulating the heat system dynamics in (5) maintains the computational tractability while providing a more accurate representation of the network. Besides, the acceptably accurate representation and excellent solvability of the VTCF model underly the motivation for using it in the proposed coordinated energy hub analysis. For the CHP operations, we adopt the convex polyhedron to model CHP heat and power curves as discussed in our previous work [30], which retains model convexity and is omitted here.

C. Energy Storage Operations
As the proposed framework considers the power ESS (PESS) in transmission and distribution systems and the gas ESS (GESS) in distribution systems, we discuss the ESS model with a general formulation, i.e., p which includes ESS power equations (6a), charge/discharge limits (6b)-(6c), state-of-charge (SOC) equations (6d), SOC limits (6e), and SOC consistency requirements (6f). Note that it is not necessary to include the binary exclusiveness for battery charging/discharging since it could be achieved by the efficiency implementation and ESS operation price [31].

D. Ambiguity Set Construction
The uncertainty studied in this work includes the variable renewable energy, while it could seamlessly extend to other uncertainties such as load and price. We employ the Wasserstein distance for the ambiguity set construction. Note that it works for both transmission and distribution systems. Consider the uncertain renewable outputP r,h has a series of historical observations {P 1 r,h ,P 2 r,h , . . .P E r,h }. Given the compact support space Ξ = {P r,h ∈ R + }, whose σ-algebra contains the true and empirical renewable generation probability distributions, i.e., P , P E ∈ P(Ξ), we define the Wasserstein metric as where Π denotes the joint probability distribution. Then we construct the ambiguity set as where φ(E) is a tunable Wasserstein ball radius and P E serves as the centroid of A. It is crucial that φ(E) controls the conservativeness of the ambiguity set and we follow an improved criterion for its selection [32], i.e., where κ is the confidence level, for which we typically use 95% in this study, and C is a scalar obtained by using the bisection search method [32] to solve where M (P E r,h ) denotes the mean value of observations.

III. T&D DECENTRALIZATION AND SOLUTION STRATEGY
The proposed distributionally robust T&D problem indicates a multi-agent and multi-objective optimization, which is hard to solve as each agent presents a two-stage DRO problem with independent ambiguity sets for uncertainties. In this section, we delineate how we solve the individual subproblem via the tailored C&CG method and decentralize the overall coordination problem using the AAL algorithm.

A. DRO Subproblem Reformulation and Solution
As linear programming is a special case of SOCP, we give a compact form for both (1) and (2) with uncertainty ω and without T&D superscripts as follows.
∀T, ∀D : min subject to where (7b) and (7d) constrain the first-and second-stage problems. We first reformulate the Wasserstein metric as where P (dω) denotes the probability distribution when ω falls in the centroid of A, given the conditional distribution Π(dω, dω E ) = E{P (dω)}. Based on (8), we rewrite the second stage with the ambiguity set as a semi-infinite program As problem (9) is always bounded with finite observations of ω, it is equivalent to dualize (9a) and replace the second stage problem, which yields subject to Q(y, ω) = min y d y , where the auxiliary variable δ ω is introduced to bound the dual norm [33]. We could efficiently solve problem (10) via the C&CG algorithm tailored in [34]. Note that we place the SOCP constraints in the first stage to generate an exact lower-bounding affine cut from the second stage, which accelerates the optimal searching of the algorithm.

B. Decentralized T&D Information Exchange
After each DRO-based subproblem is solved, we present a more general form for the revisited master problem in the C&CG method upon convergence as ∀T, ∀D : min subject to where x C denotes the common variables shared in the T&D coordination, e.g., power injection at the boundary node. Consider a centralized model including all agents' operations, subject to where Ψ describes the T&D feasible region and (12b) unifies the common variables. It is straightforward to write the individual Lagrangian functions for both agents as where λ D defines the Lagrangian multipliers for the common variable shared by each distribution network and ν is the secondorder penalty step size. Then following Fig. 1, we perform the AAL algorithm and solve each agent's problem in parallel. The merits of AAL include the double updates in one iteration and the acceleration brought by using ι in both updates for primal variables and Lagrangian multipliers, which outperforms other optimality-conditioned methods [24]. Similar to [13], the subproblem's solution process and the AAL are independent. However, the AAL's convergence depends on the convergence of the inner C&CG process, which is guaranteed by model convexity and linear recourse problems with finite support. This feature protects each agent's privacy as the solution processes do not cross impact with only limited information, e.g., boundary power injection, being exchanged.

IV. NUMERICAL EXPERIMENTS
In this section, we validate the effectiveness of the proposed market paradigm and solution strategy via quantitative analyses on a standard test system. We carry out case studies on the IEEE 6-bus transmission system connecting to two multi-energy hubs based on the IEEE 13-bus distribution system, whose data  could be retrieved from [35], [36], [37]. Fig. 2 illustrates the test system's topology. We solve all of the SOCP subproblems by Gurobi 9.1 on a Windows PC with quad-core Intel i7-6700 CPU and 8 GB of RAM.
For the renewable uncertainty considered in the proposed model, we normalize the historical data from a utility-level wind turbine [38], and two DER-level photovoltaic panels [39]. Then we construct the ambiguity set following the procedure described in Section II. D. for each distribution network. Note that for different distribution systems, the ambiguity set is independent since the training datasets from the historical observations are i.i.d. For the AAL algorithm, we set ε 1 , ε 2 = 0.01 and select ι = 0.25 and ν = 2.2, fine-tuned by computational trials.

A. Cost-Effectiveness of Considering Multi-Energy Hubs
We first validate the cost-effectiveness of considering multienergy hubs by comparing the following cases.
r Case 1: T&D coordination with multi-energy hubs r Case 2: T&D coordination with only power distribution r Case 3: Disconnected T&D operations For Case 2, we do not model the gas and heat distribution networks, while gas-fired generators ag 1 , ag 2 use the normalized gas fuel prices, and CHP is assumed as a linear-cost generator. For Case 3, we perform the DRO-based multi-period ED for both transmission and distribution systems independently using fixed power exchange obtained by results in Case 2. Table II tabulates the results of the three cases regarding the economic performance of using different coordination strategies. We use the finalized locational marginal prices in the boundary node from the transmission side upon convergence to calculate the power trade revenue, which means the downstream operator could benefit from exporting power to the upstream operator and vice versa. We also record the average power mismatch as the active power difference in the boundary node upon convergence.
As we compare Case 3 with Case 1 and Case 2, it is straightforward to find the disconnected T&D operation is far less costeffective than the T&D coordination. The renewable curtailment  also raises if the two agents perform individual optimization because the excess power from uncertain outputs in different regions could not support each other. If we compare Case 1 with Case 2, considering the multi-energy hub in the T&D coordination yields higher profits since the explicit network modeling reduces the potentially inaccurate cost estimation of other energy sectors. Moreover, taking GESSs into account contributes to better scheduling of the gas units. For the network loss, the reason why Case 1 outperforms other cases lies in the higher capacity factor of gas units and lower power imports from the transmission system. They both reduce the total electrical distance from power sources to load buses in the distribution network.
Figs. 3 and 4 depict the detailed dispatch results for Case 1 and Case 2. The gas units' dispatch is highly sensitive to the gas retail price. For example, at Hour 3-8, the gas price is the lowest in one day, making the gas units receive sufficient fuels and generate more power to charge the PESS while exporting power to the transmission grid. At noon, the PV generation increases and helps the ESS store more energy to feed the peak load at night. By comparing Case 1 and Case 2, the lack of explicit gas network modeling lowers the capability of the gas generation because no GESS would help with the gas arbitrage. Moreover, the gas units could not sufficiently capture the price fluctuation without the compressor model. Hence, the weak gas units' generation capability forces the downstream operator to reduce the power export and reap lower revenues. As the CHP is the only heat source in the heat network, it needs to maintain the power level to ensure sufficient heat supplies, which is common in residential districts [36].

B. Comparison Between Uncertainty Models
To further demonstrate the advantage of using the DRO-based uncertainty model, we test cases with other uncertainty models, including Case 4 (deterministic), Case 5 (SP), and Case 6 (RO). For Case 4, we use the K-means centroid value of historical observations to represent the uncertainty. For Case 5, we apply the Monte-Carlo sampling to generate 20 equiprobable scenarios and reduce (10) to a single-stage stochastic conic program. For Case 6, we replace the ambiguity set with the polyhedral uncertainty set, [34] which is constructed using the maximum and minimum values of the renewable forecast, i.e., . 5 illustrates the comparison. Case 4 serves as the base case since the problem is solved with a widely applied but inaccurate characterization of uncertainties. We could find the SP in Case 5 yields the most cost-effective result among the three uncertainty models, which is also closest to the base case, but with more iterations required for the AAL convergence. In each AAL iteration, the Lagrangian multiplier is updated as an expected value of all scenarios' outcomes, which hinders the optimal search. RO in Case 6 achieves the fastest convergence rate but with highly conservative solutions, as it always draws the worst-case scenario in which the renewable generation is limited. The DRO-based model in Case 1, comparatively, becomes the most competitive choice as it balances the tradeoff between computational speed and solution conservatism. Albeit solving subproblems with different models requires different times, the AAL algorithm parallelizes the subproblem solution procedures, making the iteration counter comparable.

C. SOCP Exactness Evaluation
Our work adopts the SOCP relaxation in both the power distribution and gas distribution systems as the branch flow model and the relaxed Weymouth equation. Note that the SOCP-based branch flow model has enjoyed wide recognition for its exactness in representing the actual AC power flow [40]. Similarly, the Weymouth equation can also be exactly relaxed via the conic formulation, as reported in [41]. To further validate the exactness of the SOCP approximation, Table III shows the maximum SOCP gaps for both the power and gas distribution networks.
As shown in Table III, the SOCP gap values are sufficiently small for both networks, which proves the proposed methods' exactness for approximating the DRO model's energy flows. It is also worth mentioning that this approximation cannot completely substitute actual energy flow analyses but could serve as a satisfactory starting setpoint when running energy flow calculations with improved computational speed.

D. Decentralization Convergence Analysis
Though the AAL algorithm converges very quickly, its convergence pattern is slightly different from other centralized or decentralized methods. We are interested in investigating the AAL convergence by comparing it with a centralized method, Benders decomposition (BD), and another decentralized method, ADMM. To simplify the BD, we reduce the DRO-based subproblem to the deterministic equivalent form, and hence (12) becomes a typical two-stage SP problem. Fig. 6 reports the evolution of one boundary power injection variable when h = 5 under the three algorithms. We could observe that AAL crosses the optimal solution several times because it utilizes two primal and Lagrangian updates in each iteration to facilitate the solution, which asks for more penalizing measures. While ADMM also has this pattern, its convergence is slower than AAL. The centralized approach of BD never crosses the optimal solution since the bounding cuts are optimally conditioned, which guarantees the global optimality albeit a low convergence rate.

E. Scalability Test
In this subsection, we provide an additional test on a largescale system to further validate the scalability of the proposed method. The test system is a T118E20 system consisting of one IEEE 118-bus system with high renewable penetration [42] and 20 duplicates of the multi-energy hub from Fig. 2. We place the  20 energy hubs in the original load buses and keep the renewable generators with the same dataset in Section IV to generate the ambiguity set. Fig. 7 presents the performance of the proposed method in the large-scale test case. Since the proposed method is a distributed algorithm, parallel computing is a natural way to maximize resource usage for solving the problem. Each transmission and distribution operator can solve their scheduling problem locally and then communicate the Lagrangian information. Hence, we record the highest computational time per iteration that one agent needs to take to solve the local problem. We also report the objective value deviation from the optimally converged result to demonstrate the convergence progress. From Fig. 7, we observe that by leveraging parallel computing, solving large-scale T&D scheduling problems is efficient at the minute level due to the fast convergence rate of the AAL method.

V. CONCLUDING REMARKS
This paper proposes a novel scheduling framework as a market paradigm for the T&D coordinative operation. The proliferating penetration of renewable generation and the increasingly closer connection with other energy sectors create more uncertainties to the distribution side, propagating system-level power mismatch to the transmission side and leading to suboptimal decision-making. We consider the multi-energy hub in the distribution system to fill this gap with state-of-the-art modeling of uncertainties in both upstream and downstream operators. Furthermore, the nature of T&D coordination requires privacy-preserving operations between multiple agents, while the AAL-based workflow shows its efficacy for handling such a problem. The numerical experiments confirm the effectiveness of the proposed modeling and solution strategies.
Future studies include considering a more realistic and comprehensive market structure, such as the UC and balancing market in the decentralized T&D scheme. It is also interesting to further investigate whether the multi-energy hub as a distribution system could help the transmission operator mitigate network congestions and contribute to carbon neutrality.