Risk-Averse Model Predictive Control Design for Battery Energy Storage Systems

When batteries supply behind-the-meter services such as arbitrage or peak load management, an optimal controller can be designed to minimize the total electric bill. The limitations of the batteries, such as on voltage or state-of-charge, are represented in the model used to forecast the system’s state dynamics. Control model inaccuracy can lead to an optimistic shortfall, where the achievable schedule will be costlier than the schedule derived using the model. To improve control performance and avoid optimistic shortfall, we develop a novel methodology for high performance, risk-averse battery energy storage controller design. Our method is based on two contributions. First, the application of a more accurate, but non-convex, battery system model is enabled by calculating upper and lower bounds on the globally optimal control solution. Second, the battery model is then modified to consistently underestimate capacity by a statistically selected margin, thereby hedging its control decisions against normal variations in battery system performance. The proposed model predictive controller, developed using this methodology, performs better and is more robust than the state-of-the-art approach, achieving lower bills for energy customers and being less susceptible to optimistic shortfall.


Symbol
Decision variable description p -ac power provided to the BESS p dc -dc power provided to the battery i bat -dc current provided to the battery v bat -battery terminal voltage v s -battery terminal slack voltage v 1 -dynamic battery voltage v oc -open-circuit-voltage ς -state-of-charge ("sigma") τ -peak net electrical load ("tau") Vector, vector equation, and matrix equation notation -If z is a decision variable, z ∈ R n is a column vector of n decision variables. This is normally used to indicate decision variables at discrete times (e.g., p at each time-step in a control horizon becomes p) -The expression z [1:3] produces a column vector with the elements of z indicated by the index(es) (in this case, the first three elements) -The expression z + y produces a column vector that is the element-wise addition of the vectors z and y -For the scalar value b, the expression b [1] denotes the multiplication of the constant (b), times a vector of ones [1], that produces a column vector populated with b -The vector equation z+y = b [1] denotes n equations, each with indexed variables (a.k.a, z [1] + y [1] = b, z [2] + y [2] = b, etc.) -The matrix equation A[z, y] T ≤ b [1] 1×n , where A ∈ R m×2 and b ∈ R m×1 , denotes m × n equations, each with indexed variables (a.k.a, A [1,1] z [1] + A [1,2] y [1] ≤ b [1] , A [2,1] z [1] + A [2,2] y [1] ≤ b [2] , A [1,1] z [2] + A [1,2] y [2] ≤ b [1] , etc.)

