A discussion on electricity prices, or the two sides of the coin

A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 Abstract. We examine how different pricing frameworks deal with nonconvex features typical of day-ahead energy prices when the power system is hydrodominated, like in Brazil. For the system operator, requirements of minimum generation translate into feasibility issues that are fundamental to carry the generated power through the network. When utilities are remunerated at a price depending on Lagrange multipliers computed for a system with fixed commitment, the corresponding values sometimes fail to capture a signal that recovers costs. Keeping in mind recent discussions for the Brazilian power system, we analyze mechanisms that provide a compromise between the needs of the generators and those of the system operator. After characterizing when a price supports a generation plan, we explain in simple terms dual prices and related concepts, such as minimal uplifts and bi-dual problems. We present a new pricing mechanism that guarantees cost recovery to all agents, without over compensations. Instead of using Lagrange multipliers, the price is defined as the solution to an optimization problem. The behavior of the new rule is compared to two other proposals in the literature on illustrative examples, including a small, yet representative, hydro-thermal system. 1. Price signals for power systems Setting prices at appropriate levels is a key driver for success in any business. This is particularly true for the energy sector and day-ahead prices, because of nonconvex features present in unit-commitment (UC) problems. As pointed out in [1], the generation of price signals became more complex nowadays, because distributed and renewable sources increased significantly the scale and complexity of UC problems. In this work we examine the issue under the light of the following through line: Pricing mechanisms in energy systems are like coins, and have two sides Signals for energy prices have the Independent System Operator (ISO) on the heads, and generation agents on the tail side. The latter are in charge of providing energy, the former is responsible for dispatching the generators in a manner that is both reliable and sufficient, so that the electricity flows through the network to meet the demand. Assessing a signal as a mechanism that provides a “good” price will naturally depend on from which side of the coin the appraisal is being made. Generators are concerned about making their business profitable, and expect the price to be high enough to cover the cost of all the “ingredients” necessary to produce electricity, including fixed costs for starting-up and shutting-down a utility. Rather than focusing on individual costs of the generated energy, the ISO has a Key words and phrases. Energy prices, unit-commitment problems, cost recovery, marginal cost. 1 2 JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 global view that may sometimes deviate the dispatch from what is perceived as a “merit” order. The ISO’s main interest is to ensure the generated electricity is carried through the network and reaches the end consumers reliably, day after day. These concerns of profit and feasibility, respectively of the generators and the ISO, define the two sides of our figurative coin. Our work is motivated by discussions held in Brazil before putting in place an official tool to determine the dispatch and prices for the next day in the whole country. Until 2020, both the dispatch and commercialization were done on a weekly basis. The accumulation of a series of changes that had been introduced in the power mix along the years made it necessary to reduce the time granularity significantly, passing to hourly decision-making. Brazil has approximately 130 thermal and 200 hydraulic utilities, and a fast growing segment of wind power, that in 2019 represented 9.7% of the total generated electricity. Additionally, to mitigate the environmental disruption of building new dams and reservoirs, particularly in the Amazon basin, the construction of run-of-river hydro-plants was prioritized in the last decades. The consequent loss of regularization capacity of the system was estimated at 10% in [2, Table 2]. Together with the wind uncertainty, this decreased the system’s flexibility and storage capability, having a direct impact on the role of thermal power in the Brazilian mix. Such considerations were crucial for the implementation of the day-ahead model DESSEM [3], developed by the Brazilian Electric Energy Research Center. Since January 2020 and January 2021, the Brazilian computational UC model is run by two separate national entities, the ISO and the Chamber of Commercialization of Energy (CCE), respectively. To take into account the needs of each agency, the model was formulated into two variants, an electrical one to be used by the ISO, and an energy oriented one, for the CCE. In its electrical version for the ISO, the model makes a DC representation of the network to define the optimal dispatch for each transmission bus. The energy variant, employed by the CCE, adapts the electrical model to reflect the pricing conditions in the country. Energy trading in Brazil is carried out in four regions, through regulated and free contracts that are cleared by the CCE. The second version of the computational tool aggregates the network into sub-markets that can exchange energy within some capacity bounds, with few electrical security constraints. The model provides marginal costs of operation (MCO) for each sub-market that the CCE uses as a basis for the price to be paid for differences between contracted and consumed or generated energy. The goal of reformulating the electrical unit-commitment model into an energy-based counterpart is to device a uniform pricing system, by considering the commercialized energy as equally available in all the consumption points of a given sub-market. For the highly complex hydro-thermal Brazilian system, DESSEM makes an individual representation of each reservoir and generating unit, modeling combined cycle thermal plants, pumping stations, reservoirs connected by water channels, water level constraints at some river sections, and many other features. The resulting UC problem with security constraints yields a mixed-integer linear program (MILP) with several hundreds of thousands variables (25% binary) and a similar number of constraints. In spite of its large scale, the problem is tackled with commercial solvers, combining parallel calculations with several ad-hoc iterative procedures that reduce the computational burden, see [3] and references therein. In what follows, we address from a mathematical perspective concerns raised by A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN 3 agents in the Brazilian system, particularly regarding cost recovery of the pricing system. Rather than dealing with the large-scale official formulation, we focus the discussion on simplified models that are sufficient to illustrate challenges that need to be resolved. A very fundamental question related to marginal costs is what can possibly be considered as being an incremental variation, when turning a unit on or off induces an instantaneous fixed cost. In a system without indivisibilities, a term coined by Scarf in [4], prices are set at the margin, quantifying the value of the next unit to be produced. But a UC problem presents no indivisibilities only when its feasible set satisfies the so-called integrality property [5, Appendix D]. In turn, this amounts to the MILP being equivalent to a convex optimization problem, a model of little practical value for the day-head problem, that typically involves relations such as minimum generation requirements, that make the problem nonconvex. The “IP prices” [6] are currently adopted in Brazil to define the hourly energy prices for the next day. These are multipliers of a dispatch problem that is defined using the output of the UC problem, freezing the commitment status of the utilities. Since no unit can be turned on or off, in some situations, IP prices fall short of the generator’s expectations. Counterintuitive situations that arise with IP prices are illustrated with a simple example in Appendix A. To address this issue we propose a new rule that maintains separated the two sides of the coin. Instead of using Lagrange multipliers, prices are computed directly as the solution of an optimization problem whose feasible set guarantees cost recovery to all the agents. To stabilize the output and avoid volatility observed with the IP prices, the objective function maintains the price close to some reference. The presentation is organized as follows. The two sides of the coin are formalized in mathematical terms by the primal and dual formulations of the ISO’s problem, given respectively in Sections 2 and Section 3 below. As explained in the latter, the difference between those two views, the thickness of the coin, is the duality gap. Furthermore, we show that the primal problem that corresponds to the generator’s preference, bi-dual to the UC problem, can be interpreted as a relaxed version of the ISO’s problem, with randomized decisions. Section 4 starts with the approach in [6], and succinctly reviews the literature. The discussion then focuses on the new cost-recovery rule. The work continues by comparing the merits of some pricing mechanisms for a hydro-thermal stylized system over 24h, and ends with concluding remarks. 2. Unit Commitment problems The UC optimization problem is formulated by the ISO to determine how to dispatch in the short term the energy produced by the generators. Solving the UC problem provides the ISO with a dispatch, an output of primal nature from the optimization point of view. The procedure uses the dispatch to output also an MCO as an important dual indicator. Whether the power system is run under market premises, like in New Zealand or the US, or with centralized dispatch as in Brazil, the MCO gives a price signal, akin to the shadow prices in Linear Programming. We now illustrate the main issues that arise when defining prices in a day-ahead setting, using an idealized UC model, that we named unitoy. The modeling of thermal power plants follows the official tool used in Brazil [3], see also [7]. The simplified systems below are purely thermal, as this suffices to study indivisibilities 4 JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 such as the minimum generation constraints (1). An instance including hydrogeneration is considered in the final section. 2.1. Formulation of the UC problem. The optimization horizon covers T time steps for a system with generation units gathered in the set M. The UC variables are the energy pti ∈ R generated by i at time t, and the commitment uti ∈ {0, 1} that indicates whether unit i at time t is on (uti = 1) or it is off (uti = 0). The overall generation and commitment of the ith unit are the vectors pi = (p1i , . . . , pTi ) and ui = (u1i , . . . , uTi ). Technological constraints are written abstractly as (pi , ui ) ∈ Pi ⊂ R2T . Typical relations in this set are the capacity and ramp constraints, given below:  t pi,min uti ≤ pti ≤ pti,max uti t = 1, . . . , T (1) (pi , ui ) ∈ Pi contains |pti − pt−1 | ≤ ∆pi t = 1, . . . , T , i where p0i is the initial generation level for the ith unit. Binary relations on the commitment are considered separately, by letting  Qi := (pi , ui ) ∈ Pi : ui ∈ {0, 1}T . Each unit has variable generation cost Ci (pti ) as well as fixed operational costs Fi+ and Fi− , the latter being incurred whenever the unit is turned on or off, respectively. The total operational cost for unit i is given by GCosti (pi , ui ) = T  X t=1 +  Ci (pti ) + Fi+ [uti − ut−1 ]+ + Fi− [ut−1 − uti ]+ , i i where [·] = max(·, 0) denotes the positive-part function. Given a system demand D = (D1 , . . . , DT ) ∈ RT , the UC problem is  X  GCosti (pi , ui ) min  (p ,u    i i )i∈M i∈M X (2) pti = Dt t = 1, . . . , T s.t.    i∈M   (pi , ui ) ∈ Qi , i ∈ M. The UC model minimizes the total generation cost while satisfying the demand D, respectively represented by the objective function and the first constraint in (2). This goal must be achieved respecting the second constraint, with rules for the particular features of each source of energy. The technological conditions therein are to be satisfied for one unit i, independently of the behavior of other units in the system. For a hydro-thermal system, if there is more than one hydro-plant along a river, water released by units uphill has an impact on the generation of units that are downhill. In this case the set Pi is defined for the ith cascade as a whole, and the commitment uti is a vector whose components indicate the status of each one of the units along the sequence of hydro-plants. 2.2. The value function. The computational model [3] employed by the ISO and the CCE is responsible for two tasks. The first one is to determine the commitment and levels of dispatch for each unit, by solving a UC problem with a formulation similar to the one in (2). The second one is to define the price of the generated energy. As shown by the brief review in Section 4.1, this second task is not simple, there are several proposals to try to balance the two sides of the coin mentioned A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN 5 in the introduction. The following sections give a mathematical flavor to the main difficulties encountered when addressing the issue for (2). 2.2.1. Theoretical considerations. Let (p∗ , u∗ ) denote an optimal dispatch, that is a generation plan set by the ISO after solving (2). To determine the compatibility between the ISO’s dispatch and the commercialization prices, Brazilian generators examine their payments considering the given price π = (π1 , . . . , πt ) ∈ RT+ as follows. The margin of the ith generator is given by Mi (π, p∗i , u∗i ) = π ⊤ p∗i − GCosti (p∗i , u∗i ) , (3) P where x⊤ y = j xj yj stands for the inner product of two column vectors x and y. The Brazilian cost-minimizer ISO can be interpreted as maximizing the margin for all the agents (the “welfare”). To see this relation, it suffices to add the margins and recall that any feasible dispatch satisfies the demand (assumed to be non-responsive to price): X X GCosti (p∗i , u∗i ) . Mi (π, p∗i , u∗i ) = π ⊤ D − (4) i∈M i∈M The individual view of one utility does not necessarily agree with the global view of welfare adopted by the operator. To start with, the margin (3) may not represent an actual profit for the generator, it might be negative if the price π is too low and insufficient to cover the generation cost. Price-taker generators examine the situation from their side of the coin, and compare the ISO’s dispatch (p∗i , u∗i ) with the level of production that maximizes their margin. For the ith agent this perception amounts to computing the following self-dispatch n π ⊤ pi − GCosti (pi , ui ) . (5) (Pi∗ , Ui∗ ) solving Mmax (π) := (pimax i ,ui )∈Qi The price π is said to support the generation plan when (p∗i , u∗i ) solves (5) for all i ∈ M, that is, when for all the agents there is an agreement between the ISO’s dispatch and the ideal self-dispatch above. While the plan (p∗i , u∗i ) is the view of the ISO, the generation obtained by the unit when solving (5) reflects the individual interest of the generator (individual views are myopic, there is no reason for the demand to be satisfied when adding the decentralized decisions Pi∗ ). Likewise, Mmax (π) is the best possible margin for i the generator and, since the plan (p∗ , u∗ ) is feasible for the maximization problem (5), the relation (6) Mmax (π) − Mi (π, p∗ , u∗ ) ≥ 0 i always holds. The term uplift in the literature refers to adjustments that were proposed as side payments to compensate for a positive difference between these two measures; see Section 44.1. The difference (6), considered a lost opportunity by the generator, provides an estimation of the thickness of our figurative coin. Both sides agree when (2) is convex, a property that is rare for UC problems. For general UC problems our next result characterizes a price π that supports a generation plan through the important inequality (8) defined below. 6 JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 Theorem 2.1 (Characterization of prices supporting generation plans). For the UC problem (2) consider the following value function v : RT → (−∞, +∞]  X  GCosti (pi , ui ) min   (pi ,ui )∈Qi ,i∈M i∈M X (7) v(D) =  pti = Dt , t = 1, . . . , T , s.t.   i∈M and let {(p∗i , u∗i )}i∈M denote a solution to (7). A price π supports (p∗ , u∗ ) for the demand D if and only if (8) v(D′ ) ≥ v(D) + π ⊤ (D′ − D) for all D′ ∈ RT . Proof. Suppose first the price π supports the generation plan, so that all the differences in (6) are null. By definition of the margins, (p∗ , u∗ ) satisfies (4) and, therefore, X X Mmax (π) . Mi (π, p∗i , u∗i ) = − v(D) − π ⊤ D = − i i∈M i∈M Each term of the rightmost summation can be bounded from above by noting that, for any (pi , ui ) ∈ Qi , the maximization in (5) implies that Mmax (π) ≥ π ⊤ pi − i GCosti (pi , ui ). Therefore, X X pi . GCosti (pi , ui ) − π ⊤ v(D) − π ⊤ D ≤ i∈M i∈M P In particular, for all (pi , ui ) ∈ Qi such that i∈M pi = D′ , the inequality holds and, hence, X GCosti (pi , ui ) − π ⊤ D′ . v(D) − π ⊤ D ≤ i∈M Taking the minimum for all such (pi , ui ) yields v(D) − π ⊤ D ≤ v(D′ ) − π ⊤ D′ , and (8) holds. Reciprocally, if π satisfies (8), this means that D is a global minimizer of the shifted function W(D′ ) := v(D′ ) − π ⊤ D′ . Thus, v(D) − π ⊤ D = W (D) = minD′ W (D′ ) = minD′ v(D′ ) − π ⊤ D′ . Expanding the inner product, we see that  X GCosti (pi , ui ) − π ⊤ D′ min v(D) − π ⊤ D = min ′ D P (pi ,ui )∈Qi ,i∈M i∈M t ′ t i∈M = = = pi =(D ) ,1≤t≤T min X GCosti (pi , ui ) − X GCosti (pi , ui ) − (pi ,ui )∈Qi ,i∈M min (pi ,ui )∈Qi i∈M X min i∈M =− (pi ,ui )∈Qi X t=1 i∈M T X t=1 i∈M  T X ⊤ GCosti (pi , ui ) − π pi Mmax (π) , i  πt X i∈M π t pti pti   i∈M where the last equality follows from (3). Combined with (4) the equality results in X (Mmax (π) − Mi (π, p∗i , u∗i )) = 0 , 0≤ i i∈M A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN 7 which, together with (6), ensures that all the terms are zero and therefore shows that π supports the generation plan, as claimed.  Some remarks are in order. First, the result is a characterization, the set of supporting prices for a given demand D is exactly the set in (8), of vectors defining a supporting hyperplane for the value function v at that demand D. When v is convex, such set coincides with the subdifferential  ∂v(D) = π ∈ RT | v(D′ ) ≥ v(D) + π ⊤ (D′ − D) for all D′ ∈ RT . Under very mild assumptions, the subdifferential is nonempty. Hence, whenever v is convex and (7) has a solution, there exists a supporting price. Moreover, subgradients are the Lagrange multipliers associated with the demand constraint of the linear relaxation of (7) (all 0-1 variables relaxed to the interval [0, 1]). Thanks to the characterization given by Theorem 2.1, such Lagrange multipliers are legitimate MCO, or shadow prices. Note however that, unless the value function is differentiable and the subdifferential shrinks to just the singleton derivative, there is an ambiguity on what “the” MCO should be and on which element in the set supports the generation plan (any element in the subdifferential set qualifies as such). This issue goes beyond a mathematical fanciness, see our comments in Remark 2(2.2.3). Second, convexity of the value function is not automatic. The simple examples considered in the next section illustrate nonconvex and even discontinuous situations, related to Scarf’s indivisibilities. For nonconvex v, existence of supporting prices is rare, and this is a consequence of the characterization in Theorem 2.1. Specifically, for π to support the generation plan at D, the affine function v(D) + π ⊤ (· − D) must stay below the value function. If v is locally concave, as on the left in Figure 2 below, the condition fails and no price can support the generation plan. Finally, note that since the feasible set in (7) is bounded, the function values are always finite unless the problem is infeasible, in which case v(D) = +∞ (for example if D is negative, which has not much practical meaning, but is mathematically possible). 2.2.2. Some illustrative instances. Our analysis in Figures 1 and 2 below is done by exploring the dispatch of two units when demand parses the interval [5000, 5300] and plotting the corresponding value function. We include a third unit b in the base, with continuous generation until a capacity of 5000, to eliminate border effects that modify the price if the demand is close to the minimum value in the interval. For the Brazilian system, the unit b represents Angra, the only nuclear power plant in the country. Table 1 has the information for the system that is used as a basis for our illustrations. Table 1. A very simple power system with one unit in the base and two generating units off at departure (Cit (p) = Cit p is linear and T = 1 in (2)) unit Cit 1 5 2 12 b 4 Fi+ = Fi− 0.0 0.0 0.0 pti,min 0.0 0.0 5000 pti,max 150 150 5000 u0i in (1) 0 0 1 In our plots, the abscissa parses the demand in an interval that is beyond the maximum capacity of the base unit (D ≥ 5000). Figure 1 reports with a full line 8 JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 the value function, which is convex with the given data. Shaded areas therein represent the generation level of units 1 and 2. Figure 1. Convex value function for a system as in Table 1 and T = 1 in (7). The change of slope reflects the more expensive variable cost of unit 2, that enters into operation when the capacity of the other units is insufficient to satisfy the demand (D = 5150). By Theorem 2.1, at that point any price between the two slopes supports the generation plan, and there is a full set of MCO. Having positive start-up and shut-down costs and minimum generation requirements in (1) is not uncommon. This results in value functions as those represented in the left and right plots in Figure 2, leading to situations that go against the expected behavior of prices in the energy business (and require adjustments for dispatched generators not to incur losses). The plots were obtained by changing the data from Table 1 to F1+ = F1− = 500 (left) and to F1+ = F1− = 1000, p2,min = 500 (right). On the left, approximately in the range D ∈ (5050, 5150), the change in the slopes indicates a reduction in the MCO, when unit 1 is generating instead of unit 2 (unit 1 has cheaper variable cost but high fixed cost). This drop occurs in spite of an increase in demand, a behavior that goes against the market expectation and could harm the effectiveness of demand-response programs. As remarked after Theorem 2.1, at D = 5075 no price can supports that demand, because v is locally concave. On the right, near D = 5100, unit 2 (with minimal generation but no fixed cost) starts generating instead of unit 1. The discontinuity in the value function results into a MCO taking any value in a set that is unbounded. This issue can produce an infinite signal that translates into unduly large prices in practice. 2.2.3. Remark on ambiguity of prices. Figure 1 represents the easiest configuration of a power system, with a convex value function for which all subgradients support the generation plan. In particular, at D = 5150 any price in [5, 12] will satisfy the property. Even in this most favorable situation, an ambiguity needs to be A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN 9 Figure 2. Positive values for fixed costs and minimal generation introduce nonconvexities and discontinuities in the value function. On the left, no price can support the dispatch at D = 5075. addressed, or at least acknowledged, to avoid bad surprises for the generators. Specifically, suppose that more than one value is a good candidate for a price, as in the figure with D = 5150. Suppose also the value in question is computed as a Lagrange multiplier of some optimization problem, as often in practice. In this case, the output can be different when running the same code with different solvers or on different machines, and all the obtained prices will be legitimate ones. A price that is not uniquely defined can lead generators to make wrong business plans; for instance if the price computed by the CCE is 5 and the one obtained with the same code in the generator’s computer was 12. The ambiguity persists also if the same code is run in one machine and, to accelerate the process, concurrent Linear Programming algorithms are used in parallel by the solver. A simplex method would provide a vertex with value 5 or 12, while with an interior point method the price can be any value in the interval. 3. Bi-dual considerations Generators are interested in a problem that from the optimization point of view is dual to (7). Starting from this UC problem, the ISO’s side of the coin, we now derive a certain dual function. The dual function is separable, and its evaluation amounts to computing the self-dispatch margins (5) for all the generators. These issues are discussed below, as well as examining what happens if the game is played in the opposite direction, computing a certain bi-dual problem for the ISO. 3.1. Generators have a dual point of view. In the UC problem (7), requiring satisfaction of demand couples the decisions of all the generators. Decentralization is achieved by associating a multiplier λ ∈ RT to the demand constraint. The corresponding Lagrangian function  X GCosti (pi , ui ) − λ⊤ pi L(p, u, λ) := λ⊤ D + i∈M X Li (pi , ui , λ) = λ⊤ D + i∈M JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 10 is separable into partial Lagrangians Li (pi , ui , λ) := GCosti (pi ) − λ⊤ pi , depending only on the variables of the ith generator. A problem dual to (7), given by maxi-minimizing the Lagrangian, inherits the separable structure present in the Lagrangian. Specifically, the dual problem is defined as  (9) max θ(λ) : λ ∈ RT , where θ(λ) (10) := min {L(p, Xu, λ) : (pi , ui ) ∈ Qi , i ∈ M} min Li (pi , ui , λ) . = λ⊤ D + i∈M The partial dual functions θi (λ) := (pi ,ui )∈Qi min (pi ,ui )∈Qi Li (pi , ui , λ) in the summation are the negative of the margin maximization problems (5), written with π = λ: ( max π ⊤ pi − GCosti (pi , ui ) (pi ,ui ) −θi (π) = max −Li (pi , ui , π) = = Mmax (π) . i (pi ,ui )∈Qi s.t. (pi , ui ) ∈ Qi Accordingly, a dual price also solves the minimization problem below ) ( ) ( X X max max ⊤ ⊤ Mi (λ) − λ D . Mi (λ) = − min max θ(λ) = max λ D − λ λ λ i∈M i∈M This characterization of the dual price as the one minimizing the discrepancy between the expectation of best margin from the agents and the total payment they receive. In turn, the payment is related to the margins seen by the ISO in (4), so the dual problem is in fact minimizing the differences (6). These relations justify the naming of minimal uplift that is associated with the dual prices in [8], also called convex-hull prices in the electric engineering literature. When there is no discrepancy, the minimal uplifts are all null and in the optimization jargon it is said that there is no duality gap. For the duality gap to be zero, convexity is again of paramount importance. If the primal problem (2) is convex, so is the value function (7). In this case, we claim that a dual price π ∗ supports the optimal dispatch. This can be seen by combining the various definitions with the convexity assumption, to write v(D) = min (pi ,ui )∈Qi ,i∈M max L(p, u, λ) = max λ λ min (pi ,ui )∈Qi ,i∈M L(p, u, λ) = θ(π ∗ ) . By linearity of the Lagrangian with respect to D, this means that, for any D′ ∈ RT , v(D′ ) = max {θ(λ) + λ⊤ (D′ − D)} ≥ θ(π ∗ ) + π ∗ ⊤ (D′ − D) = v(D) + π ∗ ⊤ (D′ − D) . λ Since the last inequality is (8), by Theorem 2.1, the price π ∗ supports the dispatch, as claimed. The dual prices, computed solving (9) for our three examples, are reported in Figure 3, together with the IP prices from [6] presented in Section 4. 3.2. On solution procedures. For systems with many heterogeneous units, the primal UC problem (7) may seem at first glance too difficult to solve. The ISO must decide every day the dispatch for the next day, and calculations need to be done in relatively short time reliably. This is usually addressed following two different approaches. A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN 11 One possibility is to use Lagrangian relaxation as in (9) and solve the dual of the UC described in Subsection 3.1. This is advocated in [9, 10], for example. The main idea is that the dual problem has a nice structure, being amenable to be solved using nonsmooth convex optimization methods. Another major advantage of this approach is that the computational effort scales moderately with the problem size [11]. This property may become increasingly important as the time step length in UC problem decreases, to better deal with renewable uncertainty. Additionally, as mentioned, the computed dual solution gives a price signal with minimal uplifts. In fact, Electricité de France has a long tradition of using Lagrangian relaxation, that goes back to the seminal work [12]. The approach continues to be used nowadays, adapted to the market setting, to compute half-hour dispatch and prices [13]. Due to the constraints on the operational process in France, the problem needs to be solved in less than 10 minutes, with very strong requirements both on optimality and feasibility (all schedules must be feasible and a 1% gap adds up into several millions of euros per year). In this setting, the separable structure of the dual function in (9) can be exploited by an iterative procedure based on decomposition, that maximizes the dual function. Because the dual function is concave but nonsmooth, special techniques must be put in place, that guarantee accuracy and robustness. Bundle algorithms for nonsmooth minimization are the methods of choice in this case, [14, Part II]. We do not enter into further details here. For other success stories of decomposition methods applied to energy optimization, we refer to [15] and the many references therein. Nothing guarantees that the generation associated with the dual price will be feasible for the UC problem. When the dual approach ends, a second procedure must be put in place, to find a feasible dispatch. The process of obtaining a primal feasible point after solving the dual problem is called in the literature Phase 2, feasibility recovery, or primal recovery. There are many proposals, among which we can mention [16, 17, 18, 19, 20], and also the very recent work [21]. In general, primal output associated with the dual price can be used as a starting point for this second phase. On the other hand, as well explained in [22], in the last years commercial solvers took a major step forward in their ability to solve real-world MILPs. The progress was so significant that it is now possible to solve many UC problems directly, using off-the shelf software for a variety of power systems. Notwithstanding, for large-scale systems like the Brazilian one, the security constrained UC problem [3] can be solved in reasonable times only after enhancing the commercial code with smart heuristics to reduce the computational complexity. The effectiveness of heuristics depends largely on exploiting the MILP’s combinatorial structure, for example to eliminate symmetries. Every time the UC problem is enriched with some new feature (say, to model demand response or carbon capture programs), the combinatorial structure changes and a successful heuristics can fail for the new structure. Another possible drawback is that the solution time required by MILP formulations may not scale well with problem size, see [10, 11]. There is always a trade-off to be found when solving hard complex problems, and the choice of the best methodology is often driven by the available know-how. 3.3. What primal problem corresponds to the generators’ preference? Solving the dual problem (9) gives a price deemed acceptable from the generator’s point of view, but not a feasible dispatch, the interest of the ISO. Nevertheless, 12 JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 the dual process outputs some primal pairs (p̄, ū), defined below, that solve a problem that is dual to the dual, called bi-dual. This problem convexifies the original problem by looking at all possible convex combinations. Since the relaxed constraint has T components, by Carathéodory theorem, it suffices to consider convex combinations of at most T + 1 elements. We let the unit simplex ) ( T +1 X ∆ := α ∈ [0, 1]T +1 : αk = 1 , k=1 and set K := {1, . . . , T + 1}, so that the bi-dual problem to (7) is  X X  GCosti (pi (k), ui (k)) min αk   α,(u,p)   i∈M k∈K   s.t. α ∈ ∆ (11) (p ui (k)) ∈ Qi , i ∈ M,k ∈ K  i (k),  X X   t t  α p (k) = D , t = 1, . . . , T .  k i  i∈M k∈K The bi-dual coincides with the original problem if the later has no 0-1 variables and the objective and constraints functions are convex (the optimality conditions are necessary and sufficient to characterize a minimizer). In the ideal convex setting, (11) is the dual of (9) and coincides with our initial UC problem (7): any side of the coin is dual to the other one, there is no duality gap, and the minimal uplifts are all null. By contrast, when there are indivisibilities, the duality gap may be positive and the bi-dual differs from (7). Notwithstanding, an interesting interpretation for the bi-dual is given in [23]. More precisely, when solving (9) with a bundle method [14], the output is an optimal dual price π ∗ and T + 1 points n  o p̄(k), ū(k) and ᾱk : k ∈ K that are feasible for (7), barring the demand constraint Furthermore, the corresponding convex combinations X X p̂i (ᾱ, p̄) := ᾱk p̄i (k) , ûi (ᾱ, ū) := ᾱk ūi (k) k∈K̄ k∈K̄ solve the bi-dual problem (11), see [23, Section 6]. This convex combination, called pseudo-planning in [24], satisfies the demand constraint, but is not feasible for the UC problem (2). Now the commitment component only satisfies ûi ∈ [0, 1]T , as in Table 2 below for D = 5149 and 5151. For operational purposes, an additional step to recover primal feasibility becomes necessary. The proposal [24] finds a feasible point that is closest in some sense to the pseudo-planning, a methodology that reveals superior to the one in [25], that linearizes an augmented Lagrangian. A third alternative is to consider the coefficients αk in (11) as the probability of the ISO taking a randomized decision (p(k), u(k)), satisfying the technological constraints in (7), and the relaxed version of the demand constraint in (11). The objective in (11) represents the expected value of the operational cost; see [26, Section 3.4] and [27]. With this interpretation, the bundle output with largest ᾱk can be a good starting point for recovering feasibility. As an illustration, we applied the bundle method in [28] to solve the dual problem of the system yielding the discontinuous value function on the right in Figure 2. We considered three different demand instances, D ∈ {5149, 5150, 5151}, that capture the behavior in one of the regions where the value function jumps. The respective A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN 13 prices are π = {11.67, 12, 12}, the duality gap is 0.03%, 0.00%, and 0.00% (and indication that the last two problems satisfy the integrality property). The bundle output and corresponding pseudo-plannings are also reported in Table 2. Table 2. Solution to UC problem (2) and output of the bundle method when solving (9) for the system with fixed costs and minimal generation (p∗ , u∗ ) solving (2) pseudo-planning (p̂, û) solving (11) bundle output k = 1 : ᾱ(k) D = 5149 b (5000,1) unit 1 (149,1) unit 2 (0,0) b (5000,1) 2 ) unit 1 (149,1- 300 unit 2 (0, 0) 1− 2 300 (p̄i (k), ūi (k)) b (5000,1) unit 1 (150,1) unit 2 (0,0) k = 2 : ᾱ(k) 2 300 (p̄i (k), ūi (k)) b (5000,1) unit 1 (0,0) unit 2 (0,0) D = 5150 b (5000,1) unit 1 (150,1) unit 2 (0,0) D = 5151 b (5000,1) unit 1 (51,1) unit 2 (100,1) b (5000,1) unit 1 (150,1) unit 2 (0,0) b (5000,1) unit 1 (150,1) unit 2 (1, 0.01) 0 0.99 b (5000,1) unit 1 (150,1) unit 2 (0,0) 1 b (5000,1) unit 1 (150,1) unit 2 (0,0) b (5000,1) unit 1 (150,1) unit 2 (0, 0) 0.01 b (5000,1) unit 1 (150,1) unit 2 (100, 1) Interpreting the coefficients ᾱ as probabilities, the dual results obtained when D = 5149 suggest that, with probability close to 1, the ISO should dispatch units b and 1 and leave unit 2 off, paying a price π ∗ = 11.67. This coincides with the optimal dispatch, and the situation is similar for D = 5150. When D = 5151, there is a large mismatch between the pseudo-planning and the UC solution (p∗ , u∗ ), when unit 2 enters into operation with its minimal generation (p2min = 100). In this case, the pseudo-planning assigns a low probability to turn the unit on to generate 1, instead of the minimum 100, therefore making necessary a second phase, to recover primal feasibility. 4. Comparing some pricing proposals We start explaining the mechanism from [6], implemented in Brazil, and review briefly the literature. Afterwards, we focus on a new proposal that covers costs of dispatched generators, and compare some pricing rules for a hydro-thermal system over a 24h horizon. 4.1. The IP prices and other prices in the literature. Solving (7) yields the optimal dispatch and commitment p∗ and u∗ , and this information is used to define a price. If the value function v is convex, there is an ambiguity pointed out in Remark 2(2.2.3), but at least Theorem 2.1 ensures that any subgradient is a price that supports the dispatch and yields a genuine MCO. Nonconvex configurations require some adjustments to provide proper incentives for both operation and investment. The innovative idea in [6] was to measure the rate of change of the 14 JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 integrality nature associated with the commitment by solving the following linear program called post-UC problem in [29]:  X  GCosti (pi , ui ) min   T   (pi ,ui )∈Pi ,ui ∈[0,1] ,i∈M i∈M ∗ (12) i∈M s.t. uX i = ui ,   t t  p = D , t = 1, . . . , T .  i  i∈M The proposal [6] was the first one to split the payment to generators into two parts, one variable, in line with the MCO, is the Lagrange multiplier of the demand constraint and is called IP price in the literature. The second part refers to adjustments to be paid aside of the energy price, called uplifts, to compensate differently each generator for incurred fixed costs. Such payments are measured by the multiplier associated with the constraint u = u∗ . Some drawbacks of this rule are shown in Figure 3, with the dual and IP prices for our three illustrative systems. Figure 3. Dual and IP prices (light and dark lines), obtained with (9) and (12) for the three situations in Section 2.2. For the convex instance, when D = 5150, any value in [5, 12] is a legitimate price. The dual price always increases together with the demand. In a manner consistent with the intuition, start-up and shutdown costs (middle) and minimal generation and fixed costs (right) yielded higher prices than for the convex case (left). By contrast, the IP price exhibits a decrease for the nonconvex case (middle), even if the demand increases. This is explained by the cheaper variable cost of the unit 1, not dispatched for values of D below 5050 because of its high start-up cost. A similar phenomenon is observed in the discontinuous case (right), with the minimal generation requirement for unit 1 inducing another change. Continuing with our analysis of the rule [6], in the post UC-problem (12), the demand constraint multiplier gives an MCO for a frozen configuration, no unit that is on/off is turned off/on. Since (12) is obtained from (7) by adding extra constraints for which the UC solution remains feasible, (p∗ , u∗ ) is also optimal for (12). Therefore, Theorem 2.1 states that the IP price supports the generation plan (p∗ , u∗ ), but only when the commitment variable is fixed. The ISO determines if a given unit is turned on and, moreover, generators can not choose their production level at the given IP price. As a result, if for a dispatched unit i′ the IP price falls short A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN 15 in covering generation costs, the unit is still forced to deliver the energy p∗i′ at a loss. Uplifts are meant to address such undesirable situation. Notice that in Figure 3 the IP prices, contrary to the dual scheme, do not follow all the increases in demand, a feature that can make less effective demandresponse programs. The IP prices are in general lower than the dual ones. This phenomenon was observed in [8]: the dual scheme (9) increases prices until achieving the minimum necessary uplifts that compensate generators for being dispatched with the centralized solution, instead of self-scheduling. According to [6], IP prices are employed by the New York ISO and in the Pennsylvania-New Jersey-Maryland Interconnection market. The rule was also put in practice in Brazil since January 2021, to define the commercialization prices (but not the compensations). Regarding compensations, a multiplier of the commitment constraint with negative sign would result in a surplus that penalizes efficiency. A fix suggested by [6] is to compensate only if the multiplier is positive. The sign of the compensation will depend on how the constraint was written in (12) (u−u∗ = 0 or u∗ −u = 0). Furthermore, different formulations of the technological sets Pi , even if equivalent on paper, do change the multiplier set and, therefore, also produce a different output. A simple example in Appendix A illustrates this phenomenon. Similarly to the concerns in Remark 2(2.2.3), delicate issues related to degeneracy can have a significant impact on the business, if not suitably taken into consideration. Among other selected proposals that followed [6], the pricing rule in [30] has some flavor of parametric integer programming, using Benders’ cuts to estimate derivatives for value functions like the ones in Figure 2. Although the idea is nicely illustrated for a variation of Scarf’s example, it is not clear how the approach can be put in practice for general UC problems. The proposal [8] amounts to computing prices as in Section 3. In the literature three denominations are used to refer to the same mechanism: dual prices, minimal uplift prices and convex-hull pricing, where the latter wording refers to the value function of a bi-dual problem close to (11). As explained, since dual prices minimize uplifts, they are generally higher than IP prices. The Semi-Lagrangian price [31] is computed by solving a dual problem with the demand relaxed as in (10), but maintaining the constraint as an inequality. From the nonsmooth optimization point of view, the dual update is similar to a subgradient method, with well-known issues regarding convergence properties, see [14, part II]. Furthermore, keeping the constraint makes the dual function nonseparable, which can be a drawback for large systems. The report [32] compares advantages and drawbacks of the US and European pricing systems, respectively close by the IP and dual approaches. In [32, Section 6] the author proposes a new model that reduces uplifts for the former and increase welfare for the latter. Consequences of trade-offs adopted when defining compensations were analyzed in economic terms by [33], and of inexact MILP solution in terms of surpluses and prices in [34, 35]. Those works point out that the solution of the post-UC problem varies sharply with the generation plan. Two suboptimal (p∗ , u∗ ), apparently close, give very different prices and compensations in (12). Since the output depends on the ISO’s decision of stopping the MILP solver before optimality, such abrupt changes are perceived as arbitrary redistributions of profit and surplus, and are controversial. In [36] it is argued that dual prices minimize an upper estimate of the total redistribution that is done when solutions are suboptimal. The redistribution 16 JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 is measured in terms of the duality and optimality gaps and computable estimates are derived in that work. An interesting analysis in [36] explains that, even if the MILP is fully solved and the optimality gap is null, a nonconvex UC problem will always have redistribution payments, because of positive duality gaps. This motivates the study of different UC formulations, that reduce the duality gap. For a geometric study inspecting how different, apparently equivalent, MILP formulations result in different duality gaps, we refer to [37]. We conclude this brief revision of the literature by referring to the discussion in the review [38, Section 5], where a large number of other pricing designs is examined thoroughly. A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN 17 4.2. Ensuring cost recovery with limited compensations. The increased recent research on the impact of nearly optimal solutions on energy prices confirms our comments in Remark 2(2.2.3) and the interest of stabilizing financial outcomes in the business. The new pricing rule, introduced in problem (13) below, addresses this issue by including a target price in the objective function that stabilizes the optimization process. It is stated in [38] that prices computed from optimal Lagrange multipliers of (12) suffer from volatility, which sometimes results in too large or too small compensations. The work also comments on the importance of guaranteeing that all the dispatched units recover their costs. Our new scheme keeps these observations in mind by resorting to a natural interpretation of duality. Namely, given a linear production problem under inventory constraints, the dual problem can be interpreted as searching prices for the items in the inventory in such a way that they compensate any possible production level. This goal is achievable with Linear Programming models because there is no duality gap and any real production level is allowed. With the UC problem (7), however, the situation is different, especially if the duality gap is large due to the integrality constraints, the minimal generation levels, and significant fixed costs. Indeed, when the value function is discontinuous and exhibits a jump, as in the right plot in Figure 2, only an “infinite” price could capture the instantaneous change in the optimal cost. In this setting, the price must be accompanied by some compensation. Prices and compensations computed following [6], as in Section 44.1, are multipliers of the post-UC problem (12), a relaxation of the ISO’s problem. Our proposal is to define prices considering instead an economic problem that represents directly the interest of the generators. This is consistent with the separate concerns of feasibility and profit, respectively of the ISO and the generators, mentioned in the introduction. Instead of being Lagrange multipliers of a problem having a generation plan as decision variable, in our approach prices are the decision variables of an optimization problem. Additionally, the model defines the variable Ei ≥ 0 below as the total compensation to be paid to unit i, over the whole horizon. After talking to several agents in the business, it appears that generators care more for the value of the global compensation than for knowing individual values for each time step. The economic pricing problem is formulated following some leading principles. First, by constraint [onlyF] in (13) below, uplifts cover fixed costs only, defined by FCosti (u∗i ) = T  X t=1 ∗  Fi+ [(u∗i )t − (u∗i )t−1 ]+ + Fi− [(u∗i )t−1 − (u∗i )t ]+ , where a solution (p , u∗ ) to the UC problem (2) is given. Second, by constraint [Non-conf], supply is non-confiscatory [39] and, third, uplifts avoid over compensations through constraint [NoOverComp]. Finally, constraint [Lim] ensures that total payments in uplifts do not exceed a proportion β ∈ (0, 1) of the market profit. These requirements define the feasible set in (13), to address concerns of the generators. The ISO’s view is represented in the objective function, which makes use of a target π ∗ to drive the selection. The target can be the dual price from Section 3 that minimizes the uplifts, or the IP price, or some other reference deemed valid as marginal price. The choice of what is an appropriate reference, with desirable properties, entails an agreement between all the agents. Having a target price gives a 18 JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 systematic way to output the same price. The role of the target is to stabilize the selection mechanism. This feature addresses both the issues raised in Remark 2(2.2.3) and variations resulting from suboptimal solutions. In problem (13), a difference of close, yet different, generation plans (p∗ , u∗ ) is seen through its effect on the costs only. Given a large constant M > 0 whose value can be computed from the market configuration, the mathematical formulation of our pricing problem is (13)  1  kπ − π ∗ k22 min   t  2 (π ,Ei ,si )    s.t. Ei ≤ FCosti (u∗i )si , i ∈ M [OnlyF]     T  X  ∗ ∗   π t p∗t [Non-conf]  i + Ei ≥ GCosti (pi , ui ), i ∈ M   t=1                   T X ∗ ∗ π t p∗t i + Ei ≤ GCosti (pi , ui ) + M (1 − si ), i ∈ M [NoOverComp] t=1 ! T X X X X ∗ ∗ t ∗t GCosti (pi , ui ) [Lim] π pi − Ei ≤ β i∈M π ≥ 0, t=1 i∈M E ≥ 0, i∈M si ∈ {0, 1} , i ∈ M . Closeness to the target is defined using the squared 2-norm, so (13) is a mixedinteger quadratic program. Other measures are of course possible, including the 1-norm, that yields a MILP instead. The binary variable si models whether generator i will be paid a compensation. As stated above, constraint [Lim] is used to limit the total compensation to β times overall profit of the entire market. The rationale is related to consider compensations as an internal charge to participate of the market profit. This is somewhat reminiscent of the zero-sum prices described in Section 4 of [38]. A suitable value of β can be determined by the regulating agency in agreement with the producers. The final constraint states that prices and compensations can not be negative and that variables si are all binary. Since (13) considers the overall compensation for each unit Ei , instead of a compensation for each period and unit, the problem is not large-scale, having as many integer variables as units in the system. If a compensation for each period is needed, for example to ensure that some generator recovers fixed costs in a given window of time, the variable Ei can be split into a sum of compensations for each period, as PT in Ei = t=1 eti , adding constraints to model the desired properties for the individual compensations eti (keeping in mind the problem can become ill-conditioning because many different values of eti result in the same global compensation Ei ). Our rule (13) is close in spirit to [29], but formulated for more than one period and allowing different targets than the IP price. As such, (13) also bears some resemblance with the Dual Pricing Algorithm (DPA) in Chapter 2 in the PhD dissertation [39], defined for a market with demand that is responsive to price. A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN 19 Adapted to our setting, with a constant demand, the DPA boils down to the following problem, depending on a parameter a ≥ 0: (DPA)     min    (π t ,Eit )    s.t.          X Ei + a |π t − (π IP )t |/(π IP )t t=1 i∈M T X T X ∗ ∗ π t p∗t i + Ei ≥ GCosti (pi , ui ), i ∈ M t=1 π ≥ 0, E ≥ 0. These two pricing mechanisms compare as follows. Our rule (13) separates the ISO’s concern in the objective, from the generator’s perception, formulated as constraints. The objective in (DPA) uses the scaling parameter a to try to balance proximity to the ISO’s target price with the generator’s expectations of having their cost covered by direct payments only, without compensations. The feasible set in (DPA) ensures cost recovery, the minimal acceptable revenue structure for the generators. Using the dual and IP prices as target π ∗ , we computed prices with these two rules for the system with discontinuous value function in Figure 2, letting the demand D ∈ [5000, 5300]. We let β = 0.01 and M = 109 in (13) and a = 1 in (DPA). At the critical value D = 5150, both rules set the price 5 if the reference in (13) is the IP price, the same used in (DPA). If we set the target in (13) as the dual price, the selected value is 12, as expected. 4.3. Dynamics of prices for a hydro-thermal system. The storage capability of energy in Brazil is still very large (291 GWmonth in 2020), but hydro-power has become more vulnerable to inflow variability in the last years. There is a marked seasonality in the supply, with run-of-river hydro-plants producing half of the total generation in the rainy season. Utilities like Belo Monte, Jiraú and Santo Antônio, with respective installed capacity of 11,300 MW, 3,750 MW, and 3,568 MW, contribute to preserving storage capability in other regions, but also bring down what other hydro-plants foresee as maximum consumption in dry periods. To reduce risk, [2] explains that the ISO’s resorts to some operational procedures to achieve a storage target level that is considered secure at the end of each year. For the short-term, the ISO takes operational measures that protect the system from different degrees of severity of the dry seasons especially regarding the starting levels of reservoirs at the end of the rainy season (source: energy summary released by the ISO in 2020, www.ons.org). To recreate those features, we consider a power system with three ideal thermal units, with increasingly faster ramping dynamics and decreasing capacity. All units have start-up and shut-down costs and minimal generation requirements. There is also hydro-power associated with a reservoir of volume Vht at time t. Given an initial volume Vh0 , deterministic inflows I t , and a factor ηh , converting turbine outflow into energy, the water balance equation is Vht − Vht−1 + ηh pth = I t , for t = 1, . . . , T . 20 JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 All hydro-variables are continuous and, in the technological set Ph , for given reservoir bounds Vhmin < Vhmax the volumes satisfy the relations   Vhmin safe := max Vhmin , 0.20Vhmax ≤ Vht ≤ Vhmax . The hydro-generation cost is a piecewise affine function derived from the “futurecost-function” of the mid-term planning tool that depends on the final volume VhT , as in [3, eq.(3)]. To mimic, respectively, “dry” and “wet” initial conditions for the reservoir, we solved the UC problem varying the initial condition for the reservoir, by setting Vh0dry = 0.9Vhmin safe + 0.1Vhmax and Vh0wet = 0.1Vhmin safe + 0.9Vhmax , The hourly demand follows the profile of one summer weekday in the South East sub-market in Brazil, scaled to the small system. A comment regarding the UC model is in order. Typically, there is a difference of magnitude between the futurecost-function, when compared to the generating cost of the thermal power plants. For the results that follow, (7) includes a constraint fixing the final volume of the reservoir to its starting value: VhT = Vh0 . Because of this constraint, the operational cost for hydro-power in (13) is null: GCosth (p∗ , u∗ ) = 0, which eliminates scaling issues that could distort the price in our stylized configuration. Prices IP, DPA, and Econ(IP), our rule with IP target, are reported in Figure 4 for this base case. Figure 4. Base case, solving (12) to compute the IP price π IP ; the DPA price solves (DPA) with a = 1 and the Econ(IP) price corresponds to solving (13) with π ∗ = π IP , β = 0.01 and M = 109 . The yellow line shows the demand profile. In the base case the values of inflows I t represent a favorable hydrological condition. The graph corresponds to the output obtained with Vh0 = Vh0dry , noting that all prices rules were similar with initial condition Vh0 = Vh0wet . In the base case all the rules exhibit a similar pattern for both starting volumes. The only significant difference occurs with the DPA price, that at time 13 jumps more 50% over the IP and Econ(IP) prices. Near that time, the most expensive thermal unit is turned on. When the hydrological situation is less favorable, it is expected that the pricing rules are more sensitive to the starting volume. We run a variation of the base case in which inflows are reduced in a 25% (I t = 0.75I tbase ), taking the same system and demand as before. When the starting volume is close to the maximum level Vhmax , prices exhibit a pattern similar to the one in Figure 4. By contrast, under a A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN 21 dry scenario, when Vh0 is close to the lower bound Vhmin safe , the temporal structure of prices changes, as shown by Figure 5. Figure 5. When inflows are reduced a 25%, prices are sensitive to the initial volume of reservoirs. When starting with a high volume (Vh0wet ), prices remain low until the peak hours, as in Figure 4. But if the initial condition is dry that is, when Vh0 = Vh0dry , prices increase in the first hours too, coinciding with the time when the reservoir volume remains close to its minimum value. We note in Figure 5 that again DPA augments the price significantly at a single hour. We observed the same phenomenon in many runs, indicating that DPA avoids compensations by increasing substantially the IP price at only one time. The behavior of Econ(IP) is more smooth, the rule increases the IP prices between 5% and 10% until time 19, approximately. Regarding compensations, there were none for DPA, and the IP compensations are the largest (IP prices are lowest in general). In spite of the large side payment, the IP prices never cover the cost of the fast generator, even after adding compensations. By contrast, Econ(IP) computes a compensation that always zeroes the profit of the fast generator. The margin for each unit defined as Mi (π, p∗i , u∗i ) in (3) using the three prices in Figure 5, can be seen in Figure 6, where the output was scaled with respect to the IP price, to show all the units in a similar scale. 5. Conclusion We have seen that there are many proposals to remunerate generators in the presence of nonconvexities. The positive duality gap associated to that situation measures the discrepancy between the ISO’s and the generator’s expectations. To define which trade-offs to put in place, the pricing system must be assessed from different angles, among which determining if prices follow the demand can be important in Brazil, to ensure that demand-response programs are effective. Ideally, if the system implements compensations, they should be sufficient to ensure cost recovery, without over compensating generators. The simple model derived from Table 1, when particularized in its three instances, illustrates well issues that need to be resolved to reach an agreement on those trade-offs. An additional layer of complexity arises with a hydro-thermal system, if hydro-operational costs are indirectly valuated through the final volume of reservoirs. The different pricing mechanisms explored have all pros and cons. The proposal in Section 4, by [6], keeps prices low but, as shown in Figure 6 does not recover the cost for the fast thermal unit, without significant uplift payments. In some situations, IP prices over compensate, [32]. The new approach in Subsection 44.2 22 JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 Figure 6. Margins (3), scaled to the one obtained with the IP price for the configuration in Figure 5, with one hydro and three thermal plants, of slow, middle and fast dynamics. If only the generated energy is paid, without any compensation, the margins represent the total payment received by the generator. The corresponding value is a profit if positive, as with the slow and middle thermal and the hydro units. Without compensations paid on the side, the most expensive thermal power plant, with fast dynamics, incurs losses with IP and Econ(IP), but not DPA. The rule DPA achieves this by the jump of 50% in the price of time 16. Compensations are null with DPA, as expected. Rule Econ(IP) computes a compensation that zeroes the loss of the fast thermal unit. Finally, the IP compensation, even though it is more than four times larger than Econ(IP)’s, still leaves the fast unit in a loss. presents some appeal in this sense, as it makes explicit the generator’s search for cost recovery without over compensations. The optimization problem, written directly with the prices and compensations as variables, is not too large and can be solved easily with commercial solvers. Furthermore, thanks to its objective function, the output is stabilized and possible ambiguities and degeneracies such as those pointed out in Remark 2(2.2.3) and in Appendix A, are eliminated. Finally, the formulation (13) can easily include constraints that have some extra economical or financial meaning. An interest possibility would be to add to the feasible set in (13) a bound on the payment redistribution that is inherent to MILPs, if the UC problem is inexactly solved. According to the theory in [36], such constraints, which depend on estimates of the optimality and duality gaps, could be combined the dual price as target in the objective function of (13) to help avoiding abrupt, and somewhat arbitrary, oscillations in the pricing mechanism when the UC problem is not solved until optimality. A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN 23 Acknowledgments. The corresponding author thanks the Isaac Newton Institute for Mathematical Sciences for support and hospitality during the programme The mathematics of energy systems [EPSRC grant number EP/R014604/1]. The authors acknowledge the support of TACEMM, the Transatlantic Consortium of Energy Markets and Modeling. JPL is funded by PRONEX–Optimization, CS and PS by CEPID CeMEAI [FAPESP 2013/07375-0]; additionally, CS and PS are partly funded by [CNPq 306089/2019-0 and CNPq 306089/2019-0 and FAPESP 2018/24293-0], respectively. This study was kicked off thanks to the series of workshops “Mathematical solutions to industrial problems” organized by CeMEAI in 2018 and 2019, with the support of Engie Brazil and Cepel Eletrobras. Several issues in our paper were addressed in response to the intensive work and remarks of graduate students who actively participated in both editions of the workshop, namely of J. Chavez from Unicamp, and B. Colonetti, M. Cordova and K. Vinente from UFSC, the Federal University of Santa Catarina, Brazil. References [1] A. Bublitz, D. Keles, F. Zimmermann, Ch. Fraunholz, and W. Fichtner. A survey on electricity market design: Insights from theory and real-world implementations of capacity remuneration mechanisms. Energy Economics, 80:1059 – 1078, 2019. [2] B. Bezerra, L. A. Barroso, M. Brito, F. Porrua, B. Flach, and M. V. Pereira. Measuring the hydroelectric regularization capacity of the brazilian hydrothermal system. In IEEE PES General Meeting, pages 1–7, 2010. [3] T.N. Santos, A.L. Diniz, C.H. Saboia, R.N. Cabral, and L.F. Cerqueira. Hourly pricing and day-ahead dispatch setting in Brazil: The ”DESSEM” model. Electric Power Systems Research, 189:106709, dec 2020. Part of special issue: Proceedings of the 21st Power Systems Computation Conference (PSCC 2020), edited by D. Molzahn. [4] H. E. Scarf. The allocation of resources in the presence of indivisibilities. Journal of Economic Perspectives, 8(4):111–128, nov 1994. [5] L. .V Snyder and Z.-J. M. Shen. Fundamentals of Supply Chain Theory. Wiley, jun 2019. [6] R. P. O’Neill, P. M. Sotkiewicz, B. F. Hobbs, M. H. Rothkopf, and W. R. Stewart. Efficient market-clearing prices in markets with nonconvexities. European Journal of Operational Research, 164(1):269–285, July 2005. [7] M. Cordova, W. de Oliveira, and C. Sagastizábal. Revisiting Augmented Lagrangian Duals. http://www.optimization-online.org/DB_HTML/2020/03/ 7709.html, 2020. [8] P. R. Gribik, W. W. Hogan, and S. L. Pope. Market-Clearing Electricity Prices and Energy Uplift. Technical report, John F. Kennedy School of Government, Harvard University, 2007. [9] R. Baldick. The generalized unit commitment problem. IEEE Transactions on Power Systems, 10(1):465–475, 1995. [10] X. Guan, Q. Zhai, and A. Papalexopoulos. Optimization based methods for unit commitment: Lagrangian relaxation versus general mixed integer programming. In 2003 IEEE Power Engineering Society General Meeting, volume 2, pages 1095–1100 Vol. 2, 2003. 24 JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 [11] T. Li and M. Shahidehpour. Price-based unit commitment: A case of Lagrangian relaxation versus mixed integer programming. IEEE Transactions on Power Systems, 20(4):2015–2025, 2005. [12] A. Merlin and P. Sandrin. A new method for unit commitment at electricite de france. IEEE Transactions on Power Apparatus and Systems, PAS102(5):1218–1225, 1983. [13] R. Taktak and C. D’Ambrosio. An overview on mathematical programming approaches for the deterministic unit commitment problem in hydro valleys. Energy Systems, 8(1):57–79, jan 2016. [14] J.F. Bonnans, J.Ch. Gilbert, C. Lemaréchal, and C. Sagastizábal. Numerical Optimization. Theoretical and Practical Aspects. Universitext. Springer-Verlag, Berlin, 2006. Second edition, xiv+490 pp. [15] C. Sagastizábal. Divide to conquer: decomposition methods for energy optimization. Mathematical Programming, 134:187–222, 2012. [16] F. Zhuang and F.D. Galiana. Towards a more rigorous and practical unit commitment by Lagrangian relaxation. IEEE Transactions on Power Systems, 3(2):763–773, 1988. cited By 325. [17] C. Lemaréchal, C. Sagastizábal, F. Pellegrino, and A. Renaud. Bundle methods applied to the unit-commitment problem. In System Modelling and Optimization, pages 395–402. Springer US, 1996. [18] N. J. Redondo and A. J. Conejo. Short-term hydro-thermal coordination by lagrangian relaxation: solution of the dual problem. IEEE Transactions on Power Systems, 14(1):89–95, 1999. [19] A. Borghetti, A. Frangioni, F. Lacalandra, and C. A. Nucci. Lagrangian heuristics based on disaggregated bundle methods for hydrothermal unit commitment. IEEE Power Engineering Review, 22(12):60–60, 2002. [20] A. Belloni, A.L. Diniz Souto Lima, M.E. Piñeiro Maceira, and C.A. Sagastizábal. Bundle relaxation and primal recovery in unit commitment problems. the Brazilian case. Annals of Operations Research, 120(1-4):21–44, 2003. cited By 37. [21] M. R. Scuzziato, E. C. Finardi, and A. Frangioni. Solving stochastic hydrothermal unit commitment with a new primal recovery technique based on lagrangian solutions. International Journal of Electrical Power & Energy Systems, 127:106661, may 2021. [22] R. E. Bixby, M. Fenelon, Z. Gu, E. Rothberg, and R. Wunderling. 18. mixedinteger programming: A progress report. In The Sharpest Cut, pages 309–325. Society for Industrial and Applied Mathematics, jan 2004. [23] S. Feltenmark and K. C. Kiwiel. Dual Applications of Proximal Bundle Methods, Including Lagrangian Relaxation of Nonconvex Problems. SIAM Journal on Optimization, 10(3):697–721, 2000. [24] L. Dubost, R. Gonzalez, and C. Lemaréchal. A primal-proximal heuristic applied to the French unit-commitment problem. Mathematical Programming, 104:129–151, 2005. [25] P. Carpentier, G. Gohen, J.-C. Culioli, and A. Renaud. Stochastic optimization of unit commitment: a new decomposition framework. IEEE Transactions on Power Systems, 11(2):1067–1073, may 1996. [26] L. Lasdon. Optimization theory for large systems. Macmillan Series in Operations Research, 1970. A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN 25 [27] D. Bertsekas, G. Lauer, N. Sandell, and T. Posbergh. Optimal short-term scheduling of large-scale power systems. IEEE Transactions on Automatic Control, 28(1):1–11, 1983. [28] C. Lemaréchal and C. Sagastizábal. Variable metric bundle methods: from conceptual to implementable forms. Mathematical Programming, 76:393–410, 1997. [29] R. P. O'Neill, A. Castillo, B. Eldridge, and R. B. Hytowitz. Dual pricing algorithm in ISO markets. IEEE Transactions on Power Systems, 32(4):3308– 3310, jul 2017. [30] M. Bjørndal and K. Jörnsten. Equilibrium prices supported by dual price functions in markets with non-convexities. European Journal of Operational Research, 190(3):768 – 789, 2008. [31] V. Araoz and K. Jörnsten. Semi-lagrangean approach for price discovery in markets with non-convexities. European Journal of Operational Research, 214(2):411–417, 2011. [32] M. Van Vyve. Linear prices for non-convex electricity markets: models and algorithms. Core discussion papers, Université catholique de Louvain, Center for Operations Research and Econometrics (CORE), 2011. [33] R. Sioshansi, R. O'Neill, and S.S. Oren. Economic consequences of alternative solution methods for centralized unit commitment in day-ahead electricity markets. IEEE Transactions on Power Systems, 23(2):344–352, may 2008. [34] R. Sioshansi and A. Tignor. Do centrally committed electricity markets provide useful price signals? The Energy Journal, 33(4), oct 2012. [35] R. Sioshansi. Pricing in centrally committed electricity markets. Utilities Policy, 31:143–145, dec 2014. [36] B. Eldridge, R. O'Neill, and B. F. Hobbs. Near-optimal scheduling in dayahead markets: Pricing models and payment redistribution bounds. IEEE Transactions on Power Systems, 35(3):1684–1694, may 2020. [37] C. Lemaréchal and A. Renaud. A geometric study of duality gaps, with applications. Mathematical Programming, 90(3):399–427, 2001. [38] G. Liberopoulos and P. Andrianesis. Critical Review of Pricing Schemes in Markets with Non-Convex Costs. Operations Research, 64(1):17–31, January 2016. [39] R. S. B. Hytowitz. Optimization Modeling to Address the Impacts of Electric Power Market Design on Operations and Policy. PhD thesis, Whiting School of Engineering, Johns Hopkins University, 2018. Appendix A. Example illustrating odd responses of IP prices A simple instance of (7), with two units and no minimal generation requirement, but fixed costs, suffices to show some awkward output of the post-UC problem (12), used by rule [6] to determine prices and compensations. 26 JUAN PABLO LUNA1 , CLAUDIA SAGASTIZÁBAL2 AND PAULO J. S. SILVA3 A.1. Unique IP price, nonunique compensations. With a demand D = 10, the UC problem is  2 X    ipi + 10iui min   (pi ,ui )i=1,2   i=1   2 X (14) pi = 10 s.t.     i=1    0 ≤ pi ≤ 20ui i = 1, 2   ui ∈ {0, 1} i = 1, 2 . Since all units have capacity larger than D and the first unit is cheaper, the solution is (p∗1 , u∗1 ) = (10, 1) and (p∗2 , u∗2 ) = (0, 0) . The post-UC problem (12), yielding IP prices and compensations, is  2 X    ipi + 10iui min   (pi ,ui )i=1,2   i=1    2  X   s.t. pi = D (λ) (15) i=1    0 ≤ pi ≤ 20ui i = 1, 2 (ηi )     0 ≤ ui ≤ 1 i = 1, 2 (ζi )     u1 = 1 (w1 )   u2 = 0 (w2 ) , where multipliers are denoted by variables in greek letters between parentheses. The optimality conditions, equivalent to (15) write down as follows: 0 ≤ i − λ + ηi ⊥ pi ≥ 0, i = 1, 2 0 ≤ 10i − 20ηi + ζi − wi ⊥ ui ≥ 0, D − p1 − p2 = 0 (16) i = 1, 2 0 ≤ 20ui − pi ⊥ ηi ≥ 0 i = 1, 2 0 ≤ 1 − ui ⊥ ζi ≥ 0 i = 1, 2 1 − u1 = 0 0 − u2 = 0 . The dual variables that solve this KKT system when D = 10 are not unique. This can be checked by plugging in (16) any element in the set    λ = 1 , η1 = 0 , η2 = 0 L := (λ, η, w) : , w1 ≥ 10 , w2 ≤ 20 and checking satisfaction of the the optimality conditions. The IP price, equal to λ, is unique and coincides with the variable cost of the cheaper unit 1, the only one dispatched. This value is aligned with the intuition, but compensations present a counter intuitive pattern. More precisely, the compensation for unit 1 can be any value w1 ∈ [10, +∞) while unit 2’s compensation can be any w2 ∈ (−∞, 20]. This means that the dispatched unit 1 can be charged any value w1 larger than its fixed cost 10, while the second unit that is not dispatched, could be charged or paid, according to the sign output by the solver for w2 . Even with the fix proposed in [6], that pays compensations only to dispatched generators, A DISCUSSION ON ELECTRICITY PRICES, OR THE TWO SIDES OF THE COIN 27 for unit 1, the only unit generating in the system, this means that at best it can get a zero profit, incurring in losses for any value of w1 > 10. A.2. Wild variations for equivalent formulations. The situation can be even more awkward, if (14) is reformulated in an equivalent manner, changing the capacity constraints therein from 20 to 10. Since the optimal dispatch remains the same, the corresponding post-UC problem is  2 X    ipi + 10iui min    (pi ,ui )i=1,2 i=1     2  X   s.t. pi = 10 (λ) (17) i=1    0 ≤ pi ≤ 10ui i = 1, 2 (ηi )     0 ≤ ui ≤ 1 i = 1, 2 (ζi )     u1 = 1 (w1 )   u2 = 0 (w2 ) . Both problems (15) and (17) have the same feasible set, so they are equivalent. But the corresponding dual solutions are very different. More precisely, instead of the set L, the dual variable associated with (17) is now unique, equal to λ = 101 , η1 = 100 , η2 = 91 , w1 = −1000 , w2 = 890 . Compensations are exceedingly large and, most importantly, an apparently innocuous difference in the post-UC problem made the former IP price of 1 jump to 101! A.3. Nonunique IP prices. Finally, if the demand in (14) is D = 20, the optimal generation plan is (p∗1 , u∗1 ) = (20, 1) and (p∗2 , u∗2 ) = (0, 0) . If in both the post-UC problem (15) and the optimality conditions (16) we let D = 20. we obtain a set of optimal multipliers that is parameterized by λ:    λ ≥ 1 , η1 = λ − 1 , η2 = max(0, λ − 2) L := (λ, η, w) : . w1 ≥ 30 − 20λ , w2 ≤ 20 − max(0, 20λ − 40) As commented in Remark 2(2.2.3) and, moreover, the IP price can be any value larger than 1, the generation cost of unit 1. 1 Universidade Federal do Rio de Janeiro, Brazil, 2,3 IMECC-Unicamp, Brazil