A Novel Strategy to Reduce Computational Burden of the Stochastic Security Constrained Unit Commitment Problem

: The uncertainty related to the massive integration of intermittent energy sources (e.g., wind and solar generation) is one of the biggest challenges for the economic, safe and reliable operation of current power systems. One way to tackle this challenge is through a stochastic security constraint unit commitment (SSCUC) model. However, the SSCUC is a mixed-integer linear programming problem with high computational and dimensional complexity in large-scale power systems. This feature hinders the reaction times required for decision making to ensure a proper operation of the system. As an alternative, this paper presents a joint strategy to efﬁciently solve a SSCUC model. The solution strategy combines the use of linear sensitivity factors (LSF) to compute power ﬂows in a quick and reliable way and a method, which dynamically identiﬁes and adds as user cuts those active security constraints N − 1 that establish the feasible region of the model. These two components are embedded within a progressive hedging algorithm (PHA), which breaks down the SSCUC problem into computationally more tractable subproblems by relaxing the coupling constraints between scenarios. The numerical results on the IEEE RTS-96 system show that the proposed strategy provides high quality solutions, up to 50 times faster compared to the extensive formulation (EF) of the SSCUC. Additionally, the solution strategy identiﬁes the most affected (overloaded) lines before contingencies, as well as the most critical contingencies in the system. Two metrics that provide valuable information for decision making during transmission system expansion are studied. and J.E.S.-A.; Formal analysis, C.C.M.-C., J.E.S.-A., J.M.L.-L., Á.J.-D. and J.G.V.; Funding acquisition, Á.J.-D., J.M.L.-L. and J.G.V.; Investigation, C.C.M.-C. and J.E.S.-A.; Methodology, C.C.M.-C. and J.E.S.-A.; Project administration, Á.J.-D., J.M.L.-L. and J.G.V.; Resources,C.C.M.-C., J.E.S.-A., Á.J.-D., J.M.L.-L. and J.G.V.; Software, C.C.M.-C. and J.E.S.-A.; Supervision, Á.J.-D., J.M.L.-L. and J.G.V.; Validation, C.C.M.-C. and J.E.S.-A.; Visualization, C.C.M.-C. and J.E.S.-A.; Writing—original draft, C.C.M.-C. and J.E.S.-A.; and Writing—review and editing, Á.J.-D., J.M.L.-L. and J.G.V.


