Compact Optimization Learning for AC Optimal Power Flow

Model 1 presents the AC-OPF formulation [32]. The power network can be represented as a graph $(\mathcal{N},\mathcal{E})$ where $\mathcal{N}$ denotes the set of bus indices containing generators and load units, and $\mathcal{E}$ is the set of transmission line indices between two buses. $(lij)\in\mathcal{E}$ where $l$ is a branch index connected from node $i$ to node $j$. The set $\mathcal{E}^{R}$ captures the reversed orientation of $\mathcal{E}$, i.e., $(lji)\in\mathcal{E}^{R},\>\forall(lij)\in\mathcal{E}$.

A generator output $S^{g}=p^{g}+jq^{g}$ is a complex number, where the real part $p^{g}$ is an active power generation (injection) and the imaginary part $q^{g}$ is a reactive power generation. An AC voltage $V=v\angle\theta$ is represented by a voltage magnitude $v$ and a voltage angle $\theta$. The objective function  (1) minimizes the sum of quadratic cost functions $c_{i}(\cdot)$ with respect to the active power generations $p_{i}^{g},i\in\mathcal{G}$. At the reference bus $r$, the voltage angle is set to zero as defined in constraint (2). Constraints (3a), (3b), and (4) capture the bounds on variables $p^{g}$, $q^{g}$, and $v$, respectively. Constraints (5a), (5b) represent the complex power flow for each transmission line, which is governed by Ohm’s law. Here, $Y_{l}$, $b_{l}^{c}$, and $T_{l}$ are the series admittance, line charging susceptance, and transformer parameter at each branch $l$, respectively. Constraints (6) ensures that, at each bus, the active and reactive power balance are satisfied. Here, $Y^{s}$ is the bus shunt admittance and $\mathcal{G}_{i}$ represents the set of generator indices attached to the bus $i$. Constraints (7) capture the thermal limits and ensure that the apparent power does not exceed its limit $\bar{s}_{ij}$ for every transmission line.

For the sake of simplicity, in what follows, the input configuration parameters and the optimal solution to the AC-OPF are denoted by $x$ and $y^{*}$, respectively.

The end-to-end learning approach in this paper consists in using supervised learning to find a mapping from an input $x$ to an optimal solution $y^{*}$ to AC-OPF. This work assumes that the set of commitment decisions and generator bids are predetermined. This setting is common in prior end-to-end learning studies for AC-OPF (e.g., [2, 15, 10]). Existing supervised learning schemes exploit the instance data $\{x_{i},y^{*}_{i}\}$. However, it often suffers from the high dimensionality of the output when dealing with power networks of industrial sizes. Table I reports the specifications of the nine power networks from PGLib [32] and a realistic version of the French system (denoted by France_2018) used in the experiments. The French system¹¹1Refer to [33] for more details on this test case. uses the realistic grid topology of the French transmission system captured in 2018 with annual time series data fot the renewable generation capacity and load demand.

For the PGLib test cases, each instance has different active and reactive load demands $p^{d}$ and $q^{d}$, i.e., $x:=\{p^{d},q^{d}\}$. For France_2018, the generation capacity of the renewable generators is also varying in addition to the load demands. As such, for the PGLib cases, the dimensions of $x$ is $\dim\left(x\right)=2|\mathcal{L}|$, and for France_2018, the input of the mapping is increased to $x:=\{p^{d},q^{d},\{\overline{p}_{i}^{g}\}_{i\in\mathcal{G}^{r}}\}$, where $\overline{p}$ represents the upper bounds of the active generation, and $\mathcal{G}^{r}$ is the set of renewable generator indices. Note that among 1890 generators in this system, there are 1609 renewable non-dispatchable generators including hydro, solar, and wind generators.

From Table I, observe that the dimension of the output $y$ is much higher than that of $x$: this implies that the DNN for the OPF will necessitate an excessive number of trainable parameters. It is the goal of this paper to propose a scalable approach that mitigates this curse of dimensionality.

