Computationally Efficient Market Simulation Tool for Future Grid Scenario Analysis
Abstract
This paper proposes a computationally efficient electricity market simulation tool (MST) suitable for future grid scenario analysis. The market model is based on a unit commitment (UC) problem and takes into account the uptake of emerging technologies, like demand response, battery storage, concentrated solar thermal generation, and HVDC transmission lines. To allow for a subsequent stability assessment, the MST requires an explicit representation of the number of online generation units, which affects powers system inertia and reactive power support capability. These requirements render a fullfledged UC model computationally intractable, so we propose unit clustering, a rolling horizon approach, and constraint clipping to increase the computational efficiency. To showcase the capability of the proposed tool, we use a simplified model of the Australian National Electricity Market with different penetrations of renewable generation. The results are verified by comparison to a more expressive and computationallyintensive binary UC, which confirm the validity of the approach for long term future grid studies.
Electricity market, future grid, electricity market simulation tool, optimization, scenario analysis, unit commitment, stability assessment, inertia, loadability.
[A01]Set of consumers . \nomenclature[A02]Set of generators , . \nomenclature[A03]Set of synchronous generators, . \nomenclature[A04]Set of renewable generators, . \nomenclature[A05]Set of concentrated solar thermal generators, . \nomenclature[A06]Set of synchronous generators in region . \nomenclature[A07]Set of subhorizons . \nomenclature[A08]Set of power lines , . \nomenclature[A09]Set of AC power lines, . \nomenclature[A10]Set of HVDC power lines, . \nomenclature[A11]Set of nodes . \nomenclature[A12]Set of nodes in region . \nomenclature[A13]Set of prosumers . \nomenclature[A14]Set of regions . \nomenclature[A15]Set of storage plants . \nomenclature[A16]Set of time slots .
[D01]Number of online units of generator , in BUC and in MST. \nomenclature[D02]Integer startup status variable of a unit of generator , in BUC and in MST. \nomenclature[D03]Integer shutdown status variable of a unit of generator , in BUC and in MST. \nomenclature[D04]Voltage angle at node . \nomenclature[D05]Power flow on line . \nomenclature[D08]Power loss on line . \nomenclature[D09]Power dispatch of generator . \nomenclature[D10]Grid/feedin power of prosumer . \nomenclature[D11]Power flow of storage plant . \nomenclature[D12]Battery power flow of prosumer . \nomenclature[D13]Thermal energy stored in TES of generator . \nomenclature[D14]Energy stored in storage plant . \nomenclature[D15]Battery charge state of prosumer .
[P]Fix/variable cost of a unit of generator . \nomenclature[P]Startup/shutdown cost of a unit of generator . \nomenclature[P]Load demand of consumer . \nomenclature[P]Load demand of prosumer . \nomenclature[P]Power reserve requirement of node . \nomenclature[P]Minimum/maximum limit of variable . \nomenclature[P]Total number of identical units of generator . \nomenclature[P]Rampup/down rate of a unit of generator . \nomenclature[P]Minimum up/down time of a unit of generator . \nomenclature[P]Time slot offset index. \nomenclature[P]Time resolution. \nomenclature[P]Susceptance of line . \nomenclature[P]Max. output power of renewable generator . \nomenclature[P]Max. thermal power capture by generator . \nomenclature[P]Inertia of a unit of generator . \nomenclature[P]MVA rating of a unit of generator . \nomenclature[P]Minimum synchronous inertia requirement of node . \nomenclature[P]Efficiency of component . \nomenclature[P]Aggregated PV power of prosumer . \nomenclature[P]Feedin price ratio.
[I]Number of online units of generator at start of horizon. \nomenclature[I]Power dispatch of generator at start of horizon. \nomenclature[I]Minimum number of units of generator required to remain online for time . \nomenclature[I]Minimum number of units of generator required to remain offline for time . \nomenclature[I]Energy stored in TES of at start of horizon. \nomenclature[I]Energy stored in storage plant at start of horizon. \nomenclature[I]Battery state of charge for prosumer at start of horizon.
1 Introduction
Power systems worldwide are moving away from domination by largescale synchronous generation and passive consumers.
Instead, in future grids
Future grid planning thus requires a major departure from conventional power system planning, where only a handful of the most critical scenarios are analyzed. To account for a wide range of possible future evolutions, scenario analysis has been proposed in many industries, e.g. in finance and economics [2], and in energy [3, 4]. In contradistinction to power system planning, where the aim is to find an optimal transmission and/or generation expansion plan, the aim of scenario analysis is to analyze possible evolution pathways to inform power system planning and policy making. Given the uncertainty associated with longterm projections, the focus of future grid scenario analysis is limited only to the analysis of what is technically possible, although it might also consider an explicit costing [5]. In more detail, existing future grid feasibility studies have shown that the balance between demand and supply can be maintained even with high penetration of RESs by using largescale storage, flexible generation, and diverse RES technologies [6, 7, 8, 9, 10]. However, they only focus on balancing and use simplified transmission network models (either copper plate or network flow; a notable exception is the Greenpeace panEuropean study [11] that uses a DC load flow model). This ignores network related issues, which limits these models’ applicability for stability assessment.
To the best of our knowledge, the Future Grid Research Program, funded by the Australian Commonwealth Scientific and Industrial Research Organisation (CSIRO) is the first to propose a comprehensive modeling framework for future grid scenario analysis that also includes stability assessment. The aim of the project is to explore possible future pathways for the evolution of the Australian grid out to 2050 by looking beyond simple balancing. To this end, a simulation platform has been proposed in [12] that consists of a market model, power flow analysis, and stability assessment, Fig. 1. The platform has been used, with additional improvements, to study fast stability scanning [13], inertia [14], modeling of prosumers for market simulation [15, 16], impact of prosumers on voltage stability [17], and power system flexibility using CST [18] and battery storage [19]. In order to capture the interseasonal variations in the renewable generation, computationally intensive timeseries analysis needs to be used. A major computational bottleneck of the framework is the market simulation.
Within this context, the contribution of this paper is to propose a unified generic market simulation tool (MST) based on a unit commitment (UC) problem suitable for future grid scenario analysis, including stability assessment. The tool incorporates the following key features:

