Benders' decomposition of the unit commitment problem with semidefinite relaxation of AC power flow constraints

TL;DR

This paper introduces a novel approach combining Benders' decomposition with semidefinite relaxation for the AC power flow in the unit commitment problem, improving convergence and solution quality.

Contribution

It presents a new formulation integrating semidefinite programming relaxation into Benders' decomposition for AC power flow constraints in unit commitment.

Findings

Improved convergence speed in numerical tests.

Enhanced feasibility and optimality cuts.

Effective application to 6-bus and 118-bus networks.

Abstract

In this paper we present a formulation of the unit commitment problem with AC power flow constraints. It is solved by a Benders decomposition in which the unit commitment master problem is formulated as a mixed-integer problem with linearization of the power generation constraints for improved convergence. Semidefinite programming relaxation of the rectangular AC optimal power flow is used in the subproblem, providing somewhat conservative cuts. Numerical case studies, including a 6-bus and the IEEE 118-bus network, are provided to test the effectiveness of our proposal. We show in our numerical experiments that the use of such strategy improves the quality of feasibility and optimality cuts generated by the solution of the convex relaxation of the subproblem, therefore reducing the number of iterations required for algorithm convergence.

Tables

Table 1. Table 1: Generator Data for the 3-GEN case.

Parameter	G1	G2	G3
$\alpha$ $\left[\$/{\mathrm{MW}}^2\right]$	$4\cdot {10}^{-4}$	$1\cdot {10}^{-3}$	$5\cdot {10}^{-3}$
$\beta$ $\left[\$/\mathrm{MW}\right]$	$13.7$	$40$	$17.7$
$\gamma$ $[\$]$	$177$	$130$	$137$
$\overline{P}$ $\left[\mathrm{MW}\right]$	$210$	$100$	$70$
$\underset{¯}{P}$ $\left[\mathrm{MW}\right]$	$100$	$10$	$10$
$\overline{Q}$ $\left[\mathrm{MVAR}\right]$	$210$	$100$	$70$
$\underset{¯}{Q}$ $\left[\mathrm{MVAR}\right]$	$-210$	$-100$	$-70$
$SU$ $[\$]$	$100$	$200$	$50$
$SD$ $[\$]$	$50$	$100$	$50$
$p_0$ $\left[\mathrm{MW}\right]$	$150$	$0$	$15$
$UR$ $\left[\mathrm{MW}\right]$	$55$	$50$	$15$
$DR$ $\left[\mathrm{MW}\right]$	$55$	$50$	$15$
$T_0$ $\left[h\right]$	$2$	$-1$	$1$
$T_{on}$ $\left[h\right]$	$4$	$3$	$2$
$T_{off}$ $\left[h\right]$	$4$	$2$	$2$

Table 2. Table 2: Iterations until convergence with and without modified master problem.

Case study	$M{M}^k$	$M^k$
3-GEN	6	13
IEEE-118	6	$\infty$

Table 3. Table 3: Convergence of upper and lower bounds for the IEEE-118 case.

Iteration	LB	UB
$1$	$1,615,160.46$	$4,884,906.44$
$2$	$1,615,209.08$	$4,765,981.95$
$3$	$1,615,191.01$	$3,545,654.32$
$4$	$1,615,755.13$	$1,700,015.57$
$5$	$1,699,320.68$	$1,699,482.33$
$6$	$1,699,482.33$	$1,699,482.33$

Table 4. Table 4: General data for the 3-GEN case.

$T\left[h\right]$	$Sbase\left[\mathrm{MVA}\right]$	$i_0$	$V_0[p.u.]$	$\overline{V}[p.u.]$	$\underset{¯}{V}[p.u.]$
$24$	$100$	$1$	$1$	$1.05$	$0.95$

Table 5. Table 5: Transmission line data.

Line	$\overline{F}\left[\mathrm{MW}\right]$	$r[p.u.]$	$x[p.u.]$
1-2	110	0.02	0.09
2-3	65	0.004	0.012
1-4	110	0.03	0.09
2-4	60	0.02	0.06
4-5	60	0.004	0.012
5-6	80	0.02	0.06
3-6	50	0.006	0.002

Table 6. Table 6: Load data (in MW megawatt \mathrm{MW} and M mega \mathrm{M} VAR respectively). Load on buses 1,2 and 6 are zero.

$T\left[h\right]$	dp3	dq3	dp4	dq4	dp5	dq5
1	$34.05$	$9.01$	$68.1$	$18.01$	$68.1$	$18.01$
2	$31.92$	$8.45$	$63.84$	$16.88$	$63.84$	$16.88$
3	$30.86$	$8.17$	$61.71$	$16.32$	$61.71$	$16.32$
4	$29.8$	$7.88$	$59.58$	$15.76$	$59.58$	$15.76$
5	$29.8$	$7.88$	$59.58$	$15.76$	$59.58$	$15.76$
6	$30.86$	$8.17$	$61.71$	$16.32$	$61.71$	$16.32$
7	$34.05$	$9.01$	$68.1$	$18.01$	$68.1$	$18.01$
8	$40.44$	$10.7$	$80.86$	$21.38$	$80.86$	$21.38$
9	$46.29$	$12.25$	$92.57$	$24.48$	$92.57$	$24.48$
10	$50.54$	$13.37$	$101.08$	$26.73$	$101.08$	$26.73$
11	$52.69$	$13.94$	$105.34$	$27.86$	$105.34$	$27.86$
12	$53.2$	$14.08$	$106.4$	$28.14$	$106.4$	$28.14$
13	$52.69$	$13.94$	$105.34$	$27.86$	$105.34$	$27.86$
14	$53.2$	$14.08$	$106.4$	$28.14$	$106.4$	$28.14$
15	$53.2$	$14.08$	$106.4$	$28.14$	$106.4$	$28.14$
16	$51.61$	$13.65$	$103.21$	$27.35$	$103.21$	$27.35$
17	$51.08$	$13.52$	$102.14$	$27.01$	$102.14$	$27.01$
18	$51.08$	$13.52$	$102.14$	$27.01$	$102.14$	$27.01$
19	$49.48$	$13.09$	$98.95$	$26.17$	$98.95$	$26.17$
20	$48.95$	$12.95$	$97.89$	$25.88$	$97.89$	$25.88$
21	$48.95$	$12.95$	$97.89$	$25.9$	$97.89$	$25.9$
22	$49.48$	$13.09$	$98.95$	$26.17$	$98.95$	$26.17$
23	$46.25$	$12.24$	$92.57$	$24.48$	$92.57$	$24.48$
24	$38.31$	$10.14$	$76.61$	$20.26$	$76.61$	$20.26$

Equations

$$\underset{\begin{array}{c}\mathbf{x},\mathbf{y},\mathbf{z},\\ {\mathbf{P}}_{t,g},{\mathbf{V}}_t,\mathbf{q}\end{array}}{\min }\sum _{t=1}^T\sum _{g=1}^{N_G}{\mathbf{C}}_g\bullet {\mathbf{P}}_{t,g}+{\mathbf{c}}^{⊺}\mathbf{x}+{\mathbf{u}}^{⊺}\mathbf{y}+{\mathbf{h}}^{⊺}\mathbf{z}$$