Introduction
The security constrained unit commitment (SCUC) problem is the optimization model that schedules a set of power generation units to meet an expected energy demand over a typical time horizon of 24 h while meeting several system constraints (e.g., operational limits of the generation units and transmission limits) under normal operating conditions and also under N − 1 contingencies. The SCUC problem is widely used in day ahead and real time markets to obtain an economically, reliably, and safely operation of the power system [1]. The SCUC is a complex mixed integer programming (MIP) problem since it belongs to the class of NP-hard problems [2]. This feature has challenged researchers to approach the complexity of the SCUC problem in different ways when dealing with large-scale power systems. Examples of these solution approaches include Lagrangian relaxation (LR) and Benders decomposition (BD) techniques ,which have been implemented to solve the SCUC model in efficient computing times [3][4][5][6][7]. Dealing with N − 1 security constraints posses a challenge when solving the SCUC problem. Therefore, several authors have proposed different alternatives to tackle them. In [8], the authors defined necessary and sufficient conditions to identify and eliminate inactive N − 1 security constraints. Similarly, H. Wu et al. [9] presented analytical feasibility conditions to define only the subset of security constraints of the SCUC problem feasible region. Likewise, in [10], the authors mathematically ruled out unnecessary security constraints using again necessary and sufficient conditions. However, the application of these methods could be limited for large-scale practical problems, and therefore decomposition techniques such as BD should be required in a complementary manner [9,10]. On the other hand, some authors have addressed this challenge without relying to decomposition techniques or any of the previously mentioned methods. This is the case of Capitanescu et al. [11] who developed a preventive security-constrained optimal power flow (PSCOPF) defining two constraint filtering techniques under the concepts of dominated and non-dominated constraints. Including only the latter in the PSCOPF model accelerates its solution time.
Currently, a key factor to reduce the computing times when solving the SCUC problem is the use of linear sensitivity factors (LSFs) to model the power grid. The Power Transfer Distribution Factor (PTDF) and the Line Outage Distribution Factor (LODF) estimate power flows changes under power injection variations in buses and line outages (i.e., N − 1 contingencies), respectively [12]. Nonetheless, Eslami et al. [13] proposed a new method for reducing SCUC variables and constraints under a conventional formulation of the power grid for post-contingency constraints. Several iterative methodologies based on LSFs have been developed to increase the computational efficiency of SCUC models for large-scale power systems. In [2], the authors added only binding N − 1 security constraints to the SCUC problem using the Gurobi deferred constraints function (i.e., lazy constraints). Similarly, Marín-Cano et al. [14] proposed the implementation of user cuts to improve the computational performance of the SCUC and presented indices that provide expansion signals for the transmission network. Xavier et al. [15] added effective security constraints based on the concept of priority queue, adding at most only one N − 1 security constraint per line. All these methodologies aim to generate iteratively the feasible region of the SCUC problem reducing effectively the computation times needed to solve this problem. However, it is noteworthy that none of the aforementioned methodologies has been applied in the analysis of the safe operation of electrical systems under uncertainty [16].
Recently, the evolution of power systems has imposed major challenges in the solution of the SCUC due to the uncertainty caused by the integration of large amounts of intermittent renewable resources. One way to address these challenges is through the SSCUC problem [17]. In [18], the authors evaluated security through a multi-stage stochastic programming model where random system disturbances (generation unit and transmission line failures) are modeled through scenario trees. To solve the SSCUC, they decomposed the problem using LR. The authors of [19] used BD to solve a SSCUC considering the intermittency and volatility of wind generation through scenarios. Two different procedures to tackle the SSCUC are compared in [17]. The first one is based on scenarios, whereas the second one is based on optimization intervals. The SSCUC model presented in [20] considers wind power production scenarios and non-spinning reserves. The authors relied on chance-constrained programming to solve this variant of the SSCUC. Similarly, the authors of [21] studied an N − 1 security and chance-constrained unit commitment (SCCUC) model. They formulated this SCCUC problem as a mixed-integer second-order cone program with two types of reserves to consider wind fluctuations and single component outages (i.e., generators or transmission lines outages). In [22], the authors relied on BD to solve a SSCUC model considering different types of pumped-storage hydropower plants. Finally, Park et al. [23] modeled a SSCUC using two-stage stochastic programming, considering dynamic line rating to improve the use of transmission lines during N − 1 conditions. The above references highlight the importance of the SSCUC models as a way to hedge against uncertainty in the safe operation of power systems. However, it should be pointed out that high computational speed is of great importance for operators at the moment of real-time decision-making, since the mathematical models for UC under uncertainty are computationally much more complex than their deterministic counterparts [24]. To the best of our knowledge, there are few studies that address the computational complexity of the SSCUC problem and focus on solving the problem quickly and effectively. Among them, Wang and Fu [1] proposed a fully parallel stochastic SCUC approach to obtain efficient and fast solutions in large-scale electrical systems and N. Yang et al. [25] proposed a method that aims to improve the efficiency of the SSCUC problem under constrained ordinal optimization (COO). This approach has two advantages: The first one is a recognition of discrete variables to reduce the solution space of the approximate COO rough model. The second one is the inclusion of a strategy to reduce non-effective security constraints, as presented in [26], which in turn is based on the analytical conditions proposed in [8].
Following the authors of [1,25], in this work, we propose a novel and easy-to-implement strategy to efficiently solve a SSCUC problem under the uncertainty of large penetration of intermittent energy sources. This strategy has the following features: (i) it uses sensitivity linear factors to compute power flows to reduce the number of variables and constraints when compared to conventional network modeling; (ii) it considers user cuts to dynamically add active N − 1 security constraints, as proposed in [14]; and (iii) it uses the PHA decomposition algorithm to relax the SSCUC model. As shown in [27], the PHA is able to achieve high-quality solutions for stochastic optimization in large-scale power systems. As an additional result from the proposed strategy, we define several indices that allow identifying the most critical and vulnerable lines in terms of the overloads they cause or suffer under single contingencies. Such information is useful to the system planner for short-term transmission expansion analyses.
The rest of the paper is organized as follows. Section 2 introduces some preliminary background on stochastic optimization and presents the SSCUC model studied in this paper. Section 3 presents the elements of the solution strategy proposed in the paper as well as the general framework that embeds N − 1 security cuts into the progressive hedging algorithm. Section 4 includes the results of the evaluation of the proposed strategy in the IEEE RTS-96 test system. Section 5, introduces and exemplifies transmission expansion indices calculated with the results of the proposed strategy. Finally, Section 6 concludes the paper and provides additional research extensions.

