Contingencies-Based Distributionally Robust Co-Risk Operation for Combined Electricity and Heat System

Sustainable development of combined electricity and heat system (CEHS) can effectively facilitate the energy transition from fossil energy consumption to comprehensive utilization of multiple energy. The coupled energy system, featured with the risk of components failure and the integrated uncertain renewable energy sources (RESs), can face a safety and reliability operation challenge. Therefore, achieving the tradeoff between the economy and risk aversion for the CEHS operation is an urgent issue to be settled. This paper proposes a contingencies-based distributionally robust co-risk operation model (DRCROM) for the CEHS, which provides the failure risk sample of components applied by the non-sequenced Monte Carlo simulation (MCS) method and the penalty of energy spillage (including wind spilling and load shedding) under contingencies. Moreover, the uncertain power output risk of RESs can be addressed by using the distributionally robust individual chance-constrained (DRICC) methodology with Wasserstein metric ambiguity set. As a result, incorporating the uncertain component state and power output of RESs into the risk operation can be considered in this paper. Case studies can be conducted to give risk analyses that the average violation probability <inline-formula> <tex-math notation="LaTeX">$V_{p}$ </tex-math></inline-formula> and the total co-risk dispatch cost <inline-formula> <tex-math notation="LaTeX">$C^{\text {Obj}}$ </tex-math></inline-formula> are mainly affected by <inline-formula> <tex-math notation="LaTeX">$\varepsilon $ </tex-math></inline-formula>. Contingencies can also increase risk indices and influence the dispatch energy and energy spillage.


I. INTRODUCTION
The cross-sectoral integration of energy systems is a current research area in the energy transition, which has manifested great potential in enhancing energy performance with economy, efficiency, carbon reduction, sustainability, and reliability [1]- [3]. Considering the increasing heat demand in some urban with long winters, the CEHS can promote the heat source transform from fossil fuel to electrical energy while enhancing the interaction of different energy systems and the integration of renewable energy sources (RESs), which is served as the important trend of the future energy system. However, the uncertain RESs integration and coupled multi-energy network with numerous components may bring The associate editor coordinating the review of this manuscript and approving it for publication was Wencong Su . challenges of co-risk operation [4]- [6]. Hence, the co-risk operation model for the CEHS should be developed to tackle the uncertain RESs integration and components outage.
The rapid development of the CEHS can accommodate higher penetration of renewable energy, which is feathered as variability and uncertainty in nature [7]. However, the seamless integration of renewable energy may produce the operation risk due to the uncertain power output, enhancing the risk propagation from the electricity subsystem to the heating subsystem. The uncertain RESs can be addressed by the existing approach, such as stochastic program (SP), robust optimization (RO), and distributionally robust optimization (DRO). The SP assumes a probability distribution followed by the uncertain parameter, which is hard to acquire the accurate distribution from the given historical data in practice [8]- [10]. RO mainly focuses on the realization of VOLUME 9, 2021 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ the worst-case uncertainties without probability distributions under consideration. Nevertheless, it may yield overly conservative solutions [11]- [13]. The DRO can be modeled by constructing the moment-based ambiguity set with the underlying distribution of uncertainty [14]- [16]. Then, the DRO can be solved by incorporating a tractable approximation method into the model. However, the moment-based DRO approach may just apply the moment information obtained from the history data, which cannot be used efficiently. Chance constrained programming (CCP) can guarantee a low probability of constraint violation [17]- [21]. However, the CCP emerges the poor conservation without the capture of the magnitude of violations. The DRCC incorporates the chance constraints into the DRO model to hedge against the uncertainty of probability distributions and mitigate the violation risk of overlimit chance constraints, which provides the risk constraints holding with a specified risk probability level [22]- [25]. Nevertheless, the DRCC approach cannot provide specific risk indices. The components of the CEHS would face the challenge of outage risk considering its lifetime. The analysis of components failure can be studied from two aspects, deterministic component failure [26]- [31], stochastic component failure [5], [32]. The operation state for the deterministic component failure is predefined. By investigating the contingencies, risk aversion strategies can be developed by the decision-maker in advance. Alternatively, the risk dispatch strategies are formed when the contingencies have arisen. In [31], considering the N-1 deterministic contingency of crucial devices, a planning model for an electricity and hydrogen integrated energy system (EH-IES) is proposed considering the system reliability. The uncertain component failure is investigated by mainly using the MCS approach. In [5], a novel risk assessment for the integrated electricity and heat system is proposed under contingency, in which the post-contingency is simulated by using the MCS method.
In summary, the researches above have presented more solutions for the optimal economical operation incorporating the uncertain power output and component outage into the CEHS. However, considering the operation risks brought by both the uncertain RESs integration and component failure is not focused on. Therefore, the CEHSCROM can be developed in this paper. The goal of the proposed model is to hedge against contingencies with the impacts of components failure and the risk of the uncertain power output of RESs. Compared with the existing work, the main contribution of this work can be summarized as follows: (1) The proposed CEHSCROM is built upon the contingencies-based distributionally robust co-risk optimization within the Wasserstein metric framework, which the metric-based Wasserstein ambiguity set can guarantee a comparatively stronger out-of-sample performance. The CEHS faces the challenge of operation risk, which originates from two aspects, the uncertain RESs, and the component failure. The DRCC can be employed to mitigate the operational risk caused by uncertain RESs, achieving the tradeoff between the economy and risk. The decision-maker can adjust the risk coefficient to satisfy the risk constraints. Meanwhile, the MCS is used to simulate the CEHS component failure, which is incorporated into the risk constraints.
(2) Given that the chance constraints are hard to be simultaneously violated for a larger time scale, the risk violation coefficient set in the different risk chance constraints should be diverse. Therefore, the uncertain RESs addressed can employ the distributionally robust individual chanceconstrained (DRICC) to balance the economy and risk. Since the proposed model is considered as the known NP-hard problem, the tractable conservative approximation is derived. The individual worst-case chance constraints can be approximated by the worst-case condition value at risk (WC-CVaR), which the proposed model can be reformulated into the tractable safe approximation with better conservation.
(3) To quantify the risk level of the CEHS, the risk indices are present to describe the operation risk of the CEHS. The violation probability V p is proposed to describe the violation probability of chance constraints. The risk indice C Obj can quantify the co-risk cost for the CEHS. The determinants of risk indices can be given in the risk analyses.