Parameters
All parameters can be found in Tables I-V. I. INTRODUCTION B ATTERY Energy Storage Systems (BESS) are becoming an integral part of a resilient and efficient electrical system. Distributed energy resources (DER) such as BESS are able to support the grid through advanced control and functionality [1]. In addition to responding to local conditions of voltage and frequency, energy management systems can forecast future conditions of variables such as price and load to optimally schedule BESS operation. The time-ofuse (ToU) and peak demand charge billing mechanisms are driving the adoption of BESS in commercial applications in many areas [2]. Both TOU and demand charge rate structures incentivize costumers to, in aggregate, reduce system peak demand, allowing a utility to defer or avoid costly capacity upgrades [3]. A primary concern in both applications is how to make control decisions that maximize the value of the BESS to This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ Fig. 1. Time-of-use price schedule (top) [4], and electrical load (bottom) [5]. the owner. Optimal control of BESS is a challenging problem due to their complexity and the uncertainty of the underlying chemical processes. The objective of this paper is to develop and demonstrate a methodology to design an advanced, energy management level controller for BESS. We accomplish this by first reducing model uncertainty, through model improvements, and then by shaping model uncertainty to impose a risk-averse bias on control decisions.
We introduce the following simple case study to provide a basis for an optimal control objective function. Consider a hypothetical commercial electrical customer billed for power under both ToU and demand charges. This customer decides to purchase and install a battery to reduce their electrical bill. The energy contract charges 9¢/kWh during off-peak, 11¢/kWh during partial-peak hours, and 15¢/kWh during peak [4] according to the schedule shown in Fig. 1 (top). The utility then charges $50/kW service fee (d) according to the peak net load measured during the billing period. For simplicity we will assume a weekly billing period for this example, although it would be more typical to use a monthly or longer period to set the demand charge. This price is consistent with demand charges in specific localities in California and New York [6]. The customer's week ahead load forecast is shown in Fig. 1 (bottom). These load data are adapted from the EPRI test circuit 'Ckt5' loadshape, July 20th through 26th, normalized to a 5 kW peak [5]. We will assume that the load and price are known a priori. Without the battery, total bill would be calculated according to (1).
where f bill is the total electric bill, l is the n-length vector of load, c is the n-length vector of ToU prices, d is the service fee in $/kW for peak net load measured during the billing period, and • T denotes a vector's transpose. Unless otherwise stated, we use a time-step t = 15 minutes (0.25 hours), and n = 672 (1 week). With the addition of a BESS that can supply (−), or absorb (+), power p, the customer's cost can be modified to (2).
where p is the battery system power that element wise subtracts from l when the battery system is discharging. The problem is thus formulated: design a control algorithm to optimally calculate a vector of battery system power p that minimizes the customer's bill f bill without exceeding the battery's limits, including those on charge rate and depth-ofdischarge. Further, the quality of the solution will be judged based both on how much it reduces the customer's bill and how robust it is to uncertainty. A common approach to controller design is called model predictive control (MPC). MPC is a real-time, state-feedback, optimal control approach that involves solving a finite-horizon online optimization problem at each time step that results in a sequence of future control actions as well as predictions of the future states [7]- [9]. In designing MPC, the choice of what model to use can be critical. The simplest BESS model assumes that changes in SoC are proportional to the energy charged or discharged from ac point of interconnection. This approach to optimal control represents the state-of-the-art and has been used for improving wind farm dispatch in Australian electricity markets [10], and achieving distribution feeder dispatchability [11]. Another common approach, based on the need for improved accuracy, is to use a BESS model that assumes that changes in SoC are proportional to the charge, in amp-hours, supplied or absorbed by the battery itself. While it has been used [12]- [14], this approach can be difficult because the feasible subspace it defines is fundamentally nonconvex. Historically, the only way to apply the more accurate model to calculate optimal control schedules was to either approximate the model using pseudospectral methods [12], or to use dynamic programming [13], [14]. Further, given the precision of this type of model, the performance of optimal controllers that rely on it can be sensitive to variations in battery performance. This paper makes two fundamental contributions to the state-of-the-art: 1) formulation of an optimal controller for a residential lithium-ion battery system, based on the more accurate charge reservoir model (Section III) with upper and lower bounds to check the viability of solutions found through gradient based methods (Section V), and 2) a method to modify the controller to be risk-averse to variations in battery performance (Section VI). In Section VII we demonstrate the improved controller performance, due to using the more accurate model, and demonstrate how the risk-averse modification makes it more robust to model uncertainty. Together, these contributions make up an advanced methodology for designing BESS controllers that perform better and are more robust than those designed through traditional methods. Section VIII concludes the paper with a summary of the results and a discussion of the broad applicability of the proposed control design approach.

II. ENERGY RESERVOIR MODEL
Energy Reservoir Model (ERM) refers to a class of models that calculates SoC as a linear function of energy into and out of the BESS. The ERM is widely used in battery energy storage control problems [10], [11], [15]- [17] and has the advantage of being linear in charge and discharge power. This allows for convex, and therefore computationally efficient, formulations of the optimal control problem. The ERM formulation used here is shown in (3). Definitions for parameters are given in Table II.
An 2 norm regularization is applied in (3a) and scaled by the constant 1 to even out peak battery power when it is not needed. The constraint (3b) ensures that control decisions are made based on the current estimated SoC (ς 0 ). The constraint (3c) represents the intuitive assumption that the BESS will continue to operate after the end of the current control horizon and that the next period will be similar to this one. While not necessary in closed-loop implementation, (3c) makes simulation results easier to interpret and compare.

III. CHARGE RESERVOIR MODEL
Charge Reservoir Model (CRM) refers to a class of models that calculates SoC as a function of charge (current integrated over time) into and out of the battery itself. The CRM is also used in battery energy storage control problems and has the advantage of being more accurate over a longer time horizon or over larger range of SoC [12]. The disadvantage of the CRM is that the subspace of feasible solutions is fundamentally nonconvex. Therefore, it is more computationally complex and difficult integrate into an on-board controller [13]. The CRM formulation used here is shown in (5). The parameters for this model are listed in Table III.
, p dc ∈ R n is the dc electrical power provided to battery, v bat ∈ R n is the battery terminal voltage, v s ∈ R n + is the slack voltage used in calculation of an upper bound, v oc ∈ R n+1 is the battery open-circuit voltage, and τ ∈ R is the peak power demand. The CRM objective includes a 2 norm power regularization and an 1 norm slack voltage cost, weighted by the constant 2 . The weight 2 is chosen, using a simple trial and error sweep, to be the smallest value that is still large enough to drive the slack voltage to zero under normal operation. The CRM includes constraints on inverter conversion efficiency (5b), Ohm's law relating dc power, voltage and current (5c), the battery equivalent circuit model (5d), and the open-circuit voltage curve (5e). Note that the inverter conversion efficiency (5b) is a convex inequality that collapses to equality as long as energy prices in the objective are positive. This model uses the big cell method discussed in [18]. These additional parameters and constrains allow the CRM to more accurately represent the physical dynamics of battery systems.