market structure agnostic modeling framework,

integration of various types and penetrations of RES and emerging demandside technologies,

generic demand model considering the impact of prosumers,

explicit network representation, including HVDC lines, using a DC power flow model,

explicit representation of the number of online synchronous generators,

explicit representation of system inertia and reactive power support capability of synchronous generators,

computational efficiency with sufficient accuracy.
The presented model builds on our existing research [16, 15, 19, 17, 14, 18] and combines all these in a single coherent formulation.
In more detail, to reduce the computational burden, the following techniques are used building on the methods proposed in [20, 21]:

unit clustering,

rolling horizon approach,

constraint clipping.
The computational advantages of our proposed model are shown on a simplified 14generator model of the Australian National Energy Market (NEM) as a test grid [22]. Four cases for different RES penetration are run for one to seven days horizon length, and computational metrics are reported. To reflect the accuracy of the proposed MST, system inertia and voltage stability margins are used as a benchmark. In simulations, RES and load traces are taken from the National Transmission Network Developed Plan (NTNDP) data, provided by the Australian Energy Market Operator (AEMO) [23].
The remainder of the paper is organized as follows: Literature review and related work are discussed in Section II, while Section III details the MST. A detailed description of the simulation setup is given in Section IV. In Section V results are analyzed and discussed in detail. Finally, Section VI concludes the paper.
2 Related Work
In order to better explain the functional requirements of the proposed MST, we first describe the canonical UC formulation. An interested reader can find a comprehensive literature survey in [24].
2.1 Canonical Unit Commitment Formulation
The UC problem is an umbrella term for a large class of problems in power system operation and planning whose objective is to schedule and dispatch power generation at minimum cost to meet the anticipated demand, while meeting a set of systemwide constraints. In smart grids, problems with a similar structure arise in the area of energy management, and they are sometimes also called UC [25]. Before deregulation, UC was used in vertically integrated utilities for generation scheduling to minimize production costs. After deregulation, UC has been used by system operators to maximize social welfare, but the underlying optimization model is essentially the same.
Mathematically, UC is a largescale, nonlinear, mixedinteger optimization problem under uncertainty. With some abuse of notation, the UC optimization problem can be represented in the following compact formulation [26]:
(1)  
(2)  
(3)  
(4)  
Due to the timecouplings, the UC problem needs to be solved over a sufficiently long horizon. The decision vector for each time interval consist of continuous and binary variables. The continuous variables, , include generation dispatch levels, load levels, transmission power flows, storage levels, and transmission voltage magnitudes and phase angles. The binary variables, , includes scheduling decisions for generation and storage, and logical decisions that ensure consistency of the solution. The objective (1) captures the total production cost, including fuel costs, startup costs and shutdown costs. The constraints include, respectively: dispatch related constraints such as energy balance, reserve requirements, transmission limits, and ramping constraints (2); commitment variables, including minimum up and down, and startup/shutdown constraints (3); and constraints coupling commitment and dispatch decisions, including minimum and maximum generation capacity constraints (4).
The complexity of the problem stems from the following: (i) certain generation technologies (e.g. coalfired steam units) require long startup and shutdown times, which requires a sufficiently long solution horizon; (ii) generators are interconnected, which introduces couplings through the power flow constraints; (iii) on/off decisions introduce a combinatorial structure; (iv) some constraints (e.g. AC load flow constraints) and parameters (e.g. production costs) are nonconvex; and (v) the increasing penetration of variable renewable generation and the emergence of demandside technologies introduce uncertainty. As a result, a complete UC formulation is computationally intractable, so many approximations and heuristics have been proposed to strike a balance between computational complexity and functional requirements. For example, power flow constraints can be neglected altogether (a copper plate model), can be replaced with simple network flow constraints to represent critical interconnectors, or, instead of (nonconvex) AC, a simplified (linear) DC load flow is used.
2.2 UC Formulations in Existing Future Grid Studies
In operational studies: the nonlinear constraints, e.g. ramping, minimum up/down time (MUDT) and thermal limits are typically linearized; startup and shutdown exponential costs are discretized, and; nonconvex and nondifferentiable variable cost functions are expressed as piecewise linear function [27, 20]. In planning studies, due to long horizon lengths, the UC model is simplified even further. For example: combinatorial structure is avoided by aggregating all the units installed at one location [21, 28, 29]; piecewise linear cost functions and constraints are represented by one segment only; some costs (e.g. startup, shutdown and fix costs) are ignored; a deterministic UC with perfect foresight is used, and; noncritical binding constraints are omitted [30, 31]
In contrast to operation and planning studies, the computational burden of future grid scenario analysis is even bigger, due to a sheer number of scenarios that need to be analyzed, which requires further simplifications. For example, the Greenpeace study [11] uses an optimal power flow for generation dispatch and thus ignores UC decisions. Unlike the Greenpeace study, the Irish All Island Grid Study [34] and the European project eHighway2050 [35] ignore load flow constraints altogether, however they do use a rolling horizon UC, with simplifications. The Irish study, for example doesn’t put any restriction on the minimum number of online synchronous generators to avoid RES spillage, and the eHighway2050 study uses a heuristics to include DR. The authors of the eHighway2050 study, however, acknowledge the size and the complexity of the optimization framework in long term planning, and plan to develop new tools with a simplified network representation [35].
In summary, a UC formulation depends on the scope of the study. Future grid studies that explicitly include stability assessment bring about some specific requirements that are routinely neglected in the existing UC formulations, as discussed next.
3 Market Simulation Tool
3.1 Functional Requirements
The focus of our work is stability assessment of future grid scenarios. Thus, MST must produce dispatch decisions that accurately capture the kinetic energy stored in rotating masses (inertia), active power reserves and reactive power support capability of synchronous generators, which all depend upon the number of online units and the respective dispatch levels.
For the sake of illustration, consider a generation plant consisting of three identical (synchronous) thermal units, with the following characteristics: (i) constant terminal voltage of \SI1pu; (ii) minimum technical limit ; (iii) power factor of ; (iv) maximum excitation limit ; and (v) normalized inertia constant . We further assume that in the overexcited region, the excitation limit is the binding constraint, as shown in Fig. 2. Observe that the maximum reactive power capability depends on the active power generated, and varies between at and at . We consider three cases defined by the total active power generation of the plant: (i) \SI0.8pu, (ii) \SI1.2pu, and (iii) \SI1.6pu. The three scenarios correspond to the rows in Fig. 3, which shows the active power dispatch level , reactive power support capability , online active power reserves , and generator inertia . The three columns show feasible solutions for three different UC formulations: all three units are aggregated into one equivalent unit (AGG), standard binary UC (BUC) when each unit is modeled individually, and the proposed market simulation tool (MST). A detailed comparison of the three formulations is given in Section V.
Although the results are selfexplanatory, a few things are worth emphasizing. In case (i), aggregating the units into one equivalent unit (AGG) results in the unit being shut down due to the minimum technical limit. The individual unit representation (BUC), on the other hand, does allow the dispatch of one or two units, but with significantly different operational characteristics. In cases (ii) and (iii), the total inertia in the AGG formulation is much higher, which has important implications for frequency stability. A similar observation can be made for the reactive power support capability, which affects voltage stability. Also, dispatching power from all three units results in a significantly higher active power reserve. And last, a higher reactive power generation due to a lower reduces the internal machine angle, which improves transient stability.
In conclusion, a faithful representation of the number of online synchronous machines is of vital importance for stability assessment. An individual unit representation, however, is computationally expensive, so the computational burden should be reduced, as discussed in the following section. Next, an explicit network representation is required. An AC load flow formulation, however, is nonlinear (and nonconvex), which results in an intractable mixedinteger nonlinear problem. Therefore, we use a DC load flow representation with a sufficiently small voltage angle difference on transmission lines. Our experience shows that an angle difference of \SI30\degree results in a manageable small number of infeasible operating conditions that can be dealt with separately.
3.2 Computational Speedup
The MST is based on the UC formulation using constant fixed, startup, shutdown and production costs. To improve its computational efficiency, the dimensionality of the optimization problem is reduced employing: (i) unit clustering [21] to reduce the number of variables needed to represent a multiunit generation plant; (ii) a rolling horizon approach [30, 25, 36] to reduce the time dimension; and (iii) constraint clipping to remove most nonbinding constraints.
Unit Clustering
Linearized UC models are computationally efficient for horizons of up to a few days, which makes them extremely useful for operational studies. For planning studies, however, where horizon lengths can be up to a year, or more, these models are still computationally too expensive. Our work builds on the clustering approach proposed in [21], where identical units at each generation plant are aggregated by replacing binary variables with fewer integer variables. The status of online units, startup/shutdown decisions and dispatched power are tracked by three integer variables and one continuous variable per plant per period, as opposed to three binary and one continuous variable per unit per period. Further clustering proposed in [21] is not possible in our formulation because of the explicit network representation required in the MST.
Rolling Horizon
Solving the UC as one block, especially for long horizons, is computationally too expensive. This can be overcome by breaking the problem into several smaller intervals called subhorizons [30, 25, 36]. To ensure accuracy and consistency of the solution, a proper overlap between subhorizons is maintained and the terminating state of the previous subhorizon is used as the initial condition of the next subhorizon. The minimum subhorizon length depends on the time constants associated with the decision variables. While these might be in the order of hours for thermal power plants, they can be significantly longer for energy storage. Largescale hydro dams, for example, require horizon lengths of several weeks, or even months. In our research, however, the subhorizon length is up to a few days to cater for thermal energy storage (TES) of CST plants and battery storage. The optimization of hydro dams is not explicitly considered, however it can be taken into account heuristically, if needed.
Constraint Clipping
The size of the problem can be reduced by removing nonbinding constraints, which doesn’t affect the feasible region. For instance, an MUDT constraint on a unit with an MUDT less than the time interval is redundant
3.3 MST UC Formulation
Objective function
The objective of the proposed MST is to minimize total generation cost for all subhorizons :
(5) 
where are the decision variables of the problem, and , , , and are fixed, startup, shutdown and variable cost, respectively. As typically done in planning studies [21], [33], the costs are assumed constant to reduce the computational complexity. The framework, however, also admits a piecewise linear approximation proposed in [20].
System Constraints
System Constraints
Power balance: Power generated at node must be equal to the node power demand plus the net power flow on transmission lines connected to the node:
(6) 
where represent respectively the set of generators, consumers, prosumers
Power reserves: To cater for uncertainties, active power reserves provided by synchronous generation are maintained in each region :
(7) 
For synchronous generators other than concentrated solar thermal (CST), reserves are defined as the difference between the online capacity and the current operating point. For CST, reserves can either be limited by their online capacity or energy level of their thermal energy system (TES). Variable in (7) represents the total number of online units at each generation plant, and and represent the sets of generators and nodes in region , respectively.
Minimum synchronous inertia requirement: To ensure frequency stability, a minimum level of inertia provided by synchronous generation must be maintained at all times (more details are available in [14]) in each region :
(8) 
Network constraints
Network constraints include DC power flow constraints and thermal line limits for AC lines, and active power limits for HVDC lines.
Line power constraints:
A DC load flow model is used for computational simplicity for AC transmission lines
(9) 
where the variables and represent voltage angles at nodes and , respectively.
Thermal line limits: Power flows on all transmission lines are limited by the respective thermal limits of line :
(10) 
where represents the thermal limit of line .
Generation constraints
Generation constraints include physical limits of individual generation units. For the binary unit commitment (BUC), we adopted a UC formulation requiring three binary variables per time slot (on/off status, startup, shutdown) to model an individual unit. In the MST, identical units of a plant are clustered into one individual unit [21]. This requires three integer variables (on/of status, startup, and shutdown) per generation plant per time slot as opposed to three binary variables per generation unit per time slot in the BUC, as discussed in Section III.B of A Computationally Efficient Market Model for Future Grid Scenario Studies.
Generation limits: Dispatch levels of a synchronous generator are limited by the respective stable operating limits:
(11) 
The power of RES
(12) 
Unit on/off constraints: A unit can only be turned on if and only if it is in off state and vice versa:
(13) 
In a rolling horizon approach, consistency between adjacent time slots is ensured by:
(14) 
where is the initial number of online units of generator . Equations (13) and (14) also implicitly determine the upper bound of and in terms of changes in .
Number of online units: Unlike the BUC, the MST requires an explicit upper bound on status variables:
(15) 
Rampup and rampdown limits: Ramp rates of synchronous generation should be kept within the respective rampup (16), (17) and rampdown limits (18), (19):
(16)  
(17)  
(18)  
(19) 
In the MST, a ramp limit of a power plant is defined as a product of the ramp limit of an individual unit and the number of online units in a power plant . If is binary, these ramp constraints are mathematically identical to ramp constraints of the BUC. If a ramp rate multiplied by the length of the time resolution is less than the rated power, the rate limit has no effect on the dispatch, so the corresponding constraint can be eliminated. Constraints explicitly defined for are used to join two adjacent subhorizons in the rollinghorizon approach.
Minimum up and down times: Steam generators must remain on for a period of time once turned on (minimum up time):
(20)  
(21) 
Similarly, they must not be turned on for a period of time once turned off (minimum down time):
(22)  
(23) 
Similar to the rate limits, if the minimum up and down times are smaller than the time resolution , the corresponding constraints can be eliminated. Due to integer nature of discrete variables in the MST, the definition of the MUDT constraints in the RH approach requires the number of online units for the last time interval to establish the relationship between the adjacent subhorizons. If the is smaller than time resolution , then these constraints can be eliminated.
CST constraints:
CST constraints include TES energy balance and storage limits.
TES state of charge (SOC) determines the TES energy balance subject to the accumulated energy in the previous time slot, thermal losses, thermal power provided by the solar farm and electrical power dispatched from the CST plant:
(24)  
(25) 
where, is the thermal power collected by the solar field of generator .
TES limits: Energy stored is limited by the capacity of a storage tank:
(26) 
Utility storage constraints
Utilityscale storage constraints include energy balance, storage capacity limits and power flow constraints. The formulation is generic and can capture a wide range of storage technologies.
Utility storage SOC limits determine the energy balance of storage plant :
(27)  
(28) 
Utility storage capacity limits: Energy stored is limited by the capacity of storage plant :
(29) 
Charge/discharge rates limit the charge and discharge powers of storage plant :
(30) 
where and represent the maximum power discharge and charge rates of a storage plant, respectively.
Prosumer subproblem
The prosumer subproblem captures the aggregated effect of prosumers. It is modeled using a bilevel framework in which the upperlevel unit commitment problem described above minimizes the total generation cost, and the lowerlevel problem maximizes prosumers’ selfconsumption. The coupling is through the prosumers’ demand, not through the electricity price, which renders the proposed model market structure agnostic. As such, it implicitly assumes a mechanism for demand response aggregation. The KarushKuhnTucker optimality conditions of the lowerlevel problem are added as the constraints to the upperlevel problem, which reduces the problem to a single mixed integer linear program.
The model makes the following assumptions: (i) the loads are modeled as price anticipators; (ii) the demand model representing an aggregator consists of a large population of prosumers connected to an unconstrained distribution network who collectively maximize selfconsumption; (iii) aggregators do not alter the underlying power consumption of the prosumers; and (iv) prosumers have smart meters equipped with home energy management systems for scheduling of the PVbattery systems, and, a communication infrastructure is assumed that allows a twoway communication between the grid, the aggregator and the prosumers. More details can be found in [15].
Prosumer Objective function: Prosumers aim to minimize electricity expenditure:
(31) 
where is the applicable feedin price ratio. In our research, we assumed , which corresponds to maximization of selfconsumption.
The prosumer subproblem is subject to the following constraints:
Prosumer power balance: Electrical consumption of prosumer , consisting of grid feedin power, , underlying consumption, , and battery charging power, , is equal to the power taken from the grid, , plus the power generated by the PV system, :
(32) 
Battery charge/discharge limits: Battery power should not exceed the charge/discharge limits:
(33) 
where and represent the maximum power discharge and charge rates of the prosumer’s battery, respectively.
Battery storage capacity limits: Energy stored in a battery of prosumer should always be less than its capacity:
(34) 
Battery SOC limits: Battery SOC is the sum of the power inflow and the SOC in the previous period:
(35)  
(36) 
where represents the initial SOC and is used to establish the connection between adjacent subhorizons.
4 Simulation Setup
The case studies provided in this section compare the computational efficiency of the proposed MST with alternative formulations. For detailed studies on the impact of different technologies on future grids, an interested reader can refer to our previous work [19, 17, 15, 14, 16, 18].
4.1 Test System
We use a modified 14generator IEEE test system that was initially proposed in [22] as a test bed for smallsignal analysis. The system is loosely based on the Australian National Electricity Market (NEM), the interconnection on the Australian eastern seaboard. The network is stringy, with large transmission distances and loads concentrated in a few load centres. Generation, demand and the transmission network were modified to meet future load requirements. The modified model consists of 79 buses grouped into four regions, 101 units installed at 14 generation plants and 810 transmission lines.
4.2 Test Cases
To expose the limitations of the different UC formulations, we have selected a typical week with sufficiently varying operating conditions. Four diverse test cases with different RES penetrations are considered. First, RES0 considers only conventional generation, including hydro, black coal, brown coal, combined cycle gas and open cycle gas. The generation mix consists of \SI2.31\giga\watt hydro, \SI39.35\giga\watt of coal and \SI5.16\giga\watt of gas, with the peak load of \SI36.5\giga\watt. To cater for demand and generation variations, \SI10\percent reserves are maintained at all times. The generators are assumed to bid at their respective short run marginal costs, based on regional fuel prices [37].
Cases RES30, RES50, RES75 consider, respectively, \SI30\percent, \SI50\percent and \SI75\percent annual energy RES penetration, supplied by wind, PV and CST. Normalized power traces for PV, CST and wind farms (WFs) for the 16zones of the NEM are taken from the AEMO’s planning document [23]. The locations of RESs are loosely based on the AEMO’s 100% RES study [10].
4.3 Modeling Assumptions
Power traces of all PV modules and wind turbines at one plant are aggregated and represented by a single generator. This is a reasonable assumption given that PV and WF don’t provide active power reserves, and are not limited by ramp rates, MUDT, and startup and shutdown costs, which renders the information on the number of online units unnecessary.
Also worth mentioning is that RES can be modeled as negative demand, which can lead to an infeasible solution. Modeling RES (wind and solar PV) as negative demand is namely identical to preventing RES from spilling energy. Given the high RES penetration in future grids, we model RES explicitly as individual generators. Unlike solar PV and wind, CST requires a different modeling approach. Given that CST is synchronous generation it also contributes to spinning reserves and system inertia. Therefore, the number of online units in a CST plant needs to be modeled explicitly.
An optimality gap of \SI1% was used for all test cases. Simulation were run on Dell OPTIPLEX 9020 desktop computer with Intel(R) Core(TM) i74770 CPU with \SI3.40\giga\hertz clock speed and \SI16\gigaB RAM.
5 Results and Discussion
To showcase the computational efficiency of the proposed MST, we first benchmark its performance for different horizon lengths against the BUC formulation employing three binary variables per unit per time slot and the AGG formulation where identical units at each plant are aggregated into a single unit, which requires three binary variables per plant per time slot. We pay particular attention to the techniques used for computational speedup, namely unit clustering, rolling horizon, and constraint clipping. Last, we compare the results of the proposed MST with BUC and AGG formulations for voltage and frequency stability studies.
5.1 Binary Unit Commitment (BUC)
We first run the BUC for horizon lengths varying from one to seven days, Fig. 4 (top).
As expected, with the increase in the horizon length, the solution time increases exponentially. For a sevenday horizon, the solution time is as high as \SI25000\second (\SI7\hour). Observe how the computational burden is highly dependent on the RES penetration. The variability of the RES results in an increased cycling of the conventional thermal fleet, which increases the number of on/off decisions and, consequently the computational burden. In addition to that, a higher RES penetration involves an increased operation of CST. This poses an additional computational burden due to the decision variables associated with TES that span several time slots. In summary, the computational burden of the BUC renders it inappropriate for scenario analysis involving extended horizons.
5.2 Aggregated Formulation (AGG)
Aggregating identical units at a power plant into a single unit results in a smaller number of binary variables, which should in principle reduce the computational complexity. Fig. 4 confirms that this is mostly true, however, for RES50HL7 the computation time is higher than in the BUC formulation. The reason for that is that, in this particular case, the BUC formulation has a tighter relaxation than the AGG formulation and, consequently, a smaller root node gap. Compared to the MST formulation, with a similar number of variables than the AGG formulation, the MST has considerably shorter computation time due to a smaller root node gap.
In terms of accuracy, the AGG formulation works well for balancing studies [19, 18]. On the other hand, the number of online synchronous generators in the dispatch differs significantly from the BUC, which negatively affects the accuracy of voltage and frequency stability analysis, as shown later. Due to a large number of online units in a particular scenario, a direct comparison of dispatch levels and reserves from each generator is difficult. Therefore, we compare the total number of online synchronous generators, which serves as a proxy to the available system inertia. Fig. 5 shows the number of online generators of four different RES penetration levels for a horizon length of seven days. For most of the hours there is a significant difference between the number of online units obtained from the BUC and the AGG formulation.
In conclusion, despite its computational advantages, the AGG formulation is not appropriate for stability studies due to large variations in the number of online synchronous units in the dispatch results. In addition to that, the computational time is comparable to the BUC in some cases.
We now evaluate the effectiveness of the techniques for the computational speedup.
Unit Clustering
In unit clustering, binary variables associated with the generation unit constraints are replaced with a smaller number of integer variables, which allows aggregating several identical units into one equivalent unit, but with the number of online units retained. This results in a significant reduction in the number of variables and, consequently, in the computational speedup. Compared to the BUC, the number of variables in the MST with this technique alone reduces from \SI24649 to \SI5990 for RES75 with a horizon length of seven days. Therefore, the solution time for RES75HL7 reduces from \SI25000\second in the BUC to \SI450\second in MST with unit clustering alone.
Rolling Horizon Approach
A rolling horizon approach splits the UC problem into shorter horizons. Given the exponential relationship between the computational burden and the horizon length, as discussed in Section 5.1, solving the problem in a number of smaller chunks instead of in one block results in a significant computational speedup. The accuracy and the consistency of the solution are maintained by having an appropriate overlap between the adjacent horizons. However, the overlap depends on the time constants of the problem. Long term storage, for example, might require longer solution horizons. The solution times for different RES penetrations are shown in Table 1. Observe that in the RES75 case, the effect of rolling horizon is much more pronounced, which confirms the validity of the approach for studies with high RES penetration.
RES0  RES30  RES50  RES75  
(minutes)  (minutes)  (minutes)  (minutes)  
BUC MB 7  6.92  12.95  37.11  415.25 
AGG MB 7  6.81  27.08  154.27  27.37 
MST MB 7  2.12  3.34  4.73  5.32 
BUC RH 2+1  2.38  4.03  24.18  74.70 
AGG RH 2+1  0.15  0.20  0.27  0.25 
MST RH 2+1  0.35  0.71  0.60  0.76 
Constraint Clipping
Eliminating non binding constraints can speedup the computation even further. Table 2 shows the number of constraints for different scenarios with and without constraint clipping. Observe that the number of redundant constraints is higher in scenarios with a higher RES penetration. The reason is that a higher RES penetration requires more flexible gas generation with ramp rates shorter than the time resolution (one hour in our case). Note that the benefit of constraint clipping with a shorter time resolution will be smaller.
Case  Number of Constraints  

