Abstract
This work presents a general theory for the construction of a polyhedral outer approximation of the reachable set (“polyhedral bounds”) of a dynamic system subject to time-varying inputs and uncertain initial conditions. This theory is inspired by the efficient methods for the construction of interval bounds based on comparison theorems. A numerically implementable instance of this theory leads to an auxiliary system of differential equations which can be solved with standard numerical integration methods. Meanwhile, the use of polyhedra provides greater flexibility in defining tight enclosures on the reachable set. These advantages are demonstrated with a few examples, which show that tight bounds can be efficiently computed for general, nonlinear systems. Further, it is demonstrated that the ability to use polyhedra provides a means to meaningfully distinguish between time-varying and constant, but uncertain, inputs.
Similar content being viewed by others
References
Althoff M, Stursberg O, Buss M (2008) Reachability analysis of nonlinear systems with uncertain parameters using conservative linearization. In: Proceedings of the 47th IEEE conference on decision and control, pp 4042–4048
Bemporad A, Morari M (1999) Robust model predictive control: a survey. In: Garulli A, Tesi A (eds) Robustness in identification and control, lecture notes in control and information sciences. Springer, London, pp 207–226
Ben-Israel A, Greville TNE (2003) Generalized inverses: theory and applications, 2nd edn. Springer, New York
Bertsimas D, Tsitsiklis JN (1997) Introduction to linear optimization. Athena Scientific, Belmont
Chachuat B, Villanueva M (2012) Bounding the solutions of parametric ODEs: when Taylor models meet differential inequalities. In: Bogle IDL, Fairweather M (eds) Proceedings of the 22nd European symposium on computer aided process engineering, pp 17–20
Chutinan A, Krogh BH (1998) Computing polyhedral approximations to flow pipes for dynamic systems. In: Proceedings of the 37th IEEE conference on decision and control, vol 2, pp 2089–2094
Chutinan A, Krogh BH (2003) Computational techniques for hybrid system verification. IEEE Trans Autom Control 48(1):64–75
Dang T, Maler O (1998) Reachability analysis via face lifting. In: Henzinger T, Sastry S (eds) Hybrid systems: computation and control. Springer, Berlin, pp 96–109
Dontchev AL, Lempio F (1992) Difference methods for differential inclusions: a survey. SIAM Rev 34(2):263–294
Filippov AF (1988) Differential equations with discontinuous righthand sides. Kluwer Academic, Boston
Gomez JA, Höffner K, Barton PI (2014) DFBAlab: a fast and reliable MATLAB code for dynamic flux balance analysis. BMC Bioinform 15(409):1–10
Greenstreet MR, Mitchell I (1998) Integrating projections. In: Henzinger T, Sastry S (eds) Hybrid systems: computation and control. Springer, Berlin, pp 159–174
Greenstreet MR, Mitchell I (1999) Reachability analysis using polygonal projections. In: Vaandrager FW, van Schuppen JH (eds) Hybrid systems: computation and control, lecture notes in computer science. Springer, Berlin, pp 103–116
Harrison GW (1977) Dynamic models with uncertain parameters. In: Avula XJR (ed) Proceedings of the first international conference on mathematical modeling, vol 1, pp 295–304. University of Missouri Rolla
Harwood SM, Höffner K, Barton PI (2015) Efficient solution of ordinary differential equations with a parametric lexicographic linear program embedded. Numerische Mathematik. doi:10.1007/s00211-015-0760-3
Harwood SM, Scott JK, Barton PI (2013) Bounds on reachable sets using ordinary differential equations with linear programs embedded. In: Tarbouriech S, Krstic M (eds) Nonlinear control systems, vol 9, pp 62–67
Harwood SM, Scott JK, Barton PI (2015) Bounds on reachable sets using ordinary differential equations with linear programs embedded. IMA J Math Control Inf. doi:10.1093/imamci/dnu054
Hindmarsh AC, Brown PN, Grant KE, Lee SL, Serban R, Shumaker DE, Woodward CS (2005) SUNDIALS: suite of nonlinear and differential/algebraic equation solvers. ACM Trans Math Soft 31(3):363–396
Hwang I, Stipanovic DM, Tomlin CJ (2003) Applications of polytopic approximations of reachable sets to linear dynamic games and a class of nonlinear systems. In: Proceedings of the American control conference, vol 6, pp 4613–4619
IBM: IBM ILOG CPLEX: high performance mathematical programming solver. http://www-01.ibm.com/software/commerce/optimization/cplex-optimizer/ (2014)
Kapela T, Zgliczyński P (2009) A Lohner-type algorithm for control systems and ordinary differential inclusions. Discrete Contin Dyn Syst Ser B 11(2):365–385
Kieffer M, Walter E (2004) Guaranteed nonlinear state estimator for cooperative systems. Numer Algorithms 37(1–4):187–198
Lambert JD (1991) Numerical methods for ordinary differential systems: the initial value problem. Wiley, New York
Le Guernic C (2009) Reachability analysis of hybrid systems with linear continuous dynamics. Ph.D. thesis, Université Joseph-Fourier
Lin Y, Stadtherr MA (2007) Validated solutions of initial value problems for parametric ODEs. Appl Numer Math 57(10):1145–1162
Lin Y, Stadtherr MA (2008) Fault detection in nonlinear continuous-time systems with uncertain parameters. AIChE J 54(9):2335–2345
Makino K, Berz M (1996) Remainder differential algebras and their applications. In: Berz M, Bischof C, Corliss G, Griewank A (eds) Computational differentiation: techniques, applications, and tools, Chap 5, pp 63–75. SIAM
Mangasarian OL, Shiau TH (1987) Lipschitz continuity of solutions of linear inequalities, programs, and complementarity problems. SIAM J Control Optim 25(3):583–595
Mattheij R, Molenaar J (2002) Ordinary differential equations in theory and practice. SIAM, Philadelphia
Mitchell I, Tomlin CJ (2000) Level set methods for computation in hybrid systems. In: Lynch N, Krogh BH (eds) Hybrid systems: computation and control. Springer, Berlin, pp 310–323
Mitchell IM, Bayen AM, Tomlin CJ (2005) A time-dependent Hamilton–Jacobi formulation of reachable sets for continuous dynamic games. IEEE Trans Autom Control 50(7):947–957
Mitsos A, Chachuat B, Barton PI (2009) McCormick-based relaxations of algorithms. SIAM J Optim 20(2):573–601
Nedialkov NS (1999) Computing rigorous bounds on the solution of an initial value problem for an ordinary differential equation. Ph.D. thesis, University of Toronto
Neumaier A (1990) Interval methods for systems of equations. Cambridge University Press, New York
Ramdani N, Meslem N, Candau Y (2009) A hybrid bounding method for computing an over-approximation for the reachable set of uncertain nonlinear systems. IEEE Trans Autom Control 54(10):2352–2364
Ramdani N, Nedialkov NS (2011) Computing reachable sets for uncertain nonlinear hybrid systems using interval constraint-propagation techniques. Nonlinear Anal Hybrid Syst 5:149–162
Scott JK (2012) Reachability analysis and deterministic global optimization of differential-algebraic systems. Ph.D. thesis, Massachusetts Institute of Technology
Scott JK, Barton PI (2010) Tight, efficient bounds on the solutions of chemical kinetics models. Comput Chem Eng 34(5):717–731
Scott JK, Barton PI (2013) Bounds on the reachable sets of nonlinear control systems. Automatica 49(1):93–100
Scott JK, Barton PI (2013) Improved relaxations for the parametric solutions of ODEs using differential inequalities. J Glob Optim 57:143–176
Scott JK, Chachuat B, Barton PI (2013) Nonlinear convex and concave relaxations for the solutions of parametric ODEs. Optimal Control Appl Methods 34(2):145–163
Singer AB, Barton PI (2006) Bounding the solutions of parameter dependent nonlinear ordinary differential equations. SIAM J Sci Comput 27(6):2167–2182
Singer AB, Barton PI (2006) Global optimization with nonlinear ordinary differential equations. J Glob Optim 34(2):159–190
Srinivasan B, Amrhein M, Bonvin D (1998) Reaction and flow variants/invariants in chemical reaction systems with inlet and outlet streams. AIChE J 44(8):1858–1867
Szarski J (1965) Differential inequalities. Polish Scientific, Warszawa
Wets RJB (1985) On the continuity of the value of a linear program and of related polyhedral-valued multifunctions. Math Program Study 24:14–29
Acknowledgments
The authors thank Novartis Pharmaceuticals as part of the Novartis-MIT Center for Continuous Manufacturing for funding this research, and Garrett Dowdy for some suggestions on the technical results.
Author information
Authors and Affiliations
Corresponding author
Appendix: Parametric affine relaxations
Appendix: Parametric affine relaxations
This section discusses a method for constructing affine relaxations of the dynamics which satisfy Assumption 4. The method relies on Proposition 5.1 in [17], which is repeated below. It is a general result, and to implement it, one needs specific relaxations of various arithmetic operations. Such operations are listed in Table 1 below.
Proposition 2
(Proposition 5.1 in [17]) Let \(m, n \in \mathbb {N}\). Let \(Z \subset \mathbb {R}^n\) be a nonempty open set, and let \(Y \subset \mathbb {R}^m\). Let \(\mathbf {g} : Z \rightarrow \mathbb {R}^m\) and \(f : Y \rightarrow \mathbb {R}\). Define \({Z}^{\mathbb {I}} = \{ (\mathbf {v},\mathbf {w}) \in \mathbb {R}^n \times \mathbb {R}^n : \mathbf {v}\le \mathbf {w}, [\mathbf {v}, \mathbf {w}] \subset Z \}\), and similarly define \(Y^{\mathbb {I}}\). For \(i \in \{ 1, \dots , m\}\), let \(\mathbf {g}_i^{al}\) and \(\mathbf {g}_i^{au}\) be locally Lipschitz continuous mappings \({Z}^{\mathbb {I}} \rightarrow \mathbb {R}^n\) and \(g_i^{bl}, g_i^{bu}, g_i^L, g_i^U\) be locally Lipschitz continuous mappings \({Z}^{\mathbb {I}} \rightarrow \mathbb {R}\) which satisfy
for all \((\mathbf {v},\mathbf {w}) \in {Z}^{\mathbb {I}}\).
Let \(\mathbf {f}^{al}\) and \(\mathbf {f}^{au}\) be locally Lipschitz continuous mappings \({Y}^{\mathbb {I}} \rightarrow \mathbb {R}^m,\) and \(f^{bl}\) and \(f^{bu}\) be locally Lipschitz continuous mappings \({Y}^{\mathbb {I}} \rightarrow \mathbb {R}\) which satisfy
for all \((\mathbf {v}',\mathbf {w}') \in {Y}^{\mathbb {I}}\).
Let \(h : Z \rightarrow \mathbb {R}\) be defined by \(h(\mathbf {z}) = f( \mathbf {g}(\mathbf {z}) )\). For \(i \in \{1, \dots , m\}\), let
Let \(\mathbf {h}^{al}, \mathbf {h}^{au} : {Z}^{\mathbb {I}} \rightarrow \mathbb {R}^n\) and \(h^{bl}, h^{bu} : {Z}^{\mathbb {I}} \rightarrow \mathbb {R}\) be defined by
Then, \(\mathbf {h}^{al}, \mathbf {h}^{au}, h^{bl}, h^{bu}\) are locally Lipschitz continuous mappings on \({Z}^{\mathbb {I}}\) which satisfy
for all \((\mathbf {v},\mathbf {w}) \in {Z}^{\mathbb {I}}\).
Constructing relaxations which satisfy Assumption 4 proceed by recursively applying Proposition 2, illustrated by the following example.
Example 1
This example demonstrates how to construct affine relaxations for a simple function (in fact, this function is part of the dynamics of the Lotka–Volterra problem in Sect. 5.1). Let \(f : \mathbb {R} \times \mathbb {R}^2 \ni (p,\mathbf {z}) \mapsto p z_1 (1 - z_2)\). As in the Lotka–Volterra problem, p is an uncertain parameter. The goal is to construct, for any \(p \in [p^L,p^U]\) such that \(p^L > 0\), affine relaxations of \(f(p,\cdot )\) on any interval \([\mathbf {z}^L,\mathbf {z}^U]\) such that \(z_1^L > 0\) and \(z_2^L > 1\) (in the context of Assumption 4, obtaining a nontrivial \(\widetilde{\mathbf {c}}_i^u\) is beneficial only when U is not an interval). Furthermore, one desires that these relaxations are locally Lipschitz continuous with respect to \((\mathbf {z}^L,\mathbf {z}^U)\).
The process resembles the construction of an interval enclosure of the range of f via interval arithmetic, and indeed part of the method involves interval arithmetic. Evaluation of f is broken down into a sequence of auxiliary variables called “factors,” which can be expressed as simple arithmetic operations on previously computed factors. An interval enclosure and affine relaxation of each factor can also be computed, and following the rules in Proposition 2 and Table 1, the affine relaxations will also be locally Lipschitz continuous in the manner desired. See Table 2 for the factored expression. Note that factor \(v_3\), corresponding to the parameter p, is initialized with the trivial affine relaxations \(\mathbf {0} ^{ T }\mathbf {z}+ p^L \le p \le \mathbf {0} ^{ T }\mathbf {z}+ p^U\). This ensures that the final relaxations obtained are valid for all \(p \in [p^L,p^U]\). Also, note that the restrictions \(z_1^L > 0, z_2^L > 1\), and \(p^L > 0\), simplify the evaluation and preclude the need to consider the different “branches” when constructing the affine relaxations for factors \(v_5\) and \(v_6\), as indicated in Proposition 2 (for example, this implies that \({1/2}(v_4^L + v_4^U) < 0\)). Although in general, the different cases must be taken into account.
The final factor, \(v_6\), gives the value of f, and thus one also has
for all \((p,\mathbf {z}) \in [p^L,p^U] \times [\mathbf {z}^L,\mathbf {z}^U]\). However, by virtue of Proposition 2, \(\mathbf {v}_6^{al}, \mathbf {v}_6^{au}, v_6^{bl}, v_6^{bu}\) can be considered locally Lipschitz continuous functions with respect to \((\mathbf {z}^L,\mathbf {z}^U)\).
More generally, Assumption 4 requires parametric relaxations which are continuous on \(T \times D_x^{\mathbb {I}}\). To include dependence on t, one could construct relaxations with respect to t as well, over the degenerate interval [t, t]. Then, it is clear that the final relaxations would be locally Lipschitz continuous on \((T \times D_x)^{\mathbb {I}}\). The following lemmata shows that this yields the desired properties.
Lemma 6
Assume \(m, n, p \in \mathbb {N}\). Let \(C \subset \mathbb {R}^m\) be nonempty and compact and \(D \subset \mathbb {R}^n\) be nonempty. Let \(\mathbf {g} : C \times D \rightarrow \mathbb {R}^p\) be locally Lipschitz continuous. Then for all \(\mathbf {z}\in D\), there exists a neighborhood \(N(\mathbf {z})\) and \(L > 0\) such that for all \((\mathbf {y},\mathbf {z}_1), (\mathbf {y},\mathbf {z}_2) \in C \times N(\mathbf {z}) \cap D\)
Proof
Choose \(\mathbf {z}\in D\). For each \(\mathbf {y}\in C\), let \(N(\mathbf {y},\mathbf {z})\) be a neighborhood of \((\mathbf {y},\mathbf {z})\) such that \(\mathbf {g}\) is Lipschitz continuous on \(N(\mathbf {y},\mathbf {z}) \cap (C \times D)\), with corresponding Lipschitz constant \(L(\mathbf {y})\). However, this collection of open sets forms an open cover of \(C \times \{ \mathbf {z}\}\), which is compact, and thus we can choose a finite number of these neighborhoods \(\{ N(\mathbf {y}_i, \mathbf {z}) : 1 \le i \le k\}\), such that their union, \(\widetilde{N}\), contains \(C \times \{ \mathbf {z}\}\). Let L be the (finite) maximum of the corresponding Lipschitz constants (i.e. \(L = \max \{ L(\mathbf {y}_i) : 1 \le i \le k \}\)). Note that \(\widetilde{N}\) is an open set, and \(\mathbf {g}\) is Lipschitz continuous on \(\widetilde{N} \cap (C \times D)\) with Lipschitz constant L.
We claim that there exists a \(\delta > 0\) such that \(C \times N_{\delta }(\mathbf {z}) \subset \widetilde{N}\) (where \(N_{\delta }(\mathbf {z})\) is viewed as a subset of \(\mathbb {R}^n\)). This follows from, for instance, Lemma 1 in Section 5 of [10]. The argument is that the complement of \(\widetilde{N}, \widetilde{N}^C\), is closed and disjoint from \(C \times \{\mathbf {z}\}\), and so there exists a \(\delta > 0\) such that the distance between any point in \(C \times \{ \mathbf {z}\}\) and any point in \(\widetilde{N}^C\) is greater than \(\delta \). This implies that \(C \times N_{\delta }(\mathbf {z})\) is disjoint from \(\widetilde{N}^C\), which in turn implies \(C \times N_{\delta }(\mathbf {z}) \subset \widetilde{N}\). The result follows from Lipschitz continuity on \((C \times N_{\delta }(\mathbf {z}) ) \cap ( C \times D) = C \times N_{\delta }(\mathbf {z}) \cap D\). \(\square \)
Lemma 7
Assume \(n, p \in \mathbb {N}\). Let \(T \subset \mathbb {R}\) be nonempty and compact and \(D \subset \mathbb {R}^n\) be nonempty. Define \((T \times D)^{\mathbb {I}} = \{ (s,\mathbf {v},t,\mathbf {w}) : s\le t, \mathbf {v}\le \mathbf {w}, [s,t] \times [\mathbf {v},\mathbf {w}] \subset T \times D \}\) and define \(D^{\mathbb {I}}\) similarly. Let \(\widetilde{\mathbf {g}} : (T \times D)^{\mathbb {I}} \rightarrow \mathbb {R}^p\) be locally Lipschitz continuous. Define \(\mathbf {g} : T \times D^{\mathbb {I}} \rightarrow \mathbb {R}^p\) by \(\mathbf {g}(t,\mathbf {v},\mathbf {w}) = \widetilde{\mathbf {g}}(t,\mathbf {v},t,\mathbf {w})\). Then \(\mathbf {g}\) is continuous, and for all \((\mathbf {v},\mathbf {w}) \in D^{\mathbb {I}}\) there exists a neighborhood \(N(\mathbf {v},\mathbf {w})\) and \(L > 0\) such that for all \((t,\mathbf {v}_1,\mathbf {w}_1), (t,\mathbf {v}_2,\mathbf {w}_2) \in T \times N(\mathbf {v},\mathbf {w}) \cap D^{\mathbb {I}}\)
Proof
The mapping \(\mathbf {h} : (t,\mathbf {v},\mathbf {w}) \mapsto (t, \mathbf {v}, t, \mathbf {w})\) is Lipschitz continuous, and so \(\mathbf {g}\), as the composition of locally Lipschitz \(\widetilde{\mathbf {g}}\) and \(\mathbf {h}\) is locally Lipschitz continuous on \(T \times D^{\mathbb {I}}\). Applying Lemma 6 we obtain the desired result. \(\square \)
Rights and permissions
About this article
Cite this article
Harwood, S.M., Barton, P.I. Efficient polyhedral enclosures for the reachable set of nonlinear control systems. Math. Control Signals Syst. 28, 8 (2016). https://doi.org/10.1007/s00498-015-0153-2
Received:
Accepted:
Published:
DOI: https://doi.org/10.1007/s00498-015-0153-2