II. CONTINGENCIES-BASED STOCHASTIC CO-RISK OPERATION TOWARDS CEHS A. MCS-BASED COMPONENT STATE SAMPLING METHOD
The CEHS consists of numerous devices, which can bring operation risk challenges for stochastic operation outrage. Therefore, the risk of component stochastic failure operation should be considered. In this section, the impact of contingencies caused by the component outrage is mainly focused on. Two state models, normal operation and outrage operation is depicted as follows, Assume that the CEHS is made up of m components in total. S c represents the components state. S c = 1 indicates that the c− th component is a normal operation state while S c = 0 indicates that the c− th component is an outrage operation state. θ c represents the random number of c components under the uniform distribution. ϑ c depicts failure probability of c− th component.
In this work, the outrage of each power branch, the power source, the heat source, and the heat network should be considered. Moreover, the state of the CEHS is the integration of the state of each component. The non-sequence Monte Carlo Simulation is employed to analyze the impact of component failure under contingencies.

B. STOCHASTIC RISK OPERATION MODEL FOR CEHS
The CEHS comprises the power system (PS) and heating system (HS) [33], [34]. The PS consists of B u buses, B r branches, Z wind power generators, and De end-user electric demand. The HS is made up of P pipelines, N nodes, and D h heating demand. E energy hubs (EH) that comprise the CHP and HP are introduced to the CEHS as the coupled devices. They can produce electricity and heat energy by consuming natural gas, as shown in Fig. 1. The electric power source P G ∈ R E+Z is produced by the wind power generators P W (ξ ) ∈ R Z and the CHP P CHP ∈ R E while the heat source H is supplied by the CHP H CHP ∈ R E and heat pump H HP ∈ R E . In this work, it is assumed that P L is given. Given the uncertain nature of wind power output, the wind power can be modeled by P W (ξ ) = W (µ + ξ ), where µ ∈ R Z is defined as the mean power production and ξ ∈ R Z is served as the random variable with the mean µ 0 ∈ R Z and covariance matrix ν 0 ∈ R Z ×Z . W ∈ R Z ×Z depicts the diagonal matrix of the rated capacity of stochastic wind power generators. The upward and downward reserve margins R CHP,U and R CHP,D are chosen for each CHP. The stochastic risk operation for the CEHS is stated as the following linear program model. It is noted that the adjustment power is described as the expected value while the risk aversion of wind spilling and load shedding is expressed as the expected penalties value, Here, the objection function (2) minimizes the total risk operation cost, which is the sum of the CHP operation cost, the worst expected adjustment operation cost, and the weighted worst expected risk cost of the wind spilling and loading. C CHP , C CHP , and C CHP is the generation cost, upward/downward reserve cost of CHP, respectively. Power flow with radial electricity network topologies described by linearized DistFlow branch model has been employed in distribution network [15], [16]. Equations (3)-(7) represent the power system constraints. Constraints (3)-(5) enforce the nodal balance of active and reactive power, respectively. Constraint (6) describes the relationship between the voltage and the branch power flows. Constraint (7) provides the voltage bound. Equations (8)-(12) represent the mathematical model constraints for the district heating network. Constraints (8) depict the supply and return pipelines for the nodal heat power, respectively. Constraint (9) states that the inevitable heat loss leads to the temperature drop along the pipeline, constraint (10) defines the mixture temperature at confluence nodes, and constraint (11) enforces the same temperature for the mixed fluid at the confluence node. Equations (13)- (14) represent the CHP capacity limits. Equations (16)-(17) enforce the electric and heat power balance in the EH. Constraints (15), (19) are chance constraints of the adjustment power of CHP CHP , H HP for the EH. Constraint (18) defines the electric power capacity bound for the EH. Other variables can be refered to the paper [15], [16].
Given that the exponential term in the temperature loss constraint shown in Equation (9) is non-convex, which can bring a great challenge for the solution of the proposed model. Thus, the temperature loss constraint should be approximately reformulated to solve the optimization model tractably. Utilizing the approximation that e −x = 1 − x and substitute the expression (9) into the heat loss equation (6), the heat loss can be approximately expressed as, The aforementioned stochastic risk operation model is intractable since it includes an infinite-dimensional optimization problem. To mitigate the complexity, the functional form of decision rules can be restricted to be linear. In detail, the functional form of recourse decisions y (ξ ) can be approximated by the linear decision rules of the form y (ξ ) = Y ξ [35], [36]. Taking constraint (3) with probability 1 for example, two deterministic constraints can be obtained by matching the zero-and first-order coefficients of ξ on both sides of (3). Thus, the above semi-infinite problem with the decision rule approximation can be re-expressed as the following form, The outrage of the heat pipeline will restrict the heat delivery in the pipeline and will limit the heat production from the heat source. Thus, both contingencies in heating and power distribution networks are considered in the risk operation. The modeling of contingencies in heating and power distribution networks can be presented as follows.

