Reachability Analysis Using Message Passing over Tree Decompositions

In this paper, we study efficient approaches to reachability analysis for discrete-time nonlinear dynamical systems when the dependencies among the variables of the system have low treewidth. Reachability analysis over nonlinear dynamical systems asks if a given set of target states can be reached, starting from an initial set of states. This is solved by computing conservative over approximations of the reachable set using abstract domains to represent these approximations. However, most approaches must tradeoff the level of conservatism against the cost of performing analysis, especially when the number of system variables increases. This makes reachability analysis challenging for nonlinear systems with a large number of state variables. Our approach works by constructing a dependency graph among the variables of the system. The tree decomposition of this graph builds a tree wherein each node of the tree is labeled with subsets of the state variables of the system. Furthermore, the tree decomposition satisfies important structural properties. Using the tree decomposition, our approach abstracts a set of states of the high dimensional system into a tree of sets of lower dimensional projections of this state. We derive various properties of this abstract domain, including conditions under which the original high dimensional set can be fully recovered from its low dimensional projections. Next, we use ideas from message passing developed originally for belief propagation over Bayesian networks to perform reachability analysis over the full state space in an efficient manner. We illustrate our approach on some interesting nonlinear systems with low treewidth to demonstrate the advantages of our approach.


Introduction
Reachability analysis asks whether a target set of states is reachable over a finite or infinite time horizon, starting from an initial set for a dynamical system. This problem is fundamental to the verification of systems, and is known to be challenging for a wide variety of models. This includes cyber-physical systems, physical and biological processes. In this paper, we study reachability analysis algorithms for nonlinear, discrete-time dynamical systems. The key challenge in analyzing such systems arises from the difficulty of representing the reachable sets of these systems. As a result, we resort to over-approximations of reachable sets using tractable set representations such as intervals [16], ellipsoids, polyhedra [19], and low degree semi-algebraic sets [2]. Whereas these representations are useful for reachability analysis, they also trade off the degree of overapproximation in representing various sets against the complexity of performing operations such as intersections, unions, projections and image computations over these sets. The theory of abstract interpretation allows us to design various abstract domains that serve as representations for sets of states in order explore these tradeoffs [17,18,34]. However, for nonlinear dynamical systems, these representations often become too conservative or too expensive as the number of state variables grow.
In this paper, we study reachability analysis using the idea of tree decompositions over the dependency graph of a dynamical system. Tree decompositions are a well-known idea from graph theory [37], used to study properties of various types of graphs. The treewidth of a graph is an intrinsic property of a graph that relates to how "far away" a given graph is from a tree. For instance, trees are defined to have a treewidth of 1. Many commonly occurring families of graphs such as series-parallel graphs have treewidth 2 and so on. Formally, a tree decomposition of a graph is a tree whose nodes are associated with subsets of vertices of the original graph along with some key conditions that will be described in Sect. 2. We use tree decompositions to build an abstract domain. The abstraction operation projects a set of states in the full system state space along each of the nodes of the tree, yielding various projections of this set. The concretization combines projections back into the high dimensional set. We study various properties of this abstract domain. First, we characterize abstract elements that can potentially be generated by projecting some concrete elements along the nodes of the tree (so called canonical elements, Definition 10). Next we characterize those sets which can be abstracted along the tree decomposition and reconstructed without any loss in information (tree decomposable sets, Definition 11). In this process, we also derive a message passing approach wherein nodes of the tree can exchange information to help refine sets of states in a sound manner. However, as we will demonstrate, the abstraction is "lossy" in general since projections of tree decomposable sets are not necessarily tree decomposable. We discuss some interesting ways in which precision can be regained by carefully analyzing this situation.
We combine these ideas together into an approach for reachability analysis of nonlinear systems using a grid domain that represents complex non convex sets as a union of fixed size cells using a gridding of the state-space. Although such a domain would be prohibitively expensive, we show that the tree decomposition abstract domain can drastically cut down on the complexity of computing reachable set overapproximations in this domain, yielding precise reachable set estimation for some nonlinear systems with low treewidth. We demonstrate our approach using a prototype implementation to show that for a restricted class of systems whose dependency graphs have low treewidth, our approach can be quite efficient and precise at the same time. Although some interesting systems have low treewidth property, it is easy to see that many systems will have treewidths that are too high for our approach. Our future work will consider how systems whose dependency graphs do not have sufficiently low treewidth can still be tackled in a conservative manner using some ideas from this paper.