$$\sum _{\forall g\in {\mathcal{G}}_i}\left(\mathbf{A}\bullet {\mathbf{P}}_{t,g}+{\underline{p}}_g\cdot x_{t,g}\right)+{\mathbf{Y}}_i\bullet {\mathbf{V}}_t$$

$$\sum _{\forall g\in {\mathcal{G}}_i}\left({\stackrel{\Delta }{q}}_{t,g}+{\underline{q}}_g\cdot x_{t,g}\right)+{\tilde{\mathbf{Y}}}_i\bullet {\mathbf{V}}_t$$

$$\mathbf{A}\bullet {\mathbf{P}}_{t,g}-\left({\overline{p}}_g-{\underline{p}}_g\right)x_{t,g}$$

$${\stackrel{\Delta }{q}}_{t,g}-\left({\overline{q}}_g-{\underline{q}}_g\right)x_{t,g}$$

$$x_{t,g}-{\text{1}}_{(t>1)}\cdot x_{t-1,g}-{\text{1}}_{(t=1)}\cdot x_{0,g}$$

$${\text{1}}_{(t>1)}\cdot x_{t-1,g}+{\text{1}}_{(t=1)}\cdot x_{0,g}-x_{t,g}$$

$${\omega }_{(t,T_{on})}\left(x_t-{\text{1}}_{(t>1)}\cdot x_{t-1}\right)-\sum _{j=1}^{{\omega }_{(t,T_{on})}}{x}_{j+t-1,g}$$

$${\omega }_{(t,T_{off})}\left({\text{1}}_{(t>1)}\cdot x_{t-1}-x_t\right)-\sum _{j=1}^{{\omega }_{(t,T_{off})}}\left(1-x_{j+t-1,g}\right)$$

$$\sum _{g=1}^{N_G}\left(\mathbf{A}\bullet {\mathbf{P}}_{t,g}-\left({\overline{p}}_g-{\underline{p}}_g\right)\cdot x_{t,g}\right)$$

$$\left(\mathbf{A}\bullet {\mathbf{P}}_{t,g}+{\underline{p}}_g\cdot x_{t,g}\right)-{\text{1}}_{(t=1)}\cdot p_{0,g}-{\text{1}}_{(t>1)}\left(\mathbf{A}\bullet {\mathbf{P}}_{t-1,g}+{\underline{p}}_g\cdot x_{t-1,g}\right)$$

$${\text{1}}_{(t=1)}\cdot p_{0,g}-\left(\mathbf{A}\bullet {\mathbf{P}}_{t,g}+{\underline{p}}_g\cdot x_{t,g}\right)+{\text{1}}_{(t>1)}\left(\mathbf{A}\bullet {\mathbf{P}}_{t-1,g}+{\underline{p}}_g\cdot x_{t-1,g}\right)$$

$$-{\overline{F}}_{i,j}⩽{\mathbf{Y}}_{i,j}\bullet {\mathbf{V}}_t$$

$${\underline{V}}_i^2⩽{\mathbf{E}}_i\bullet {\mathbf{V}}_t$$

$${\mathbf{P}}_{t,g},{\mathbf{V}}_t$$

$$\mathbf{q}$$

$$\text{rank}({\mathbf{V}}_t)$$

$$\mathbf{x},\mathbf{y},\mathbf{z}$$

$$\mathbf{x}$$

$$\mathbf{z}$$

$${\stackrel{\Delta }{p}}_{t,g}\doteq p_{t,g}-{\underline{p}}_g$$

$$\mathbf{A}\doteq\left[\begin{array}[]{cc}0&1/2\\ 1/2&0\end{array}\right],\textrm{ and}\;\mathbf{P}_{t,g}\doteq\left[\begin{array}[]{cc}\overset{\Delta}{p}_{t,g}^{2}&\overset{\Delta}{p}_{t,g}\\ \overset{\Delta}{p}_{t,g}&1\end{array}\right],$$

$${\stackrel{\Delta }{p}}_{t,g}=\mathbf{A}\bullet {\mathbf{P}}_{t,g}=\frac{1}{2}{\stackrel{\Delta }{p}}_{t,g}+\frac{1}{2}{\stackrel{\Delta }{p}}_{t,g}.$$

$${\stackrel{\Delta }{q}}_{t,g}\doteq q_{t,g}-{\underline{q}}_g$$

$$\mathbf{q}\doteq {\left[{\stackrel{\Delta }{q}}_{1,1},\dots ,{\stackrel{\Delta }{q}}_{1,N_G},\dots ,{\stackrel{\Delta }{q}}_{T,1},\dots ,{\stackrel{\Delta }{q}}_{T,N_G}\right]}^{⊺}.$$

$${\mathbf{v}}_t\doteq {\left[e_{t,1},\dots ,e_{t,N_B},f_{t,1},\dots ,f_{t,N_B}\right]}^{⊺}\in {\mathbb{R}}^{2\cdot N_B},$$

$${\mathbf{V}}_t={\mathbf{v}}_t{\mathbf{v}}_t^{⊺}.$$

$${\alpha }_g\cdot p_{t,g}^2+{\beta }_g\cdot p_{t,g}+{\gamma }_g,$$

$${\alpha }_g\cdot {\left({\stackrel{\Delta }{p}}_{t,g}+{\underline{p}}_{t,g}\right)}^2+{\beta }_g\cdot \left({\stackrel{\Delta }{p}}_{t,g}+{\underline{p}}_{t,g}\right)+{\gamma }_g.$$

$$\mathbf{C}_{g}\doteq\left[\begin{array}[]{cc}\alpha_{g}&\alpha_{g}\cdot\underline{p}_{g}+\beta_{g}/2\\ \alpha_{g}\cdot\underline{p}_{g}+\beta_{g}/2&0\end{array}\right]$$

$$\mathbf{c}\doteq {\left[c_{1,1},\dots ,c_{1,N_G},\dots ,c_{T,1},\dots ,c_{T,N_G}\right]}^{⊺}\in {\mathbb{R}}_{+}^{T\cdot N_G}$$

$$c_{t,g}\doteq {\alpha }_g\cdot {\underline{p}}_g^2+{\beta }_g\cdot {\underline{p}}_g+{\gamma }_g,$$

$${\mathbf{Y}}_i\doteq -\sum _{\forall j\in {\Omega }_i}{\mathbf{Y}}_{i,j}$$

$$\mathbf{Y}_{i,j}\doteq\frac{1}{2}\left[\begin{array}[]{cc}G_{i,j}+G_{i,j}^{\intercal}&B_{i,j}^{\intercal}-B_{i,j}\\ B_{i,j}-B_{i,j}^{\intercal}&G_{i,j}+G_{i,j}^{\intercal}\end{array}\right]\;\in\;\mathbb{R}^{(2\cdot N_{B})\times(2\cdot N_{B})}$$