1) CONTINGENCIES FOR THE HS
To incorporate contingencies into the proposed model, the component state indices S i can be introduced to describe the device state of the district heating network. In this work, the state indices S i can be obtained by using the MCS method. The pipeline constraint can be reformulated as, where M i for the i-th heat pipeline constraints is defined as the ''big M'' value that is large enough to make the constraint nonbinding. When S PL n = 0 stands for the outrage state of the heat pipeline, the constraints (17) (26) where N MCS represents the sampling size.

2) CONTINGENCIES FOR THE PS
Similarly, a component state S i can be defined for the power distribution network. The branch constraints can be rewritten as, Meanwhile, the outrage of CHP and EHs is considered. Therefore, the constraints of the outrage of the components in PS are expressed as follows,

D. RISK INDICES FOR CEHS
In this section, risk indices are presented to describe the risk assessment of CEHS. The total risk cost C Obj for the CEHS can be depicted as the risk indice, which is expressed as follows, where y (p CHP , P W ,S , P L,S ), T N means the number of sampling, C obj j indicates the j− th cost of the CEHS co-risk operation cost. The quality of the solutions is measured by the empirical out of sample cost.
To quantitatively reflect the risk level, the average violation probability V p is defined as follows, where T N means the number of sampling, A j (y) ξ > b j (x (y)) indicates violating chance constraints. The risk indice V p can quantify the risk level.
The solution framework of the proposed approach for the CEHS can be depicted in Fig. 2, which provides the detailed process of solving the risk optimization model.