Related Work
As mentioned earlier, the concept of tree decompositions and treewidth originated in graph theory [37]. The concept of treewidth gained popularity when it was shown that many NP-complete problems on graphs such as graph coloring could be solved efficiently for graphs with small treewidths [5]. Courcelle showed that the problem of checking if a given graph satisfies a formula in the monadic second order logic of graphs can be solved in linear time on graphs with bounded treewidth [15]. Several NP-complete problems such as 3-coloring can be expressed in this logic. Tree decompositions are also used to solve inference problems over Bayesian networks leading to representations of the Bayesian networks such as junction trees that share many of the properties of a tree decomposition [29]. In fact, belief propagation over junction trees is performed by passing messages that marginalize the probability distributions at various nodes of the tree. This is analogous to the message passing approach described here.
Tree decomposition techniques have been applied to model checking problems over finite state systems. For instance, Obdržálek show that the μ-calculus model checking problem can be solved in linear time in the size of a finite-state system whose graph has a bounded treewidth [35]. However, as Ferrara et al. point out, requiring the state graph of a system to have a bounded treewidth is often restrictive [24]. Instead, they study concurrent finite state systems wherein the communication graph has a bounded tree width. However, they conclude that while it is more reasonable to assume that the communication graph has a bounded tree width, it does not confer much advantages to verification problems. For instance, they show that the unrolling of these systems over time potentially results in unbounded treewidth. In this paper, we consider a different approach wherein we study the treewidth of dependency graphs of the system. We find that many systems have small treewidth and exploit this property. At the same time, we note that some of the benchmarks studied have "sparse" dependency graphs but treewidths that are too large for our approach.
Tree decomposition techniques have also been studied in static analysis of programs. The control and data flow graphs of structured programs without goto-statements or exceptional control flow are known to have small treewidth that can be exploited to perform compiler optimizations such as register allocation quite efficiently [38]. Chatterjee et al. have shown how to exploit small treewidth property of the control flow graphs of procedures in programs to perform interprocedural dataflow analysis by modeling the execution of programs with procedures as recursive state machines [11]. However, this approach seems restricted to control dominated properties such as sequence of function calls. In a followup work, they study control and data flow analysis problems for concurrent systems, wherein each component has constant treewidth [10]. In contrast, our approach studies dynamical system and consider tree decompositions of the data dependency graph.
The use of message passing in this paper closely resembles past work by Gulwani and Jojic [27]. Therein, a program verification problem involving the verification pre/post and intermediate assertions in a program is solved by passing messages that can propagate information between assertions along program paths in a randomized fashion. The approach is shown to be similar to loopy belief propagation used in Bayesian inference. The key differences are (a) we use data dependencies and tree decompositions rather than control flow paths to pass information along; and (b) we formally prove properties of the message passing algorithm.
Our approach is conceptually related to a well-known idea of speeding up static analysis of large programs using "packing" of program variables [4,28]. This approach was used successfully in the Astreé static analyzer [3,4,21]. Therein, clusters of variables representing small sets of dependent local and global are extracted. The remaining program variables are abstracted away and the abstract interpretation process is carried out over just these variables. The usefulness of this approach has borne out in other abstract interpretation efforts, including Varvel [28]. The key idea in this paper can be seen as a formalization of the rather informal "clustering" approach using tree decompositions. We demonstrate theoretical properties as well as the ability to pass messages to improve the results of the abstract interpretation.
The use of the dependency graph structure to speed up reachability analysis approaches has been explored in the past for speeding up Hamilton-Jacobi-based approaches by Mo Chen et al. [12] as well as flowpipe based approaches by Xin Chen et al. [13]. Both approaches consider the directed dependency graph wherein x i is connected to x j if the former appears in the dynamical update equation of the latter variable. The approaches perform a strongly connected component (SCC) decomposition and analyze each SCC in a topological sorted order. However, this approach breaks as soon as the system has large SCCs, which is common. As a result, Xin Chen et al. show how SCCs can themselves be broken into numerous subsets at the cost of a more conservative solution. In contrast, the tree decomposition approach can be applied to exploit sparsity even when the entire dependency graph is a single SCC.

Preliminaries
In this section, we will describe the system model under analysis, the dependency graph structure and the basics of tree decompositions. Let X : {x 1 , . . . , x n } be a set of system variables and x : X → R represent a valuation to these system variables. Let D be the domain of all valuations of X, that describes the state space of the system. For convenience let x i denote x(x i ). Also, let W : {w 1 , . . . , w m } represent disturbance variables and w : W → R represent a vector of m ≥ 0 external disturbance inputs that take values in some compact disturbance space W.