$$G_{i,j}\doteq g_{i,j}\left({\mathbf{\xi }}_i{\mathbf{\xi }}_i^{⊺}-{\mathbf{\xi }}_i{\mathbf{\xi }}_j^{⊺}\right),B_{i,j}\doteq b_{i,j}\left({\mathbf{\xi }}_i{\mathbf{\xi }}_i^{⊺}-{\mathbf{\xi }}_i{\mathbf{\xi }}_j^{⊺}\right),$$

$${\tilde{\mathbf{Y}}}_i\doteq -\sum _{\forall j\in {\Omega }_i}{\tilde{\mathbf{Y}}}_{i,j}$$

$$\widetilde{\mathbf{Y}}_{i,j}\doteq\frac{1}{2}\left[\begin{array}[]{cc}B_{i,j}+B_{i,j}^{\intercal}&G_{i,j}-G^{\intercal}_{i,j}\\ G_{i,j}^{\intercal}-G_{i,j}&B_{i,j}+B_{i,j}^{\intercal}\end{array}\right]\;\in\;\mathbb{R}^{(2\cdot N_{B})\times(2\cdot N_{B})}$$

$$\omega_{(t,\tau)}=\left\{\begin{array}[]{lr}\min\left(\tau-T_{0},\tau\right)&:t=1,\\ \min\left(\tau,T-t+1\right)&:t>1.\end{array}\right.$$

$$\mathbf{E}_{i}\doteq\left[\begin{array}[]{cl}\bm{\xi}_{i}\bm{\xi}_{i}^{\intercal}&\mathbf{0}\\ \mathbf{0}&\bm{\xi}_{i}\bm{\xi}_{i}^{\intercal}\end{array}\right].$$

Full text

††journal: Electric Power Systems Research. ©2020 CC-BY-NC-ND 4.0

1 Introduction

Unit commitment (UC) plays a central role in the operation of electric power systems and markets as the problem of finding the most economical on/off status and output dispatch of power generators connected to the network over the course of (most commonly) the next day [1], such that physical and engineering constraints associated with power generation, transmission, and consumption components of the system are observed. It is a combinatorial optimization problem whose computational complexity scales exponentially with the size of the system under consideration. Moreover, nonlinearity of power flow physics makes this a challenging problem despite its long history in the academic literature [2], more so in light of the increasing trend in the penetration of non-dispatchable renewable energy sources (RES).

The taxonomy of unit commitment is broadly described in terms of its optimization constraints and solution methods, regardless of the context of its application, which usually defines whether the objective function reflects the economic welfare or costs associated with fuel, and unit startup and shutdown. With respect to its solution, given the inherent non-convex nature of the problem, algorithms proposed in the literature can be roughly classified, as cleverly noted by Fattahi et al. [3], into “methods based on a single convex [approximation] model, …a series of convex models, [or] heuristics and local-search …” in which the key driver to their choice should be the set of constraints considered in the formulation. The basic set will have included unit ramping, minimum up- and down-time, limits on generation, and power balance equations in discrete time. Other constraint sets may also be formulated, such as reserve, security and, more generally, power flow constraints in either DC or (rectangular) AC formulations. More sophisticated models may also consider the stochasticity of RES.

1.1 Related work

A properly comprehensive bibliographical review of solution methods proposed for the UC problem and its many variants is out of our intended scope for the paper. Therefore, this section traces back past work related to a lineage of algorithms of interest to AC transmission-constrained problem formulations, especially those employing convex relaxation. Since its inception, several techniques to solve the UC problem have been reported [4], with dynamic programming (DP) having been explored in early works. Ayoub and Patton [5] presented a fully linearized formulation framework [5], for which a so-called security function is proposed such that cold starts are penalized. This function has also been used to assess system security in an hourly, probabilistic basis [6]. Snyder et al. [7] and Kumar and Palanisamy [8] have considered “a priori information” for the removal of infeasible paths to allow for computational tractability. Because of its combinatorial nature, the use of heuristics has then been shown to be suitable for large instances of the problem [9], e.g. Sasaki et al. [10] proposed the use of greedy algorithms ordered by the average operation cost–these solutions are then evaluated for full supply of demand, and tested for feasibility by means of an optimal power flow (OPF).

In Lagrangian relaxation (LR) applications to UC it was common to relax reserve and demand coupling constraints in order to create individual subproblems for every unit. Other common approaches include the representation of individual unit subproblems as mixed-integer problems [11] with techniques to select the Lagrange multipliers that maximize lower bounds produced by the relaxation, as well as the identification of identical solutions [12] that cause dual solutions to be far from the optimum by means of successive subproblems.

Other efficient methods to solve linear UC formulations were based on decomposition, such as Benders’ [13]. One common approach of Benders’ decomposition-based methods is to use piecewise-linear approximations of quadratic costs [14]. Fu et al. [15], in one of the seminal works presenting AC transmission constraints, proposed to approximate the solution to the security-constrained UC master problem using LR and DP, with transmission constraints being checked in the subproblem in a rectangular AC power flow model, resulting in an optimization-simulation decomposition. Constraint violations are added back to the master problem as feasibility cuts in addition to cuts representing network contingencies. Nasri et al. [16] considered the uncertainty of wind power generation in real-time operation as the second of a two-stage stochastic program, with DC network-constrained day-ahead market clearing formulated in the first stage. In their proposed Benders’ decomposition, the UC master problem is followed by a “sufficiently convexified” AC network-constrained subproblem whose “asymptotic convexification” depends on the number of scenarios. Castillo et al. [17] proposed a constraint-set linearization-based solution to the mixed-integer non-linear (MINLP) formulation of the UC problem with AC power flow constraints based on outer approximation (later extended by Liu et al. [18] in a decomposition framework), where an analysis of the economic and operation advantages of AC network-constrained UC over its DC counterpart is provided.

Conic relaxation of polynomial optimization problems has driven interest in different applications in power systems, whether in the form of second-order conic (SOCP) or semidefinite (SDP) programs. Quadratic relaxations of binary variables [19] and rectangular formulations of AC power flow [20, 21, 22] can be formulated as non-convex quadratically-constrained quadratic problems (QCQP) suitable for Shor’s scheme-based SDP relaxations of MINLP formulations of the UC problem [23] with local resolution of relaxed infeasible solutions by means of approximation heuristics. Jabr [24] proposes a rank-1 constraint imposed on the SDP relaxation with an iterative convex rank reduction procedure. Local resolution of the relaxation is done by a greedy heuristic based on average fuel costs at full load in order to yield feasible schedules. Mhanna and Jabr [25] propose the use of a heuristic based on a pruning algorithm, with complementary repair mechanisms that exploit the nature of problem constraints, and the SDP formulation itself. Another single convex model worth citing is that proposed in the work of Fattahi et al. [3] for the DC network-constrained UC problem, in which reformulation-linearization is used to introduce relaxed non-convex quadratic inequalities in order to tighten the SDP relaxation seeking to avoid local feasibility resolution. More recently, Zohrizadeh et al. [26] propose a power loss penalization method based on a series of SOCP relaxations that require an initialization “sufficiently close to the feasible set” for convergence.

1.2 Contributions

This paper contributes to the long line of methods for the solution of the UC problem building upon the recent developments in AC optimal power flow, and conic relaxation, with a QCQP formulation of the transmission constraints projected onto the space of semidefinite cones, for which the discrete unit commitment decisions are solved in a Benders’ decomposition approach in order to overcome solutions with integrality gap due to a relaxation, while taking advantage of a series of convex SDP relaxations of the network. Furthermore, and in specific terms, the following contributions can be highlighted:

Not any sort of initialization with $\epsilon$-feasibility is required.

• 2.

A linearization of the power generation constraints is included in the master problem as a means to approximately represent the optimal power flow subproblem, and thus further constrain the search for integer solutions.

• 3.

In the proposed formulation, the rank of the SDP matrices associated with active power are guaranteed to be 1 at the end of every iteration of the algorithm.

• 4.

Consequently, a rank reduction procedure (RRP) is iteratively performed upon voltage variables, the only set of matrix variables for which rank-1 is not guaranteed.

Finally, two numerical case studies based on commonly available data sets are presented as evidence to illustrate the effectiveness of the proposed methodology.

2 Unit Commitment Problem Formulation

We formulate the UC problem as a mixed-integer SDP problem of finding the unit commitment and active power generation dispatches that minimize costs subject to AC power flow constraints:

[TABLE]

subject to:

[TABLE]

The UC problem (1)–(18) is defined over a variable space consisting of active and reactive power, bus voltages, and discrete decisions describing unit status, and startup/shutdown events. Besides the constraints associated with AC power flow, the constraint set includes minimum up and down times, ramp, and spinning reserve inequalities. Its formulation enables the application of a Benders’ decomposition technique [27] by means of a convex relaxation of the rectangular AC OPF problem based on [28], that seeks to provide optimal solutions to active power generation and UC decisions subject to AC power flow constraints over the voltage variable space.

2.1 Problem dimensions

The UC problem is defined over a set of $g=1,2,\ldots,N_{G}$ generators, $i=1,2,\ldots,N_{B}$ buses, and discrete time stages $t=1,2,\ldots,T$, most commonly in hourly steps, and lines $l=1,2,\ldots,N_{L}$ for every undirected pair $(i,j)$ of buses $1\leqslant i\leqslant N_{B}$, and $1\leqslant j\leqslant N_{B}$, $i\neq j$ for which there is a connecting transmission line.

2.2 Problem variables

2.2.1 Startup/shutdown

Binary vectors $\mathbf{x}$, $\mathbf{y}$, and $\mathbf{z}$ represent up/down statuses, and startup and shutdown events, respectively, for a given generator $g$ at every $t$. If generator $g$ is up-and-running during $t$ then $x_{t,g}$ equals $1$, otherwise it is [math]. On the other hand, if generator $g$ is brought up at time $t$, then $y_{t,g}$ equals $1$, or [math] otherwise. Analogously, if generator $g$ is shutdown at time $t$, then $z_{t,g}$ equals $1$, or [math] otherwise. The indicator function $\textbf{1}_{(t>1)}$ is equal to $1$ for time steps $t>1$, or [math] otherwise. The value of $x_{0,g}$ is assumed to be known. Therefore:

[TABLE]

2.2.2 Power generation

Active power generation $p_{t,g}$ is subject to minimum $\underline{p}_{g}$ and maximum $\overline{p}_{g}$ limits. Let:

[TABLE]

define the incremental active power generation if generator $g$ is up-and-running. Then $\overset{\Delta}{p}_{t,g}$ can be alternatively written as a Frobenius product $\mathbf{A}\bullet\mathbf{P}_{t,g}$, where:

[TABLE]

and the Frobenius product of two matrices is defined as the sum of element-wise products, such that:

[TABLE]

Similarly, reactive power generation is represented as the difference with the minimum reactive power so:

[TABLE]

[TABLE]

2.2.3 Voltage

Real $e_{t,i}$ and imaginary $f_{t,i}$ parts of voltage at bus $1\leqslant i\leqslant N_{B}$ define the voltage vector as follows:

[TABLE]

such that the voltage variable $\mathbf{V}_{t}\in\mathbb{R}^{2\cdot N_{B}\times 2\cdot N_{B}}$ is given by the outer product of the voltage vector:

[TABLE]

2.3 Objective function

The objective function in (1) represents costs associated with active power generation, as well as generator startup and shutdown events. Let the generation costs of an up-and-running generator $g$ during step $t$ be approximated by a quadratic function of the form:

[TABLE]

then it can alternatively be expressed in terms of incremental power generation:

[TABLE]

In (1) the incremental costs are expressed as a Frobenius product $\mathbf{C}_{g}\bullet\mathbf{P}_{t,g}$, such that:

[TABLE]

and the minimum generation costs are given by the inner product $\mathbf{c}^{\intercal}\mathbf{x}$, such that:

[TABLE]

where:

[TABLE]

is constant over time, and $\alpha_{g}$, $\beta_{g}$, and $\gamma_{g}$ are the coefficients of the quadratic active power generation cost function for generator $g$.

Finally, vectors $\mathbf{u},\mathbf{h}\in\mathbb{R}_{+}^{T\cdot N_{G}}$ represent cost coefficients of startup and shutdown events for every generator $g$ and time $t$, respectively.

2.4 Power balance constraints

2.4.1 Active power balance

In Eq. (2) active power balance constraints are expressed as the difference between active load demand $dp_{t,i}$ and the sum of total active power generation and net active power flow to bus $i=1,2,\ldots,N_{B}$ in every $t=1,2,\ldots,T$. Active power generation is expressed as the sum of incremental $\mathbf{A}\bullet\mathbf{P}_{t,g}$ and minimum generation $\underline{p}_{g}\cdot x_{t,g}$ at all generators $g$ connected to bus $i$, as represented by set $\mathcal{G}_{i}$. Net active power flow injection to bus $i$ is given by $\mathbf{Y}_{i}\bullet\mathbf{V}_{t}$, such that:

[TABLE]

where $\Omega_{i}$ is the set of bus indexes of all buses directly connected to bus $i$, so that the pair $(i,j)$ exists, and:

[TABLE]

is the admittance matrix, and:

[TABLE]

where $g_{i,j}$ and $b_{i,j}$ represent line conductance and susceptance, respectively, and $\bm{\xi}_{i}\in\mathbb{R}^{N_{B}}$ is the standard basis vector.

2.4.2 Reactive power balance

In Eq. (3) reactive power balance constraints are expressed as the difference between reactive load demand $dq_{t,i}$ and the sum of total active power generation and net reactive power flow to bus $i=1,2,\ldots,N_{B}$ in every $t=1,2,\ldots,T$. Reactive power generation is expressed as the sum of incremental $\overset{\Delta}{q}_{t,g}$ and minimum generation $\underline{q}_{g}\cdot x_{t,g}$ at all generators $g$ connected to bus $i$. Net reactive power flow injection to bus $i$ is given by $\widetilde{\mathbf{Y}}_{i}\bullet\mathbf{V}_{t}$, such that:

[TABLE]

and:

[TABLE]

2.5 Active power constraints

Maximum active and reactive power generation constraints are expressed in Eq. (4) and Eq. (5) respectively, for $g=1,2,\ldots,N_{G}$, and $t=1,2,\ldots,T$ in terms of incremental power generation.

2.6 Startup and shutdown constraints

In Eq. (6) startup constraints are represented for every $g=1,2,\ldots,N_{G}$, and $t=1,2,\ldots,T$. Analogously, shutdown constraints are represented in Eq. (7) along the same dimensions.

2.7 Minimum up and down time constraints

In Eq. (8), and Eq. (9) minimum up $T_{on}$ and down $T_{off}$ time constraints are represented, respectively, for every $g=1,2,\ldots,N_{G}$, and $t=1,2,\ldots,T$, where the sum of time steps (second term) in up (down) time after a startup (shutdown) event occurs (first term) is equal to the minimum number of time steps $\omega_{(t,\tau)}$ required for a unit to be up (down). This formulation is an adaptation of that presented in [23]. Constant value $T_{0}$ represents a unit status immediately before the beginning of the considered time horizon—it is positive if unit was up, or negative otherwise, where:

[TABLE]

2.8 Spinning reserve constraints

In Eq. (10) constraints on residual active power generation availability $SR_{t}$, also known as spinning reserve, are imposed for every $t=1,2,\ldots,T$.

2.9 Maximum ramp up and down rate constraints

Active power generation is subject to maximum ramp up $RU_{g}$ and down $RD_{g}$ rates between consecutive time steps, as expressed by Eq. (11), and Eq. (12), respectively, for every $g=1,2,\ldots,N_{G}$, and $t=1,2,\ldots,T$. Where $p_{0,g}$ represents active power dispatch and its current value is assumed to be known.

2.10 Power flow constraints

In Eq. (13) maximum power flow constraints $\overline{F}_{i,j}$ are imposed on power flow $\mathbf{Y}_{i,j}\bullet\mathbf{V}_{t}$ in both directions of every line $l=1,2,\ldots,N_{L}$ such that $(i,j)$ exists.

Remark 1**.**

This constraint represents a limit on active power flow in the network, and does not cover limits on reactive power flow.

2.11 Voltage constraints

The module of the voltage is subject to minimum $\underline{V}_{i}$ and maximum $\overline{V}_{i}$ limits in every $i=1,2,\ldots,N_{B}$, as expressed in Eq. (14), where:

[TABLE]

2.12 Rank constraints

In Eq. (17) a rank-1 constraint is imposed on $\mathbf{V}_{t}$ as a consequence of the transformation in Eq. (19), and the relaxation $\mathbf{V}_{t}-\mathbf{v}\mathbf{v}^{\intercal}\succeq 0$ implicit to Eq. (15).

On the other hand, if matrix $\mathbf{P}_{t,g}\succeq 0$, and therefore, its terms $P_{t,g,1,2}=\overset{\Delta}{p}_{t,g}$ and $P_{t,g,2,2}$ equal $1$, then the following holds for the determinant of a semidefinite matrix:

[TABLE]

Remark 2**.**

If $\alpha>0$, which can be reasonably assumed since thermoelectric power generation can be satisfactorily approximated by convex functions, we can deduce that the relationship $\mathbf{C}_{g}\bullet\mathbf{P}_{t,g}\geqslant\alpha_{t,g}\cdot p_{t,g}^{2}+\beta_{t,g}\cdot p_{t,g}$ is implicit to $\mathbf{P}_{t,g}\succeq 0$, and, in a minimization algorithm, this relationship will be active since $P_{t,g,1,1}$ does not have any other lower bound. Furthermore, from such observation we can also deduce that in the optimal solution matrix $\mathbf{P}_{t,g}$ is rank-$1$.

3 Solution Methodology

In order to solve the mixed-integer SDP relaxation of the UC problem defined in (1)–(18), we propose the use of a generalized Benders’-like decomposition scheme. Initially, a feasible UC schedule is obtained with a relaxed bus voltage active power dispatch. In [28], the OPF problem is numerically perturbed in order to yield rank-2 solutions, whereas in this paper we choose to use a rank reduction procedure (RRP) that avoids such problem-dependent perturbation, and also provides low-rank solutions to the voltage variables with optimal active power generation and UC dispatches.

In the algorithm presented in the following sections, a decomposition scheme for polynomial optimization problems [29] is adapted according to an application of the generalized Benders’ technique for mixed-integer variables [27] to the UC problem, therefore rendering a Benders’-like decomposition method for mixed-integer semidefinite programming (MI-SDP) problems. The method is summarized in Fig. 1.

3.1 Benders’-like decomposition

We can write the UC problem (1)–(18) in matrix notation in order to simplify the development of our Benders’ decomposition scheme. With that goal in mind, let $\mathbf{C}\in\mathbb{R}^{(2\cdot T\cdot N_{G})\times(2\cdot T\cdot N_{G})}$ represent the incremental cost matrix for all generators and time steps, such that:

[TABLE]

Consequently, we need to represent incremental active power generation variables $\mathbf{P}\in\mathbb{R}^{(2\cdot T\cdot N_{G})\times(2\cdot T\cdot N_{G})}$ as a block-diagonal matrix along the same dimensions, with $\mathbf{P}_{t,g}\in\mathbb{R}^{2\times 2}$ in its diagonal, such that:

[TABLE]

A voltage matrix $\mathbf{V}\in\mathbb{R}^{(2\cdot T\cdot N_{B})\times(2\cdot T\cdot N_{B})}$ is constructed in analogous terms:

[TABLE]

Then the resulting equivalent UC problem in matrix notation with relaxation of the rank constraints is expressed as follows:

[TABLE]

[TABLE]

[TABLE]

Let $\bm{\varrho}_{(t,g)}\doteq\bm{\xi}_{N_{G}\cdot\left(t-1\right)+g}$ define a standard basis vector of appropriate dimension.

Active power balance constraints (2) are expressed in Eq. (21) by linear mappings $\mathcal{A}\mathbf{1}(\cdot)$ and $\mathcal{Y}\mathbf{1}\left(\cdot\right)$, representing incremental active power generation and net power flow injection, respectively. They are composed of matrices $\mathbf{A1}_{r}\in\mathbb{R}^{(2\cdot T\cdot N_{G})\times(2\cdot T\cdot N_{G})}$ and $\mathbf{Y1}_{r}\in\mathbb{R}^{(2\cdot N_{B})\times(2\cdot N_{B})}$, respectively:

[TABLE]

where $r=N_{B}\cdot(t-1)+i$, and $\otimes$ is the Kronecker product. In addition, matrix $\underline{\mathbf{P1}}\in\mathbb{R}^{(T\cdot N_{B})\times(T\cdot N_{G})}$ represents minimum active power generation, such that:

[TABLE]

and, finally:

[TABLE]

Reactive power balance constraints (3) are expressed in Eq. (22) by vectors $\mathbf{q}$, $\underline{\mathbf{q}}$ representing incremental and minimum reactive power generation respectively. Vector $\underline{\mathbf{q}}$ is composed by:

[TABLE]

and linear mapping $\mathcal{Y}\mathbf{2}(\cdot)$ representing minimum reactive power generation and net power flow injection, is composed by matrices $\mathbf{Y2}_{r}$:

[TABLE]

where $r=N_{B}\cdot(t-1)+i$, also the reactive load demand is:

[TABLE]

Active power constraints (4) are expressed in Eq. (23) in terms of linear mapping $\mathcal{A}\mathbf{2}(\cdot)$, which is composed of matrices $\mathbf{A2}_{r}\in\mathbb{R}^{(2\cdot T\cdot N_{G})\times(2\cdot T\cdot N_{G})}$, and $\overline{\mathbf{P1}}\in\mathbb{R}^{(T\cdot N_{B})\times(T\cdot N_{G})}$:

[TABLE]

Reactive power constraints (5) are expressed in Eq. (24) with upper bound vector:

[TABLE]

Startup (6) and shutdown (7) contraints are expressed in Eqs. (25) and (26), respectively, in terms of $\mathbf{U}\in\mathbb{R}^{(T\cdot N_{G})\times(T\cdot N_{G})}$ and $\mathbf{a1}\in\mathbb{R}^{T\cdot N_{G}}$:

[TABLE]

Minimum up (8) and down (9) time constraints are expressed in Eqs. (27) and (28), respectively, in terms of matrices $\mathbf{M1},\mathbf{M2}\in\mathbb{R}^{(T\cdot N_{G})\times(T\cdot N_{G})}$, such that:

[TABLE]

and:

[TABLE]

where $\mathbf{a2}\in\mathbb{R}^{T\cdot N_{G}}$ is defined as:

[TABLE]

Spinning reserve constraints (10) are expressed in Eq. (29) in terms of $\mathcal{A}\mathbf{3}(\cdot)$, which is composed of time matrices $\mathbf{A3}_{t}\in\mathbb{R}^{(2\cdot T\cdot N_{G})\times(2\cdot T\cdot N_{G})}$, $\overline{\mathbf{P2}}\in\mathbb{R}^{T\times(T\cdot N_{G})}$, and $\mathbf{r}\in\mathbb{R}^{T}$:

[TABLE]

Maximum ramp up (11) and down (12) constraints are expressed in Eqs. (30) and (31), respectively, in terms of linear mapping $\mathcal{A}\mathbf{4}(\cdot)$, which is composed of time matrices $\mathbf{A4}_{N_{G}(t-1)+g}\in\mathbb{R}^{(2\cdot T\cdot N_{G})\times(2\cdot T\cdot N_{G})}$:

[TABLE]

and $\overline{\mathbf{P3}}\in\mathbb{R}^{(T\cdot N_{B})\times(T\cdot N_{G})}$:

[TABLE]

as well as $\mathbf{g1}$ and $\mathbf{g2}\in\mathbb{R}^{T}$:

[TABLE]

for $t=1,\ldots,T$ and $g=1,\ldots,N_{G}$. Where $p_{0,g}$ is active power dispatch at the previous period assumed known.

Power flow constraints (13) are expressed in Eq. (32) in terms of linear mapping $\mathcal{Y}\mathbf{3}\left(\cdot\right)$, which is composed of time matrices $\mathbf{Y3}_{r}\in\mathbb{R}^{(2\cdot N_{B})\times(2\cdot N_{B})}$, and $\overline{\mathbf{f}}\in\mathbb{R}^{2\cdot T\cdot|\mathcal{L}|}$:

[TABLE]

where $r=N_{L}(t-1)+l$, for all $t=1,2,\ldots,T$, and $(i,j)\in\mathcal{L}$.

Voltage constraints (14) are expressed in Eq. (33) in terms of linear mapping $\mathcal{Y}\mathbf{4}\left(\cdot\right)$, which is composed of time matrices $\mathbf{Y4}_{r}\in\mathbb{R}^{2\cdot N_{B}\times 2\cdot N_{B}}$, and $\mathbf{b}\in\mathbb{R}^{2\cdot T\cdot N_{L}+2\cdot N_{B}}$:

[TABLE]

where $r=\left(N_{B}-1\right)t+i$, for all $t=1,2,\ldots,T$, and $i=\{1,\ldots,N_{B}\}\backslash i_{slack}$.

3.1.1 Unit Commitment Master Problem

Problem (20)–(37) is decomposed into a master problem ($M^{(k)}$) and a subproblem ($S^{(k)}$) at every $k$-th iteration, where the latter receives solutions from the former, returning either feasibility or optimality cuts to be added as constraints to the UC master problem, formulated as follows:

[TABLE]

3.1.2 Optimal Power Flow Subproblem

Subproblem $S^{(k)}$ is parameterized by the solution $\widetilde{\mathbf{x}}^{(k)}$ to master problem $M^{(k)}$, being modified accordingly to return feasibility and optimality cuts. This is accomplished by adding a slack variable $s$ penalized in the objective function by a factor $\sigma\gg 0$. Good numerical stability in the optimization process is obtained in our computational experiments by choosing

[TABLE]

which represents both an upper bound to objective function (1) and a sufficiently big number. If either $s>0$ or its corresponding constraint in the dual problem ($SD^{(k)}$) is active in the optimal solution, then ($S^{(k)}$) is infeasible, and the solution ($\mu_{1}^{(k)},\ldots,\mu_{7}^{(k)}$) to ($SD^{(k)}$) is used to construct a feasibility cut (41) to be added to ($M^{(k)}$) as a constraint. If otherwise $s=0$ in the optimal solution, then the OPF subproblem ($S^{(k)}$) is feasible and the solution ($\lambda_{1},\ldots,\lambda_{7}$) is used to construct an optimality cut (42) to be added to ($M^{(k)}$) as a constraint:

[TABLE]

where $\widehat{\mathbf{e}}_{1},\ldots,\widehat{\mathbf{e}}_{9}$ are vectors of ones with appropriate dimension in every equation. In order to obtain a convex subproblem the rank-1 constraint over $\mathbf{V}$ in Eq. (36) is relaxed. The corresponding dual of ($S^{(k)}$) is formulated as follows:

[TABLE]

s.t.

[TABLE]

Feasibility cuts are defined below:

[TABLE]

as well as optimality cuts:

[TABLE]

Depending on the nature of the cut generated, it is added to either the set of feasibility cuts:

[TABLE]

or the set of optimality cuts:

[TABLE]

In OPF problems, zero duality gaps of SDP relaxations are not generally guaranteed, as extensively recorded in the recent literature. In [30] the authors have shown that low order SDP relaxations are sufficient to obtain the global solution of a wide variety of OPF problems. There exist, however, some cases for which SDP relaxations do not guarantee zero duality gaps, as shown in [31] and [32]. Different techniques for dealing with non-zero duality gaps in such problems have been subject of open research. In [28] the authors have shown that minor modifications to line resistances in IEEE systems yield rank-2 solutions that enable further rank-1 solution recovery. It has been observed in [32] that non-convexities in the transmission system can be mitigated by the penalization of apparent losses in order to reduce the difference of angles in the corresponding lines. An analogous approach is presented in [30], where the authors propose higher order relaxation of the insulting lines. Other alternative is to use SDP relaxations in a branch-and-bound scheme [33]. In this work, a rank reduction procedure is proposed as means to deal with non-zero duality gap of SDP relaxations, as described in Section 3.2.

3.1.3 Modified Master Problem

A linearization of ($S^{(k)}$) is performed into the modified master problem ($MM^{(k)}$) as a means to represent the subproblem into the master problem, thus constraining integer solutions, and consequently improving convergence of the overall problem, especially in the initial iterations, where, if we don’t use this linearization in the master problem, it will result in many infeasible unit commitment solutions to the subproblems. At some point of the iterations it’s expected that this linearization will be less tight than the Benders’ cuts and could also be neglected. Therefore, variables $p_{t,g}$ are introduced in the master problem, as well as a corresponding relaxed model representing active power generation constraints from the OPF subproblem, consisting of relaxed power balance equations:

[TABLE]

[TABLE]

for $t=1,2,\ldots,T$, where $l_{t}^{(k-1)}$, representing power losses updated from ($S^{k-1}$), is assumed as a ratio of total load at the beginning of each iteration. At optimality, Eq. (43) and Eq. (44) are active. Also, spinning reserve constraints

[TABLE]

the equivalent for reactive generation:

[TABLE]

for $t=1,2,\ldots,T$, minimum and maximum active and reactive power generation

[TABLE]

[TABLE]

ramp up

[TABLE]

and, finally, ramp down constraints

[TABLE]

for $t=1,2,\ldots,T$ and $g=1,2,\ldots,N_{G}$. Additionally, objective function (1) is approximated by a linear function as illustrated in Fig. 2. This linear approximation represents a lower bound to the quadratic fuel costs, and is calculated by using the mean value theorem of convex functions:

[TABLE]

Furthermore, without any loss of generality, we choose:

[TABLE]

and finally obtain a linear lower bound of the cost function represented by $w+\mathbf{c}^{\intercal}\mathbf{x}$:

[TABLE]

where:

[TABLE]

and

[TABLE]

for $t=1,2,\ldots,T$, and $g=1,2,\ldots,N_{G}$.

Finally, the modified master problem is formulated below:

[TABLE]

where $\mathbf{H}\mathbf{p}+\mathbf{G}\mathbf{x}\leq\bm{\pi}$ represents constraints (43), (45), (47), (49), and (50), $\widetilde{\mathbf{H}}\mathbf{q}+\widetilde{\mathbf{G}}\mathbf{x}\leqslant\widetilde{\bm{\pi}}$ represents constraints (44), (46), (48), and feasibility $\mathcal{F}_{i}$ and optimality $\mathcal{O}_{i}$ cuts are added iteratively. Note that $w\geqslant 0$ is incurred from the cost function lower bound in Eq. (51).

3.2 Rank Reduction Procedure

As demonstrated in Section 2.12, the rank of $\mathbf{P}^{*}$ is guaranteed to be one at every iteration. A rank reduction procedure [34] is necessary, however, such that a perturbation matrix $\mathbf{D}_{r}$ in the rank reduction iteration $r$ is added to the voltage matrix $\mathbf{V}^{*}$ that solves problem (20)–(37) at the optimum of the Bender’s algorithm. The rank reduction procedure listed in Algorithm 1 maintains feasibility as explained in the remainder of this section. Moreover, because $\mathbf{V}$ does not appear in the objective function, new solutions resulting from the rank reduction procedure are still optimal. A central idea to this procedure is to calculate $\mathbf{D}_{r}$ such that at every iteration, $\mathbf{V}^{*}_{r}$ does not violate any constraints, has its rank reduced, and be positive semidefinite.

An upper bound on the rank of the final optimal solution to subproblem ($S^{(k)}$) can be expressed as a function of the number of problem constraints, as shown in [35]:

[TABLE]

such that $\mathbf{V}^{*}_{r}$ is updated as follows:

[TABLE]

where $\omega_{r}$ is a scalar parameter that guarantees that every update to $\mathbf{V}^{*}_{r}$ remains positive semidefinite. It follows that the perturbation matrix $\mathbf{D}_{r}$ must be subject to the following constraints:

[TABLE]

where $\mathbf{D}_{r}$ has the following form:

[TABLE]

where $\mathbf{R}_{r}$ is obtained from the Cholesky decomposition of $\mathbf{V}^{*}_{r}=\mathbf{R}_{r}\mathbf{R}_{r}^{\intercal}$. By the distributive property of the Frobenious product we obtain the following relationship:

[TABLE]

In order to satisfy (58), $\mathbf{Z}_{r}$ is constructed as a [math]-$1$ matrix such that there is a one-valued element for every zero-valued element in the corresponding coordinate of matrices $\left(\mathbf{R}_{r}^{\intercal}\mathbf{Y}_{i}\mathbf{R}_{r}\right)$, $\left(\mathbf{R}_{r}^{\intercal}\mathbf{Y}_{i,j}\mathbf{R}_{r}\right)$, and $\left(\mathbf{R}_{r}^{\intercal}\mathbf{E}_{l}\mathbf{R}_{r}\right)$. This guarantees that (54), (55), and (56) are satisfied.

Substituting (57) in (53), we obtain:

[TABLE]

From (59) we conclude that the necessary condition for positive semidefiniteness of the new solution is given by:

[TABLE]

If $\mathbf{Z}_{r}\preceq\mathbf{0}$, then a multiplication by $-1$ follows. This is characterized by the value assumed by parameter $\phi_{z}$. If $\mathbf{Z}_{r}$ is indefinite, then the algorithm is not able to reduce its rank at the $r$-th iteration. However it is still possible to reduce its rank in the following iterations since $\mathbf{Z}_{r}\neq\mathbf{Z}_{r+1}$ because of the added perturbations, and also because the SDP formulation of the network constraints guarantees that at least one rank-1 solution exists. In that sense, additional $\varsigma$ iterations are provisioned. Parameter $\omega_{r}$, on the other hand, is a scale of perturbation matrix $\mathbf{D}_{r}$ that induces rank reduction of $\mathbf{V}^{*}_{r}$ in at least one degree by reducing the rank of $\left(\mathbf{I}-\omega_{r}\mathbf{Z}_{r}\right)$. This is achieved by:

[TABLE]

where $\lambda_{z}$ is an eigenvalue of $\mathbf{Z}_{r}$.

Remark 3**.**

If $\text{rank}(\mathbf{V}_{r}^{*})>1$ then feasibility of the solution must be recovered. This can be achieved, without any guarantees of global optimality, by means of a local resolution, either approximately by, e.g. a convex-concave procedure [36], or locally, using a method to solve the AC OPF with strong local convergence properties, e.g. an infeasible primal-dual interior-point method with line search filter [37].

4 Numerical Results

The Benders’-like decomposition algorithm described in Section 3.1 was implemented in MATLAB® 8.1 with solution of the mixed-integer master problem by the IBM® CPLEX® 12.6 solver using the branch-and-cut algorithm, and solution of the OPF subproblem by SDP relaxation using SDPA 7.3.8.

In order to validate the methodology we have constructed two numerical cases: 3-GEN and IEEE-118. The purpose of the small-scale 3-GEN case is to didactically serve to showcase the proposed solution to the UC problem. It has a 6-bus network with limits on power flow, and 3 generators with different cost functions, as listed on Table 1. In this case generator G1 is less costly and has higher capacity, whereas generators G2 and G3 present equivalent costs depending on the range of operation. On the other hand, the purpose of the case study based on the IEEE 118-bus dataset111Problem data can be obtained from the following file available online: http://motor.ece.iit.edu/Data/SCUC_118.xls and also by request to the authors. is to test the effectiveness of the proposed methodology in a larger-scale problem. In order to speed up the convergence of the algorithm, we adopted a strategy that assumes initial estimated power losses to be between 5% and 10% of total system load.

In Fig. 3 we can observe the evolution of both the upper and lower bounds during the algorithm iterations while solving the 3-GEN case. As expected, the gap between the two bounds is very large in the first iteration. However, in only 6 iterations upper and lower bounds converge to an infinitesimal gap with 3 feasibility cuts and 3 optimality cuts.

A characteristic of Benders’-like decomposition methods is that the solutions of the integer master problem are evaluated as they serve as parameters of the subproblem. This interaction between the master and subproblem explains the behaviour observed in Fig. 4, which illustrates how the unit commitment schedule of generator G2 is modified over the iterations as it adapts to the inclusion of power system constraints—e.g. in this case, feasibility cuts associated with spinning reserve and ramp rate constraints.

Unit commitment and power generation results are illustrated in Fig. 5a and Fig. 6a, respectively. It can be observed in Table 1 that generator G1 is the least costly, resulting in it being dispatched at all times as illustrated in Fig. 5, whereas more costly generators G2 and G3 complement generation to meet the load demand and spinning reserve requirements. This observed behaviour in the final solution is a numerical evidence of consistency towards solution optimality and feasibility. In order to illustrate the effectiveness of the proposed methodology, it is compared with a simple heuristic deemed herein as MILP-AC. In this heuristic, the unit commitment problem is solved for linearized network constraints with DC power flow, and piecewise-linear costs as a mixed-integer linear program whose unit commitment results are shown in Fig. 5b. Upon optimality, a power factor of .85 is used to calculate active and reactive power generation set points from the optimal DC power flow. These set points are used to initialize the AC power flow for the final generation dispatches as shown in Fig. 6b. The MILP-AC heuristic was implemented in the Python programming language using PyPSA [38] and Pyomo [39] with IBM® CPLEX®. Pyomo was used to define the piecewise-linear cost functions, and spinning reserve and initial ramping constraints not readily available in PyPSA version 0.17 at the time of writing. In the 3-GEN case, both the Benders’ $MM^{k}$ decomposition with SDP relaxation of AC power flow constraints and the MILP-AC heuristic yielded UC and active power feasible solutions, with a cost difference of 0,52% lower for MILP-AC at $93,404.52.

In the IEEE-118 case, generation unit configurations similar to those of the 3-GEN case were used. Convergence of the larger-scale case is also achieved in only six iterations as illustrated in Table 3. In order to test the proposed methodology for its effectiveness in solving the AC optimal power flow, we reduced the limits of some power lines. This additional restraining is dealt with effectively by the convex SDP relaxation. As a comparison, in this modified IEEE 118-bus case the MILP-AC solution is infeasible for active power generation after the ensuing AC power flow with a 14% cost increase at $1,958,123.15.

The low number of iterations required for convergence of our Benders’-like decomposition algorithm is a result of the use of the modified master problem ($MM^{(k)}$), as evidenced by additional numerical experiments. We verified that convergence was reached in 13 iterations for the same 3-GEN case study if modification of the master problem is not present as described in Section 3.1.3, as illustrated in Table 2. It can be observed that convergence of the IEEE-118 case is not even possible without the use of ($MM^{(k)}$). It is important to note that in these numerical experiments with the use of ($M^{(k)}$), we did not consider estimates for power losses, which can constitute an additional difficulty for convergence.

Additionally, we have tested our methodology without RRP in order to assess its effectiveness in the presented case studies. In the 3-GEN case, all voltage matrices were rank-1. In the IEEE-118 case, on the other hand, the maximum rank of voltage matrices were rank-13 at the final of the iterations.

5 Conclusion

In this paper we presented a formulation of the UC problem with AC power flow constraints over the voltage variable space solved by SDP relaxation in a Benders’-like decomposition approach. Our proposed methodology incorporates linearized formulations of constraints of the OPF subproblem into a modified master problem, thus improving the quality of generated cuts. This modification of the UC master problem has allowed convergence of the algorithm in very few iterations, with low rank voltage matrices resulted from our proposed RRP, in both the three-generator small-scale and the IEEE 118-bus-based larger-scale case studies, as illustrated in the numerical experiments. Our Benders’-like decomposition of the problem with convex relaxation of AC OPF constraints enables the incorporation of contingency analysis to account for scenarios of system failure analogously to the work of Fu et al. [40] in future developments. Moreover, because of the nature of the cuts generated from a SDP relaxation, it should be expected such cuts to be more conservative as compared to traditional linear formulations in the literature. Likewise, it is also expected to increase efficiency of the algorithm by adapting strategies proposed in the literature to improve cut selection in linear formulations [41, 42, 43, 44]. Future work also includes an assessment of the feasibility and optimality of the solutions found by our method in comparison to other solvers, as well as the study of the effects of OPF rank-1 solution recovery approaches on the solution to the UC problem.

Appendix A Additional 3-GEN Case Data

Acknowledgements

This research was funded in part by the São Paulo Research Foundation (FAPESP) under grant numbers 13/07570-7 and 14/08972-4. The authors would like to thank Professor Mohammad Shahidehpour of Illinois Institute of Technology for his invaluable insights in the early versions of this manuscript.