III. WASSERSTEIN DISTRIBUTIONALLY ROBUST SOLUTION METHODOLOGY
For notation convenience, a general form of the model problem can be compactly reformulated as, where x P CHP , R CHP,U , R CHP,D , V , T S n , T R n , P HP i , p gas , P in , P EH , H EH , c C CHP , C CHP,U , C CHP,D , d C W , ωC P , is defined. can represent the set of all x, y , which can match the constraints (3)-(14) [15], [16].

A. DR-BASED MODEL FORMULATION
Considering the uncertain wind power output, it is hard to obtain the exact probability distribution estimated from the historical wind power data. Thus, the model problem will lack the fundamental input. To achieve the solution, a remedy that replaces the unknown P with the discrete empirical distribution P N can be applied. Therefore, the Wasserstein distributionally robust methodology in this paper can be used to tackle the infinite dimensionality of the problem, which can hedge against all distributions in the neighborhood of P N by incorporating the Wasserstein metric. The type-1 Wasserstein distance between P and P N on R Z can be defined via, for all distributions P, P N ∈M ( ) while is defined by the polytope := ξ ∈R Z : H ξ ≤ h where represents the joint distribution of ξ 1 and ξ 2 with marginals P and P N , which can be quantified as an optimal transportation cost that minimizes the cost from one distribution described by P to another on described by P N . The norm · encodes the transportation cost. Then, the ambiguity set can be defined by using the Wasserstein metric, which can be viewed as the Wasserstein ball of radius ρ centered at the empirical distribution P N on the training dataset N . The Wasserstein distributionally robust form of the proposed model problem (31)-(32) can be recast as which minimizes the worst-case expected risk operation cost and satisfies the constraints for all distributions in the ambiguityset . ε j reflects the risk coefficient of individual chance constraints.

B. WORST-CASE TRACTABLE REFORMULATION
The reformulation of the worst-case expected loss of the proposed model problem is considered [35], [36]. The objective function (35) can be derived as follows, where · stands for the dual norm. λ represents the nonnegative scalar variable used in the reformulation of the supremum of the expected term. γ and s are auxiliary variables for ease of notation. The proposed optimization problem can be reformulated into a tractable program. Moreover, two conservative approximations can be given for the feasible set, of a generic joint chance constraint of the form (36). The joint chance constraint (36) can be divided into K linear inequalities, and the matrix A( ) and the vector b(x) can be decomposed as The set of individual violation tolerances k ≥0, k ≤ K is given and K k=1 k = is set. Thus, the Bonferroni's inequality can be exploited to split joint chance constraint up into a family of K simpler but more conservative individual chance constraints, which is approximated as follows [35], B can be conservatively approximated by using the worstcase CVaR method, which is illustrated as follows, In the last, the conic reformulation for the set BC is admitted in Appendix [35].

IV. CASE STUDY A. CASE MODEL DESCRIPTION
The case study can be conducted using the proposed CEHS model that composes the power and heat systems. The PS is based on the modified IEEE 9-bus power system, while the HS is modified from the 32-node Barry Island heating system [15], [38], as illustrated in Fig. 3. Three EHs are installed in the CEHS to connect the PS and HS, which serve as the coupled components to achieve energy exchange between different energy systems. The PS is a 12.66kV radial network with nine buses, five power sources, and 32 branches. The power sources contain two wind generators integrated into the PS and three CHPs installed in the EHs. Moreover, the heat sources in HS are composed of three CHPs and HPs, which merge the EHs. The RESs power can be based on the hourly wind power measurements provided in the 2020 European Energy Market (EEM20) forecasting competition [39]. The wind spilling cost is 1000 $/MWh, and the load shedding cost is 500 $/MWh.

B. COMPARISON WITH SAA AND RO METHODOLOGY
In thissection, the co-risk cost analyses are implemented to investigate the dispatch performance using the proposed approach, compared with the RO and Sample Average Approximation (SAA) approaches. Firstly, the appropriate violation probability ε = 0.01 and the Wasserstein radius    ρ = 0.0001 are set. The dispatch energy cost C D, E , the energy spillage cost C S , and the total cost C obj can be obtained by using the three approaches. The computational results for the three approaches are summarized in Table 2.
As seen from Table 2, the DRCC significantly outperforms the SAA method in the energy spillage cost C S (including the wind spilling and load shedding), achieving a 32.53% reduction. Similarly, the total cost C obj can reduce 26.85% using the DRCC method, while the dispatch energy cost C D, S reduces a little. The main reason is the incorporated distribution information of wind power output, and the less conservative SAA causes the poor dispatch result. Meanwhile, compared with the RO method, the 22.96% energy spillage cost C S and the 18.9% total cost C obj using the DRCC can be reduced, respectively. It owes to the over-conservative nature of the RO method, which leads to a broader range of the stochastic variable.