Definition 1 (Dynamical Model).
A model Π is a tuple X, W, D, W, f, X 0 , U , wherein X, W, D, W are as defined above, f is an arithmetic expression over variables in X, W describing the dynamics, X 0 is a set of possible initial valuations (states) and U is a designated set of unsafe states.
The dynamics are given by x(t + 1) = eval(f, x, w), wherein eval evaluates a given an expression f , a set of valuations to the system variables x ∈ D and disturbances w ∈ W, and returns a new set of valuations for each variable in X, denoted by x(t + 1).
For simplicity, we write f (x, w) to denote eval(f, x, w) for a function expression f . A state of the system is a valuation x : X → R such that x ∈ D. Given a finite sequence of disturbance inputs w(0), . . . , w(T ), for some T ≥ 0 and w(i) ∈ W for all i ∈ [0, T ], an execution of the system is a sequence of states x(0), . . . , x(T + 1), such that x(0) ∈ X 0 , x(t) ∈ D for t ∈ [0, T + 1] and x(t+1) = f (x(t), w(t)) for all t ∈ [0, T ]. According to these semantics, the system may fail to have an execution for a given disturbance sequence w(t), t ∈ [0, T ] and initial state x(0) if for some state x(t), we have f (x(t), w(t)) ∈ D.
A state x(t) is reachable (at time t) if there is an execution of the form x(0), . . . , x(t), satisfying the constraints above. We say that the unsafe state U is reachable iff some state x ∈ U is reachable. Furthermore, we say that U is reachable within a finite time horizon T , iff some state x ∈ U is reachable at time t ∈ [0, T ]. Example 1. Consider a nonlinear example of a dynamical model Π with state space x : (x 1 , x 2 , x 3 ) and w : (w 1 ). The dynamics can be written as parallel assignments to the state variables: The assignments are all evaluated in parallel to update the current state x(t) to a new state x(t + 1). The domain D is x i ∈ [−3, 3] for i = 1, 2, 3 and the disturbance We will now define the dependency (hyper)graph of the system Π. For convenience, we write the update function (expression) f of a system Π in terms of individual updates (f 1 , . . . , f n ), wherein x j = f j (x, w). We say that system variable x i (or disturbance variable w j ) is a proper input to the expression f k if x i (or w j ) occurs as a subterm in f k . Let inps(f k ) denote the set of all proper input variables to the function (expression) f k .
As an example, consider X = {x 1 , . . . , x 4 } and W = {w 1 , w 2 } and the expression f : x 1 x 4 −w 1 . The proper inputs to f are {x 1 , x 4 , w 1 }. We exclude cases such as g : sin 2 (x1) + cos 2 (x1) sin 2 (x2) + cos 2 (x2) that has {x 1 , x 2 } as proper inputs. However a simplification using elementary trigonometric rules can eliminate them. We will assume that all expressions are simplified to involve the least number of variables.

Definition 2 (Dependency Hypergraph).
A dependency hypergraph of a system Π has vertices V : X ∪ W , given by the union of the system and disturbance variables with hyperedge set E ⊆ 2 V given by E = {e 1 , . . . , e n }, wherein for each update x k := f k (x, w) (k = 1, . . . , n), we have the hyperedge e k : {x k } ∪ inps(f k ). In other words, each update x k := f k (x, w) yields an edge that includes x k along with all the system/disturbance variables that are proper inputs to f k . Example 2. The dependency hypergraph for the system from Example 1 has the vertices V : {x 1 , x 2 , x 3 , w 1 } and the edges {e 1 : {x 1 , x 2 }, e 2 : {x 2 , w 1 } and e 3 : {x 2 , x 3 }}.

Tree Decomposition
We will now discuss tree decompositions and the associated concept of treewidth of a hypergraph G : (V, E). The tree decomposition will be applied to the dependency hypergraphs (Definition 2) for systems Π (Definition 1).

Definition 3 (Tree Decomposition and Treewidth). Given a hypergraph
G : (V, E), a tree decomposition is a tree T : (N, C) and a mapping verts : N → 2 V , wherein N is the set of tree nodes, C is the set of tree edges and verts(·) associates each node u ∈ N with a set of graph vertices verts(n) ⊆ V . The tree decomposition satisfies the following conditions: 3. For each vertex v, for any two nodes n 1 , n 2 such that v ∈ verts(n 1 ) and v ∈ verts(n 2 ), then v ∈ verts(n) for each node n along the unique path between n 1 and n 2 in the tree. Stated another way, the subset of nodes N v : The width of a tree decomposition is given by max{|verts(n)| | n ∈ N } − 1. In other words, we find the node n in the tree whose associated set of vertices has the largest cardinality. We subtract one from this maximal cardinality to obtain the treewidth. A tree decomposition is optimal for a graph G if no other tree decomposition exists with a strictly smaller width. The treewidth of a hypergraph G is given by width of an optimal tree decomposition.
It is easy to show that if the graph G is a tree, it has treewidth 1. Likewise, a cycle has tree width 2.
Example 3. The tree decomposition of the hypergraph G from Example 2 has three nodes {n 1 , n 2 , n 3 } with edges (n 1 , n 2 ) and (n 2 , n 3 ). The nodes along with the associated vertex sets are as follows: Although the tree decomposition is not a rooted tree, we often designate an arbitrary node r ∈ N as the root node, and consider the tree T as a rooted tree with root r.
Finding a Tree Decomposition: Interestingly, the problem of finding the treewidth of a graph is itself a NP-hard problem. However, many practical approaches exist for graphs with small treewidths. For instance, Bodlaender presents an algorithm that runs in time O(k O(k 3 ) ) to construct a tree decomposition of width at most k or conclude that the treewidth of the graph is at least k + 1 [6]. Such an approach can be quite useful if a given graph is suspected to have a small tree width in the first place. Besides this, many efficient algorithms exist to approximate the treewidth of a graph to some constant factor. A detailed survey of these results is available elsewhere [7,8]. Open-source packages such as HTD can compute treewidth for graphs with thousands of nodes [1]. Finally, we note that if a tree decomposition of width k can be found, then one can be found with at most |V | nodes.

Lemma 1.
Let T be a tree decomposition for a (multi)graph G with vertices V and treewidth k. There exists a tree decompositionT of G with the same treewidth k, and at most |V | nodes.
A proof is provided in the extended version of the paper.

Abstract Domains Using Tree Decompositions
In this section, we will define abstract domains using tree decompositions of the dependency hypergraph of the system under analysis. Let Π be a transition system over system variables X. The concrete states are given by x ∈ D, wherein x : X → R maps each state variable x j ∈ X to its value x(x j ) (denoted x j ).