IV. EXTENDED CRM FOR SIMULATION
To perform a pseudo-empirical analysis of the optimal schedules calculated from each model we simulate how the battery system would respond to each control signal using an extended CRM that incorporates additional constraints and parameters to improve its accuracy. The simulation model uses  slightly different functions and parameters, enabling an analysis of the effects of model and parameter uncertainty on controller performance. The modified constraints are shown in (6). The parameters for these modified constraints are shown in Table IV.  with a piecewise cubic-spline fit, as has been shown to be highly accurate [19].
The simulation timestep is 1 second, meaning that it is executed 900 times between controller time steps (with t = 15 minutes). The extended CRM is implemented in simulation using the Battery-Inverter fleet model discussed in [20]. The resulting schedules are distinguished by the tags 'calculated', which stands for the optimal schedules calculated using the ERM or CRM, and 'achieved', which stands for the results of simulating the calculated schedule using the extended CRM. The discrepancy between 'calculated' and 'achieved' schedules is a result of inaccurate parameters and unrepresented battery system characteristics in the ERM and CRM models.

V. BOUNDING THE GLOBAL MINIMUM
The nonlinear CRM optimization problem shown in (5) is non-convex. Further, it can be shown that the Lagrangian of this problem is not pseudoconvex as defined in [21]. If it had either of these properties then we would know that any minimum found would be in the set of global minima but as it is, we cannot make this guarantee. Because of this some argue that gradient based methods such as Newton-Raphson are not viable for CRM optimization due to local minima in the solution space [14]; however, we find this not to be the case. Our contribution to the state-of-the-art is to bound the global minimum of this problem such that if we find a local minimum inside this range, we can be confident that it is, or is close to, the globally optimal solution. An upper bound to a minimization problem can be found by restricting the feasible set (adding additional constraints) while a lower bound can be calculated by expanding the feasible set (relaxing or removing constraints) [22].
To calculate a convex lower bound we relax the non-convex constraints to their convex hulls. First, the constraint (5b) is modified to the include positive and negative dc power (7a). We then relax the ohm's power law constraint (5c) to a convex space bounded by eight affine surfaces as shown in Fig. 3 and represented in (7b) and (7c). To do this while maintaining feasibility we split the dc power into separate positive and negative decision variables. Finally, we relax the open-circuitvoltage constraint (5e) to the convex hull shown in Fig. 2  and represented in (7d). The resulting convex problem in (7) provides a lower bound on the global minimum of (5).
. . . (5d) and (5f) through (5n) unchanged where p + dc and p − dc are the charge and discharge dc powers respectively.
To calculate an upper bound we restrict battery terminal voltage to a constant (v bat = v ocmin ). The slack voltage allows the battery equivalent circuit (8b) to become nonbinding, thereby making the fixed voltage restriction feasible. To efficiently calculate the global minimum of the upper bound we use a piecewise linear approximation of the open circuit-voltage function, as shown in (8c) through (8j). This approximation makes the upper bound a mixed-integer nonlinear program (MINLP) with a convex objective and convex constraint functions, for which there exist effective exact solution algorithms [22]. Given that the approximation embodied in the piecewise linearization is fairly accurate, we expect that the resulting solution will be a useful upper bound to the original, non-convex problem.
With these bounds and properties established, we can use an off-the-shelf primal-dual, interior-point method to solve the optimal control problem using the CRM. The freely downloadable modeling language Pyomo [23], [24] and nonlinear solver Ipopt [25] are used to implement this algorithm efficiently and the code and data we used are available as supplemental material attached to the digital version this article. The nonlinear solver Gurobi is used for calculation of the upper bound as it is able to efficiently work with integer variables [26].

VI. REDUCING CONTROL SENSITIVITY TO UNCERTAINTY
Parameters such as capacity are functions of many physical mechanisms we do not consider in the model. To consider this uncertainty, we break capacity 1 into it's mean value and a random component as shown in (9).
whereC cap ∼ N (μ = 0, σ = 2.6Ah) is the random component of the battery's capacity, assumed to be a zero-mean, normal distribution. Price arbitrage has symmetric risk, meaning that overestimating capacity is just as bad as underestimating capacity. Peak demand charge management, in contrast, has asymmetric risk, in that the down-side of overestimating capacity is larger in magnitude than the down-side of underestimating capacity. Because of asymmetric risk, we expect a risk-neutral controller to have a skewed performance distribution as shown in Fig. 4. We modify the proposed controller to consider asymmetric risk. The value-at-risk constraint needed for this modification is shown in (10). (10) whereĈ cap is the value-at-risk capacity, and P is the probability function. Because of our assumption that the capacity has a normal distribution, calculating the value-at-risk is trivial in that we can use lookup tables to determine how many standard deviations from the mean will yield a risk of 0.13% (−3σ from the table supplied in [27]). By usingĈ cap = 135.2 Ah−3×2.6 Ah = 127.4 Ah, we design the controller to consistently underestimate the battery's capacity, thereby making its control decisions robust to fluctuations in capacity. Due to this modification, we expect that a risk-averse controller will have performance distribution with reduced downside risk, in exchange for slightly reduced average performance, as illustrated in Fig. 4.
As described above some services have asymmetric risk of overestimating or underestimating capacity. The degree of this asymmetry corresponds to the potential advantage of riskaverse control. For example: backup power for a critical load has an extreme imbalance in cost between not having sufficient energy and having more than enough energy. Because of this, the intuitive risk-averse control solution is to maintain full charge at all times.
Additionally, the accuracy and precision of the BESS model is important to consider as illustrated in Fig. 5. Model accuracy refers to how close the mean estimated capacity parameter 1 The uncertainty of coulombic efficiency and internal resistance parameters are also listed in Table III. The effects of these uncertainties is negligible for the example application when compared to the capacity. Relative benefits of risk-averse control based on model accuracy/precision, given an asymmetric risk application.
is to the capacity expressed during the control horizon (in our case, the extended CRM used in simulation). Precision is a measure of how consistent the BESS capacity is in this application. If the model consistently underestimates or overestimates available energy, low accuracy with high precision, then the marginal benefit of risk-averse control will be negligible. If the model has high accuracy and precision then there is only a small margin for improvement. For controllers with low accuracy and precision, performance consistency comes at a large cost in average performance. It is the case where the BESS model has high accuracy, and relatively low precision where risk-averse control is most useful because it is able to hedge decisions for uncertainties in performance. Such controllers can achieve consistency with only a small sacrifice in average performance.
To assess the sensitivity of the proposed controller to off normal circumstances we use the optimistic shortfall. Optimistic shortfall is the difference between expected controller performance and achieved controller performance, which in our case is the total bill achieved minus the calculated optimal bill. In Fig. 4, the risk-neutral controller would have significant optimistic shortfall in the 'extreme case' (C cap = −3σ ), whereas the risk-averse controller would have low optimistic shortfall throughout the normal range of battery capacity.

VII. RESULTS
This section explains the results of simulated control action calculated using the ERM, CRM, and Risk-Averse (RA) CRM as discussed above, applied to the extended CRM model of the BESS as a pseudo-empirical analysis of controller performance and sensitivity to parameter uncertainty. Openloop control is assessed first to provide a baseline controller performance and a clear picture of the effect that uncertainty has on control. The closed-loop, model predictive controllers are then assessed for their ability to reduce the effects of uncertainty under normal operations and reduce optimistic shortfall. The effects of extreme case parameters are then analyzed to illustrate that risk-averse control design is needed to make the controller robust to the normal variations in system performance. Last, the performance of the proposed risk-averse controller is analyzed. A summary of the customer bill, % savings, and optimistic shortfall under each simulation scenario is shown in Table VI.

A. Open-Loop Control
The optimal 'calculated' schedules, along with the 'achieved' schedules, for the customer using the ERM and CRM in open-loop are shown in Fig. 6. The resulting net load curves for their control schedules are shown in Fig. 7 and 8. While the ERM is clearly more computationally efficient than the CRM, optimal schedules can be calculated using either model in just a few seconds on a mid-range laptop (hardware used for this study: i7-7600U CPU at 2.8 GHz) meaning either approach could be used for on-board control.
For the customer introduced in Section I, the baseline cost of electrical service is $311 ($61 energy, $250 peak demand). The schedule calculated using the ERM reduces this by 11.6% to $275. The schedule calculated using the CRM reduces the cost of electrical service by 13.3% to $270. These bill reductions come primarily from the BESS reducing the peak electrical load by 14.38% (ERM) and 16.0% (CRM) respectively. As this bill falls between the calculated upper and lower  bounds on the global minimum, we are confident that the minimum calculated is, or is close to, the global. While a $5 improvement in savings over the ERM does not sound significant in absolute terms, it is important to remember the scale of power systems. With approximately 5 million commercial customers in the U.S. currently eligible for tariffs with a demand charge rate of at least $15/kW [6], a 14.7% improvement in cost savings, over the ERM, from a simple change in software would have a significant impact. Note that a IEEE 1547-2018 compliant inverter would be able support local grid voltage with volt-var [1] while applying this control schedule.
While the CRM more accurately forecasts SoC, in this case, improved accuracy makes the CRM more susceptible than the ERM to overestimating future SoC and hence not being able to supply sufficient energy during the critical peak. This phenomena is illustrated in Fig. 8 where the achieved net load, derived by simulating the extended CRM using the calculated power schedule, has a peak significantly higher than the calculated net load. The gap between calculated and achieved net load schedules comes from the BESS being unable to supply sufficient energy to shave the complete peak, needing to curtail its discharge prematurely. This gap creates a large optimistic shortfall, where the achieved bill is $5.43 higher than the calculated bill. We demonstrate in the next section that this optimistic shortfall can be mostly eliminated with closed loop control.

B. Closed-Loop Control
Closed-loop control recalculates the optimal schedule at each time step. The net loads achieved by both closed-loop ERM and CRM based controllers are shown in Fig. 9. When implemented on the ERM, closed-loop control generates a small negative optimistic shortfall (optimistic surplus). This is because, as it starts to shave the peak load at a level based   Fig. 4). Probability refers to the likelihood that the capacity is less than the signified value. When total bills resulting from different capacity values are very close (within 25¢), their probabilities are added, and the total bill is averaged. on its underestimation of capacity, the SoC is updated and the controller has more energy to work with than expected. It then applies this excess energy to a discharge during the window of peak ToU price that is coincident with peak load. The CRM based model predictive controller reduces the optimistic shortfall from $5.43 (open-loop, see Fig. 8) to $0.00 (closed-loop). This is a result of the open-loop controller not supplying sufficient charge to reach ς max before the beginning of the peak. The closed-loop controller is able to adjust for the insufficient charge and have enough energy to shave the peak completely.

C. Risk-Averse Closed-Loop Control
The physical parameters of a BESS vary under normal operation and these variations can have a large impact on the optimistic shortfall of a controller. Fig. 10 shows the sensitivity of the total bill achieved by the CRM and RA CRM due to variations in capacity. When the battery's capacity is at its mean value and above (μ, +1σ , +2σ , and +3σ ), the riskneutral CRM has a slight performance advantage. However, when the battery's capacity is below expectations (−1σ , −2σ , and −3σ ), the risk-neutral CRM's performance drops off, producing an optimistic shortfall up to $22.98 in the 'extreme case' atC cap = −3σ , while the performance of the risk-averse controller does not decline. In terms of peak net load reduction, the risk-neutral controller is only 50% confident it will reduce the peak by 16.0%, but it has a roughly one in six chance it will reduce the peak less than 14%. In contrast, the risk-averse controller is 99.87% confidante that it can reduce the peak by 15.4%. This achieves the goal of making the controller more robust to battery model uncertainty.

VIII. CONCLUSION
In this paper we develop and demonstrate an advanced methodology for designing BESS controllers under ToU price arbitrage and peak demand charge management applications. A state-of-the-art ERM is used as the baseline for control performance comparison. The proposed CRM based model predictive controller outperforms the ERM based controller by achieving a lower total electric bill when pseudo-empirically applied in an example scenario. Because peak load management has asymmetric risk for overestimating available energy, we then shape the uncertainty of the CRM to consistently underestimate capacity. This risk-averse CRM yields better controller performance than the ERM and is more robust to variations in BESS performance than the CRM. This methodology for designing BESS controllers can be applied in a broad range of energy storage applications, wherever the risk profile of a scheduled service is asymmetric. Incremental improvements in controller performance can reduce the cost of deploying storage to make the grid more efficient and resilient.