C. OUT OF SAMPLE PERFORMANCE
In this section, the out of sample co-risk analyses are implemented to investigate the sample performance of the proposed approach. The component outrage sample can be obtained by using the MCS. The proposed method can be employed to solve the CEHS co-risk operation problem based on the training data. The component failure sample can be accurately obtained by increasing the sample size. The objective cost C obj can be shown in Fig. 4. Fig. 4 visualizes the total co-risk cost affected by the sample T N and violation ratio ε. It is noted that the total risk cost can slightly grow when increasing the sample size. The total cost C obj is $242.09 with the sample T N = 50 while the total cost C obj is $252.49 with T N = 1000. Meanwhile, the Wasserstein radius ρ increases from 0.0001 to 0.0025, the total cost C obj can achieve a 25.27% increase (e.g., ρ = 2.5× 10 −3 ). It is noted that the Wasserstein radius ρ can affect the conservation of the solution. It shows that the determinant of the total cost is ρ.

D. CO-RISK OPERATION ANALYSIS BASED ON CONTINGENCIES
The robustness of the proposed approach under the Wasserstein ambiguity set can be weighted by the radius ρ describing the Wasserstein distance between P and P N . Meanwhile, the risk assessment of violating chance constraints can be quantified by adjusting the violation ratio ε. Therefore, picking the appropriate Wasserstein radius ρ and violation ratio ε can achieve the tradeoff between robustness and risk-averse. The individual chance constraints are incorporated into the CEHS model considering the difference of each risk chance constraint. Thus, different chance constraints can be defined by picking different violation ratios ε. The appropriate Wasserstein radius ρ describing the most robust performance can be chosen from the result, while the violation ratio ε can be determined by the risk tolerance for the decision-maker.
In this section, four scenarios for the CEHS are designed to investigate the impact on co-risk operation under different contingencies. The detailed co-risk operation analysis can be given as follows.

1) SCENARIO 1: THE CEHS CO-RISK OPERATION WITHOUT CONTINGENCIES
In scenario 1, the risk caused by the contingencies for the CEHS is not included. The challenge of operation risk originates from the uncertain RESs integration. Namely, overlimit risk caused by uncertain power output can be weighted by the 156216 VOLUME 9, 2021  average violation probability V p . Adjusting the Wasserstein radius ρ and violation ratio ε can obtain the optimal corisk dispatch energy. In this scenario, the co-risk analysis can be implemented from two risk indices, including the average violation probability V p and the total co-risk dispatch cost C Obj , which is shown in Fig. 5 and Fig. 6. The dispatch energy for the CEHS can be shown in Fig. 7. The energy spillage can be depicted in Fig. 8. Fig. 5 illustrates the relationship between the average violation probability V p , violation ratio ε, and Wasserstein radius ρ. Note that as ρ increases, V p slightly increases. It is shown that ρ can affect the conservation of the solutions, which implies that the larger Wasserstein radius can satisfy more chance constraints. However, it is little effect on violation of chance constraints. Meanwhile, the average violation probability V p marginally declines as the violation ratio ε increases from 0.01 to 0.15. Nevertheless, when ε increases from 0.15 to 0.45, V p significantly decreases. The reason is that ε can affect the violation of chance constraints. Namely, the larger ε produces a significant influence on V p . Fig. 6 visualizes the average total cost C Obj affected by Wasserstein radius ρ and violation ratio ε. When 0.01 ≤ ε ≤ 0.05, the C Obj can show a marked increase with ρ increasing. It is noted that as the Wasserstein radius ρ increases, C Obj slightly increases. When ρ increases from 0.0005 to 0.0025(ε = 0.45), C Obj can only increase $0.034. The reason is that the Wasserstein radius ρ can affect the conservation of the solution without changing the dispatch result. As ε increases, the total risk cost C Obj decreases. When ε satisfies 0 ≤ ε ≤ 0.15, C Obj decreases quickly. However, the relative decline for the total cost C Obj is slight (ε > 0.15).    7 shows that the power output of the CHP and HP is affected by the Wasserstein radius ρ and the violation ratio ε. As seen from Fig. 7, the Wasserstein radius ρ and the violation ratio ε have little impact on the power output of CHP, which achieves the upper bound owing to the heavier load. The HP is seen as the critical adjustable flexible resource device. Increasing Wasserstein radius ρ can slightly reduce the power output of HP. However, the power output of HP increases from 0.52 MW to 0.61 MW when the violation ratio ε grows from 0.01 to 0.45. Therefore, the power output of HP is the leading dispatch energy, which is affected by ε. It shows that violation chance constraints can significantly influence the power output of HP. Fig. 8 presents energy spillage (including wind spilling and load shedding) with different ε and ρ. As seen from Fig. 8, the wind spilling is more than load shedding. Moreover, compared with ε, ρ has an essential impact on energy spillage. The main reason is that ρ can influence the conversation of the solution, which can produce an optimal dispatch solution.