Definition 5 (Extensions). Let R be a set of states involving just the variables in the set
In other words, the extension of a set embeds each element in the larger dimensional space defined by J 2 allowing "all possible values" for the dimensions in J 2 \ J 1 . We will use the notation ext(S) to denote the set ext X (S), i.e, its extension to the entire set of state variables X. For a state x S , we will use ext(x S ) denote ext({x S }).

Definition 6 (Product (Join) of Sets). Let
Let T : (N, C) be a tree decomposition of the dependency hypergraph of the system. Recall that for each node n ∈ N we associate a set of system/disturbance variables denoted by verts(n). Let verts X (n) denote the set of system variables: verts(n)∩X. We say that an update function x k := f k (x, w) is associated with a node n in the tree iff {x k } ∪ inps(f k ) ⊆ verts(n).

Lemma 2. For every system variable
Proof. This follows from those of a tree decomposition that states that every hyperedge in the dependency hypergraph must belong to verts(n) for at least one node n ∈ N .

Abstraction and Concretization
We consider subsets of the concrete states for the system Π, i.e, the set 2 D , ordered by set inclusion as our concrete domain. Given a tree decomposition, T , we define an abstract domain through projection of a concrete set along verts(n) for each node n of T .

Definition 7 (Abstract Domain). Each element s of the abstract domain
A T is a mapping that associates each node n ∈ N with a set s(n) ⊆ proj(D, verts X (n)). For We will use the notation proj(S, n) for a node n ∈ N to denote proj(S, verts X (n)).

Definition 8 (Abstraction Map).
Given a tree decomposition T , the abstraction map α T takes a set of states S ⊆ D and produces a mapping that associates tree node n ∈ N to a projection of S along the variables verts X (n). Formally, Thus, an abstract state s is a map that associates each node n of the tree to a set s(n) ⊆ D n . We now define the concretization map γ T .

Definition 9 (Concretization Map).
The concretization γ T (s) of an abstract state is defined as γ T (s) : n∈N ext(s(n)). In other words, we take s(n) for every node n ∈ N , extend it to the full dimensional space of all system variables and intersect the result over all nodes n ∈ N . Example 4. Consider a simple tree decomposition T with 2 nodes n 1 , n 2 and a single edge (n 1 , n 2 ). Let verts(n 1 ) : {x 1 , x 2 } and verts(n 2 ) : {x 2 , x 3 }. Let the domain D be the set x i ∈ {1, 2, 3} for i = 1, 2, 3. We use the notation ( to denote a state x that maps x 1 to the value v 1 , x 2 to the value v 2 and so on. 3)}. We have that s : α(S) is the mapping that projects S onto the dimensions (x 1 , x 2 ) for node n 1 and (x 2 , x 3 ) for node n 2 : Likewise, we verify that the concretization map γ(s) will yields us: For convenience, if the tree T is clear from the context, we will drop the subscripts to simply write α and γ for the abstraction and concretization map, respectively. Proof. Let S, s be such that α(S) s. Therefore, proj(S, n) ⊆ s(n) ∀n ∈ N by the definition of . Pick any, x ∈ S. First, proj(x, n) ∈ proj(S, n) and therefore, proj(x, n) ∈ s(n) for all n ∈ N . Thus, x ∈ ext(s(n)) for each node n ∈ N . Therefore, x ∈ n∈N ext(s(n)), and hence, x ∈ γ(s), by defn. of γ. Therefore, S ⊆ γ(s).
The meet operation is defined as s 1 s 2 : λn. s 1 (n) ∩ s 2 (n), and likewise, the join is defined as s 1 s 2 : λn. s 1 (n) ∪ s 2 (n). We recall two key facts that follow from Galois connection between α and γ.
1. For any set S ⊆ D, we have S ⊆ γ(α(S)). Abstracting a concrete set and concretizing it back again "loses information". To see why, we start from α(S) α(S) and apply the Galois connection to derive S ⊆ γ(α(S)). 2. Likewise, for any abstract domain object s ∈ A, we have α(γ(s)) s. I.e, for any element s, taking its concretization and abstracting it "gains information". To prove this, we start from γ(s) ⊆ γ(s) and conclude that α(γ(s)) s.

Canonical Elements and Message Passing
In the tree decomposition, various nodes share information about the subsets of vertices associated with each node. Since the subsets have elements in common, it is possible that a node n 1 has information about a variable x 2 that is also present in some other node n 2 of the tree. We will now see how to take an abstract element s and refine each s(n) by exchanging information between nodes in a systematic manner.

Definition 10 (Canonical Elements
). An abstract element s is said to be canonical if and only if for each edge (n 1 , n 2 ) ∈ C in the tree: proj(s(n 1 ), CV X (n 1 , n 2 )) = proj(s(n 2 ), CV X (n 1 , n 2 )) .
In other words, if we took the common variables verts X (n 1 ) ∩ verts X (n 2 ), the set s(n 1 ) projected along these common variables is equal to the projection of s(n 2 ) along the common variables. Example 6. Consider the abstract element s 1 from Example 5: 1}. Therefore, s 1 fails to be canonical.
The key theorem of tree decomposition is that a canonical element in the abstract domain can be seen as the projection of a concrete set S along verts X (n) for each node n of the tree. To prove that we will first establish a useful property of a canonical element s.