without CC  with CC  % reduction  
RES0  64555  61332  4.99 
RES30  62520  56617  9.44 
RES50  62500  56777  9.15 
RES75  62740  57017  9.12 
5.3 MST Computation Time and Accuracy
The proposed MST outperforms the BUC and AGG in terms of the computational time by several orders of magnitude, as shown in Fig. 4 (bottom). The difference is more pronounced at higher RES penetration levels. For RES75, the MST is more than \SI500 times faster than the BUC. In terms of the accuracy, the MST results are almost indistinguishable from the BUC results, as evident from Fig. 5 that shows the number of online synchronous units for different RES penetration levels. Minor differences in the results stem from the nature of the optimization problem. Due to its mixedinteger structure, the problem is nonconvex and has therefore several local optima. Given that the BUC and the MST are mathematically not equivalent, the respective solutions might not be exactly the same. The results are nevertheless very close, which confirms the validity of the approach for the purpose of scenario analysis. The loadability and inertia results presented later further support this conclusion.
5.4 Stability Assessment
To showcase the applicability of the MST for stability assessment, we analyze system inertia and loadability that serve as a proxy to frequency and voltage stability, respectively. More detailed stability studies are covered in our previous work, including smallsignal stability [13], frequency stability [14], and voltage stability [17].
System inertia
Fig. 6 (bottom) shows the system inertia for the BUC, AGG and the proposed MST, respectively, for RES0. Given that the inertia is the dominant factor in the frequency response of a system after a major disturbance, the minuscule difference between the BUC and the MST observed in Fig. 6 validates the suitability of the MST for frequency stability assessment. The inertia captured by the AGG, on the other hand, is either over or under estimated and so does not provide a reliable basis for frequency stability assessment.
Loadability Analysis
The dispatch results from the MST are used to calculate power flows, which are then used in loadability analysis
6 Conclusion
This paper has proposed a computationally efficient electricity market simulation tool based on a UC problem suitable for future grid scenario analysis. The proposed UC formulation includes an explicit network representation and accounts for the uptake of emerging demand side technologies in a unified generic framework while allowing for a subsequent stability assessment. We have shown that unit aggregation, used in conventional planningtype UC formulations to achieve computational speedup, fails to properly capture the system inertia and reactive power support capability, which is crucial for stability assessment. To address this shortcoming, we have proposed a UC formulation that models the number of online generation units explicitly and is amenable to a computationally expensive timeseries analysis required in future grid scenario analysis. To achieve further speedup, we use a rolling horizon approach and constraint clipping.
The effectiveness of the computational speedup techniques depends on the problem structure and the technologies involved so the results cannot be readily generalized. The computational speedup varies between \SI20 to more than \SI500 times, for a zero and \SI75% RES penetration, respectively, which can be explained by a more frequent cycling of the conventional thermal units in the highRES case. The simulation results have shown that the computational speedup doesn’t jeopardize the accuracy. Both the number of online units that serves as a proxy for the system inertia and the loadability results are in close agreement with more detailed UC formulations, which confirms the validity of the approach for long term future grid studies, where one is more interested in finding weak points in the system rather than in a detailed analysis of an individual operating condition.
Footnotes
 We interpret a future grid to mean the study of national grid type structures with the transformational changes over the longterm out to 2050.
 For the sake of brevity, by RES we mean “unconventional” renewables like wind and solar, but excluding conventional RES, like hydro, and dispatchable unconventional renewables, like concentrated solar thermal.
 An interested reader can refer to [32] for a discussion on binding constraints elimination for generation planning.
 This is especially the case when the time resolution is coarse. In our studies, the time step is one hour. In operational studies, where the resolution can be as short as five minutes, constraint clipping is less useful.
 All the constraints must be satisfied in all time slots , however, for sake of notational brevity, this is not explicitly mentioned.
 Priceresponsive users equipped with smallscale PVbattery systems.
 A sufficiently small () voltage angle difference over a transmission line is used to reduce the number of nonconvergent AC power flow cases.
 For the sake of brevity, by RES we mean “unconventional” renewables like wind and solar, but excluding conventional RES, like hydro, and dispatchable unconventional renewables, like concentrated solar thermal.
 The loadability analysis is performed by uniformly increasing the load in the system until the load flow fails to converge. The loadability margin is calculated as the difference between the base system load and the load in the last convergent load flow iteration.