2) SCENARIO 2: THE CEHS CO-RISK OPERATION WITH PS CONTINGENCIES
In scenario 2, the risk source comes from the uncertain RESs integration and components failure of PS mainly caused by the power transmission outrage. The DRCC approach can be employed to describe the risk level associated with the overlimit risk. Meanwhile, the power transmission failure risk can be weighted by using the MCS sample. Consequently, the risk assessment for the CEHS is quantified by incorporating both the overlimit risk and power transmission failure risk into the model. The Wasserstein radius ρ and violation ratio ε can affect the overlimit risk assessment. Moreover, simulating the power transmission failure can obtain the component failure risk sample. In this scenario, the average violation probability V p is served as the co-risk assessment indice, and the determinants of the risk assessment can be depicted in Fig. 9. The total co-risk dispatch cost C Obj can be considered the quantitated risk indice, as described in Fig. 10. The co-risk dispatch energy and the energy spillage for the CEHS can be shown in Fig. 11 and Fig. 12. Fig. 9 visualizes the average violation probability V p affected by Wasserstein radius ρ and violation ratio ε. Note that the average violation probability V p slightly declines as the violation ratio ε increases from 0.01 to 0.15, and the meshes in Fig. 9 shows fluctuations. Moreover, ε increases from 0.15 to 0.45, V p barely shifts. Meanwhile, Wasserstein radius ρ is little impact on V p . Compared with scenario 1, the average violation probability V p is higher in scenario 2 due to the involved component failure risk in PS.
As seen from Fig. 10, the total co-risk dispatch cost C Obj , which is served as one of the critical risk indices, is affected by ρ and ε. Increasing Wasserstein radius ρ can significantly increase the C Obj (ε = 0.01). When 0 ≤ ε ≤ 0.15, C Obj decreases quickly. C Obj is almost no reduction with ε over 0.15. With the comparison of scenario1, the total co-risk dispatch cost C Obj is higher in scenario 2 due to the energy spillage in PS considering the component failure risk.
Similarly, Fig. 11 illustrates that the dispatch energy is affected by the Wasserstein radius ρ and the violation ratio ε, including the power output of the CHP and HP. As seen from Fig. 11, ρ and ε are little impact on the power output of the CHP. The power output of HP is mainly affected by ε. Moreover, the dispatch energy is more than that in scenario one due to wind power output reduction considering the power transmission line failure. Fig. 12 presents the energy spillage affected by the Wasserstein radius ρ and the violation ratio ε. The two risk   parameters are little impact on the wind spilling. However, the load shedding is mainly affected by ρ.