Stochastic SCUC Formulation
In this section, we model the SSCUC problem as a two-stage stochastic program. We first set the scene for this by introducing some important background on stochastic optimization, and then we introduce the detailed SSCUC model. For easy reference, Appendix A presents the nomenclature used in the modeling of the SSCUC.

Background
Stochastic models to address the SCUC problem have been widely used by system operators, especially when wind generation is considered [1]. From a mathematical point of view, the SSCUC problem can be formulated as a two-stage stochastic MIP problem based on scenarios. In the first stage, the generation units are committed, while, in the second (uncertainty) stage, the power output level of each generator to meet the demand is decided, taking into account the N − 1 criterion. In sum, the SSCUC problem aims to minimize the expected operation cost (Equation (1)), while satisfying the generation and transmission constraints in all scenarios (ξ ∈ Ξ). The SSCUC problem represented through its extensive form (EF) follows: The decision variable x of the first stage represents the on/off state of a generator (i.e., its commitment), whch is implicitly independent of the realization of the scenarios (ξ) [28]. The decision variable of the second stage y(ξ) represents the power dispatched by each generation unit in each scenario (ξ) (with probability p ξ for ξ ∈ Ξ). Equations (2) and (3) represent the constraints associated with the first and second stage, respectively. A is the matrix of coefficients associated with variable x. The coefficients c and d are the costs associated with the binary variable x and the continuous variable y, respectively.
An alternative representation of the SSCUC problem is through the Extensive Form of the Scenario (EFS) as follows: Unlike Equations (1)-(3), in Equations (4)-(7) the variables of the first stage are calculated for each realization of the scenario (ξ). Additionally, it includes (7), known as the non-anticipativity constraint, which requires that the variable x does not depend on the scenarios for all feasible solutions [28].

Mathematical Model of the SSCUC Problem
The first and second stages of the SSCUC are based on the formulation and modeling in [29,30], respectively. Constraints of the first stage include binary variable logical implications, start-up and shut down costs and minimum down and up time. On the other hand, constraints of the second stage consider minimum and maximum generator outputs as well as ramping and transmission limits. Contrary to H. Pandžić et al. [29], the formulation of transmission constraints uses LSF, PTDF for normal operation and LODF for operation under N − 1 line contingencies. Details of the LSF calculation are provided in [2,14].
The objective function of the SSCUC minimizes the expected value of the operating costs of the system, consisting on the costs of the generation of thermal plants C i (t, ξ), the costs associated with the variable of involuntary load shedding L sh s (t, ξ) and the costs associated with wind power curtailment Q w (t, ξ). This objective function is subject to the constraints given by (9)-(30), as shown below. Equations (9)-(24) are well-known constraints of the UC formulation that can be found in classical references of this problem (e.g., [29]). We add to this formulation net-power balance equations (constraints (25) and (26)), involuntary load shedding (constraint (29)) and wind power curtailment (constraint (30)). Remarkably, the constraints in Equation (27) and (28) resort to the PTDF and LODF to reduce the complexity of modelling the pre-and post-contingency power flows on lines.
This objective function is subject to the constraints given by (9)-(30).
• Binary variable logic: • Operating costs of thermal plants: • Total power output of the thermal generator i, expressed as the total sum of the generation level in each segment b of the cost curve: • Minimum generator output constraint: • Maximum generator output constraint: • Initial on-off status of generator i at t = 0 (15), minimum up (16) and down time constraints (17): • Ramping constraints: • Generator off counter set-up constraints. In this case, symbols | and ∧ indicate the logical conditions IF and AND, respectively: • Net power balance constraints: • Power flow constraints under normal (27) and post-contingency (28) operation conditions. In this case, PTDF l,s is the power transfer distribution factor of line l for a power injection at bus s, and LODF l,k is the line outage distribution factor for line l when line k is out of service.

