Abstract 1 Introduction 2 Preliminaries 3 The Decomposition Method 4 Proof of Theorem 1 5 Other applications of o-minimality References

Verification of Linear Dynamical Systems via O-Minimality of the Real Numbers

Toghrul Karimov ORCID Max Planck Institute for Software Systems, Saarland Informatics Campus, Saarbrücken, Germany
Abstract

A discrete-time linear dynamical system (LDS) is given by an update matrix Md×d, and has the trajectories s,Ms,M2s, for sd. Reachability-type decision problems of linear dynamical systems, most notably the Skolem Problem, lie at the forefront of decidability: typically, sound and complete algorithms are known only in low dimensions, and these rely on sophisticated tools from number theory and Diophantine approximation. Recently, however, o-minimality has emerged as a counterpoint to these number-theoretic tools that allows us to decide certain modifications of the classical problems of LDS without any dimension restrictions. In this paper, we first introduce the Decomposition Method, a framework that captures all applications of o-minimality to decision problems of LDS that are currently known to us. We then use the Decomposition Method to show decidability of the Robust Safety Problem (restricted to bounded initial sets) in arbitrary dimension: given a matrix M, a bounded semialgebraic set S of initial points, and a semialgebraic set T of unsafe points, it is decidable whether there exists ε>0 such that all orbits that begin in the ε-ball around S avoid T.