3) SCENARIO 3: THE CEHS CO-RISK OPERATION WITH HS CONTINGENCIES
In scenario 3, the uncertain RESs integration and components failure of HS mainly caused by the heat pipe outrage can bring about the operation risk of the CEHS. As the analysis mentioned above, the overlimit risk is solved using the DRCC method, while the heat pipe failure risk is simulated using the MCS. In this scenario, the risk level for the CEHS is quantified by incorporating both the overlimit risk and heat   pipe failure risk into the model. In this scenario, the co-risk analysis can be implemented from two risk indices, including V p and C Obj , which is shown in Fig. 13 and Fig. 14. The corisk dispatch energy and energy spillage for the CEHS can be shown in Fig. 15 and Fig. 16, respectively. Fig. 13 and Fig. 14 show influence factors for the two risk indices. The determinants of the average violation probability V p is similar to that in scenario 2. In contrast to the risk of power transmission line failure, the risk caused by heat pipe failure is less than that in scenario 2 due to the slow change of temperature in the pipe. Fig. 15 gives the dispatch energy, and Fig. 16 presents the energy spillage. The determinations of dispatch energy and energy spillage are similar to those in scenario 2. Moreover, the power output of HP is mainly influenced by ε. The two risk parameters hardly affect the wind spilling. The load shedding is primarily affected by ρ. Both dispatch energy and energy spillage are less than that in scenario 2.

4) SCENARIO 4: THE CEHS CO-RISK OPERATION CONSIDERING BOTH PS AND HS CONTINGENCIES
In scenario 4, the operation challenge is brought by the uncertain RESs integration and components failure caused by the power transmission outrage for the PS and heat pipe outrage for the HS. As the analysis mentioned above, the corisk operation problem can incorporate the overlimit, power transmission, and heat pipe failure risks into the model. In this scenario, the co-risk analysis can be implemented from two risk indices, including V p and the C Obj , which is shown in Fig. 17 and Fig. 18. The co-risk dispatch energy and energy spillage for the CEHS can be shown in Fig. 19 and Fig. 20, respectively.     Fig. 18 depict influence factors for the two risk indices. The determinants of the average violation probability V p is similar to those in scenario 2 and scenario 3. In contrast to scenario 2 and scenario 3, the risk level is higher in scenario 4 due to both power transmission line failure and heat pipe failure. Fig. 19 gives the dispatch energy, and Fig. 20 presents the energy spillage in this scenario. The influence factors of dispatch energy and energy spillage are similar to those in scenarios 2 and 3. Moreover, the power output of HP is mainly influenced by ε. The two risk parameters hardly affect the wind spilling. The load shedding is primarily affected by ρ. Both dispatch energy and energy spillage are less than those in scenario 2 and scenario 3.
As seen from the analysis mentioned above, the risk indices V p and C Obj are mainly affected by ε. It shows that the violation of chance constraints can produce more risk and accelerate the total risk cost. Contingencies can also increase risk levels and the entire risk cost. The risk level and the risk cost in scenario 2, scenario 3, and scenario 4 are higher than those in scenario 1. Contingencies can also influence dispatch energy and energy spillage. Moreover, contingencies in PS can cause more risk than in other single scenarios.

V. CONCLUSION
Considering the high percentage of uncertain renewable integration and component failure, a contingencies-based DRCROM for the CEHS is present and provides co-risk operation analysis of the numerical case and the quantified risk level. The main conclusions are summarized as follows: (1) Considering the higher share of uncertain renewable energy and uncertain component failure, a contingenciesbased DRCROM is formulated to hedge against the risk caused by the uncertainty. Since the chance constraints violated may be different, the individual chance constraints can be incorporated into the risk constraints, and ε of chance constraints is different. Moreover, ε can be adjusted to enforce the risk level within the tolerance. Therefore, the DRICC methodology can be employed to handle the uncertain power output, which can lower the operation risk level. The WC-CVaR approximation can be applied to reformulate the original model into the tractable problem. The MCS is employed to simulate the component failure and obtain the sample, which is also incorporated into the risk constraints.
(2) To validate the proposed model approach, a modified Barry Island combined electricity and heat system can be considered. The numerical result shows that the risk indices V p and C Obj are mainly affected by ε. Contingencies can increase risk levels and the total risk cost in contrast to scenario 1. ρ and ε impacts on dispatch energy and energy spillage. Moreover, the power output of HP is mainly influenced by ε.
In further work, the component outage risk for the CEHS should be further investigated. The risk induced by the uncertain outage for the CEHS component can be discussed, considering the component failure distribution information disobeying a certain probability. Moreover, the extreme event caused by the carbon emission can be incorporated into the risk analysis, affecting the normal operation for the CEHS and bringing about more risk.

APPENDIX
The conic reformulation for the set BC is admitted as follows,