Solution Strategy
The proposed solution strategy embeds the dynamic generation of N − 1 security constraints into the PHA. We first describe these two elements independently and then present the embedding procedure.

Method for Adding Binding N−1 Security Constraints
LODFs allow the modeling of N − 1 contingencies with low computational cost [2,12,14]. From the standpoint of computational efficiency, it makes no sense to consider all possible N − 1 contingencies, since not all of them generate violations on the maximum thermal transmission limit of the rest of system lines. Therefore, to take into account only the binding N − 1 security constraints that define the feasible region of the SCUC problem, we follow the method proposed by Marín-Cano et al. [14]. Such approach is based on the concept of user cuts. These cuts are linear constraints strategically and dynamically added by the user, based on the implicit information of the problem [31], with the objective of tightening model (8)- (30).
Algorithm 1 shows the method proposed by Marín-Cano et al. [14] to add (only) sufficient and necessary N − 1 constraints to the model through user cuts. In Step 1, the SCR l,k,t (ξ) parameter is set to zero to discard the N − 1 security constraints (28), since the objective is to solve only the base case (normal operation) of the system during the first iteration of the process. In the main loop (Steps 2-14), the model defined by Equations (8)-(30) is solved iteratively (Step 3). This process, Steps 4-9, estimates useful parameters that allow finding those N − 1 security constraints directly linked to the feasible region of the SSCUC problem. In Step 4, parameter OP l,k,t (ξ) is reset to zero, with the aim of storing the overload value of only new N − 1 security constraints. In Step 5, the parameterF l,k (t, ξ) estimates the post-contingency power flows of the lines based on the optimal flows in the normal operating state of the line l and the contingency k, which are F * l (t, ξ) and F * k (t, ξ), respectively. These values were previously calculated in Step 3. OnceF l,k (t, ξ) is estimated, Step 6 compares it with the maximum transmission capacity for each line F max l . If this capacity is surpassed, Step 8 calculates and stores the values of the overloads generated at the lines in parameter OP l,k,t (ξ). Furthermore, the sum of this parameter allows calculating the total system overload TO in Step 9. This value is useful to define the stopping condition of the algorithm (Step 12). In Step 10, the addition of N − 1 security constraints (28) is performed, for combinations of l, k, t, ξ where SCR l,k,t (ξ) is equal 1; this last parameter allows to record and count all the overloads in each iteration. In Step 12, if there have been no overloads since the first iteration of the algorithm, then OP l,k,t (ξ) is zero, thus TO = 0.
Step 14 is the end of the iterative process when TO is zero or less than a tolerance value tol; in this case, we set tol = 1 −7 MW. For a detailed description of this method, the interested reader is referred to [14].