Lemma 3.
For every canonical element s ∈ A, node n ∈ N and element x n ∈ s(n), we have that ext(x n ) ∩ γ(s) = ∅.
Stated another way, the lemma claims that for any canonical s, any x n ∈ s(n) can be extended to form some element of γ(s). A proof is provided in the extended version.

Theorem 2. An element s is canonical (Definition 10) if and only if s = α(S) for some concrete set S.
Ideally, in abstract interpretation, we would like to work with abstract domain objects that satisfy s = α(γ(s)). One way to ensure that is to take any given domain element s 0 and simply calculate out α(γ(s 0 )) by applying the maps. However, γ(s 0 ) in our domain takes lower dimensional projections and reconstructs a set in the full states pace. It may thus be too expensive to compute. Fortunately, canonical objects satisfy the equality s = α(γ(s)). Therefore, given any object s ∈ A that is not necessarily canonical, we would like to make it canonical: I.e, we seek an objectŝ such that γ(ŝ) = γ(s), butŝ is canonical. As mentioned earlier, directly computingŝ = α(γ(s)) can be prohibitively expensive, depending on the domain. We now describe a message passing approach.
First, we convert the tree T to a rooted tree by designating an arbitrary node r ∈ N as the root of the tree.
Message Passing along Edges: Let (n 1 , n 2 ) be an edge of the tree and s be an abstract element. A message from n 1 to n 2 is defined as the set msg(s, n 1 → n 2 ) : proj(s(n 1 ), CV(n 1 , n 2 )). In other words, we project the set s(n 1 ) along the dimensions that are common to (n 1 , n 2 ).
Once a node n 2 receives M : msg(s, n 1 → n 2 ), it processes the message by updating s(n 2 ) as s(n 2 ) := s(n 2 ) ∩ ext verts(n2) (M ). In other words, it intersects the message (extended to the dimensions in n 2 ) with the current set that is associated with n 2 .
A message msg(s, n 1 → n 2 ) is given by the set proj(s(n 1 ), This results in the new abstract object s wherein the element ( Upwards Message Passing: The upwards message passing works from leaves up to the root of the tree according to the following two rules: 1. First, each leaf of the tree n passes a message to its parent n p . The parent node n p intersects its current value s(n p ) with the message to update its current set. 2. After a node has received (and processed) a message from all its children, it passes a message up to its parent, if one exists.
The upwards message passing terminates at the root since it does not have a parent to send a message to.
Example 8. Going back to Example 7, we designate n 2 as the root and the upwards pass sends the messages msg(s, n 1 → n 2 ) and msg(s, n 3 → n 2 ). This results in the following updated element: 1. To initialize, the root sends a message to all its children. 2. After a node has received (and processed) a message from its parent, it sends a message to all its children.
The overall procedure to make a given abstract object s canonical is as follows: (a) perform an upwards message passing phase and (b) perform a downwards message passing phase.
Example 9. Going back to Example 8, the downward message passing phase sends messages from n 2 → n 1 and n 2 → n 3 . The resulting elementŝ is On the other hand, it is important to perform message passing upwards first and then downwards second. Reversing this does not yield a canonical element. For instance going back to Example 7, if we first performed a downwards pass from n 2 , the result is unchanged: Performing an upwards pass now yields the element s 2 : However this is not canonical, since the element (  3) in s 2 (n 1 ) violates the requirement over the edge (n 1 , n 2 ).
Letŝ be the resulting abstract object after the message passing procedure finishes.

Theorem 3. The result of message passingŝ is a canonical object, and it satisfies γ(ŝ) = γ(s).
Proof (Sketch). First, we note that whenever a message is passed for an abstract value s from node m to n along an edge (m, n) resulting in a new abstract value s : (P1) γ(s ) = γ(s); and (P2) the projection of s (n) along the dimensions CV(m, n) is now contained in that of s (m) along CV(m, n). Furthermore, property (P2) remains unchanged regardless of any future messages that are passed along the tree edges.
Next, it is shown that after each upwards pass, when a message is passed, property (P2) (stated above) holds for each node m and its parent node n since a message is passed from m to n. During the downwards pass, property (P2) holds for each node n and its child node m in the tree. Combining the two, we note that for each edge (m, n) in the tree, we have property (P2) in either direction guaranteeing that proj(s * (m), CV(m, n)) = proj(s * (n), CV(m, n)), for the final result s * , or in other words that s * is canonical.

Decomposable Sets and Post-conditions
We have already noted that for any concrete set over S ⊆ D, the process of abstracting it by projecting into nodes of a tree T , and re-concretizing it is "lossy": I.e, S ⊆ γ(α(S)). In this section, we study "tree decomposable" concrete sets S for which γ(α(S)) = S. Ideally, we would like to prove that if a set S is tree decomposable then so is the set post(S, Π) of next states. However, we will disprove this by showing a counterexample. Nevertheless, we will present an analysis of why this fact fails and suggest approaches that can "manage" this loss in precision.

Definition 11 (Decomposable Sets). We say that a set S is tree decomposable given a tree T iff γ(α(S)) = S.
This is in fact a "global" definition of decomposability. In fact, a nice "local" definition can be provided that is reminiscent of the notion of conditional independence in graphical models. We will defer this discussion to an extended version of this paper due to space limitations. 2)} and tree T below: 2) .}. We note that the set S is not tree decomposable. On the other hand, one can verify that the set 2)} is tree decomposable.
The following lemma will be quite useful.