Keywords and phrases:
Linear dynamical systems, reachability problems, o-minimality
Category:
Track B: Automata, Logic, Semantics, and Theory of Programming
Copyright and License:
[Uncaptioned image] © Toghrul Karimov; licensed under Creative Commons License CC-BY 4.0
2012 ACM Subject Classification:
Theory of computation Logic and verification
Acknowledgements:
The author thanks the anonymous reviewers for their many helpful comments.
Funding:
The author was supported by the DFG grant 389792660 as part of the TRR 248 (see https://perspicuous-computing.science).
Related Version:
ArXiv Version: https://arxiv.org/abs/2410.13053
Editors:
Keren Censor-Hillel, Fabrizio Grandoni, Joël Ouaknine, and Gabriele Puppis

1 Introduction

Linear dynamical systems (LDS) are mathematical models widely used in engineering and sciences to describe systems that evolve over time. A discrete-time LDS is given by an update matrix M. The orbit of a point s under M is the infinite sequence of vectors (Mns)n. In formal verification, the most prominent decision problem of linear dynamical systems is the (point-to-set) Reachability Problem: given Md×d, sd, and a semialgebraic target set T, decide whether there exists n such that MnsT. In terms of program verification, this is equivalent to the Termination Problem for linear loops: given a program fragment of the form

x1,,xds1,,sd
𝐰𝐡𝐢𝐥𝐞¬P(x)𝐝𝐨xMx

where x1,,xd are variable names, si is the initial value of xi, x=(x1,,xd), M is a linear update and P is a Boolean combination of polynomial inequalities, decide whether the loop terminates. The restriction of the Termination Problem to P defined by a single linear equality (i.e., the Reachability Problem where T is restricted to hyperplanes) is Turing-equivalent to the famously open Skolem Problem of linear recurrence sequences (LRS): given an LRS (un)n defined by a recurrence relation

un+d=a1un+d1++adun

as well as the initial values u0,,ud1, decide whether there exists n such that un=0. The Skolem Problem is known to be decidable111Decidability for d4 applies to LRS over real algebraic numbers. For LRS over algebraic numbers, decidability is known for d3. only for LRS with d4 [18]. Consequently, the Termination Problem for linear loops with a linear guard is currently open for 5 or more program variables. Similarly, the Termination Problem with P defined by a single linear inequality (i.e., the Reachability Problem with halfspace targets) is Turing-equivalent to the Positivity Problem: given an LRS (un)n, decide whether un0. The Positivity Problem subsumes the Skolem Problem [12, Chap. 2.5] and is decidable for LRS with d5, but is also known to be hard with respect to certain long-standing open problems in Diophantine approximation [20]: a resolution (in either direction) of the decidability of the Positivity Problem would entail major breakthroughs regarding approximabiliity of Lagrange constants for a large class of transcendental numbers, which are currently believed to be out of reach.

Decidability of the Reachability Problem, which is formidably difficult in full generality, has also been studied under various geometric restrictions on T. The state of the art in this direction is that the Reachability Problem is decidable in arbitrary dimension for T restricted to the class of tame targets [13], i.e. for T that can be constructed through the usual set operations from semialgebraic sets that either (i) have dimension one (i.e. are “string-like”), or (ii) are contained in a three-dimensional subspace of d. In particular, the Reachability Problem is decidable for arbitrary T in dimension d3. These decidability results rely on Baker’s theorem [25] as well its p-adic analogue [26], and are tight in the sense that allowing T to be two-dimensional or contained in a four-dimensional subspace yields a decision problem that is Diophantine-hard similarly to the Positivity Problem [12, Chap. 8].

As discussed above, the Reachability Problem is intimately related to number theory, and, mathematically speaking, becomes harder and harder as we increase the dimension d. Recently, however, a string of results have emerged that show decidability of various reachability-type problems of LDS in arbitrary dimension. First, Akshay et al. [2] showed that the robust variants of the Skolem Problem as well as the Positivity Problem are decidable: given a matrix M, an initial point s, and a hyperplane or a halfspace T, we can decide whether there exists ε>0 such that for all s in the ε-ball around s, the orbit (Mns)n avoids T. The proof is based on classical analyses of LRS. A second interesting result was given by Kelmendi in [14], who showed that given M, an initial point s, and a semialgebraic target T, the frequency of visits

μ=limn1nk=0n1𝟙(MksT)

exists and can be effectively compared against 0. The proof is variation on the fundamental result [21] that ultimate positivity is decidable for simple (also known as diagonalisable) linear recurrence sequences. Lastly, [10] showed that the Pseudo-Orbit Reachability Problem is decidable for diagonalisable M and a single starting point s: given such M,s and a semialgebraic target T, we can decide whether for every ε>0 there exists (xn)n such that (i) x0=s, (ii) xn+1Mxn<ε for all n, and (iii) xnT for some T. Such (xn)n is called an ε-pseudo-orbit of s under M. The proof of decidability uses geometry of pseudo-orbits of diagonalisable M and the classical tools of linear dynamical systems in combination with o-minimality. Briefly, o-minimality of real numbers equipped with arithmetic and exponentiation tell us that every subset of d definable using first-order logic and the aforementioned operations has finitely many connected components.

The thesis of this paper is that every single variant of the Reachability Problem that is decidable in an arbitrary dimension d can, in fact, be explained using o-minimality (Section 5). Although a concept originating in model theory, a branch of mathematical logic, o-minimality has recently had spectacular applications to Diophantine geometry, including counting rational points in algebraic varieties as well as the proofs of Mordell-Lang and André-Oort conjectures [22]. To apply o-minimality, we introduce the Decomposition Method222A similar decomposition is used in [4]., a framework in which all decidability results mentioned above (and more) can be proven. As an application of the Decomposition Method, we study the Robust Safety Problem, which is motivated by situations in which a system must remain in a safe set forever, even when an adversarial perturbation is applied to the initial state. For effectiveness reasons, let us work with real algebraic numbers333If we move from the reals to the complex numbers, i.e. consider the Robust Safety Problem in ¯d, all of our results still apply. In this setting semialgebraic sets are defined by identifying d with 2d., denoted by ¯, which subsume rationals and can be effectively represented and manipulated in computer memory. For sd and ε>0, denote by B(s,ε) the open ε-ball around s. Further write B(S,ε) for sSB(s,ε). The Robust Safety Problem is to decide, given M(¯)d×d, a semialgebraic set S of initial points, and a semialgebraic set T of unsafe points, whether the sequence (MnB(S,ε))n avoids T: that is, whether MnB(s,ε) does not intersect T for all sS and n. Our main result is the following, which, for bounded S, (i) characterises the largest ε>0 such that (MnB(S,ε))n avoids T, and (ii) shows decidability of the Robust Safety Problem.

Theorem 1.

Let M(¯)d×d, Sd be non-empty, semialgebraic and bounded, and Td be semialgebraic. Further let

μ1 =sup{ε0:MnB(S,ε) does not intersect T for all n},
μ2 =sup{ε0:MnB(S,ε) does not intersect T for all sufficiently large n}.

Then μ1¯, μ2(¯){} and is effectively computable, μ1 can be approximated (both from above and below) to arbitrary precision, and it is decidable whether μ1=0. Moreover, for every ε(0,μ2) we can effectively compute N such that (MnB(S,ε))nN avoids T.

Intuitively, μ1,μ2 are the best possible margin of safety and “asymptotic safety”, respectively. From Theorem 1 it follows that the Robust Safety Problem is decidable: we simply check whether μ1=0. For positive instances, i.e. when μ1>0, we can compute ε0 that is arbitrarily close to the best possible safety margin by approximating μ1 from below. Finally, for every positive εμ2 we can decide whether (MnB(S,ε))n avoids T: for ε>μ2, the answer is immediately negative. For ε(0,μ2), we can compute N and then check whether MnB(S,ε) intersects T for some n<N. For ε=μ2, however, we do not know how to decide whether (MnB(S,ε))n avoids T: in this case we have to contend with the hard Diophantine approximation problems that also arise in the classical Reachability Problem.

Related work

At the time of writing, to the best of our knowledge, there are only two published works that apply o-minimality to verification of (discrete-time) linear dynamical systems: [3] studies o-minimal and semialgebraic invariants of LDS (see Section 5), and [10] studies the Pseudo-Orbit Reachability Problem for LDS. In the unpublished extension [9] of [10], it is described how the decidability of the Pseudo-Orbit Reachability Problem can be adapted to study the Robust Safety Problem for singleton S. The latter problem, however, turns out to be much simpler, and can be solved more generally and directly using the Decomposition Method.

Theorem 1 is a generalisation of the aforementioned results of Akshay et al. [2] to bounded S and semialgebraic T. A similar result is [8], which shows the Pseudo-Orbit Reachability Problem is decidable for hyperplane/halfspace targets, without any restrictions on M. Both of these heavily rely on the fact that T is defined by a single linear (in)equality, and their approach does not generalise to targets defined by non-linear or multiple inequalities.

Some other connections between o-minimality and dynamical (or, more generally, cyber-physical) systems have already been established. In [15], Lafferriere et al. show that o-minimal hybrid systems always admit a finite bisimulation; these do not include linear dynamical systems. In [19], Miller studies expansions of the field of real numbers with trajectories of LDS from the perspective of definability and o-minimality.

2 Preliminaries

2.1 Notation and conventions

We denote by 𝒊 the imaginary number and by 𝕋 the unit circle in . We write 𝟎 for a vector of all zeros as well as a zero matrix. In both cases, the dimensions will be clear from the context. For xd and ε>0, we define B(x,ε)={yd:xy<ε} where is the 2-norm. We work exclusively with open 2-balls. We abbreviate B(𝟎,ε) to B(ε). For xd and Yd, we define d(x,Y)=inf{yx:yY}{}. For a set Y and ε>0, we let B(Y,ε)={x:d(x,Y)<ε}.

For sets of vectors X,Y, we write X+Y for {x+y:xX,yY}. For a matrix M and a set X of vectors, we write MX to mean {Mx:xX}. Finally, for a set of matrices and a set X of vectors we define X={Mx:M,xX}. For a (not necessarily invertible) matrix Cd×d and Xd, we define C1X to be {yd:CyX}. We write Cn for (Cn)1.

We denote by Cl(X) the closure of X in a topological space that will be either explicit or clear from the context. We will only be working with Euclidean topology and induced subset topologies. When we say that an object X is effectively computable, we mean that a representation of X in a scheme that will be clear from the context is effectively computable.

2.2 Linear algebra

For a matrix Ad×d, A=maxx𝟎Ax/x. The matrix norm is sub-multiplicative: for all A,B of matching dimensions, ABAB.

Let Xidi×di for 1ik. We write diag(X1,,Xk) for the block-diagonal matrix in d×d, where d=d1++dk, constructed from X1,,Xk respecting the order. For a,b, let Λ(a,b)=[abba], and for θ, let R(θ)=Λ(sin(θ),cos(θ)). A real Jordan block is a matrix

J=[ΛIΛIΛ]d×d (1)

where I is an identity matrix, and either (i) Λ,I1×1, or (ii) Λ=ρR(θ) with I2×2, ρ>0 and θ. A matrix J is in real Jordan form if J=diag(J1,,Jl) where each Ji is a real Jordan block. Given M(¯)d×d, we can compute P,J(¯)d×d such that J is in real Jordan form and M=P1JP [6].

2.3 Logical theories

A structure 𝕄 consists of a universe U, constants c1,,ckU, predicates P1,,Pl where each PiUμ(i) for some μ(i)1, and functions f1,,fm where each fi has the type fi:Uδ(i)U for some δ(i)1. By the language of the structure 𝕄, written 𝕄, we mean the set of all well-formed first-order formulas constructed from symbols denoting the constants c1,,ck, predicates P1,,Pl, functions f1,,fm, as well as the equality symbol =, the quantifier symbols , and the connectives ,,¬. A theory is simply a set of sentences, i.e. formulas without free variables. The theory of the structure 𝕄, written Th(𝕄), is the set of all sentences in the language of 𝕄 that are true in 𝕄. A theory 𝒯

  • is decidable if there exists an algorithm that takes a sentence φ and decides whether φ𝒯, and

  • admits quantifier elimination if for every formula φ with free variables x1,,xn, there exists a quantifier-free formula ψ with the same free variables such that the sentence

    x1,,xn:(φ(x1,,xn)ψ(x1,,xn))

    belongs to 𝒯.

We will be working with the following structures and their theories.

  • Let 0=;0,1,<,+,, which is the ordered ring of real numbers. We will denote the language of this structure by or, called the language of ordered rings. Observe that using the constants 0,1 and the addition, we can obtain any constant c. Hence every atomic formula in or with k free variables is equivalent to p(x1,,xk)0, where p is a polynomial with integer coefficients and is either > or the equality. By the Tarski-Seidenberg theorem, Th(0) admits quantifier elimination and is decidable [5].

  • Let exp=;0,1,<,+,,exp, the ring of real numbers augmented with the function xex. The theory of this structure is model-complete, meaning that quantifiers in any formula can be eliminated down to a single block of existential quantifiers, and decidable assuming Schanuel’s conjecture [23, 16].

We say that a structure 𝕊 expands 𝕄 if 𝕊 and 𝕄 have the same universe and every constant, function, and relation of 𝕄 is also present in 𝕊. We will only need structures expanding 0. A set XUd is definable in a structure 𝕄 if there exist k0, a formula φ in the language of 𝕄 with d+k free variables, and a1,,akU such that for all x1,,xdU, φ(x1,,xd,a1,,ak) holds in 𝕄 if and only if (x1,,xd)X. We say that X is definable in 𝕄 without parameters if we can take k=0 above. Similarly, a function is definable (without parameters) in 𝕄 if its graph is definable (without parameters) in 𝕄.

A structure 𝕄 expanding 0 is o-minimal if every set definable in 𝕄 has finitely many connected components. The structures 0 and exp, as well as the expansion of exp with bounded trigonometric functions [23], are o-minimal.

2.4 Semialgebraic sets and algebraic numbers

A set Xd is semialgebraic if it is definable in 0 without parameters. By quantifier elimination, every semialgebraic set can be defined by a Boolean combination of polynomial equalities and inequalities with integer coefficients. We say that Zd is semialgebraic if Z~={(x1,y1,,xd,yd):(x1+𝒊y1,,xd+𝒊yd)Z} is a semialgebraic subset of 2d. We represent Z by a formula defining Z~. A function f:XY is semialgebraic if its graph {(x,f(x)):xX}X×Y is semialgebraic. We also identify a×b with ab, and define semialgebraic subsets of a×b accordingly.

A number z is algebraic if there exists a polynomial p[x] such that p(z)=0. The set of all algebraic numbers is denoted by ¯. Computationally, we will be working with real algebraic numbers. The number x¯ will be represented by a formula φor defining the singleton set {x}. In this representation, arithmetic on real algebraic numbers is straightforward, and first-order properties of a given number can be verified using the decision procedure for Th(0).

3 The Decomposition Method

We say that a matrix M is a scaling matrix if all eigenvalues of M belong to 0, and a rotation matrix if M is diagonalisable and all eigenvalues of M have modulus 1. A decomposition of Md×d is a pair of matrices (C,D) such that C is a scaling matrix, D is a rotation matrix, and M=CD=DC.

Lemma 2.

Every Md×d has a decomposition (C,D). If M(¯)d×d, then a decomposition satisfying C,D(¯)d×d can be effectively computed.

Proof.

Write M=P1JP, where J is in real Jordan form. If M(¯)d×d, then we additionally ensure the same for P and J. Once we construct a decomposition (CJ,DJ) of J, we have the decomposition (P1CJP,P1DJP) of M. Suppose J=diag(J1,,Jm), where each Ji is a real Jordan block. It suffices to construct decompositions (Ci,Di) of Ji for 1im: then (diag(C1,,Cm),diag(D1,,Dm)) is a decomposition of J.

If Ji has a single real eigenvalue, then we simply take D=I and C=Ji. Now suppose J is of the form (1), where Λ=ρΛ(a,b)=ρR(θ) for some ρ>0 and I is the 2×2 identity matrix. Note that the entries of Λ(a,b) as well as ρ are real algebraic in case M(¯)d×d. Let N be the nilpotent matrix satisfying Ji=diag(Λ(a,b),,Λ(a,b))+N. We compute the decomposition (Ci,Di) with Di=aI and Ci=diag(Λ(0,b),,Λ(0,b))+N. The only eigenvalue of Ci is ρ, and the eigenvalues of Di are e𝒊θ,e𝒊θ𝕋¯. Since Di is a scalar multiple of I, it commutes with Ci.

We will apply the Decomposition Method to the Robust Safety Problem as follows. Let Md×d and S,Td. We have that for every ε>0, MnB(S,ε) does not intersect T for all n if and only if

DnB(S,ε) does not intersect CnT

for all n. We will study the sequences (DnB(S,ε))n and (CnT)n separately. The sequence (DnB(S,ε))n is, in the parlance of dynamical systems theory, “uniformly recurrent”; this is proven using Kronecker’s theorem in Diophantine approximation. The sequence (CnT)n, on the other hand, converges to a limit shape in a strong sense thanks to o-minimality. We will then combine our analyses of the two sequences to prove Theorem 1.

3.1 Applications of Kronecker’s theorem

We next develop the tools necessary for analysing the sequence (DnB(S,ε))n where D is a rotation matrix and S is bounded and semialgebraic. Our main tool is Kronecker’s classical theorem in simultaneous Diophantine approximation. For x,y define [[x]]y to be the closest distance from x to a multiple of y, and write [[x]]=[[x]]1.

Theorem 3 (Kronecker, see [7, Chap. III.5]).

Let λ1,,λl and x1,,xl be real numbers such that for all u1,,ul,

u1λ1++ulλlu1x1++ulxl.

Then for every ε>0 there exist infinitely many n such that [[nλixi]]<ε for all i.

Let β1,,βl𝕋¯, which in our case will be the eigenvalues of a rotation matrix D. Write β=(β1,,βl)𝕋l and βk=(β1k,,βlk) for k. We will use Kronecker’s theorem to analyse the orbit (βn)n in 𝕋l, from which we will deduce various properties of the sequence (Dn)n in d×d. Let

G(β)={(k1,,kl)l:β1k1βlkl=1}

which is called the group of multiplicative relations of (β1,,βl). Observe that G(β) is an abelian subgroup of l and thus has a finite basis of ml elements. Such a basis can be computed using the result [17] of Masser, which places an effective upper bound M(β) on the smallest bit length of a basis of G(β) in terms of the description length of β. To compute a basis, it remains to enumerate all linearly independent subsets of l of bit length at most M(β), and select a maximal one that only contains multiplicative relations satisfied by β.

Let X(β)={z𝕋l:G(β)G(z)}𝕋l. We next prove that (βn)n is uniformly recurrent in X(β).

Lemma 4.

The set X(β) is compact and semialgebraic with an effectively computable representation. Moreover, for every open subset O of 𝕋l containing some αX(β) there exist infinitely many n such that βnO.

Proof.

To compute a representation of X(β), first compute a basis V={v1,,vm} of G(β) as described above. Write vi=(vi,1,,vi,l) for 1im. Then X(β)=i=1lXi where

Xi={(z1,,zl)𝕋l:z1vi,1zlvi,l=1}.

It remains to observe that each Xi is closed and semialgebraic.

To prove the second claim, let α=(e𝒊a1,,e𝒊al)X(β), where ai for all i. Write βi=e𝒊bi where bi for all i. For all n,

βnα=max1il|e𝒊nbie𝒊ai|max1il[[nbiai]]2π=2πmax1il[[nbi/(2π)ai/(2π)]].

We will apply Kronecker’s theorem to the right-hand side of the last equality and to show that it can be made arbitrarily small for infinitely many n. To do this, first we need to show that for all u1,,ul,

u1b12π++ulbl2πu1a12π++ulal2π.

Suppose u1b12π++ulbl2π, and observe that this is equivalent to e𝒊u1b1e𝒊ulbl=1. Since αX(β), G(β)G(α) and hence e𝒊u1a1e𝒊ulal=1. The latter, in turn, is equivalent to u1a12π++ulal2π. Applying Theorem 3 we obtain that for every ε>0 there exist infinitely many n such that [[nbi/(2π)ai/(2π)]]<ε. Hence for every ε>0 there exist infinitely many n such that βnα<ε. It remains to construct, for the given open set O, a value ε>0 such that {z𝕋l:zα<ε}O.

Corollary 5.

The set X(β) is the topological closure444This holds both in 𝕋l and l, since the former is itself closed in the latter. of {βn:n}.

Proof.

For every n, G(β)G(βn) and hence βnX(β). From the density of (βn)n in X(β) (Lemma 4) it follows that X(β) is the closure of (βn)n in 𝕋l.

3.2 Dynamics of rotation matrices

In this section, let D(¯)d×d be a rotation matrix and 𝒟d×d be the closure of {Dn:n} in d×d.

Lemma 6.

We have the following.

  1. (a)

    The set 𝒟 is compact, semialgebraic, and effectively computable.

  2. (b)

    For every A𝒟, det(A)=1. Moreover, there exists b>0 such that A,A1b for all A𝒟.

  3. (c)

    For every non-empty open subset O of 𝒟 there exist infinitely many n such that DnO.

Proof.

Let β1,,βl𝕋¯ be the eigenvalues of D. Write βi=e𝒊bi, where bi, and let βn=(β1n,,βln) for n. Define f:𝕋ld×d by

f(e𝒊θ1,,e𝒊θl)=diag(R(θ1),,R(θl),I),

where I is the identity matrix of the suitable dimension, and let U=f(𝕋l). We have that f is a continuous bijection (i.e. a homeomorphism) between 𝕋l and U. Next, write D=P1JP, where P(¯)d×d and J=f(β1,,βl)(¯)d×d. Observe that for all n, Jn=f(βn). Denote by 𝒥 the closure of {Jn:n}. Since f is a homeomorphism, 𝒥=f(X(β)), where X(β) is the closure of {βn:n} (Corollary 5). Because f is a semialgebraic function (Section 2.4) and X(β) is semialgebraic and effectively computable (Lemma 4), we conclude the same for 𝒥. To prove (a) it remains to observe that 𝒟=P1𝒥P and hence 𝒟 is semialgebraic and effectively computable. Since X(β) is compact, 𝒥 and hence 𝒟 are also compact.

Due to the structure of J, for every n and xd, det(Jn)=1 and

Jnx=Jnx=x

which implies that Jn=Jn=1. Choose b=P1P. We have that for every n, det(Dn)=det(P1)det(Jn)det(P)=1, DnP1JnPb, and similarly Dnb. Statement (b) then follows from the continuity of the determinant, the matrix norm, and the matrix inverse. Finally, by Lemma 4, for every non-empty open subset O of X(β) there exist infinitely many n such that βnO. Statement (c) then follows from the fact that f is a homeomorphism.

We can now give the main result of this section, which is about the orbit of an open set under a rotation matrix.

Lemma 7.

Let Od be open.

  1. (a)

    For all n, DnO𝒟O.

  2. (b)

    For every y𝒟O there exists ε>0 such that DnOB(y,ε) for infinitely many n.

Proof.

Statement (a) follows immediately from the fact that Dn𝒟 for all n. To prove (b), consider y𝒟O. By definition of 𝒟, there exist A𝒟 and zO such that y=Az. Let ε1,ε2>0 be such that OB(z,ε1+ε2). We choose ε=ε1/b, where b is such that X,X1b for all X𝒟 (Lemma 6 (b)). Recalling that every X𝒟 is invertible, construct δ>0 be such that for all X𝒟 with XA<δ, we have that yXB(z,ε2).

Consider Dn𝒟 such that DnA<δ. By Lemma 6 (c), there exist infinitely many such n. We have that

DnODnB(z,ε1)+DnB(ε2)y+DnB(ε1).

Since Dnb, we have that DnB(ε1)B(ε1/b)=B(ε). Therefore, DnOB(y,ε).

3.3 Real exponentiation and o-minimality

We now develop the tools necessary for studying dynamics of scaling matrices. The main results of this section are certain (effective) convergence arguments in o-minimal structures.

We say that a sequence (Zn)n of subsets of d is definable in a structure 𝕊 expanding 0 if there exist k0, a1,,am, and φ𝕊 with d+m+1 free variables such that for all n and x=(x1,,xd)d, xZn if and only if φ(x1,,xd,n,a1,,am) holds in 𝕊. The following form of convergence is also known as Kuratowski convergence, and captures the kind of set-to-set convergence with which we will work.

Definition 8.

Let (Zn)n be a sequence of subsets of d, and

L={xd:lim infnd(x,Zn)=0}.

We say that (Zn)n converges if L={xd:limnd(x,Zn)=0}, in which case L is called the limit shape of (Zn)n.

We immediately have the following.

Lemma 9.

Suppose L is the limit shape of (Zn)n. Then

  1. (a)

    L is closed,

  2. (b)

    L= if and only if for every compact X, there exists N such that ZnX= for all nN, and

  3. (c)

    for every compact Xd and ε>0, there exists N such that for all nN,

    ZnXB(LX,ε).

    In particular, if LX=, then ZnX= for all sufficiently large n.

Proof.

Statement (a) is immediate from the definition, and (b) follows from the Bolzano-Weierstraß theorem. To prove (c), fix compact X and ε>0. The set YXB(LX,ε) is compact. For every yY, since yL, there exist Ny and an open subset Byy of Y such that ZnBy= for all nNy. Hence 𝒞={By:yY} is an open cover of Y. Let  be a finite sub-cover. We can then take N=max{Ny:y}. The second part of (c) follows from the fact that B(,ε)=.

We mention that by the properties above, for (Zn)n contained in a compact set X our notion of convergence coincides with convergence with respect to the Hausdorff metric. Next, we recall that a function definable in an o-minimal structure is ultimately monotonic [24, Sec. 4.1]. The following is an immediate consequence, but we give a little more detail to illustrate the role of o-minimality.

Theorem 10.

Every (Zn)n definable in an o-minimal expansion 𝕊 of 0 converges.

Proof.

Let a1,,am and φ𝕊 be a formula with d+m+1 free variables such that

(x1,,xd)Znφ(x1,,xd,n,a1,,am) holds in 𝕊

for all (x1,,xd)d and n. Consider xk with lim infnd(x,Zn)=0 and Δ>0. Let

T={t0y=(y1,,yd):d(x,y)<Δφ(y1,,yd,t,a1,,am)}

which is definable in 𝕊. By o-minimality, T is a finite union of interval subsets of 0. Since lim infnd(x,Zn)=0, T must be unbounded. Therefore, T must enclose an interval of the form [N,). That is, d(x,Zn)<Δ for all sufficiently large n. It follows that limnd(x,Zn)=0.

In the cases we will encounter, the limit shape L will be a semialgebraic set whose representation can be computed effectively. To show this, we first need a lemma.

Lemma 11.

Let k>0, φor with k+m+1 free variables, and ρ1,,ρm>0¯. Define

A={(x1,,xk)kN.nN:φ(x1,,xk,n,ρ1n,,ρmn)}

and

B={(x1,,xk)kN.nN:¬φ(x1,,xk,n,ρ1n,,ρmn)}.

We have the following.

  1. (a)

    AB=k.

  2. (b)

    Both A and B are semialgebraic sets whose representations can be computed effectively.

  3. (c)

    Given (x1,,xk)Bk, we can effectively compute N such that for all nN, φ(x1,,xk,n,ρ1n,,ρmn) does not hold.

Proof.

Let (x1,,xk)k, and define T={t0:φ(x1,,xk,t,ρ1t,,ρmt)}. By o-minimality of exp, either T is bounded, or it encloses an interval of the form [N,). That is, (x1,,xk) belongs to either A or B. This proves (a).

Next, using quantifier elimination in Th(0), compute a formula

iIjJpi,j(x1,,xk,y0,,ym)Δi,j0 (2)

in or equivalent to φ(x1,,xk,y0,,ym), where each pi,j is a polynomial with rational coefficients and Δi,j{>,=}. Fix iI,jJ and write

pi,j(x1,,xk,n,ρ1n,,ρmn)=l=1dhl(x1,,xk)ql(n)Rln

where R1>>Rd>0, each Ri real algebraic, and hl,ql are non-zero polynomials. Without loss of generality we can further assume that ql(n)>0 for all sufficiently large n. Let

Ai,j={(x1,,xk)k:pi,j(x1,,xk,n,ρ1n,,ρmn)Δi,j0 for all sufficiently large n}.

Observe that whether (x1,,xk)Ai,j only depends on the sign of hl(x1,,xk) for 1ld. In particular, if Δi,j is the equality, then Ai,j is defined by l=1dhl(x1,,xk)=0. Otherwise, Ai,j is defined by

l=1d(hl(x1,,xk)>0s=1l1hs(x1,,xk)=0).

Hence Ai,j is semialgebraic with an effectively computable representation. It remains to observe that A=iIjJAi,j.

To prove (c), fix x=(x1,,xk)Bk. Let φor be the formula obtained by substituting the values of x1,,xk into (2), which will be of the form

ψ(y0,,ym)iIjJvi,j(y0,,ym)Δi,j0

for non-zero polynomials vi,j with rational coefficients and Δi,j{>,=}. We have to construct N such that for all nN, ψ(n,ρ1n,,ρmn) does not hold. To do this, it suffices to construct, for each iI and jJ, an integer Ni,j such that either vi,j(n,ρ1n,,ρmn)Δi,j0 holds for all nNi,j, or it does not hold for all nNi,j. We can then take N=maxi,jNi,j.

Fix i,j, and write vi,j(n,ρ1n,,ρmn)=l=1dql(n)Rln where R1>>Rd>0, and each ql is a non-zero polynomial with rational coefficients. Compute rationals M1,M2,c>0 and Ni,j such that

  • R1>M1>M2>R2,

  • for all tNi,j, |q1(t)|>c, and

  • cM1n>(d1)ql(n)M2n for all nNi,j.

The sign of vi,j(n,ρ1n,,ρmn) is stable for nNi,j and the same as the sign of q1(n).

The following is our main effective convergence result.

Lemma 12.

Let (Zn)n be a sequence of subsets of d. Suppose there exist φor and ρ1,,ρm>0¯ such that for all x=(x1,,xd)d and n,

xZnφ(x1,,xd,n,ρ1n,,ρmn). (3)

Then (Zn)n converges, and the limit shape L of (Zn)0 is semialgebraic and can be effectively computed from φ,ρ1,,ρm.

Proof.

The sequence (Zn)n is definable in exp and thus converges by Theorem 10. Let

X={(ε,x1,,xd)ε>0 and N.nN:d((x1,,xd),Zn)<ε}.

Invoking Lemma 11 with k=d+1 we conclude that X is semialgebraic and effectively computable. It remains to observe that the limit shape is

L={(x1,,xd)ε>0:(ε,x1,,xd)X}

which is semialgebraic and can be effectively computed from X.

3.4 Dynamics of scaling matrices

Now consider a scaling matrix C(¯)d×d and a semialgebraic set Td. In order to apply the Decomposition Method, one of the prerequisites is to understand the sequence of sets (CnT)n. We have the following.

Lemma 13.

Suppose the non-zero eigenvalues of C are ρ1,,ρm. We can compute φor such that for all n and (x1,,xd)d,

xCnTφ(x1,,xd,n,ρ1n,,ρmn).
Proof.

Note that ρ1,,ρm must be positive by the assumption that C is a scaling matrix. The proof follows immediately from the real Jordan form and the definition of CnT (Section 2.1).

Lemma 14.

The sequence (CnT)n converges to a limit shape L that is semialgebraic and an effectively computable.

Proof.

Let ρ1,,ρm and φ be as above, and apply Lemma 12.

4 Proof of Theorem 1

We now prove our main result. Given M(¯)d×d, non-empty and bounded semialgebraic S, and semialgebraic T, let (C,D) be a decomposition of M with C,D(¯)d×d. Recall the definitions of μ1,μ2 from the statement of Theorem 1, and that for all n and ε,

MnB(S,ε) intersects TDnB(S,ε) intersects CnT.

Let 𝒟 be the closure of {Dn:n} (Lemma 6) and L be the limit shape of (CnT)n (Lemma 14). Define μ3=sup{ε0:𝒟B(S,ε) intersects L}.

Figure 1 depicts our application of the Decomposition Method. The bounded initial set is the triangle S. We are interested in the sequences (DnB(S,ε))n for various ε0: specifically, we want to understand the time steps n at which DnB(S,ε) intersects CnT. We see that ε=μ3 is the critical value: in this case 𝒟B(S,ε), which envelopes DnB(S,ε) for all n, just about touches the limit shape L. We have that for ε(0,μ3), for all sufficiently large n, CnT is very close to L and thus is well-separated from DnB(S,ε). On the other hand, for ε>μ3, once again CnB(S,ε) stabilises around L for large values of n, and by uniform recurrence of DnB(S,ε) in 𝒟B(S,ε) (Lemma 7), DnB(S,ε) intersects CnT for infinitely many n. We next formalise these arguments.

Figure 1: Illustration of the Decomposition Method. Here 𝒟={R(θ):θ}2×2 is the group of all 2×2 matrices whose action is a rotation, the larger annulus (light blue) is 𝒟B(S,μ3), and the smaller one (blue) is 𝒟S. The sequence (CnT)n (represented by the curvy line) converges to L. The length of the dashed line segment in μ3.
Claim 0.

μ2=μ3.

Proof.

Suppose there exists ε(0,μ3). We will show that εμ2. Let X be the closure of 𝒟(S,ε), which is disjoint from L. Applying Lemma 9 (c), X(CnT)= for all sufficiently large n, which implies that DnB(S,ε)(CnT)=. Since the choice of ε was arbitrary, it follows that μ3μ2.

Now suppose ε>0 is such that ε>μ3. Take xL(𝒟B(S,ε)). Invoking Lemma 7 (with O=B(S,ε)), let ε>0 be such that DnB(S,ε)B(x,ε) for infinitely many n. Since xL, by definition of convergence, there exists N such that for all nN, CnT intersects B(x,ε). Therefore, for infinitely many nN, DnB(S,ε) intersects CnT. That is, εμ2. It follows that μ2μ3.

Observe that the time step n does not appear in the definition of μ3, which we now know to be equal to μ2. Since both 𝒟 and L are semialgebraic, we can use a decision procedure for the theory of 0 to check whether μ2 is finite. If this is the case, we can compute a formula φor that defines the set {μ2}. Thus μ2 must be algebraic.

Now consider ε(0,μ2). We have to construct N such that (MnB(S,ε))nN avoids T. Using Lemma 13 we can construct a formula φor such that for all n, φ(ε,n,ρ1n,,ρmn) holds if and only if 𝒟B(S,ε) intersects CnT, where ρ1,,ρm are the non-zero eigenvalues of C. By the choice of ε, φ(ε,n,ρ1n,,ρmn) evaluates to false for all sufficiently large n. Therefore, we can compute the desired N using Lemma 11 (c).

We now move on to μ1. For n, let

εn=sup{ε0:DnB(S,ε) does not intersect CnT}.

Observe that each εn is algebraic, μ1=infnεn, and μ2=lim infnεn. Note that each εn is non-negative, but it is possible that εn>μ2. We perform a case analysis based on whether μ2 is infinite.

Case I

Suppose μ2 is infinite, which is the case if and only if L=. By Lemma 9 (b), for every compact X, CnT does not intersect X for all sufficiently large n. Choose any positive ε, and let X=Cl(𝒟B(S,ε)). Compute N (using Lemma 13, as above) such that CnT does not intersect X for all nN. It follows that εnε for all nN and

μ1=infnεn=sup{ε0,,εN1}.

Thus we can explicitly compute μ1, and it is algebraic. We can also decide whether μ1=0.

Case II

Now suppose μ2¯, i.e. L is non-empty. Then μ1¯ since infnεn is either equal to εn for some n, or to lim infnεn. It remains to show how to decide whether μ1=0 and approximate it to arbitrary precision. Let ε(0,μ2) and X be the closure of 𝒟B(S,ε). Using Lemma 9 (c), Lemma 13, and Lemma 11 compute N such that for all nN, CnT is disjoint from 𝒟B(S,ε). Then εnε for all nN. Compute ξ=min{ε0,,εN1}. We have that μ1=0 if and only if ξ=0, which gives us the desired decision procedure. The approximation μ~=min{ξ,ε}, on the other hand, satisfies μ~μ1 and, since μ1μ2, has the absolute error of at most μ2ε. Therefore, we can obtain arbitrarily good sandwiching approximations (μ~,μ~+μ2ε) of μ1 by taking εμ2 from below.

The reason we are unable to determine μ1 exactly is as follows. Suppose we compute ε0,,εN for some large N, and observe that ε0>>εN>μ2. There are two possibilities: either εn remains above μ2 for all n, in which case μ1=μ2, or εn<μ2 for some n, in which case μ1 is an element of the sequence (εn)n. At the moment we do not whether it is possible to determine which of the above is the case by looking at ε0,,εN.

5 Other applications of o-minimality

In this section we briefly discuss how o-minimality and the Decomposition Method are related to various other problems of linear dynamical systems.555A complete account of the claims in this section, as well as details of how the Pseudo-Orbit Reachability Problem can be solved using the Decomposition Method, will be given in a future paper. We first discuss semialgebraic and o-minimal invariants. An invariant of the LDS with update matrix Md×d is a set d such that M. We say that an invariant of M separates an initial set S from a target set T if S and T=. Intuitively, such certifies that all trajectories starting at S are safe: for all sS and n, MnsT. In [3], the authors prove the following. Let M(¯)d×d, T be semialgebraic, and S={s} be a singleton with s(¯)d.

  • It is decidable whether M has an invariant definable in exp.

  • If such exists, then it can be taken to be semialgebraic.

Their arguments can be implemented in the framework of the Decomposition Method. We can moreover easily generalise from singleton to bounded S.

Proposition 15 (See [12, Chap. 9.5]).

Let Sd be bounded and semialgebraic, T be semialgebraic, and M(¯)d×d with a decomposition (C,D) over (¯)d×d. Further let 𝒟 be the closure of {Dn:n} and L be the limit shape of (CnT)n.

  1. 1.

    For all n, MnS intersects T if and only if DnS intersects CnT.

  2. 2.

    We can compute N such that either (i) for all nN, 𝒟S intersects CnT, or (ii) 𝒟S does not intersect CnT for all nN.

  3. 3.

    There exists an exp-definable invariant of M separating S from T if and only if (i) holds, in which case can be taken to be semialgebraic.

Statement (1) is immediate from the definition of a decomposition, and (2) is an application of Lemma 11: the dichotomy (i-ii), without effectiveness of N, can be shown using just o-minimality. Proving (3) requires retracing the arguments of [3]. Intuitively, Proposition 15 states that we can topologically separate MnS from T using an exp-definable set if and only if we can topologically separate 𝒟S (and hence DnS) from CnT for all sufficiently large n. Combining this with our analysis of the Robust Safety Problem yields the following.

Proposition 16.

Let M,S,T be as above. If there exists ε>0 such that MnB(S,ε) does not intersect T for all n, then there exists a semialgebraic invariant separating S from T.

The idea of the proof is that if M,S,T is a positive instance of the Robust Safety Problem, then there exists ε>0 such that for all sufficiently large n, 𝒟B(S,ε) does not intersect CnT; this is even stronger than the topological separation described in (ii) of Proposition 15 (b). The converse of Proposition 16 does not hold already in dimension d=1.

We next discuss the main result of [14] that, given M(¯)d×d, an initial point s(¯)d, and semialgebraic Td, the frequency

μ=limn|{0k<n:MksT}|n

is well-defined, can be approximated to arbitrary precision, and can be effectively compared against zero. The method of [14] is to first use lower bounds on sums of S-units [11], a deep result from algebraic number theory, to express μ in terms λ1,,λm𝕋¯. The second step is to use Weyl’s equidistribution theorem to write μ as an integral, which can be evaluated to arbitrary precision using numerical techniques. We can use o-minimality and the Decomposition Method to give a fully geometric version of the first step.

Proposition 17.

Let Td be semialgebraic, M(¯)d×d with a decomposition (C,D) over (¯)d×d, and s(¯)d. Further let 𝒟 be the closure of {Dn:n} and L be the limit shape of (CnT)n. Then

μ=limn|{0k<n:MksT}|n=limn|{0k<n:DnsCnT}|n

is exactly equal to the measure of L(𝒟s) in a suitable probability space over 𝒟s.

The aforementioned probability space is derived from the Haar measure on 𝒟. We sketch the proof of Proposition 17. Write X=𝒟s. Because X is bounded, by Lemma 9, (CnT)X converges uniformly (with respect to the Hausdorff metric) to LX as n. Thus when measuring the frequency μ, we can pass from the sequence (CnT)n to the limit shape LX. We thus need to understand the frequency with which (Dns)n hits LX, which, by Weyl’s equidistribution theorem, can expressed as an integral over X. To check whether μ>0, we simply check, using tools from semialgebraic geometry, whether LX has full dimension in X.

It is now understood that o-minimality gives us strong ergodic properties of linear dynamical systems, despite the fact that they are typically not compact. This idea is pursued further in the recent work [1], where it shown how to compute an integral representation for the mean payoff μ=1nk=0n1w(Mn), where M(¯)d×d and wd×d is a weight function definable in exp. Let (C,D) be a decomposition of M. The idea is to write

μ=1nk=0n1f(CnDn)=1nk=0n1fn(Dn)

where fn:𝒟 with fn(X)=CnX. Then the sequence of functions (fn)n, when viewed as a sequence of subsets of 𝒟×, converges pointwise to a limit function f:𝒟 (Lemma 12). The next step is to argue that, in fact, μ=limn1nk=0n1f(Dn), after which we can compute the integral representation of μ using the fact that (Dn)n is ergodic by Weyl’s equidistribution theorem. We can also approximate μ to arbitrary precision using the aforementioned integral representation. This, however, requires assuming Schanuel’s conjecture since in the most general setting all we know about f is that it is definable in exp.

References

  • [1] Rajab Aghamov, Christel Baier, Toghrul Karimov, Joël Ouaknine, and Jakob Piribauer. Linear dynamical systems with continuous weight functions. arXiv preprint 2405.06512, 2025.
  • [2] S. Akshay, Hugo Bazille, Blaise Genest, and Mihir Vahanwala. On robustness for the Skolem, positivity and ultimate positivity problems. Logical Methods in Computer Science, 20, 2024. doi:10.46298/LMCS-20(2:11)2024.
  • [3] Shaull Almagor, Dmitry Chistikov, Joël Ouaknine, and James Worrell. O-minimal invariants for discrete-time dynamical systems. ACM Transactions on Computational Logic (TOCL), 23(2):1–20, 2022. doi:10.1145/3501299.
  • [4] Shaull Almagor, Edon Kelmendi, Joël Ouaknine, and James Worrell. Invariants for continuous linear dynamical systems. In 47th International Colloquium on Automata, Languages, and Programming (ICALP 2020), pages 107–1. Schloss Dagstuhl – Leibniz-Zentrum für Informatik, 2020.
  • [5] Jacek Bochnak, Michel Coste, and Marie-Françoise Roy. Real algebraic geometry, volume 36. Springer Science & Business Media, 2013.
  • [6] Jin-Yi Cai. Computing Jordan normal forms exactly for commuting matrices in polynomial time. International Journal of Foundations of Computer Science, 5(04):293–302, 1994. doi:10.1142/S0129054194000165.
  • [7] John Cassels. An Introduction to Diophantine Approximation. Cambridge University Press, 1957.
  • [8] Julian D’Costa, Toghrul Karimov, Rupak Majumdar, Joël Ouaknine, Mahmoud Salamati, Sadegh Soudjani, and James Worrell. The pseudo-skolem problem is decidable. LIPIcs, 202, 2021. doi:10.4230/LIPICS.MFCS.2021.34.
  • [9] Julian D’Costa, Toghrul Karimov, Rupak Majumdar, Joël Ouaknine, Mahmoud Salamati, and James Worrell. The pseudo-reachability problem for diagonalisable linear dynamical systems. arXiv preprint arXiv:2204.12253, 2022.
  • [10] Julian D’Costa, Toghrul Karimov, Rupak Majumdar, Joël Ouaknine, Mahmoud Salamati, and James Worrell. The pseudo-reachability problem for diagonalisable linear dynamical systems. In 47th International Symposium on Mathematical Foundations of Computer Science, 2022.
  • [11] Jan-Hendrik Evertse. On sums of S-units and linear recurrences. Compositio Mathematica, 53(2):225–244, 1984.
  • [12] Toghrul Karimov. Algorithmic verification of linear dynamical systems. PhD thesis, Saarland University, 2024.
  • [13] Toghrul Karimov, Engel Lefaucheux, Joël Ouaknine, David Purser, Anton Varonka, Markus Whiteland, and James Worrell. What’s decidable about linear loops? Proceedings of the ACM on Programming Languages, 6(POPL):1–25, 2022. doi:10.1145/3498727.
  • [14] Edon Kelmendi. Computing the density of the positivity set for linear recurrence sequences. Logical Methods in Computer Science, 19, 2023.
  • [15] Gerardo Lafferriere, George Pappas, and Shankar Sastry. O-minimal hybrid systems. Mathematics of control, signals and systems, 13:1–21, 2000. doi:10.1007/PL00009858.
  • [16] A. Macintyre, A. Wilkie, and P. Odifreddi. On the decidability of the real exponential field. Kreisel’s Mathematics, 115:451, 1996.
  • [17] D. W. Masser. Linear relations on algebraic groups. In Alan Baker, editor, New Advances in Transcendence Theory, pages 248–262. Cambridge University Press, 1988. doi:10.1017/CBO9780511897184.016.
  • [18] M. Mignotte, T. N. Shorey, and R. Tijdeman. The distance between terms of an algebraic recurrence sequence. Journal für die reine und angewandte Mathematik, 349, 1984.
  • [19] Chris Miller. Expansions of o-minimal structures on the real field by trajectories of linear vector fields. Proceedings of the American Mathematical Society, 139(1):319–330, 2011.
  • [20] Joël Ouaknine and James Worrell. Positivity problems for low-order linear recurrence sequences. In Proceedings of the Twenty-Fifth Annual ACM-SIAM Symposium on Discrete Algorithms, December 2013. doi:10.1137/1.9781611973402.27.
  • [21] Joël Ouaknine and James Worrell. Ultimate positivity is decidable for simple linear recurrence sequences. In International Colloquium on Automata, Languages, and Programming, pages 330–341. Springer, 2014. doi:10.1007/978-3-662-43951-7_28.
  • [22] J. Pila, G. Jones, and A. Wilkie. O-minimality and Diophantine geometry. In Proceedings ICM, 2014.
  • [23] Lou van den Dries and Chris Miller. On the real exponential field with restricted analytic functions. Israel Journal of Mathematics, 85:19–56, 1994.
  • [24] Lou van den Dries and Chris Miller. Geometric categories and o-minimal structures. Duke Mathematical Journal, 84(2):497–540, 1996. doi:10.1215/S0012-7094-96-08416-1.
  • [25] G. Wüstholz and A. Baker. Logarithmic forms and group varieties. Journal für die reine und angewandte Mathematik, 442:19–62, 1993.
  • [26] Kunrui Yu. p-Adic logarithmic forms and group varieties II. Acta Arithmetica, 89(4):337–378, 1999. URL: http://eudml.org/doc/207276.