Note that one can recover the whole optimal solution estimates from the active generations and voltage magnitudes by solving a power flow problem as in [2]. This could be combined with the methods proposed herein to provide further reduction in the output space at the cost of a more costly training or inference procedure. However considering the full output space, i.e., $y:=\{p^{g},q^{g},v,\theta\}$, makes it possible to use Lagrangian Duality [3, 15] or Primal-Dual Learning [9] frameworks to improve the learning procedure further.

The idea underlying Compact Learning is to jointly learn the principal components and the mapping between the original inputs and the outputs in the subspace of the principal components. Figure 3 contrasts the overall architecture of the Compact Learning with the conventional plain learning approach. The plain approach learns a mapping $y=M_{\phi}(x)$ from an input configuration $x$ to an optimal solution $y$, where $\phi$ are the associated trainable parameters. Since $y$ is high-dimensional, it is desirable to have a large number of parameters to obtain accurate solution estimates. Compact Learning in contrast learns a mapping $z=R_{\phi}(x)$ from an input configuration $x$ into an output $z$ in the subspace of principal components before recovering an optimal solution prediction as $y=Wz$. The output dimension of $R_{\phi}$, $\dim\left(z\right)$, in Compact Learning should be substantially smaller than the dimension of $y$, making it possible to reduce the number of trainable parameters. In the upcoming sections, the way of learning $R$ and $W$ simultaneously is detailed.

Let $p$ and $d$ be the dimensions of $z$ and $y$, respectively, i.e., $p=\dim\left(z\right)$ and $d=\dim\left(y\right)$. The goal of Compact Learning is to have $p$ substantially smaller than $d$, i.e., $p\ll d$. Matrix $W\in\mathbb{R}^{d\times p}$ is unitary, i.e., $W^{\top}W=I$, and its columns are composed of the orthonormal principal components of the output space associated with $p$ largest eigenvalues. It is obviously possible to obtain $W$ through PCA, but this computation takes a substantial amount memory when there are many instances.

Instead, Compact Learning uses the Generalized Hebbian Algorithm (GHA), also called Sanger’s rule [39], to learn $W$ in a stochastic manner. GHA is specified in Algorithm 1. Using the mini-batch of optimal solutions $\{y^{*(i)}\}$, it first updates the element-wise running mean $\mu$ and variance $\sigma^{2}$ of optimal solution using the momentum parameter $\beta$. The momentum parameter $\beta$ is set to zero for the first iteration, and from the second iteration it is set to a value ranged $(0.0,1.0)$. The optimal solution $y^{*(i)}$ is normalized using the running mean and variance as shown in line 5 in Algorithm 1. GHA uses Gram-Schmidt process to find out the $\Delta W$ of the $p$ leading principal components, where $\mathcal{LT}[\cdot]$ represents the lower triangular matrix. Once $\Delta W$ is determined, $W$ is updated by adding $\gamma\Delta W$, where $\gamma$ is the learning rate that decreases as the epoch number $e$ increases, i.e.,

In the update rule, $\gamma_{\text{init}}$ and $\gamma_{\text{min}}$ are hyperparameters representing the initial and minimum learning rate, respectively.

One can simply use deterministic PCA to define $W$ instead of using GHA. However this requires significant computing time before starting the training, especially when dealing with a large number of data instances of high dimensionality. In contrast, GHA amortizes the computational effort over the training procedure leading to a better handling of the largest test cases.

The training process of Compact Learning is summarized in Algorithm 2 which jointly learns the principal components $W$ and the mapping $R$. Each iteration samples a mini-batch of size $B$ and applies the GHA Algorithm 1 (line 3 in Algorithm 2). The mapping $R$ then produces $z^{(i)}$ for each input $x^{(i)}$ and the output $y^{(i)}$ is obtained through the mapping $Wz^{(i)}$ and a denormalization step (line 5 in Algorithm 2). Finally, the parameters $\phi$ associated with $R$ are updated using backpropagation from the loss function $\mathcal{L}$, which captures the mean absolute error between the prediction $y^{(i)}$ and the optimal solution $y^{*(i)}$ (ground truth).