Lemma 4. Let S 1 , S 2 be tree decomposable sets over T . Their intersection is tree decomposable.
Let Π be a transition system over system variables in x ∈ D. For a given set S ⊆ D, us define the post-condition post(S, Π) to be the set of states reachable in one step starting from some state in S: Let us also consider a transition relation R over pairs of states (x, x ) ∈ D⊗D: The relation R can be viewed as the intersection of n relations: R : xj ∈X R j , wherein In other words, R j is a component of R that models the update of the system variable x j . Also for each x j ∈ X, let e j : inps(f j ) ∪ x j be the inputs to the update function f j and the node x j itself. Given the tree T , we define the extended tree T as having the same node set N and edge set C as T . However, verts T (n) = verts T (n) ∪ {x j |x j ∈ verts T (n)}. Note that T with the labeling verts T satisfies all the condition of a tree decomposition for a graph G save the addition of vertices x i in each node of the tree. We will write verts (n) to denote the set verts T (n).

Lemma 5. The transition relation R of a system Π is tree T decomposable.
The proof is provided in the extended version and is done by writing R as an intersection of tree decomposable relations R j , and appealing to Lemma 4.
First, we show the negative result that the image of a tree (T ) decomposable set under a tree (T ) decomposable transition relation is not tree decomposable, in general.
Example 11. Let X = {x 1 , x 2 , x 3 } and consider again the tree decomposition from Example 10. Let S be the set {( x1 * , x2 * , x3 * )}, wherein we use the wild card character as notation that can be substituted for any element in the set {1, 2}. Therefore, we take S to be a set with 8 elements. Clearly S is tree decomposable in the tree T from Example 10.
Consider the transition relation R that will be written as the intersection of three transition relations: Clearly R is tree T decomposable. We can now compute the post-condition of S under this relation. The reader can verify the post-conditionŜ : 2)}. However,Ŝ is not tree decomposable. We note thatŝ : α(Ŝ) is the setŝ(n 1 ) : As noted above, the set R is tree T decomposable. If S is tree decomposable, we can extend S to a set S : ext X (S) that is now defined over X ∪ X and is also tree decomposable. As a result S ∩ R is also tree decomposable. However, the postcondition of S is the set proj(S ∩ R, X ). Thus, the key operation that failed was the projection operation involved in computing the post-condition. This suggests a possible solution to this issue albeit an expensive one: at each step, we maintain the reachable states using both current and next state variables, thus avoiding projection. In effect, the reachable states at the i th step will be entire trajectories of the system expressed over variables X 0 ∪ X 1 ∪ · · · X i . This is clearly not practical. However, a more efficient solution is to note that some of the current state variables can be projected out without losing the tree decomposability property. Going back to Example 11, we note that we can safely project away {x 1 , x 3 }, while maintaining the new reachable set in terms of (x 2 , x 1 , x 2 , x 3 ). In this way, we may recover the lost precision back.
In conclusion, we note that tree decompositions may lose precision over postconditions. However, the loss in precision can be avoided if carefully selected "previous state variables" are maintained as the computation proceeds. The question of how to optimally maintain this information will be investigated in the future.

Grid-Based Interval Analysis
We now combine the ideas to create a disjunctive interval analysis using tree decompositions. The main idea here is to apply tree decompositions not to the concrete set of states but to an abstraction of the concrete domain by grid-based intervals.
We will now describe the interval-based abstraction of sets of states dynamical system Π in order to perform over-approximate reachability analysis. Let us fix a system Π : x, w, D, W, f, X 0 , U as defined in Definition 1. We will assume that the domain of state variables D is a hyper-rectangle given by D : for each j = 1, . . . , n. In other words, each system variable x j lies inside the interval [L(x j ), U(x j )]. Likewise, we will assume that W : We will consider a uniform cell decomposition wherein each dimension is divided into some natural number M > 0 of equal sized subintervals. The i th subinterval of variable x j is denoted as subInt(x j , i), and is given by [L(x j ) + iδ j , L(x j ) + (i + 1)δ j ] for i = 0, . . . , M − 1 and δ j : . Similarly, we will define subInt(w k , i) for disturbance variables w k whose domains are also divided into M subdivisions. The overall domain D × W is therefore divided into M m+n cells wherein each cell is indexed by a tuple of natural numbers i : i 1 , . . . , i n , i n+1 , . . . , i n+m , such that i j ∈ {0, . . . , M − 1} and the cell corresponding to i is given by:

Definition 12 (Grid-Based Abstract Domain).
The grid based abstract domain is defined by the set C : P(i ∈ {0, . . . , M} m+n ), wherein each abstract domain element is a set of grid cells. The sets are ordered simply by set inclusion ⊆ between sets of grid cells. The abstraction map α C : P(D) → C is defined as follows: The concretization map γ C is defined above in (1).