References
 Australian Energy Market Operator (AEMO), “Black System South Australia 28 September 2016. Third Preliminary Report.” Report, 2016.
 L. Fahey and R. M. Randall, Learning From the Future. Wiley, 1998.
 J. Foster, C. Froome, C. Greig, O. HoeghGuldberg, P. Meredith, L. Molyneaus, T. Saha, L. Wagner, and B. Ball, “Delivering a competitive Australian power system Part 2: The challenges, the scenarios,” Report, 2013.
 G. Sanchis, “eHighway2050: Europe’s future secure and sustainable electricity infrastructure. Project results.” Report, 2015.
 B. Elliston, J. Riesz, and I. MacGill, “What cost for more renewables? The incremental cost of renewable generation â An Australian National Electricity Market case study,” Renewable Energy, vol. 95, pp. 127 – 139, 2016.
 M. Wright and P. Hearps, “Zero Carbon Australia Stationary Energy Plan,” The University of Melbourne Energy Research Institute, Report, 2010.
 B. Elliston, I. MacGill, and M. Diesendorf, “Least cost 100% renewable electricity scenarios in the Australian National Electricity Market,” Energy Policy, vol. 59, pp. 270–282, 2013.
 C. Budischak, D. Sewell, H. Thomson, L. Mach, D. E. Veron, and W. Kempton, “Costminimized combinations of wind power, solar power and electrochemical storage, powering the grid up to 99.9% of the time,” Journal of Power Sources, vol. 225, pp. 60–74, 2013.
 I. G. Mason, S. C. Page, and A. G. Williamson, “A 100% renewable electricity generation system for New Zealand utilising hydro, wind, geothermal and biomass resources,” Energy Policy, vol. 38, pp. 3973–3984, Aug. 2010.
 AEMO, “100 Per Cent Renewable Study  Modelling Outcomes,” Report, 2013.
 Greenpeace Germany, “PowE[R] 2030  A European Grid for 3/4 Renewable Energy by 2030,” Report, 2014.
 H. Marzooghi, D. J. Hill, and G. Verbič, “Performance and stability assessment of future grid scenarios for the Australian NEM,” in 2014 Australasian Universities Power Engineering Conference (AUPEC). IEEE, Sep. 2014, pp. 1–6.
 R. Liu, G. Verbič, and J. Ma, “Fast Stability Scanning for Future Grid Scenario Analysis,” 2016. [Online]. Available: http://arxiv.org/abs/1701.03436
 A. S. Ahmadyar, S. Riaz, G. Verbič, J. Riesz, and A. Chapman, “Assessment of minimum inertia requirement for system frequency stability,” in 2016 IEEE International Conference on Power System Technology (POWERCON), Sep. 2016.
 H. Marzooghi, S. Riaz, A. C. Chapman, G. Verbič, and D. J. Hill, “Generic demand modelling considering the impact of priceresponsive users for future grid scenarios,” 2016. [Online]. Available: http://arxiv.org/abs/1605.05833
 H. Marzooghi, G. Verbič, and D. J. Hill, “Aggregated demand response modelling for future grid scenarios,” Sustainable Energy, Grids and Networks, vol. 5, pp. 94–104, 2016.
 S. Riaz, H. Marzooghi, G. Verbič, A. C. Chapman, and D. J. Hill, “Impact study of prosumers on loadability and voltage stability of future grids,” in 2016 IEEE International Conference on Power System Technology (POWERCON), Sep. 2016.
 S. Riaz, A. C. Chapman, and G. Verbič, “Evaluation of concentrated solarthermal generation for provision of power system flexibility,” in Proc. Power Syst. Comput. Conf, Conference Proceedings.
 ——, “Comparing utility and residential battery storage for increasing flexibility of power systems,” in Power Engineering Conference (AUPEC), 2015 Australasian Universities. IEEE, Conference Proceedings.
 M. CarriÃ³n and J. M. Arroyo, “A computationally efficient mixedinteger linear formulation for the thermal unit commitment problem,” Power Systems, IEEE Transactions on, vol. 21, no. 3, pp. 1371–1378, 2006.
 B. S. Palmintier and M. D. Webster, “Heterogeneous unit clustering for efficient operational flexibility modeling,” IEEE Transactions on Power Systems, vol. 29, no. 3, pp. 1089–1098, 2014.
 M. Gibbard and D. Vowles, “Simplified 14generator model of the SE Australian power system,” The University of Adelaide, Report, 2010.
 Australian Energy Market Operator (AEMO), “National Transmission Network Development Plan 2016,” Report, 2016.
 M. Tahanan, W. van Ackooij, A. Frangioni, and F. Lacalandra, “Largescale Unit Commitment under uncertainty,” 4OR  A Quarterly Journal of Operations Research, vol. 13, no. 2, pp. 115–171, Jan. 2015.
 F. Ramos, C. Cañizares, and K. Bhattacharya, “Effect of price responsive demand on the operation of microgrids,” in Proc. Power Syst. Comput. Conf, Conference Proceedings, pp. 1–7.
 D. Bertsimas, E. Litvinov, X. A. Sun, J. Zhao, and T. Zheng, “Adaptive Robust Optimization for the Security Constrained Unit Commitment Problem,” IEEE Transactions on Power Systems, vol. 28, no. 1, pp. 52–63, Feb. 2013.
 J. M. Arroyo and A. J. Conejo, “Optimal response of a thermal unit to an electricity spot market,” IEEE Transactions on Power Systems, vol. 15, no. 3, pp. 1098–1104, 2000.
 K. Hara, M. Kimura, and N. Honda, “A method for planning economic unit commitment and maintenance of thermal power systems,” IEEE Transactions on Power Apparatus and Systems, vol. PAS85, no. 5, pp. 427–436, 1966.
 N. Langrene, W. v. Ackooij, and F. Breant, “Dynamic constraints for aggregated units: Formulation and application,” IEEE Transactions on Power Systems, vol. 26, no. 3, pp. 1349–1356, 2011.
 A. Tuohy, P. Meibom, E. Denny, and M. O’Malley, “Unit commitment for systems with significant wind penetration,” Power Systems, IEEE Transactions on, vol. 24, no. 2, pp. 592–601, 2009.
 B. Palmintier and M. Webster, “Impact of unit commitment constraints on generation expansion planning with renewables,” in 2011 IEEE Power and Energy Society General Meeting, Conference Proceedings.
 B. Palmintier, “Flexibility in generation planning: Identifying key operating constraints,” in Power Systems Computation Conference (PSCC), 2014, Conference Proceedings, pp. 1–7.
 L. Zhang, T. Capuder, and P. Mancarella, “Unified unit commitment formulation and fast multiservice lp model for flexibility evaluation in sustainable power systems,” IEEE Transactions on Sustainable Energy, vol. 7, no. 2, pp. 658–671, 2016.
 EirGrid and System Operator for Northern Ireland (SONI), “All Island TSO Facilitation of Renewables Studies,” Report, 2010.
 C. Pache, P. Panciatici, S. Lumbreras, F. Echavarren, and L. Sigrist, “eHighway2050 D8.7: Recommendations about critical aspects in longterm planning methodologies,” Report, 2015.
 IRE, “All Island Grid Study,” Report, 2008.
 ACIL Tasman, “Fuel resource, new entry and generation costs in the NEM,” Report, 2009.