Progressive Hedging Algorithm
The extensive form of a stochastic optimization problem consists of a single mathematical model where all constraints are written for all possible scenarios [32]. For large-scale problems, such formulation does not guarantee a solution due to the lack of computational resources [32]. As a solution alternative, PHA is a scenario-based decomposition method initially proposed by Rockafellar and Wets [33]. In PHA, the stochastic problem is solved independently scenario by scenario, by relaxing the non-anticipativity constraint. This decomposition obtains good solutions with high computational efficiency [28]. Furthermore, PHA has proven to be much more stable than other widely used algorithms in the literature such as BD [32,34]. For these reasons, PHA has recently been used as a promising solution method for the stochastic operation of electrical power systems, such as unit commitment [27,32,[35][36][37] and the short-and medium-term hydro-thermal planning [34,38]. Algorithm 2, taken from Gade et al. [28], describes the PHA for two-stage Stochastic MIP problems. In PHA, ρ is a penalty (scalar or vector) of the same length of the non-anticipativity variable x [28,32]. In Step 1, the iteration counter ν and weight vector w 0 (ξ) are initialized. In Steps 2-4, the initial solution of the subproblems of each scenario is obtained.
Step 5 updates the iteration counter ν.

In
Step 6, the algorithm defines the current best estimate of a non-anticipativity solution (x ν ) that does not depend on the realization of scenario ξ.
Step 7 estimates multipliers w(ξ) that are used to update the non-anticipativity condition. In Steps 8-10, the multipliers w(ξ) are used together with a proximal term ( ρ 2 x −x ν 2 ) to find an optimal solution when the x values do not depend on the scenarios. The main loop (Steps 5-10) is repeated until the convergence condition is met (Step 11).
is the same for all scenarios ξ ∈ Ξ then stop and reportx ν ; 12 else go to 5;

Adjustment of ρ Parameter and Convergence Improvement
As reported in [39,40], an important factor for the performance of the PHA is its high sensitivity to the value of the penalty parameter ρ. This parameter determines the step length of the W dual price update, and its inappropriate choice could lead to non-convergence or slow convergence times [32]. Currently, there is no general rule for the selection of ρ based on theoretical analysis [40]. However, high values of ρ could give faster convergence rates in Lagrangian methods [39]. An adaptive method, as well as a comprehensive summary of the different strategies for the adjustment of ρ can be found in [40].
Different strategies for adjusting the penalty factor ρ have been used when dealing with the stochastic UC problem. The authors of [27] used a value of ρ proportional to the hedging variable (i.e., ρ is proportional to the cost of the on/off state variables). An adjustment of ρ based on locational marginal prices for different nodes and time periods was proposed by Ryan et al. [32]. Alternatively, Li et al. [35] proposed a method for using the shadow prices of the coverage variable as a penalty factor ρ. Finally, Ordoudis et al. [36] defined a new way of adjusting ρ as a fraction of the objective function proportional to the cost of the hedging variable averaged over all scenarios.
Additionally, the cyclic behavior of the PHA does not guarantee convergence while solving the Stochastic UC problem, due to the binary nature of the on/off state variables. To tackle this problem, some heuristic variable setting procedures have been proposed to improve the convergence speed of the PHA [27,35]. In this work, we follow the rounding technique proposed by Li et al. [35] to guarantee convergence at the end of the PHA. Settingx when its value is close to 0 or 1, this rounding technique has proven to be flexible, easy to implement, and it has been used successfully by other authors [36]. In this technique, two thresholds values (α, β) define three categories to round the value ofx, as follows: Although this rounding technique can speed-up the computation time of the PHA, it could also generate infeasibility [35]. To avoid this behavior, threshold values of β must be adjusted close to 0 to guarantee enough availability of powered-up units within the system [35,36]. On the other hand, α can be chosen with more flexibility since a large availability of power units guarantees feasibility.

Integrating N−1 Security Constraints Within the PHA
As a novel and effective strategy to address the computational challenge imposed by the SSCUC problem, Algorithm 3 shows the integration of the iterative method proposed by Marín-Cano et al. [14] within the PHA. Algorithm 3: Method for Adding N − 1 Constraints Embedded Within the PHA. Steps 1 and 2 initialize the necessary working parameters for the N − 1 constraint addition method and for the PHA, respectively. Step 3 updates the iteration for the main iterative process of the PHA where every sub-problem of each scenario (Step 4) is solved. Steps 5-17 account for the proposed N − 1 security constraint addition method, which is embedded within the iterative PHA process.
Here, x ν (ξ) represents the hedging variable of PHA in each iteration ν and each scenario ξ. In this case, the hedging variable is the binary variable (x ∈ {0, 1}) that assigns the on/off status of each generation unit. Once a solution for each scenario has been found through the N − 1 cut-adding method (Step 17), Steps 18-20 are the continuation of the PHA to find an optimal solution x ν (ξ) for all scenarios meeting the non-anticipativity condition (Step 19), and the penalty factor (Step 20) and weight updating procedures (Step 21), respectively. It is important to note that, to improve the PHA convergence, in Step 20, the penalty factor ρ ν is updated in every PHA iteration ν, according to the strategies presented in Section 3.2.1. The procedure stops when the stopping condition (convergence of x) is met and the valuex ν reported (Step 21).

Results
The test system used to evaluate the proposed strategy is the IEEE RTS-96 (see Figure 1) used previously whose information is fully available in [30]. This system is made up of 73 buses, 120 transmission lines, 96 thermal generators, and 19 wind farms (with an installed capacity of 6900 MW and located at buses 102, 114, 116, 118, 121, 119, 123, 202, 212, 213, 219, 220, 223, 301, 306 and 309. There are 51 loads with a maximum demand of 7539 MW and an average demand of 6258 MW. A time horizon of 24 h is considered. The SSCUC problem for this system considers 10 scenarios each with a given probability representing different generation outputs of the wind farms throughout the 24 h of the day. A detailed description of the scenarios under consideration, as well as the Excel sheets containing all input data, can be consulted in [30]. The proposed computational strategy was tested on two computers: an Intel Xeon E5 @ 2.40GHz server with 44 cores and 256 GB of RAM memory (hereafter, server) and an Intel Core i7 @ 3.4GHz desktop computer with 8 cores and 8 GB of RAM memory (hereafter, desktop). The proposed approach was implemented using the commercial algebraic modeling system GAMS version 24.8.5, with CPLEX V12.6.1 as optimization engine.

Parameters Setting for the PHA
One of the first parameters to be adjusted is the penalty factor ρ. Values of ρ based on energy production costs of the hedging variable are suggested as a good option to foster the convergence of the PHA [36]. Additionally, we selected two of the updating strategies cited in Section 3.2.1: ρ 1 is the strategy proposed by Watson and Woodruff [27], whereas ρ 2 is the one proposed by Ordoudis et al. [36].
Another important factor is the selection of the hedging variable. Some authors have carried out experiments to study which is the best hedging variable within the Stochastic UC problem [35,36]. Under different adjustment strategies of ρ, Ordoudis et al. [36] showed that the best-performing hedging variable is the on/off status variable, rather than the start and stop down variable. Accordingly, in our computational experiments, we used as hedging variable the on/off state variable (x ∈ {0, 1}).
Finally, to guarantee convergence and feasible PHA solutions in the different tests of the proposed computational strategy, we set parameters α and β in the context of the so-called rounding technique. It is advisable to assign 0 to β to guarantee whole availability of generation units to have a higher chance of reaching feasibility as mentioned in [35,36]. On the other hand, α must be set for every problem. In this case, an initial value of α is set and small variations over this value are performed.
As illustrated in [36], small increments of α can provide high quality solutions in short PHA running times. Therefore, in our experiments, we set β to zero, and the values of α were set heuristically. Beginning at a given (relatively high value) of 0.6, we ran PHA reducing the value of α iteratively in steps of 0.01 until reaching a good trade-off between solution quality and running time. In this case, we stopped the search at α = 0.53. Below this value, PHA presented convergence problems. The results presented in this section (on the server) were obtained with this value. It is noteworthy that the value of α was set to a different value (α = 0.57) for the experiments performed on the desktop.

Formulations
To validate the performance of the strategy proposed in this work, we carried out a comparative analyses between the following three formulations: Formulation A: Extensive Formulation of the SSCUC without the N − 1 cut-adding method. Formulation B: Extensive Formulation of the SSCUC with the N − 1 cut-adding method. Formulation C: This formulation corresponds to the computational strategy proposed in this work. That is, the SSCUC formulation with the N − 1 cut-adding method embedded within PHA. Additionally, we compared two methods for adjusting ρ in this formulation (ρ 1 and ρ 2 ), as defined in Section 3.2.1.
As a benchmark, we first solved the EF of the SSCUC problem (Formulation A). This gives some key metrics to compare the performance of the other formulations: the optimality gap and the objective function of the optimal solution of the problem (i.e., the best relaxed solution). As stopping criterion, the minimum relative optimality gap defined in GAMS was set to 1% for all tested formulations. The optimality gap was calculated as the difference between the solution of a formulation and the best relaxed solution, divided by the best relaxed solution. Table 1 presents the values of the objective function, the computing time , and the optimality gap of Formulation A, obtained with the server.

Analysis and Comparison of Results
Ten scenarios were used to consider the uncertainty associated to wind production within the IEEE-RTS-96 test case, which were taken from Pandzic et al. [30]. Figure 2 compares the gap (%) (blue bars) and running time (s) (green line), taken to solve the SSCUC model for each formulation on the server. As shown by Ordoudis et al. [36], PHA is a promising algorithm with great potential to efficiently address the stochastic unit commitment problem in large-scale power systems as the number of scenarios increases. Figure 2 shows that Formulation A takes a (large) computing time of 61, 222 s (approximately 17 h) to reach an optimality gap of 1.1% with an objective function of 1.1110 MUS$. Meanwhile, the performance of Formulation B is notably better; it took only 7245 s (approximately 2 h) to find a solution with an optimality gap of 0.94% (which corresponds to an objective function of 1.1109 MUS$). This is 8.5 times faster with respect to the running time of Formulation A for a solution of even better quality. This first result illustrates that the inclusion of the N − 1 cut-adding method proposed by Marín-Cano et al. [14] generates a more compact SSCUC model. Additionally, considering explicitly single contingencies for all 118 lines generates 3.3 million N − 1 security constraints (in Formulation A). By contrast, the N − 1 cut-adding method of Formulation C generated only 0.116 million binding security constraints defining the feasible region of the SSCUC problem. As a result, the cut-adding method embedded in PHA reduces the size of the model (of Formulation A) by 99.6% (by ignoring unnecessary security constraints within this model). As a consequence, the comparison of Formulations A and C (with ρ 1 and ρ 2 ) shows that the latter formulation achieves results of comparable quality in much faster running times. Even though the solutions reached by Formulation C have slightly greater optimality gaps (compared to the other formulations), it reaches acceptable sub-optimal solutions (with a gap of 1.62% for Formulation C-ρ 1 and 1.85% for Formulation C-ρ 2 ) in roughly 1200 s (approximately 20 min). Remarkably, Formulation C (with ρ 1 or ρ 2 ) is on average 50 times faster than Formulation A, and it is six times faster than Formulation B.   (25.8 min). Notably, the results of both variants of Formulation C on the server (see Figure 2) and on the desktop computers (see Figure 3) are comparable in terms of both solution quality and running time. This, once again, shows the contribution of our novel strategy when it comes to obtaining high-quality solutions in short times, even when using computing resources of limited capacity.
On the other hand, Formulation C was unable to get solutions with gaps below 1% on the server (see Figure 2), as well as to get gaps below 2%, on the desktop computer (see Figure 3). This behavior is due to the fact that the N − 1 SSCUC problem is an MIP problem with high dimensionality and PHA is a method that guarantees convergence and optimality only for convex problems. Therefore, some authors (e.g., Guo et al. [41]) have proposed integrating PHA with dual decomposition to have exact solutions for stochastic mixed-integer programming problems. However, Formulation C offers a good trade-off between the degree of suboptimality of the solution and the elapsed running time. The choice of the penalty factor updating strategies (Formulation C-ρ 1 or Formulation C-ρ 2 ) offers an additional control parameter, since the better is the quality of the solution (with Formulation C-ρ 1 ) the longer is the running time. As a reference, detailed results of the experiments presented in this section are given in Appendix B.

Expansion Transmission Indices
In addition to establishing the safe operation of the system, an added value of the N − 1 outage addition method proposed by Marín-Cano et al. [14] is the identification of the most severe line faults and vulnerable lines under contingency analysis. This additional information serves as indicative signals of needed expansion of the transmission system. These indices are computed considering the number of overloads and most severe outage of the network in scenario ξ (i.e., SCR l,k,t (ξ)).
The most severe line failures are those whose contingency generates the greatest amount of overloads on other lines l. Index L S k estimates the expected value of the number of lines l affected by contingency k, and it is obtained by rounding the weighted sum of the parameter SCR l,k,t (ξ) over the sets of time periods t, contingencies l and scenarios ξ. Briefly, L S k indicates which line contingencies are the most severe, generating the greatest number of overloads on other lines.
As an example, Figure 4 shows the most severe line faults that cause the greatest overload impact on the system, according to parameter L S k . Outages of lines L119, L120, L118 and L64 (highlighted in orange in Figure 1) would generate, respectively, 176, 176, 81 and 60 thermal overloads in the transmission lines that remain in operation on the system. Index L V l , estimates the expected value of the number of times that line l is overloaded when the contingency analysis is performed. This index show how vulnerable a line is when other lines are out of service. That is, L V l indicates the number of overloads that a line experiences when transmission contingencies are produced in the system. This parameter is calculated by rounding the weighted sum of SCR l,k,t (ξ) over the sets of time periods t, contingencies k and scenarios ξ. Figure 5 shows the most vulnerable lines of the IEEE-RTS-96 test system based on Index L V l . In this system, the lines with the highest number of overloads for different N − 1 contingencies are L11, L118, L12 and L115 (highlighted in green in Figure 1 Finally, the results of these indices suggest which lines to prioritize for the reinforcement of the transmission network through parallel circuits on some of the lines with the highest values of L V l and L S k . This could guide short-term transmission expansion analyses aimed at reducing the risks of post-contingency overloads on other elements of the network. That, in turn, increases the security, reliability and operating margins of the power system to deal with an increasing demand in the short and medium terms.

Conclusions
This work presented a novel strategy to efficiently solve the Stochastic Security Constrained Unit Commitment (SSCUC) problem (coined as Formulation C through the paper). It combines Linear Sensitivity Factors (LSF), an efficient cut-adding method and the progressive hedging algorithm (PHA). The use of the cut-adding method reduces drastically (by 99%) the size of the problem by only including active N − 1 security constraints in the model. Embedding (for the first time) this method into the PHA provides an effective and efficient solution approach, which is able to obtain solution of comparable quality but is 50 times faster than the extensive formulation (EF) of the SSCUC problem (Formulation A through the paper). Computational experiments also showed that the proposed solution strategy is stable to the choice of the updating procedure of the ρ penalty factor, a key parameter for the successful convergence of the PHA. Additionally, the fine tuning of parameters α and β depends on the problem characteristics (i.e., power system topology) and the need for high quality solutions of the SSCUC problem in short running times. The setting of these values allow the exploration of the trade-off between running time and solution quality. This characteristic of the proposed approach is valuable and can be exploited by the system operator when very short running times are needed.
Similarly, computational experiments on a (rather standard) desktop computer showed its value when it comes to obtain high quality solutions of the SSCUC problem with limited computational resources (where the EF or even the cut-adding method without PHA are not implementable). Additionally, as a side result of the proposed solution strategy, it is possible to identify the most affected (overloaded) lines before contingencies, as well as the most critical contingencies in the system. This is possible thanks to the values of two indices calculated with the information of the N − 1 cuts added during the execution of the PHA iterative process. These two indices provide valuable information for decision-making during short-to medium-term transmission system expansion studies.
Although the PHA technique is more stable than other step-wise decomposition algorithms such as Benders decomposition, it is advisable to test the performance of the N − 1 security constraint addition method with this type of technique. This suggest the integration of PHA and a promising decomposition method for Stochastic Mixed-Integer Programming: (i.e., Fenchel decomposition based on PHA), which has not yet been reported on the analysis of stochastic models for the operation of electrical power systems. Likewise, to further reduce the running time of the proposed approach, parallel versions of the solution strategy could be tested.
Finally, the method proposed in the paper can be extended to take into account other sources of uncertainty such as photovoltaic generation and demand, as well as batteries to hedge against them. For this, new constraints must be taken into account, such as the charging and discharging of batteries during the time horizon and the expected generation of the photovoltaic plants, among others.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript: