Enforcing robust control guarantees within neural network policies

Priya L. Donti, Melrose Roderick, Mahyar Fazlyab, J. Zico Kolter

Introduction

The field of robust control, dating back many decades, has been able to provide rigorous guarantees on when controllers will succeed or fail in controlling a system of interest. In particular, if the uncertainties in the underlying dynamics can be bounded in specific ways, these techniques can produce controllers that are provably robust even under worst-case conditions. However, as the resulting policies tend to be simple (i.e., often linear), this can limit their performance in typical (rather than worst-case) scenarios. In contrast, recent high-profile advances in deep reinforcement learning have yielded state-of-the-art performance on many control tasks, due to their ability to capture complex, nonlinear policies. However, due to a lack of robustness guarantees, these techniques have still found limited application in safety-critical domains where an incorrect action (either during training or at runtime) can substantially impact the controlled system.

In this paper, we propose a method that combines the guarantees of robust control with the flexibility of deep reinforcement learning (RL). Specifically, we consider the setting of nonlinear, time-varying systems with unknown dynamics, but where (as common in robust control) the uncertainty on these dynamics can be bounded in ways amenable to obtaining provable performance guarantees. Building upon specifications provided by traditional robust control methods in these settings, we construct a new class of nonlinear policies that are parameterized by neural networks, but that are nonetheless provably robust. In particular, we project the outputs of a nominal (deep neural network-based) controller onto a space of stabilizing actions characterized by the robust control specifications. The resulting nonlinear control policies are trainable using standard approaches in deep RL, yet are guaranteed to be stable under the same worst-case conditions as the original robust controller.

We describe our proposed deep nonlinear control policy class and derive efficient, differentiable projections for this class under various models of system uncertainty common in robust control. We demonstrate our approach on several different domains, including synthetic linear differential inclusion (LDI) settings, the cart-pole task, a quadrotor domain, and a microgrid domain. Although these domains are simple by modern RL standards, we show that purely RL-based methods often produce unstable policies in the presence of system disturbances, both during and after training. In contrast, we show that our method remains stable even when worst-case disturbances are present, while improving upon the performance of traditional robust control methods.

Related work

We employ techniques from robust control, (deep) RL, and differentiable optimization to learn provably robust nonlinear controllers. We discuss these areas of work in connection to our approach.

Robust control. Robust control is concerned with the design of feedback controllers for dynamical systems with modeling uncertainties and/or external disturbances (Zhou and Doyle, 1998; Başar and Bernhard, 2008), specifically controllers with guaranteed performance under worst-case conditions. Many classes of robust control problems in both the time and frequency domains can be formulated using linear matrix inequalities (LMIs) (Boyd et al., 1994; Kothare et al., 1996); for reasonably-sized problems, these LMIs can be solved using off-the-shelf numerical solvers based on interior-point or first-order (gradient-based) methods. However, providing stability guarantees often requires the use of simple (linear) controllers, which greatly limits average-case performance. Our work seeks to improve performance via nonlinear controllers that nonetheless retain the same stability guarantees.

Reinforcement learning (RL). In contrast, RL (and specifically, deep RL) is not restricted to simple controllers or problems with uncertainty bounds on the dynamics. Instead, deep RL seeks to learn an optimal control policy, represented by a neural network, by directly interacting with an unknown environment. These methods have shown impressive results in a variety of complex control tasks (e.g., Mnih et al. (2015); Akkaya et al. (2019)); see Buşoniu et al. (2018) for a survey. However, due to its lack of safety guarantees, deep RL has been predominantly applied to simulated environments or highly-controlled real-world problems, where system failures are either not costly or not possible.

Efforts to address the lack of safety and stability in RL fall into several main categories. The first tries to combine control-theoretic ideas, predominantly robust control, with the nonlinear control policy benefits of RL (e.g., Morimoto and Doya (2005); Abu-Khalaf et al. (2006); Feng et al. (2009); Liu et al. (2013); Wu and Luo (2013); Luo et al. (2014); Friedrich and Buss (2017); Pinto et al. (2017); Jin and Lavaei (2018); Chang et al. (2019); Han et al. (2019); Zhang et al. (2020)). For example, RL has been used to address stochastic stability in H∞H_{\infty} control synthesis settings by jointly learning Lyapunov functions and policies in these settings (Han et al., 2019). As another example, RL has been used to address H∞H_{\infty} control for continuous-time systems via min-max differential games, in which the controller and disturbance are the “minimizer” and “maximizer” (Morimoto and Doya, 2005). We view our approach as thematically aligned with this previous work, though our method is able to capture not only H∞H_{\infty} settings, but also a much broader class of robust control settings.

Another category of methods addressing this challenge is safe RL, which aims to learn control policies while maintaining some notion of safety during or after learning. Typically, these methods attempt to restrict the RL algorithm to a safe region of the state space by making strong assumptions about the smoothness of the underlying dynamics, e.g., that the dynamics can be modeled as a Gaussian process (GP) (Turchetta et al., 2016; Akametalu et al., 2014) or are Lipschitz continuous (Berkenkamp et al., 2017; Wachi et al., 2018). This framework is in theory more general than our approach, which requires using stringent uncertainty bounds (e.g. state-control norm bounds) from robust control. However, there are two key benefits to our approach. First, norm bounds or polytopic uncertainty can accommodate sharp discontinuities in the continuous-time dynamics. Second, convex projections (as used in our method) scale polynomially with the state-action size, whereas GPs in particular scale exponentially (and are therefore difficult to extend to high-dimensional problems).

A third category of methods uses Constrained Markov Decision Processes (C-MDPs). These methods seek to maximize a discounted reward while bounding some discounted cost function (Altman, 1999; Achiam et al., 2017; Taleghan and Dietterich, 2018; Yang et al., 2020). While these methods do not require knowledge of the cost functions a-priori, they only guarantee the cost constraints hold during test time. Additionally, using C-MDPs can yield other complications, such as optimal policies being stochastic and the constraints only holding for a subset of states.

Differentiable optimization layers. A great deal of recent work has studied differentiable optimization layers for neural networks: e.g., layers for quadratic programming (Amos and Kolter, 2017), SAT solving (Wang et al., 2019), submodular optimization (Djolonga and Krause, 2017; Tschiatschek et al., 2018), cone programs (Agrawal et al., 2019), and other classes of optimization problems (Gould et al., 2019). These layers can be used to construct neural networks with useful inductive bias for particular domains or to enforce that networks obey hard constraints dictated by the settings in which they are used. We create fast, custom differentiable optimization layers for the latter purpose, namely, to project neural network outputs into a set of certifiably stabilizing actions.

Background on LQR and robust control specifications

In this paper, our aim is to control nonlinear (continuous-time) dynamical systems of the form

for some design parameter α>0\alpha>0. (This particular condition implies exponential stability with a rate of convergence α\alpha. See, e.g., Haddad and Chellaboina (2011) for a more rigorous definition of (local and global) exponential stability. Condition (2) comes from Lyapunov’s Theorem, which characterizes various notions of stability using Lyapunov functions.) For certain classes of bounded dynamical systems, time-invariant linear control policies u(t)=Kx(t)u(t)=Kx(t), and quadratic Lyapunov functions V(x)=xTPxV(x)=x^{T}Px, it is possible to construct such guarantees using semidefinite programming. For instance, consider the class of norm-bounded LDIs (NLDIs)

2 LQR control objectives

In addition to designing for stability, it is often desirable to optimize some objective characterizing controller performance. While our method can optimize performance with respect to any arbitrary cost or reward function, to make comparisons with existing methods easier, for this paper we consider the well-known infinite-horizon “linear-quadratic regulator” (LQR) cost, defined as

Enforcing robust control guarantees within neural networks

We now present the main contribution of our paper: A class of nonlinear control policies, potentially parameterized by deep neural networks, that is guaranteed to obey the same stability conditions enforced by the robustness specifications described above. The key insight of our approach is as follows: While it is difficult to derive specifications that globally characterize the stability of a generic nonlinear controller, if we are given known robustness specifications, we can create a sufficient condition for stability by simply enforcing that our policy satisfies these specifications at all tt. For instance, given a known Lyapunov function, we can enforce exponential stability by ensuring that our policy sufficiently decreases this function (e.g., satisfies Equation (2)) at any given x(t)x(t).

In the following sections, we present our nonlinear policy class, as well as our general framework for learning provably robust policies using this policy class. We then derive the instantiation of this framework for various settings of interest. In particular, this involves constructing (custom) differentiable projections that can be used to adjust the output of a nominal neural network to satisfy desired robustness criteria. For simplicity of notation, we will often suppress the tt-dependence of xx, uu, and ww, but we note that these are continuous-time quantities as before.

Given a dynamical system of the form (1) and a quadratic Lyapunov function V(x)=xTPxV(x)=x^{T}Px, let

We note that this policy class is differentiable if the projections can be implemented in a differentiable manner (e.g., using convex optimization layers (Agrawal et al., 2019), though we construct efficient custom solvers for our purposes). Importantly, as all policies in this class satisfy the stability condition (2) for all states xx and at all times tt, these policies are certifiably robust under the same conditions as the original (linear) controller for which the Lyapunov function V(x)V(x) was constructed.

Since πθ\pi_{\theta} is differentiable, we can solve this problem via a variety of approaches, e.g., a model-based planning algorithm if the true dynamics are known, or virtually any (deep) RL algorithm if the dynamics are unknown.While this problem is infinite-horizon and continuous in time, in practice, one would optimize it in discrete time over a large finite time horizon.

This general procedure for constructing stabilizing controllers is summarized in Algorithm 1. While seemingly simple, this formulation presents a powerful paradigm: by simply transforming the output of a neural network, we can employ an expressive policy class to optimize an objective of interest while ensuring the resultant policy will stabilize the system during both training and testing.

We instantiate our framework by constructing “safe” sets C(x)\mathcal{C}(x) and their associated (differentiable) projections PC(x)\mathcal{P}_{\mathcal{C}(x)} for three settings of interest: NLDIs, polytopic linear differential inclusions (PLDIs), and H∞H_{\infty} control settings. As an example, we describe this procedure below for NLDIs, and refer readers to Appendix B for corresponding formulations for the additional settings we consider.

2 Example: NLDIs

In order to apply our framework to the NLDI setting (3), we first compute a quadratic Lyapunov function V(x)=xTPxV(x)=x^{T}Px by solving the optimization problem (6) for the given system via semidefinite programming. We then use the resultant Lyapunov function to compute the system-specific “safe” set C(x)\mathcal{C}(x), and then create a fast, custom differentiable solver to project onto this set.

Consider the NLDI system (3), some stability parameter α>0\alpha>0, and a Lyapunov function V(x)=xTPxV(x)=x^{T}Px with PP satisfying Equation (4). Assuming PP exists, define

We seek to find a set of actions such that the condition (2) is satisfied along all possible trajectories of (3). A set of actions satisfying this condition at a given xx is given by

Let S:={w:∥w∥2≤∥Cx+Du∥2}\mathcal{S}:=\{{w:\|w\|_{2}\leq\|Cx+Du\|_{2}}\}. We can then rewrite the left side of the above inequality as

by the definition of the NLDI dynamics and the closed-form minimization of a linear term over an L2L_{2} ball. Rearranging yields an inequality of the desired form. We note that by definition of the specifications (4), there is some KK corresponding to PP such that the policy u=Kxu=Kx satisfies the exponential stability condition (2); thus, Kx∈CNLDIKx\in\mathcal{C}_{\text{NLDI}}, and CNLDI\mathcal{C}_{\text{NLDI}} is non-empty. Further, as the above inequality represents a second-order cone constraint in uu, this set is convex in uu. ∎

We further consider the special case where D=0D=0, i.e., the norm bound on ww does not depend on the control action. This form of NLDI arises in many common settings (e.g., where ww characterizes linearization error in a nonlinear system but the dynamics depend only linearly on the action), and is one for which we can compute the relevant projection in closed form (as described shortly).

Consider the NLDI system (3) with D=0D=0, some stability parameter α>0\alpha>0, and Lyapunov function V(x)=xTPxV(x)=x^{T}Px with PP satisfying Equation (4). Assuming PP exists, define

The result follows by setting D=0D=0 in Theorem 1 and rearranging terms. As the above inequality represents a linear constraint in uu, this set is convex in uu. ∎

2.2 Deriving efficient, differentiable projections

For the general NLDI setting (3), we note that the relevant projection PCNLDI(x)\mathcal{P}_{\mathcal{C}_{\text{NLDI}}(x)} (see Theorem 1) represents a projection onto a second-order cone constraint. As this projection does not necessarily have a closed form, we must implement it using a differentiable optimization solver (e.g., Agrawal et al. (2019)). For computational efficiency purposes, we implement a custom solver that employs an accelerated projected dual gradient method for the forward pass, and employs implicit differentiation through the fixed point equations of this solution method to compute relevant gradients for the backward pass. Derivations and additional details are provided in Appendix C.

In the case where D=0D=0 (see Corollary 1.1), we note that the projection operation PCNLDI-0(x)\mathcal{P}_{\mathcal{C}_{\text{NLDI-0}}(x)} does have a closed form, and can in fact be implemented via a single ReLU operation. Specifically, defining ηT:=2xTPB\eta^{T}:=2x^{T}PB and ζ:=−xT(2PA+αP)x−2∥GTPx∥2∥Cx∥2\zeta:=-x^{T}(2PA+\alpha P)x-2\|G^{T}Px\|_{2}\|Cx\|_{2}, we see that

Experiments

Having instantiated our general framework, we demonstrate the power of our approach on a variety of simulated control domains.Code for all experiments is available at https://github.com/locuslab/robust-nn-control In particular, we evaluate performance on the following metrics:

Average-case performance: How well does the method optimize the performance objective (i.e., LQR cost) under average (non-worst case) dynamics?

Worst-case stability: Does the method remain stable even when subjected to adversarial (worst-case) dynamics?

In all cases, we show that our method is able to improve performance over traditional robust controllers under average conditions, while still guaranteeing stability under worst-case conditions.

We evaluate our approach on five NLDI settings: two synthetic NLDI domains, the cart-pole task, a quadrotor domain, and a microgrid domain. (Additional experiments for PLDI and H∞H_{\infty} control settings are described in Appendix I.) For each setting, we choose a time discretization based on the speed at which the system evolves, and run each episode for 200 steps over this discretization. In all cases except the microgrid setting, we use a randomly generated LQR objective where the matrices Q1/2Q^{1/2} and R1/2R^{1/2} are drawn i.i.d. from a standard normal distribution.

Synthetic NLDI settings. We generate NLDIs of the form (3) with s=5s=5, a=3a=3, and d=k=2d=k=2 by generating matrices A,B,G,CA,B,G,C and DD i.i.d. from normal distributions, and producing the disturbance w(t)w(t) using a randomly-initialized neural network (with its output scaled to satisfy the norm-bound on the disturbance). We investigate settings both where D≠0D\neq 0 and where D=0D=0. In both cases, episodes are run for 2 seconds at a discretization of 0.01 seconds.

Cart-pole. In the cart-pole task, our goal is to balance an inverted pendulum resting on top of a cart by exerting horizontal forces on the cart. For our experiments, we linearize this system as an NLDI with D≠0D\neq 0 (see Appendix D), and add a small additional randomized disturbance satisfying the NLDI bounds. Episodes are run for 10 seconds at a discretization of 0.05 seconds.

Planar quadrotor. In this setting, our goal is to stabilize a quadcopter in the two-dimensional plane by controlling the amount of force provided by the quadcopter’s right and left thrusters. We linearize this system as an NLDI with D=0D=0 (see Appendix E), and add a small disturbance as in the cart-pole setting. Episodes are run for 4 seconds at a discretization of 0.02 seconds.

Microgrid. In this final setting, we aim to stabilize a microgrid by controlling a storage device and a solar inverter. We augment the system given in Lam et al. (2016) with LQR matrices and NLDI bounds (see Appendix F). Episodes are run for 2 seconds at a discretization of 0.01 seconds.

2 Experimental setup

Robust MBP (ours): A model-based planner that assumes the true dynamics are known.

Robust PPO (ours): An RL approach based on PPO (Schulman et al., 2017) that does not assume known dynamics (beyond the bounds used to construct the robust policy class).

Robust MBP is optimized using gradient descent for 1,000 updates, where each update samples 20 roll-outs. Robust PPO is trained for 50,000 updates, where each update samples 8 roll-outs; we choose the model that performs best on a hold-out set of initial conditions during training. We note that while we use PPO for our demonstration, our approach is agnostic to the particular method of training, and can be deployed with many different (deep) RL paradigms.

We compare our robust neural network-based method against the following baselines:

Robust LQR: Robust (linear) LQR controller obtained via Equation (6).

Robust MPC: A robust model-predictive control algorithm (Kothare et al., 1996) based on state-dependent LMIs. (As the relevant LMIs are not always guaranteed to solve, our implementation temporarily reverts to the Robust LQR policy when that occurs.)

RARL: The robust adversarial reinforcement learning algorithm (Pinto et al., 2017), which trains an RL agent in the presence of an adversary. (We note that unlike the other robust methods considered here, this method is not provably robust.)

LQR: A standard non-robust (linear) LQR controller.

MBP and PPO: The non-robust neural network policy class π^θ(x)\hat{\pi}_{\theta}(x) optimized via a model-based planner and the PPO algorithm, respectively.

In order to evaluate performance, we train all methods on the dynamical settings described in Section 5.1, and evaluate them on two different variations of the dynamics:

Original dynamics: The dynamical settings described above (“average case”).

Adversarial dynamics: Modified dynamics with an adversarial test-time disturbance w(t)w(t) generated to maximize loss (“worst case”). We generate this disturbance separately for each method described above (see Appendix G for more details).

Initialization states are randomly generated for all experiments. For the synthetic NLDI and microgrid settings, these are generated from a standard normal distribution. For both cart-pole and quadrotor, because our NLDI bounds model linearization error, we must generate initial points within a region where this linearization holds. In particular, the linearization bounds only hold for a specified L∞L_{\infty} ball, BNLDIB_{\text{NLDI}}, around the equilibrium. We use a simple heuristic to construct this ball and jointly find a smaller L∞L_{\infty} ball, BinitB_{\text{init}}, such that there exists a level set LL of the Robust LQR Lyapunov function with Binit⊆L⊆BNLDIB_{\text{init}}\subseteq L\subseteq B_{\text{NLDI}} (details in Appendix H). Since Robust LQR (and by extension our methods) are guaranteed to decrease the relevant Lyapunov function, this guarantees that these methods will never leave BNLDIB_{\text{NLDI}} when initialized starting from any point inside BinitB_{\text{init}} – i.e., that our NLDI bounds will always hold throughout the trajectories produced by these methods.

3 Results

Table 1 shows the performance of the above methods. We report the integral of the quadratic loss over the prescribed time horizon on a test set of states, or indicate cases where the relevant method became unstable (i.e., the loss became orders of magnitude larger than for other approaches). (Sample trajectories for these methods are also provided in Appendix H.)

These results illustrate the basic advantage of our approach. In particular, both our Robust MBP and Robust PPO methods show improved “average-case” performance over the other provably robust methods (namely, Robust LQR and Robust MPC). As expected, however, the non-robust LQR, MBP, and PPO methods often perform better within the original nominal dynamics, as they are optimizing for expected performance but do not need to consider robustness under worst-case scenarios.

However, when we apply allowable adversarial perturbations (that still respect our disturbance bounds), the non-robust LQR, MBP, and PPO approaches diverge or perform very poorly. Similarly, the RARL agent performs well under the original dynamics, but diverges under adversarial perturbations in the generic NLDI settings. In contrast, both of our provably robust approaches (as well as Robust LQR) remain stable under even “worst-case” adversarial dynamics. (We note that the baseline Robust MPC method goes unstable in one instance, though this is due to numerical instability issues, rather than issues with theoretical guarantees.)

Figure 1 additionally shows the performance of all neural network-based methods on the test set over training epochs. While the robust and non-robust MBP and PPO approaches both converge quickly to their final performance levels, both non-robust versions become unstable under the adversarial dynamics very early in the process. The RARL method also frequently destabilizes during training. Our Robust MBP and PPO policies, on the other hand, remain stable throughout the entire optimization process, i.e., do not destabilize during either training or testing. Overall, these results show that our method is able to learn policies that are more expressive than traditional robust methods, while guaranteeing these policies will be stable under the same conditions as Robust LQR.

Conclusion

In this paper, we have presented a class of nonlinear control policies that combines the expressiveness of neural networks with the provable stability guarantees of traditional robust control. This policy class entails projecting the output of a neural network onto a set of stabilizing actions, parameterized via robustness specifications from the robust control literature, and can be optimized using a model-based planning algorithm if the dynamics are known or virtually any RL algorithm if the dynamics are unknown. We instantiate our general framework for dynamical systems characterized by several classes of linear differential inclusions that capture many common robust control settings. In particular, this entails deriving efficient, differentiable projections for each setting, via implicit differentiation techniques. We show over a variety of simulated domains that our method improves upon traditional robust LQR techniques while, unlike non-robust LQR and neural network methods, remaining stable even under worst-case allowable perturbations of the underlying dynamics.

We believe that our approach highlights the possible connections between traditional control methods and (deep) RL methods. Specifically, by enforcing more structure in the classes of deep networks we consider, it is possible to produce networks that provably satisfy many of the constraints that have typically been thought of as outside the realm of RL. We hope that this work paves the way for future approaches that can combine more structured uncertainty or robustness guarantees with RL, in order to improve performance in settings traditionally dominated by classical robust control.

This work was supported by the Department of Energy Computational Science Graduate Fellowship (DE-FG02-97ER25308), the Center for Climate and Energy Decision Making through a cooperative agreement between the National Science Foundation and Carnegie Mellon University (SES-00949710), the Computational Sustainability Network, and the Bosch Center for AI. This material is based upon work supported by the National Science Foundation Graduate Research Fellowship Program under Grant No. DGE1745016. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the author(s) and do not necessarily reflect the views of the National Science Foundation.

We thank Vaishnavh Nagarajan, Filipe de Avila Belbute Peres, Anit Sahu, Asher Trockman, Eric Wong, and anonymous reviewers for their feedback on this work.

References

Appendix A Details on robust control specifications

As described in Section 3.1, for many dynamical systems of the form (1), it is possible to specify a set of linear, time-invariant policies guaranteeing infinite-horizon exponential stability via a set of LMIs. Here, we derive the LMI (4) provided in the main text for the NLDI system (3), and additionally describe relevant LMI systems for systems characterized by polytopic linear differential inclusions (PLDIs) and for H∞H_{\infty} control settings.

Consider the general NLDI system (3). We seek to design a time-invariant control policy u(t)=Kx(t)u(t)=Kx(t) and a quadratic Lyapunov function V(x)=xTPxV(x)=x^{T}Px with P≻0P\succ 0 for this system that satisfy the exponential stability criterion V˙(x)≤−αV(x),  ∀t\dot{V}(x)\leq-\alpha V(x),\;\forall t. We derive an LMI characterizing such a controller and Lyapunov function, closely following and expanding upon the derivation provided in Boyd et al. .

Specifically, consider the NLDI system (3), reproduced below:

The time derivative of this Lyapunov function along the trajectories of the closed-loop system is

The exponential stability condition V˙(x)≤−αV(x)\dot{V}(x)\leq-\alpha V(x) is thus implied by inequality

Additionally, the norm bound on ww can be equivalently expressed as

Using the S-procedure, it follows that for some λ≥0\lambda\geq 0, the following matrix inequality is a sufficient condition for exponential stability:

Using Schur Complements, this matrix inequality is equivalent to

Left- and right-multiplying both sides by P−1P^{-1}, and making the change of variables S=P−1S=P^{-1}, Y=KSY=KS, and μ=1/λ,\mu=1/\lambda, we obtain

Using Schur Complements again on this inequality, we obtain our final system of linear matrix inequalities as

where then K=YS−1K=YS^{-1} and P=S−1P=S^{-1}. Note that the first matrix inequality is homogeneous; we can therefore assume μ=1\mu=1 (and therefore, λ=1\lambda=1), without loss of generality.

A.2 Exponential stability in PLDIs

Consider the setting of polytopic linear differential inclusions (PLDIs), where the dynamics are of the form

where then K=YS−1K=YS^{-1} and P=S−1P=S^{-1}. The derivation of this LMI follows similarly to that for exponential stability in NLDIs, and is well-described in Boyd et al. .

Consider the following H∞H_{\infty} control setting with linear time-invariant dynamics

In cases such as these with larger or more unstructured disturbances, it may not be possible to guarantee asymptotic convergence to an equilibrium. In these cases, our goal is to construct a robust controller with bounds on the extent to which disturbances affect some performance output (e.g., LQR cost), as characterized by the L2\mathcal{L}_{2} gain of the disturbance-to-output map. Specifically, we consider the stability requirement that this L2\mathcal{L}_{2} gain be bounded by some parameter γ>0\gamma>0 when disturbances are present, and that the system be exponentially stable in the disturbance-free case. This requirement can be characterized via the condition that for all tt and some σ≥0\sigma\geq 0,

We note that when E(x(t),x˙(t),u(t))≤0\mathcal{E}(x(t),\dot{x}(t),u(t))\leq 0 for all tt, both of our stability criteria are met. To see this, note that integrating both sides of (A.12) from to ∞\infty and ignoring the non-negative terms on the left hand side after integration yields

This is precisely the desired bound on the L2\mathcal{L}_{2} gain of the disturbance-to-output map (see Khalil and Grizzle ). We also note that in the disturbance-free case, substituting w=0w=0 into (A.12) yields

where the last inequality follows from the non-negativity of the LQR cost; this is precisely our condition for exponential stability.

We now seek to design a time-invariant control policy u(t)=Kx(t)u(t)=Kx(t) and quadratic Lyapunov function V(x)=xTPxV(x)=x^{T}Px with P≻0P\succ 0 that satisfies the above condition. In particular, we can write

As in Appendix A.1, we left- and right-multiply both sides by P−1P^{-1}, and make the change of variables S=P−1S=P^{-1}, Y=KSY=KS, and μ=1/σ\mu=1/\sigma to obtain

Using Schur Complements again, we obtain the LMI

where then K=YS−1K=YS^{-1}, P=S−1P=S^{-1}, and σ=1/μ\sigma=1/\mu.

Appendix B Derivation of sets of stabilizing policies and associated projections

We describe the construction of the set of actions C(x)\mathcal{C}(x), defined in Equation (7), for PLDI systems (A.9) and H∞H_{\infty} control settings (A.11). (The relevant formulations for the NLDI system (3) are described in the main text.)

For the general PLDI system (A.9), relevant sets of exponentially stabilizing actions CPLDI{\mathcal{C}_{\text{PLDI}}} are given by the following theorem.

Consider the PLDI system (A.9), some stability parameter α>0\alpha>0, and a Lyapunov function V(x)=xTPxV(x)=x^{T}Px with PP satisfying (A.10). Assuming PP exists, define

We seek to find a set of actions such that the condition (2) is satisfied along all possible trajectories of (A.9), i.e., for any allowable instantiation of (A(t),B(t))(A(t),B(t)). A set of actions satisfying this condition at a given xx is given by

by definition of the PLDI dynamics and of the convex hull. Thus, if we can ensure

then we can ensure that exponential stability holds. Rearranging this condition and writing it in matrix form yields an inequality of the desired form. We note that by definition of the specifications (A.10), there is some KK corresponding to PP such that the policy u=Kxu=Kx satisfies all of the above inequalities; thus, Kx∈CPLDI(x)Kx\in{\mathcal{C}_{\text{PLDI}}(x)}, and CPLDI(x){\mathcal{C}_{\text{PLDI}}(x)} is non-empty. Further, as the above inequality represents a linear constraint in uu, this set is convex in uu. ∎

We note that the relevant projection PCPLDI(x)\mathcal{P}_{\mathcal{C}_{\text{PLDI}}(x)} represents a projection onto an intersection of halfspaces, and can thus be implemented via differentiable quadratic programming [Amos and Kolter, 2017].

For the H∞H_{\infty} control system (A.11), relevant sets of actions satisfying the condition (A.12) are given by the following theorem.

Consider the system (A.11), some stability parameter α>0\alpha>0, and a Lyapunov function V(x)=xTPxV(x)=x^{T}Px with PP satisfying Equation (A.18). Assuming PP exists, define

We seek to find a set of actions such that the condition E(x,x˙,u)≤0\mathcal{E}(x,\dot{x},u)\leq 0 is satisfied along all possible trajectories of (A.11), where E\mathcal{E} is defined as in (A.12). A set of actions satisfying this condition at a given xx is given by

Expanding and rearranging terms, this becomes

We note that by definition of the specifications (A.18), there is some KK corresponding to PP such that the policy u=Kxu=Kx satisifies the conditions above (see (A.17)); thus, Kx∈CH∞Kx\in\mathcal{C}_{H_{\infty}}, and CH∞\mathcal{C}_{H_{\infty}} is non-empty. We note further that CH∞\mathcal{C}_{\mathcal{H}_{\infty}} is an ellipsoid in the control action space, and is thus convex in uu. ∎

Appendix C A fast, differentiable solver for second-order cone projection

In order to construct the robust policy class described in Section 4 for the general NLDI system (3) and the H∞H_{\infty} setting (A.11), we must project a nominal (neural network-based) policy onto the second-order cone constraints described in Theorem 1 and Appendix B.2, respectively. As this projection operation does not necessarily have a closed form, we implement it via a custom differentiable optimization solver.

More generally, consider a set of the form

where for brevity we define G=[AcT]G=\begin{bmatrix}A\\ c^{T}\end{bmatrix} and h=[bd]h=\begin{bmatrix}b\\ d\end{bmatrix}, and where 1F\mathbf{1}_{\mathcal{F}} denotes the indicator function for membership in the set F\mathcal{F}.

We describe our fast solution technique for computing this projection, as well as our method for obtaining gradients through the solution.

and the dual problem is given by max⁡μmin⁡x,zL(x,z,μ)\max_{\mu}\min_{x,z}\mathscr{L}(x,z,\mu). To form the dual problem, we minimize the Lagrangian with respect to xx and zz as

We note that the first term on the right side is minimized at x⋆(μ)=y+GTμx^{\star}(\mu)=y+G^{T}\mu. Thus, we see that

where the last identity follows from definition of the second-order cone F\mathcal{F}. Hence the negative dual problem becomes

where Lf=λmax⁡(GGT)L_{f}=\lambda_{\max}(GG^{T}) is the Lipschitz constant of ff, and PF\mathcal{P}_{\mathcal{F}} is the projection operator onto F\mathcal{F} (which has a closed form solution; see Bauschke ). Letting mf=λmin⁡(GGT)m_{f}=\lambda_{\min}(GG^{T}) denote the strong convexity constant of ff, the momentum parameter is then scheduled as [Nesterov, 2013]

After computing the optimal dual variable μ⋆\mu^{\star}, i.e., the fixed point of (C.10), the optimal primal variable can be recovered via the equation x⋆=y+GTμ⋆x^{\star}=y+G^{T}\mu^{\star} (as can be observed from the first-order conditions of the Lagrangian (C.4)).

C.2 Obtaining gradients

In order to incorporate the above projection into our neural network, we need to compute the gradients of all problem variables (i.e., GG, hh, and yy) through the solution x⋆x^{\star}. In particular, we note that x⋆x^{\star} has a direct dependence on both GG and yy, and an indirect dependence on all of GG, hh, and yy through μ⋆\mu^{\star}.

To compute the relevant gradients through μ⋆\mu^{\star}, we apply the implicit function theorem to the fixed point of the update equations (C.10). Specifically, as these updates imply that μ⋆=ν⋆\mu^{\star}=\nu^{\star}, their fixed point can be written as

Define M:=\frac{\partial\mathcal{P}_{\mathcal{F}}(\cdot)}{\partial(\cdot)}\big{|}_{(\cdot)=\mu^{\star}-\frac{1}{L_{f}}\nabla f(\mu^{\star})}, and note that ∇f(μ⋆)=GGTμ⋆+Gy+h\nabla f(\mu^{\star})=GG^{T}\mu^{\star}+Gy+h. The differential of the above fixed-point equation is then given by

Rearranging terms to separate the differentials of problem outputs from problem variables, we see that

where II is the identity matrix of appropriate size.

We can then form the relevant gradient terms directly as

Accounting additionally for the direct dependence of some of the problem variables on x⋆x^{\star} (recalling that x⋆=y+GTu⋆x^{\star}=y+G^{T}u^{\star}), the desired gradients are then given by

Appendix D Writing the cart-pole problem as an NLDI

D.2 Obtaining C𝐶C and D𝐷D

We solve separately for each FiF_{i} to minimize the difference between the right and left sides of Equation (D.5) (while enforcing that the right side is larger than the left side) over a discrete grid of points within x‾≤x≤xˉ\underline{x}\leq x\leq\bar{x} and u‾≤u≤uˉ\underline{u}\leq u\leq\bar{u}. By assuming that FiF_{i} is symmetric, we are able to cast this as a linear program in the upper triangular entries of FiF_{i}.

To obtain the matrices CC and DD used for the cart-pole experiments in the main paper, we let xˉ=[1.520.21.5]T\bar{x}=\begin{bmatrix}1.5&2&0.2&1.5\end{bmatrix}^{T}, uˉ=10\bar{u}=10, x‾=−xˉ\underline{x}=-\bar{x}, and u‾=−uˉ\underline{u}=-\bar{u}. As each entry-wise difference in Equation (D.5) contained exactly three variables (i.e., a total of three entries from xx and uu), we solved each entry-wise linear program over a mesh grid of 50 points per variable.

Appendix E Writing quadrotor as an NLDI

E.2 Obtaining C𝐶C and D𝐷D

We obtain the matrices CC and DD via a similar method as described in Appendix D, though in practice we only consider the linearization error with respect to xx (i.e., since the dynamics are linear with respect to uu, we have D=0D=0). We let xˉ=[110.150.60.61.3]\bar{x}=\begin{bmatrix}1&1&0.15&0.6&0.6&1.3\end{bmatrix} and x‾=−xˉ\underline{x}=-\bar{x}. As for cart-pole, each entry wise difference in the equivalent of Equation (D.5) contained exactly three variables (i.e., a total of three entries from xx and uu), and each entry-wise linear program was solved over a mesh grid of 50 points per variable.

Appendix F Details on the microgrid setting

To construct an NLDI of the form (3) for this system, we directly use the AA, BB, and GG matrices given in Lam et al. . We generate CC i.i.d. from a normal distribution and let D=0D=0, to represent the fact that the disturbance ww and the entries of the state xx are correlated, but that ww is likely not correlated with the actions uu. Finally, we let QQ and RR be diagonal matrices with 1 in the entries corresponding to quantities represented in the performance index yy, and with 0.1 in the rest of the diagonal entries, to emphasize that the variables in yy are the most important in describing the performance of the system.

Appendix G Generating an adversarial disturbance

In the NLDI settings explored in our experiments, we seek to construct an “adversarial” disturbance w(t)w(t) that obeys the relevant norm bounds ∥w(t)∥2≤∥Cx(t)+Du(t)∥2\|w(t)\|_{2}\leq\|Cx(t)+Du(t)\|_{2} while maximizing the loss. To do this, we use a model predictive control method where the actions taken are w(t)w(t). Specifically, for each policy π\pi, we model w(t)w(t) as a neural network specific to that policy. Every 10 steps of a roll-out, we optimize w(t)w(t) through gradient descent to maximize the loss over a horizon of 40 steps, subject to the constraint ∥w(t)∥2≤∥Cx(t)+Du(t)∥2\|w(t)\|_{2}\leq\|Cx(t)+Du(t)\|_{2}.

Appendix H Additional experimental details

Initial states. To pick initial states in our experiments, for the synthetic settings, we sample each attribute of the state i.i.d. from a standard Gaussian distribution. For cart-pole and planar quadrotor, we sample uniformly from bounds chosen such that the non-robust LQR algorithm (under the original dynamics) did not go unstable. For cart-pole, these bounds were chosen to be px∈p_{x}\in, φ∈[−0.1,0.1]\varphi\in[-0.1,0.1], p˙x=φ˙=0\dot{p}_{x}=\dot{\varphi}=0. For planar quadrotor, these bounds were px,pz∈p_{x},p_{z}\in, φ∈[−0.05,0.05]\varphi\in[-0.05,0.05], px˙=pz˙=φ˙=0\dot{p_{x}}=\dot{p_{z}}=\dot{\varphi}=0.

Constructing NLDI bounds. Given these initial states, for the cart-pole and quadrotor settings, we needed to construct our NLDI disturbance bounds such that they would hold over the entire trajectory of the robust policy; if not, the robustness specification (A.8) would not hold, and our agent might in fact increase the Lyapunov function. To ensure this approximately, we used a simple heuristic: we ran the (non-robust) LQR agent for a full episode with 50 different starting conditions, and constructed an L∞L_{\infty} ball around all states reached in any of these trajectories. We then used these L∞L_{\infty} balls on the states to construct the matrices CC and DD for our disturbance bounds, using the procedure described in Appendices D and E.

Computing infrastructure and runtime. All experiments were run on an XPS 13 laptop with an Intel i7 processor. The planar quadrotor and synthetic NLDI experiment with D=0D=0 took about 1 day to run (since the projections were simple half-space projections), while all the other synthetic domains and cart-pole took about 3 days to run. The majority of the run-time was in computing the adversarial disturbances for test-time evaluations.

Hyperparameter selection. For our experiments, we did not perform large parameter searches. The learning rate we chose for our model-based planner, (both robust and non-robust) remained constant for the different domains; we tried learning rates of $1\times10−31\text{\times}{10}^{-3},1\times10−41\text{\times}{10}^{-4},1\times10−51\text{\times}{10}^{-5}andfoundand found1\text{\times}{10}^{-3}workedbestforthenon−robustversionandworked best for the non-robust version and1\text{\times}{10}^{-4}$ worked best for the robust version. For our PPO hyperparameters, we simply used those used in the original PPO paper.

One parameter we had to tune for each environment was the time step. In particular, we had to pick a time step high enough that we could run episodes for a reasonable total length of time (within which the non-robust agents would go unstable), but low enough to reasonably approximate a continuous-time setting (since, for our robustness guarantees, we assume the agent’s actions evolve in continuous time). Our search space was small, however, consisting of 0.050.05, 0.020.02, 0.010.01, and 0.0050.005 seconds.

Trajectory plots. Figure H.1 shows sample trajectories of different methods in the cart-pole domain under adversarial dynamics. The non-robust LQR and model-based planning approaches both diverge and the non-robust PPO doesn’t diverge, but doesn’t clearly converge after 10 seconds. The robust methods, on the other hand, all clearly converge after 10 seconds.

Runtime comparison. Tables H.1 and H.2 show the evaluation and training time of our methods and the baselines over 50 episodes run in parallel. In the NLDI cases where D=0D=0, i.e., Generic NLDI (D=0D=0) and Quadrotor, our projection adds only a very small computational cost. In the other cases, the additional computational cost is more significant, but our method is still far less expensive than the Robust MPC method.

In addition to the NLDI settings explored in the main text, we test the performance of our method on PLDI and H∞H_{\infty} control settings. As for the experiments in the main text, we choose a time discretization based on the speed at which the system evolves, and run each episode for 200 steps over this discretization. In both cases, we use a randomly generated LQR objective where the matrices Q1/2Q^{1/2} and R1/2R^{1/2} are drawn i.i.d. from a standard normal distribution.

Synthetic PLDI setting. We generate PLDI instances (A.9) with s=5s=5, a=3a=3, and L=3L=3. Specifically, we generate convex hull matrices (A1,B1),…,(A3,B3)(A_{1},B_{1}),\ldots,(A_{3},B_{3}) i.i.d. from normal distributions, and generate (A(t),B(t))(A(t),B(t)) by using a randomly-initialized neural network with softmax output to weight the convex hull matrices. Episodes were run for 2 seconds at a discretization of 0.01 seconds.

Synthetic \boldmathH∞\boldmath{H_{\infty}} setting. We generate H∞H_{\infty} control instances (A.11) with s=5s=5, a=3a=3, and d=2d=2 by generating matrices A,BA,B and GG i.i.d. from normal distributions. The disturbance w(t)w(t) was produced using a randomly-initialized neural network, with its output scaled to satisfy the L2\mathcal{L}_{2} bound on the disturbance. Specifically, we scaled the output of the neural network to satisfy an attenuating norm-bound on the disturbance; at time tt, the norm-bound was given by 20×f(2×t/T)20\times f(2\times t/T), where TT is the time horizon and ff is the standard normal PDF function. Episodes were run for T=2T=2 seconds at a discretization of 0.01 seconds.

Results are given in Figure I.1 and Table I.1.

Appendix J Notes on linearization via PLDIs and NLDIs

While we linearize the cart-pole and quadrotor dynamics via NLDIs in our experiments, we note that these dynamics can also be characterized via PLDIs. More generally, in this section, we show how we can use the framework of PLDIs to model linearization errors arising in the analysis of nonlinear systems.

for some z=tξz=t\xi, where t∈t\in. Now, let p=s+ap=s+a. Defining the Jacobian of ff as

and recalling that f(0)=0f(0)=0, we can rewrite (J.2) as

We now seek to bound the Jacobian using polytopic bounds. To do this, note that we can write

where AκA_{\kappa}’s are the vertices of the polytope in (J.6), i.e.,

Together, Equations (J.2), (J.4), (J.7), and (J.8) characterize the original nonlinear dynamics as a PLDI.

We note, however, that we can in fact express this PLDI more concisely as an NLDI. More precisely, we would like to find matrices A,B,CA,B,C parameterizing the form of NLDI below, which is equivalent to that presented in Equation (3) (see Chapter 4 of Boyd et al. ):

It can shown that the solution to the SDP

yields the matrices AA, BB, and CC with V=CTCV=C^{T}C and W=BBTW=BB^{T}, which can be used to construct NLDI (J.9). While the NLDI here is more concise than the PLDI, the trade-off is that the NLDI norm bounds obtained via this method may be rather loose. As such, for our settings, we obtain NLDI bounds numerically (see Appendices D and E), as these are tighter than NLDI specifications obtained via the above method (though they are potentially slightly inexact). An alternative approach would be to examine how to tighten the conversion from PLDIs to NLDIs, which has been explored in other work (e.g. Kuiava et al. ).