Definition 13 (Interval Propagator). An interval propagator (IP) is a higher order function that takes in the description of a function f with k realvalued inputs and p real valued outputs, and an interval
and outputs an interval (hyperrectangle over R p ) IntvlProp(f, I) such that the following soundness guarantees hold: In practice, interval arithmetic approaches have been used to build sound interval propagators [33]. However, they suffer from issues such as the wrapping effect that make their outputs too conservative. This can be remedied by either (a) performing a finer subdivision of the inputs (i.e, increasing M ) to ensure that the intervals I being input into the IntvlProp are sufficiently small to guarantee tight error bounds; or (b) using higher order arithmetics such as affine arithmetic or Taylor polynomial arithmetic [25,32].
The interval propagator serves to define an abstract post-condition operation over sets of cellsŜ ⊆ C. Given such a set,Ŝ, we compute the post condition in the abstract domain. Informally, the post condition is given (a) by iterating over each cell in S; and (b) computing the possible next cells using IntvlProp. Formally, we define the abstract post operation as follows: Given this machinery, an abstract T -step reachability analysis is performed in the standard manner: (a) abstract the initial state; (b) compute post condition for T steps; and (c) check for intersections of the abstract states with the abstraction of the unsafe set. We can also define and use widening operators to make the sequence of iterates converge. The grid based abstract domain can offer some guarantees with respect to the quality of the abstraction. For instance, we can easily bound the Hausdorff distance between the underlying concrete set and the abstraction as a function of the discretization sizes δ j . However, the desirable properties come at a high computational cost since the number of cells grows exponentially in the number of system and disturbance variables.

Tree Decomposed Analysis
We now consider a tree-decomposed approach based on the concept of nodal abstractions. The key idea here is to perform the grid-based abstraction not on the full set of system and disturbance variables, but instead on individual nodal abstractions over a tree decomposition T .

Definition 14 (Nodal Abstractions).
A nodal abstraction Nodal Abstraction(Π, n) corresponding to a node n ∈ N is defined as follows 1. The set of system variables are given by X n : verts X (n) with domain given by D n : proj(D, X n ). 2. The initial states are given by proj(X 0 , X n ).
3. The unsafe set is given by proj(U, X n ). 4. The set of disturbance variables are Y n : verts W (n) with domain given by W n : proj(W, W n ). 5. The updates are described by a relation R(X n , X n ) that relate the possible current states X n and next states X n . The relation is constructed as a conjunction of assertions over variables w) is associated with the node n, we add the conjunct x i = f i (X n , W n ), noting that the proper inputs to f i are contained in verts(n). (b) Otherwise, x i ∈ proj(D, {x i }) that simply states that the next state value of the variable x i is some value in its domain.
Given a system Π, the nodal abstraction is a conservative abstraction, and therefore, it preserves reachability properties.

Lemma 6.
For any reachable state x of Π at time t, its projection proj(x, X n ) is a reachable state of NodalAbstraction(Π, n) at time t.
Since each nodal abstraction involves at most ω + 1 variables, the abstraction at each node can involve at most M ω+1 cells where ω is the tree width. Also, note that a tree decomposition can be found with tree width ω that has at most |X| + |W | nodes. This implies that the number of nodal abstractions can be bounded by (|X| + |W |).
Let Π(n) : NodalAbstraction(Π, n) be the nodal abstraction for tree node n ∈ N . For each node n ∈ N , we instantiate a grid based abstract domain for Π(n) ranging over the variables verts X (n). At the i th step of the reachability analysis, we maintain a map s i each node n to a set of grid cells s i (n) defined over verts(n).
The message passing is performed not over projections of concrete states but over cells belonging to the grid based abstract domain. Nevertheless, we can easily extend the soundness guarantees in Theorem 3 to conclude soundness of the composition.
Once again, we can stop this process after T steps or use widening to force convergence. We now remark on a few technicalities that arise due to the way the tree decomposition is constructed.

Intersections with Unsafe Sets:
Checking for a non-empty intersection with the unsafe sets may require constructing concrete cells over the full dimensional space if the unsafe sets are not tree decomposable for the tree T . However in many cases, the unsafe states are specified as intervals over individual variables, which yields a tree decomposable set. In such cases, we need to intersect the abstraction at each node with the unsafe set and perform message passing to make it canonical before checking for emptiness.

Handling Guards and Invariants:
We have not discussed guards and invariants. It is assumed that such guards and invariants are tree decomposable over the tree T . In this case, we can check which abstract cells have a non-empty intersection with the guard using message passing. The handling of transition systems with guards and invariants will be discussed as part of future extensions.