The predictions from Compact Learning may violate the physical and engineering constraints of the AC-OPF. Some applications are required to address these infeasibilities and this study considers two post-processing methods for this purpose: (1) solving the power flow problem; and (2) warm-starting the exact AC-OPF solver.

The performance of Compact Learning is demonstrated using nine test cases from PGLIB v21.07 and the realistic version of the French power system (denoted as France_2018) as described in Table I. A total of 52,000 instances were generated. 50,000 instances are used for training and the remaining 2,000 instances are tested for reporting the performance results. For the PGLib cases, these instances were obtained by perturbing the load demands as Eq (9) in Section IV: For France_2018, 100,000 instances for training are generated by perturbing the instances in September. Specifically, the upper bounds of active generation of the wind and hydro renewable generators are perturbed by replacing with the samples from $\mathcal{N}(\overline{p}^{g}_{i},0.2(\overline{p}^{g}_{i}-\underline{p}^{g}_{i}))$. Also, the upper bounds of solar generators and load demands are perturbed by multiplying factor sampled from multivariate Gaussian $\mathcal{N}(\mathbf{1},\Sigma)$ where $\Sigma$ is based on the correlation coefficient of $0.8$. Those perturbed upper bound of active generation is ensured to be greater than or equal to the lower bound of it. Also, 2000 realistic test instances are extracted from the instances in October for reporting the performance. As such, for this test case, the distribution of test instances is not necessarily the same as that of the training instances. Note that accurate renewable forecasting would improve the quality of the training instances, but this experiment setting for France_2018 should provide realistic and difficult circumstances where the model may be deployed for real operations. The instances were solved using Pyomo [41] and Ipopt v3.12 [34] with the HSL ma27 linear solver.

The performance of Compact Learning is compared with the plain approach that directly outputs the optimal solution. As in [10], four fully-connected layers followed by ReLU activations are used for the mapping functions for both Compact Learning and plain learning approaches. For the plain approach, two distinct models are experimented; the first model, which is named Plain-Large, has $d$ hidden nodes for each fully-connected layer (where $d$ is the output dimension). The second baseline, which is named Plain-Small, has the $p$ hidden nodes for each layer (where $p$ is the number of principal components considered in Compact Learning). The Compact Learning model has the same number of weight parameters as Plain-Small, but the last layer of the Compact Learning model is learned through GHA. Indeed, Plain-Small has an encoder-decoder structure as the number of the hidden nodes is smaller than that of the output. The ratio of $p$ to $d$, which is also called principal component ratio, is set to $5\%$ for six smaller test cases (up to 6515_rte) and to $1\%$ for four bigger cases. Mini-batch of $64$ instances is used, and the maximum epoch is set to $1,000$. The models are trained using the Adam optimizer [42] with a learning rate of $1\mathrm{e}{\textrm{-}4}$, which is decreased at $900$ epochs by $0.1$. The overall implementation used PyTorch and the models were trained on a machine with a NVIDIA Tesla V100 GPU and Intel Xeon 2.7GHz. For GHA (Algorithm 1), the momentum parameter $\beta$ is set to $0.9999$ from the second iteration. The initial and minimum learning rate ($\gamma_{\text{init}}$ and $\gamma_{\text{min}}$) are set to $1\mathrm{e}{\textrm{-}4}$, and $1\mathrm{e}{\textrm{-}8}$, respectively. The parameter $\epsilon$ to prevent ill-conditioning is set to $1\mathrm{e}{\textrm{-}8}$.

Table III reports the accuracy of the models for predicting optimal solutions. Five distinct models with randomly initialized trainable parameters per method were trained: the results report the average results and the standard deviations (in parenthesis). The table shows the averaged value of optimality gaps and maximum constraint violations on the 2,000 test instances. The optimality gap is calculated as $100\times|\frac{f(y)-f(y^{*})}{f(y^{*})}|$ where $f(y)$ is the objective function (1). It also reports the maximum constraint violations (in per unit), which is computed as