Experimental Evaluation
In this section, we describe an experimental evaluation of our approach over a set of benchmark problems. Our evaluation is based on a C++-based prototype implementation that can read in the description of a nonlinear dynamical system over a set of system and disturbance variables. The dynamics can currently include polynomials, rational functions and trigonometric functions. Our implementation uses the MPFI library to perform interval arithmetic over the grid cells [36]. We use the HTD library to compute tree decompositions [1]. The system then computes a time-bounded reachable set over the first T steps of the system's execution. Currently, we plot the results and compare the reachable set estimates against simulation data. We also compare the reachable sets computed by the tree decomposition approach against an approach without using tree decompositions. However, we note that the latter approach timed out on systems beyond 4 state variables. Table 1 presents the results over a small set of challenging nonlinear systems benchmarks along with a comparison to two other approaches (a) the approach without tree decomposition and (b) the tool SAPO [22] which computes time bounded reachable sets for polynomial systems using the technique of parallelotope bundles described by Dreossi et al. [23]. The benchmarks range in number of system variables from 3 to 20 state variables. We describe the sources for each benchmark where appropriate. Note that the SAPO tool does not handle nonpolynomial dynamics or time varying disturbances at the time of writing.
The treewidths range from 1 for the simplest system (Example 1) to 3 for the 7-state Laub Loomis oscillator example [30]. We note that the tree decomposition was constructed within 0.01 s for all the examples. We also note that systems with as many as 20 state variables are handled by our approach whereas the monolithic approach cannot handle systems beyond 4 state variables. We now compare the results of our approach to that of the monolithic approach on the two cases where the latter approach completed. System # 1: Consider again the system from Example 1 with 3 state variables and 1 disturbance. We have already noted a tree decomposition of tree width 1 for this example.

System # 2:
In this example, we consider a system over 4 state variables {x, y, z, w} and one disturbance variable w 1 .
x := 0.5x + y + 0.05xy − w 1 , y := −0.7y − 0.03x, z := z − 0.4y, w := w − 0.05xw We obtain a tree decomposition of width 2, wherein the nodes include n 1 : {x, y, w 1 }, n 2 : {y, z} and n 3 : {x, w} with the edges (n 1 , n 2 ) and (n 1 , n 3 ). Figure 1 compares the resulting reachable sets for the tree decomposed reachability analysis versus the monolithic approach. We note differences between the two reachable sets but the loss in precision is not significant.

Coordinated Vehicles:
In this example, we study nonlinear vehicle models of vehicles executing coordinated turns. Each vehicle has states (x i , y i , v x,i , v y,i , ω), representing positions, velocities and the rate of change in the yaw angle, respectively, with a disturbance w i . The dynamics are given by The vehicles are loosely coupled with ω i representing the turn rate of the i th vehicle and ω 0 that of the "lead" vehicle. The i th vehicle tries to gradually align its turn rate to that of the lead vehicle. This model represents a simple scenario of loosely coupled systems that interact using a small set of state variables. Applications including models of cardiac cells that are also loosely coupled through shared action potentials [26].  Table 1 reports results from models involving 1, 2 and 4 vehicles. Since they are loosely coupled, the treewidth of these models is 2.

Laub-Loomis Model:
The Laub-Loomis model is a molecular network that produces spontaneous oscillations for certain values of the model parameters. The model's description was taken from Dang et al. [20]. The system has 7 state variables each of which was subdivided into 100 cells yielding a large state space with 10 14 cells. We note that the tree width of the graph is 3, yielding nodes with upto 4 variables in them.
Comparison with SAPO. SAPO is a state-of-the-art tool that uses polytope bundles and Bernstein polynomials to represent and propagate reachable sets for polynomial dynamical systems [22,23]. We compare our approach directly on SAPO for identical models and initial sets. Note that SAPO does not currently handle non-polynomial models or models with time-varying disturbances. Table 1 shows that SAPO is orders of magnitude faster on all the models, with the sole exception of the 1D-Lattice-10 model. Figure 2 shows the comparison of the reachable sets computed by our approach (shaded blue region) against those computed by SAPO (black rectangles) for five different models. We note that for three of the models compared, neither reachable set is contained in the other. For the one dimensional lattice model, SAPO produces a better reachable set, whereas our approach is better for the influenza model. We also note that both for our approach the precision can be improved markedly by increasing the number of subdivisions, albeit at a large computational cost that depends on the treewidth of the model. The same is true for SAPO, where the number of directions and the template sizes have a non-trivial impact on running time.
Models with Large Treewidths. We briefly report on a few models that we attempted with large treewidths. For such models, our approach of decomposing the space into cells becomes infeasible due to the curse of dimensionality.
A model of how honeybees select between different sites [9,23] has 6 variables and its tree width is 5 with a single tree node containing all state variables. However, the large treewidth is due to two terms in the model which are replaced by disturbance variables that overapproximate their value. This brings down the treewidth to 3, making it tractable for our approach. Details of this transformation are discussed in our extended version. Treewidth reduction using abstractions is an interesting topic for future work.
We originally proposed to analyze a 2D grid lattice model taken from Vleck et al [39]. However, a 2D 10×10 lattice model has a dependency hypergraph that forms a 10×10 grid with treewidth 10. Likewise, the 17-state crazyflie benchmark for SAPO [22] could not be analyzed by our approach since its treewidth is too large.

Conclusions
We have shown how tree decompositions can define an abstract domain that projects concrete sets along the various subsets of state variables. We showed how message passing can be used to exchange information between these subsets. We analyze the completeness of our approach and show that the abstraction is lossy due to the projection operation. We show that for small tree width models, a gridding-based analysis of nonlinear system can be used whereas such approaches are too expensive when applied in a monolithic fashion. For the future, we plan to study tree decompositions for abstract domains such as disjunctions of polyhedra, parallelotope bundles and Taylor models. The process of model abstraction to reduce treewidth is another interesting future possibility. The images or other third party material in this chapter are included in the chapter's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the chapter's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.