where $g_{i}(y)$ and $h_{i}(y)$ are, respectively, the inequality and equality constraints enumerated in Model 1. The first two sets of columns represent the performance of the plain approaches. Note that, because of the high dimensionality of the output and the limited GPU memory, Plain-Large is only applicable to the four smaller test cases (up to 3022_goc): this is the limitation that motivated this study. When comparing the two plain models, Plain-Large performs better than Plain-Small as Plain-Small trades off the accuracy for scalability. Compact Learning almost always performs better than the two plain approaches. In particular, it produces predictions with significantly fewer violations (sometimes by an order of magnitude), while also delivering smaller optimality gaps on the larger test cases.

Table IV shows the number of trainable parameters of the Compact Learning model and the plain approaches. The architectures of Compact Learning and Plain-Small are exactly the same, except for the last layer: hence the number of trainable parameters of Compact Learning is smaller than those for Plain-Small by the dimension of $W$. Table IV clearly shows that Compact Learning has the smallest number of trainable parameters. Overall, Table III and Table IV indicate that Compact Learning provides an accurate and scalable approach to predict AC-OPF solutions.

One model among five trained models was randomly chosen for testing post-processing approaches and the results were evaluated on the same 2,000 test instances. Table V reports the constraint violations after applying the power flow model seeded with the predictions from the Compact Learning and Plain-Small models. The table also reports the time to solve the power flow problem. The results show that Compact Learning produces power flow solutions with the smallest constraint violations, sometimes by an order of magnitude. The results of Compact Learning are particularly impressive because the majority of the constraints are satisfied after applying the power flow. Note also that the power flow is fast enough to be used during real-time operations, opening interesting avenues for the use of learning and optimization in practice.

Table VI and Figure 4 report the results for warm-starts. The proposed warm-starting approach, WS:Compact(P+D), is compared with the following warm-starting strategies:

• •

  Flat Start: $p^{g}$, $q^{g}$, $v$ are started with their minimum values, and initial $\theta$ is set to zero. This is a default setting without warm-start.

• •

  WS:DC-OPF(P): Motivated from [20], the primal solution of the DC-OPF is used as a warm-starting point for solving AC-OPF.

• •

  WS:AC-OPF(P): The primal solution of the AC-OPF is used as a warm-starting point for solving AC-OPF again.

• •

  WS:Plain-Small(Large)(P): The primal predictions from the plain approaches are used as warm-starting points.

• •

  WS:Plain-Small(Large)(P+D): The primal and dual predictions from the plain approaches are used as warm-starting points.

• •

  WS:Compact(P): The primal predictions from Compact Learning are used as warm-starting points.

• •

  WS:Compact(P+D): The proposed warm-starting approach. The primal and dual predictions from Compact Learning are used as warm-starting points.

To predict dual solutions, an additional mapping function consisting of four fully-connected layers is trained for Compact Learning and the plain approaches. The sizes of these networks are the same as those for the primal solutions. When using both the primal and dual warm-starting, the initial barrier parameter of Ipopt is set to $1\mathrm{e}{\textrm{-}6}$ because the initial warm-starting point is closer to the optimal than the flat start. The convergence tolerance is set to $1\mathrm{e}{\textrm{-}4}$ for all cases.

Table VI reports the elapsed times and the corresponding speed-up ratio for the warm-start approaches. The first key observation is that it is critical to use both primal and dual warm-starts: only using the primal predictions is not effective in reducing the computation times of Ipopt. This is not too surprising given the implementation of interior-point methods. Primal-dual warm-starts however produce significant benefits. WS:Compact(P+D) produces the best results for all test cases. In particular, it yields a speed-up of $15.26$× for 6515_rte. Also, even for the realistic French system (France_2018), WS:Compact(P+D) gives a speed-up of $4.84$×. This is significant given the realism of this test case and highlights the potential of the combination of Compact Learning and optimization to deploy AC-OPF in real operations. Observe also that WS:Compact(P+D) strongly dominates the other approaches, including its (WS:Plain-Small(P+D) counterpart. This was anticipated by its prediction errors reported earlier. Figure 4 depicts these results visually: it plots the number of AC-OPF instances solved over time by the various warm-starting methods. The plot clearly demonstrates the benefits of Compact Learning for predicting both primal and dual solutions.